#include "itkCurvatureAnisotropicDiffusionImageFilter.h"
#include "itkGradientMagnitudeRecursiveGaussianImageFilter.h"
#include "itkSigmoidImageFilter.h"
#include "itkFastMarchingImageFilter.h"
#include "itkBinaryThresholdImageFilter.h"
#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"
#include "itkRescaleIntensityImageFilter.h"

#include "itkImage.h"
#include "itkThresholdSegmentationLevelSetImageFilter.h"
#include "itkZeroCrossingImageFilter.h"

#include "itkCommand.h"
#include "itkVTKImageExport.h"
#include "itkVTKImageImport.h"
#include "itkCastImageFilter.h"

#include "vtkSliderRepresentation.h"
#include "vtkSliderRepresentation2D.h"
#include "vtkSliderWidget.h"

#include "vtkImageImport.h"
#include "vtkImageExport.h"
#include "vtkRenderer.h"
#include "vtkRenderWindow.h"
#include "vtkRenderWindowInteractor.h"
#include "vtkActor.h"
#include "vtkPolyData.h"
#include "vtkPolyDataMapper.h"
#include "vtkContourFilter.h"
#include "vtkImageData.h"
#include "vtkDataSet.h"
#include "vtkProperty.h"
#include "vtkImagePlaneWidget.h"
#include "vtkImageViewer2.h"
#include "vtkCamera.h"
#include "vtkCellPicker.h"
#include <vtkSmartPointer.h>
#include <vtkCommand.h>

#include <vtkSeedWidget.h>
#include <vtkSeedRepresentation.h>
#include <vtkPointHandleRepresentation2D.h>
#include <vtkProperty2D.h>

//Slider event
class vtkSliderCallback2 : public vtkCommand
{
public:
  static vtkSliderCallback2 *New() 
    { return new vtkSliderCallback2; }
  void SetImageViewer(vtkImageViewer2 *viewer)
    { this->Viewer =  viewer; }
  virtual void Execute(vtkObject *caller, unsigned long , void* )
    {
      vtkSliderWidget *slider = static_cast<vtkSliderWidget *>(caller);
      vtkSliderRepresentation *sliderRepres = static_cast<vtkSliderRepresentation *>(slider->GetRepresentation());
      int pos = static_cast<int>(sliderRepres->GetValue());

    this->Viewer->SetSlice(pos);
    }
protected:
  vtkImageViewer2 *Viewer;
};



//connect vtk to itk
template <typename ITK_Exporter, typename VTK_Importer>
void ConnectPipelines(ITK_Exporter exporter, VTK_Importer* importer)
{
  importer->SetUpdateInformationCallback(exporter->GetUpdateInformationCallback());
  importer->SetPipelineModifiedCallback(exporter->GetPipelineModifiedCallback());
  importer->SetWholeExtentCallback(exporter->GetWholeExtentCallback());
  importer->SetSpacingCallback(exporter->GetSpacingCallback());
  importer->SetOriginCallback(exporter->GetOriginCallback());
  importer->SetScalarTypeCallback(exporter->GetScalarTypeCallback());
  importer->SetNumberOfComponentsCallback(exporter->GetNumberOfComponentsCallback());
  importer->SetPropagateUpdateExtentCallback(exporter->GetPropagateUpdateExtentCallback());
  importer->SetUpdateDataCallback(exporter->GetUpdateDataCallback());
  importer->SetDataExtentCallback(exporter->GetDataExtentCallback());
  importer->SetBufferPointerCallback(exporter->GetBufferPointerCallback());
  importer->SetCallbackUserData(exporter->GetCallbackUserData());
}

// number of seed points
#define N_SEEDS 13



int main( int argc, char *argv[] )
{
  // seed points
  long int seed_points[N_SEEDS][3] = {
										{	160,  242, 221	},
										{	150,  290, 198	},
										{	189,  222, 198	},
										{	287,  177, 198	},
										{	108,  223, 198	},
	                                    {	134,  342, 155	},
										{	192,  257, 155	},
										{	125,  175, 155	},
										{	252,  173, 155	},
										{	137,  309, 105	},
										{	88,  220, 105	},
										{	90,  265, 56	},
										{	89,  269, 10	}
									  };

  typedef   float           InputPixelType;
  const     unsigned int    Dimension = 3;
  typedef itk::Image< InputPixelType, Dimension >  InputImageType;

  typedef    unsigned char    OutputPixelType;
  typedef itk::Image< OutputPixelType, Dimension >  OutputImageType;

  typedef itk::BinaryThresholdImageFilter< InputImageType, 
                      InputImageType    >    ThresholdingFilterType;
  ThresholdingFilterType::Pointer thresholder1 = ThresholdingFilterType::New();


  const InputPixelType  timeThreshold = 85;
  
  thresholder1->SetLowerThreshold(           0.0  );
  thresholder1->SetUpperThreshold( timeThreshold  );
  thresholder1->SetOutsideValue(  0  );
  thresholder1->SetInsideValue(  255 );
  
  typedef  itk::ImageFileReader< InputImageType > ReaderType;
  ReaderType::Pointer reader = ReaderType::New();
  
  reader->SetFileName( "D:/Hybrid Method/Hybrid/BSP/Diffusion_filter/3D.dcm" ); // input
 

  typedef   itk::GradientMagnitudeRecursiveGaussianImageFilter< 
                               InputImageType, 
                               InputImageType >  GradientFilterType;

  typedef   itk::SigmoidImageFilter<
                               InputImageType, 
                               InputImageType >  SigmoidFilterType;

  GradientFilterType::Pointer  gradientMagnitude = GradientFilterType::New();

  SigmoidFilterType::Pointer sigmoid = SigmoidFilterType::New();

  sigmoid->SetOutputMinimum(  0.0  );        
  sigmoid->SetOutputMaximum(  1.0  );        

  typedef  itk::FastMarchingImageFilter< InputImageType, 
                              InputImageType >    FastMarchingFilterType;

  FastMarchingFilterType::Pointer  fastMarching = FastMarchingFilterType::New();
  
  gradientMagnitude->SetInput( reader->GetOutput() );
  sigmoid->SetInput( gradientMagnitude->GetOutput() );
  fastMarching->SetInput( sigmoid->GetOutput() );
  thresholder1->SetInput( fastMarching->GetOutput() );

  const double sigma = 1.3;                 

  gradientMagnitude->SetSigma(  sigma  );

  const double alpha =  -0.05;                 
  const double beta =  0.3;              

  sigmoid->SetAlpha( alpha );
  sigmoid->SetBeta(  beta  );

  typedef FastMarchingFilterType::NodeContainer           NodeContainer;
  typedef FastMarchingFilterType::NodeType                NodeType;
  NodeContainer::Pointer seeds = NodeContainer::New();

  InputImageType::IndexType  seedPosition, seedPosition2;	
  const double initialDistance = 2.0;

  NodeType node, node1, node2;
  const double seedValue = -1*initialDistance;
 	
  for (int i = 0; i<N_SEEDS; i++)	{
  	seedPosition[0] = seed_points[i][0];
  	seedPosition[1] = seed_points[i][1];
  	seedPosition[2] = seed_points[i][2];
  	node.SetValue( seedValue );
  	node.SetIndex( seedPosition );
  	seeds->InsertElement( i, node );
  	std::cout << "Seed Position " << seedPosition << std::endl;
  }

  fastMarching->SetTrialPoints(  seeds  );

  fastMarching->SetOutputSize( 
           reader->GetOutput()->GetBufferedRegion().GetSize() );

  const double stoppingTime = 100;

  fastMarching->SetStoppingValue(  stoppingTime  );
     
  ThresholdingFilterType::Pointer thresholder2 = ThresholdingFilterType::New();
  thresholder2->SetLowerThreshold( -1000.0 );
  thresholder2->SetUpperThreshold(     0.0 );

  thresholder2->SetOutsideValue(  255  );
  thresholder2->SetInsideValue(  0 );

  typedef  itk::ThresholdSegmentationLevelSetImageFilter< InputImageType, 
    InputImageType > ThresholdSegmentationLevelSetImageFilterType;
  ThresholdSegmentationLevelSetImageFilterType::Pointer thresholdSegmentation =
    ThresholdSegmentationLevelSetImageFilterType::New();

  thresholdSegmentation->SetPropagationScaling( 1.0 );                //////////////////////////////
  thresholdSegmentation->SetCurvatureScaling( 1.0 );                  //////////////////////////////

  thresholdSegmentation->SetMaximumRMSError( 0.01 );               //////////////////////////////
  thresholdSegmentation->SetNumberOfIterations( 100 );            //////////////////////////////

  thresholdSegmentation->SetUpperThreshold( 123 );                    //////////////////////////////
  thresholdSegmentation->SetLowerThreshold( 118 );                    //////////////////////////////
  thresholdSegmentation->SetIsoSurfaceValue(0.0);                     /////////////////////////////

  thresholdSegmentation->SetInput( thresholder1->GetOutput() );
  thresholdSegmentation->SetFeatureImage( reader->GetOutput() );
  thresholder2->SetInput( thresholdSegmentation->GetOutput() );
  
  reader->Update();

  typedef itk::VTKImageExport< InputImageType > ExportFilterType;

  ExportFilterType::Pointer itkExporter = ExportFilterType::New();
    
  itkExporter->SetInput( thresholder2->GetOutput() );

  vtkImageImport* vtkImporter = vtkImageImport::New();  
  ConnectPipelines(itkExporter, vtkImporter);
  
  //vtkImporter->Update();
 
  vtkImageViewer2 *imageViewer = vtkImageViewer2::New();
  imageViewer->SetInput(vtkImporter->GetOutput());
  imageViewer->SetColorLevel(127);
  imageViewer->SetColorWindow(255);

  vtkRenderWindowInteractor* iren = vtkRenderWindowInteractor::New();
  
  imageViewer->SetupInteractor(iren);
  imageViewer->GetRenderWindow()->SetMultiSamples(0);
  imageViewer->GetRenderWindow()->SetSize(600, 600);

  imageViewer->Render();
  imageViewer->GetRenderer()->ResetCamera();

  vtkSliderRepresentation2D *SliderRepres = vtkSliderRepresentation2D::New();
  int min = imageViewer->GetSliceMin();
  int max = imageViewer->GetSliceMax();
  SliderRepres->SetMinimumValue(min);
  SliderRepres->SetMaximumValue(max);
  SliderRepres->SetValue(static_cast<int>((min + max) / 2));
  
  SliderRepres->GetPoint1Coordinate()->SetCoordinateSystemToNormalizedDisplay();
  SliderRepres->GetPoint1Coordinate()->SetValue(0.92, 0.9);
  SliderRepres->GetPoint2Coordinate()->SetCoordinateSystemToNormalizedDisplay();
  SliderRepres->GetPoint2Coordinate()->SetValue(0.92, 0.1);
  SliderRepres->SetSliderLength(0.01);
  SliderRepres->SetSliderWidth(0.03);
  SliderRepres->SetEndCapLength(0.01);
  SliderRepres->SetEndCapWidth(0.05);
  SliderRepres->SetTubeWidth(0.005);
  SliderRepres->SetLabelFormat("%3.0lf");
  SliderRepres->SetTitleHeight(0.03);
  SliderRepres->SetLabelHeight(0.03);

  vtkSliderWidget *SliderWidget = vtkSliderWidget::New();
  SliderWidget->SetInteractor(iren);
  SliderWidget->SetRepresentation(SliderRepres);
  SliderWidget->KeyPressActivationOff();
  SliderWidget->SetAnimationModeToAnimate();
  SliderWidget->SetEnabled(true);
  
  vtkSliderCallback2 *SliderCb = vtkSliderCallback2::New();
  SliderCb->SetImageViewer(imageViewer);
  SliderWidget->AddObserver(vtkCommand::InteractionEvent, SliderCb);  

  imageViewer->SetSlice(static_cast<int>(SliderRepres->GetValue()));
  imageViewer->SetSliceOrientationToXY();

  SliderRepres->Delete();

  SliderCb->Delete();

  iren->Start();

  iren->Delete();

  return 0;
}