#define DEBUGing

#include "itkSWaveletTest.h"

#include <algorithm>
#include <fstream>



int main( int argc, char * argv [] )
{
  SphereMeshSourceType::Pointer  mySphereMeshSource = SphereMeshSourceType::New();

  PointType center; 
  center.Fill( 0.0 );

  VectorType scale;
  scale.Fill( 1.0 );
  
  int n = 4;

  mySphereMeshSource->SetCenter( center );
  mySphereMeshSource->SetResolution( n );
  mySphereMeshSource->SetScale( scale );

  mySphereMeshSource->Modified();

  try
    {
      mySphereMeshSource->Update();
    }
  catch( itk::ExceptionObject & excp )
    {
      std::cerr << "Error during Update() " << std::endl;
      std::cerr << excp << std::endl;
    }


  // testing task 3:
  std::cerr<<"Testing task 3:"<<std::endl; 
  std::cerr<<"    Visualize the scaling function at coarsest level on the n-th resolution subdivision."<<std::endl;
  std::cerr<<"    the n-th subdivision mesh is written in finestSubdibision.meta."<<std::endl;
  std::cerr<<"    the scalar scaling function is written in ScalingFnKWVisu.txt."<<std::endl;

  std::vector< PointType > v = mySphereMeshSource->GetVerts();
  std::vector< PointType >::const_iterator itv = v.begin();
  std::vector< PointType >::const_iterator itvEnd = v.end();

  std::vector< double > f( v.size() );
  std::vector< double >::iterator itf = f.begin();

  for ( ; itv != itvEnd; ++itv, ++itf ) {
    *itf = (*itv)[2];
  }

  mySphereMeshSource->SetScalarFunction( f );

  mySphereMeshSource->SWaveletTransform();

  double a[]={1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
  std::vector<double> va(a, a+12);
  mySphereMeshSource->SetAllCoarsestScalingCoefficients( va );
  mySphereMeshSource->SetAllWaveletCoefficientsOfAllScales( 0 );

  mySphereMeshSource->inverseSWaveletTransform();
  mySphereMeshSource->WriteReconstructedScalarFunctionForKWVisu("ScalingFnKWVisu.txt");

  return EXIT_SUCCESS;
}