]> Creatis software - creaContours.git/blobdiff - lib/Interface_ManagerContour_NDimensions/wxContourMainFrame.cxx
Added ITK Segmentation Algorithm
[creaContours.git] / lib / Interface_ManagerContour_NDimensions / wxContourMainFrame.cxx
index a3aba8df4f7abfb90b94decc93fc664c475cdd58..e81e73084b2cd0225ef3ecb515e20eb8240cd06d 100644 (file)
@@ -915,6 +915,454 @@ vtkImageData* wxContourMainFrame::getImageData(){
        return _theViewPanel->getImageData();
 }
 
+void wxContourMainFrame::onSegmentationOneSliceITK(wxString distance, wxString sigma, wxString alfa, wxString beta, wxString propagation, wxString iterations, wxString inflation)
+{
+       //JCP 20-10-08 Undo redo implementation
+       saveState();
+       //JCP 20-10-08 Undo redo implementation
+       
+       wxBusyCursor wait;
+       int                                     x                                       = _theViewPanel->GetX();
+       int                                     y                                       = _theViewPanel->GetY();
+       int                                     z                                       = _theViewPanel->GetZ();
+       SegmentationOneSliceITK( x,y,z,distance, sigma, alfa, beta, propagation, iterations, inflation);
+       RefreshInterface();
+}
+
+void wxContourMainFrame::SegmentationOneSliceITK(int x, int y, int z, wxString distanc, wxString sigm, wxString alf, wxString bet, wxString prop, wxString iter, wxString inflation)
+{
+       
+       int typeofcontour = 1;
+       //Image Data
+       vtkImageData    *imagedata      = getImageData();
+       
+       //Tipo de pixeles a utilizar internamente en ITK
+  typedef   float  InternalPixelType;
+  const     unsigned int    Dimension = 2;
+  typedef itk::Image< InternalPixelType, Dimension >  InternalImageType;
+
+  //Tipo de pixeles de salida 1
+  typedef unsigned char OutputPixelType;
+  typedef itk::Image< OutputPixelType, Dimension > OutputImageType;
+
+  //Tipo de pixeles de salida 2
+  typedef unsigned short OutputPixelType2;
+  typedef itk::Image< OutputPixelType2, Dimension > OutputImageType2;
+
+  //Definición del thresholder
+  typedef itk::BinaryThresholdImageFilter< 
+                        InternalImageType, 
+                        OutputImageType    >    ThresholdingFilterType;
+  
+  //Definición del primer filtro de conversión de pixeles
+  typedef itk::CastImageFilter<
+               OutputImageType, OutputImageType2 >  CastFilterType;
+
+  //Definición del segundo tipo de conversión de pixeles
+  typedef itk::CastImageFilter<
+               OutputImageType2, InternalImageType >  CastFilterType2;
+
+  //Tercer tipo de conversión
+  typedef itk::RescaleIntensityImageFilter< 
+                               InternalImageType, 
+                               OutputImageType >   CastFilterType3;
+
+  //Cuarto tipo de conversión
+  typedef itk::RescaleIntensityImageFilter< 
+                               OutputImageType, 
+                               OutputImageType >   CastFilterType4;
+
+  ThresholdingFilterType::Pointer thresholder = ThresholdingFilterType::New();
+                        
+  thresholder->SetLowerThreshold( 0.0 );
+  thresholder->SetUpperThreshold( 128 );
+
+  thresholder->SetOutsideValue(  255  );
+  thresholder->SetInsideValue(  0 );
+
+  //Definción de conexiónes entre VTK e ITK y el writer
+  typedef itk::VTKImageToImageFilter<OutputImageType2> ConnectorType;
+  typedef itk::ImageToVTKImageFilter<OutputImageType> ConnectorType2;
+  typedef  itk::ImageFileWriter<  OutputImageType  > WriterType;
+
+  ConnectorType::Pointer connector= ConnectorType::New();
+  ConnectorType2::Pointer connector2= ConnectorType2::New();
+  
+
+  CastFilterType::Pointer filter=CastFilterType::New();
+  CastFilterType2::Pointer filter2=CastFilterType2::New();
+
+  connector->SetInput( imagedata );
+  filter2->SetInput(connector->GetOutput());
+
+  typedef   itk::CurvatureAnisotropicDiffusionImageFilter< 
+                               InternalImageType, 
+                               InternalImageType >  SmoothingFilterType;
+
+  SmoothingFilterType::Pointer smoothing = SmoothingFilterType::New();
+
+  typedef   itk::GradientMagnitudeRecursiveGaussianImageFilter< 
+                               InternalImageType, 
+                               InternalImageType >  GradientFilterType;
+
+  typedef   itk::SigmoidImageFilter<                               
+                               InternalImageType, 
+                               InternalImageType >  SigmoidFilterType;
+
+  GradientFilterType::Pointer  gradientMagnitude = GradientFilterType::New();
+
+  SigmoidFilterType::Pointer sigmoid = SigmoidFilterType::New();
+
+  sigmoid->SetOutputMinimum(  0.0  );
+  sigmoid->SetOutputMaximum(  255.0  );
+
+
+  typedef  itk::FastMarchingImageFilter< 
+                              InternalImageType, 
+                              InternalImageType >    FastMarchingFilterType;
+
+
+  FastMarchingFilterType::Pointer  fastMarching = FastMarchingFilterType::New();
+
+  typedef  itk::GeodesicActiveContourLevelSetImageFilter< InternalImageType, 
+                InternalImageType >    GeodesicActiveContourFilterType;
+  GeodesicActiveContourFilterType::Pointer geodesicActiveContour = 
+                                     GeodesicActiveContourFilterType::New();
+
+  typedef  itk::ZeroCrossingImageFilter< 
+                              InternalImageType, 
+                              InternalImageType >    ZeroCrossingFilterType;
+ZeroCrossingFilterType::Pointer zeroCrossing =
+                                                                       ZeroCrossingFilterType::New();
+
+const double propagationScaling = atof( prop );
+
+  geodesicActiveContour->SetPropagationScaling( propagationScaling );
+  geodesicActiveContour->SetCurvatureScaling( 1.0 );
+  geodesicActiveContour->SetAdvectionScaling( 1.0 );
+
+  geodesicActiveContour->SetMaximumRMSError( 0.02 );
+  int it=atoi( iter );
+  geodesicActiveContour->SetNumberOfIterations( it );
+
+  smoothing->SetInput( filter2->GetOutput() );
+  gradientMagnitude->SetInput( smoothing->GetOutput() );
+  sigmoid->SetInput( gradientMagnitude->GetOutput() );
+  fastMarching->SetInput( sigmoid->GetOutput() );
+  geodesicActiveContour->SetInput(  fastMarching->GetOutput() );
+  geodesicActiveContour->SetFeatureImage( sigmoid->GetOutput() );
+  
+  zeroCrossing->SetInput( geodesicActiveContour->GetOutput() );
+  //thresholder->SetInput( zeroCrossing->GetOutput() );
+  thresholder->SetInput( geodesicActiveContour->GetOutput() );
+  connector2->SetInput( thresholder->GetOutput()  );
+  
+
+  smoothing->SetTimeStep( 0.125 );
+  smoothing->SetNumberOfIterations(  5 );
+  smoothing->SetConductanceParameter( 9.0 );
+
+
+  const double sigma = atof( sigm );
+  gradientMagnitude->SetSigma(  sigma  );
+
+  const double alpha =  atof( alf );
+  const double beta  =  atof( bet );
+
+  sigmoid->SetAlpha( alpha );
+  sigmoid->SetBeta(  beta  );
+  
+  typedef FastMarchingFilterType::NodeContainer  NodeContainer;
+  typedef FastMarchingFilterType::NodeType       NodeType;
+
+  NodeContainer::Pointer seeds = NodeContainer::New();
+
+  InternalImageType::IndexType  seedPosition;
+  seedPosition[0] = x;
+  seedPosition[1] = y;
+
+  const double initialDistance = atof( distanc );
+
+  NodeType node;
+
+  const double seedValue = - initialDistance;
+
+  node.SetValue( seedValue );
+  node.SetIndex( seedPosition );
+       
+  seeds->Initialize();
+  seeds->InsertElement( 0, node );
+
+  fastMarching->SetTrialPoints(  seeds  );
+
+  fastMarching->SetSpeedConstant( 1.0 );
+  
+  fastMarching->SetOutputSize( 
+           connector->GetOutput()->GetBufferedRegion().GetSize() );
+  
+  fastMarching->SetStoppingValue( 800 );
+  try
+    {
+               
+       connector2->Update();
+       vtkImageData *idata = connector2->GetOutput();
+
+       vtkMarchingContourFilter* cntVTK = vtkMarchingContourFilter::New( );
+       
+       cntVTK->SetInput( idata );
+       
+       cntVTK->SetNumberOfContours( 1 );
+       cntVTK->SetValue( 0, 255 );
+       cntVTK->Update( );
+       cntVTK->UpdateInformation();
+               
+       vtkCleanPolyData* cpd = vtkCleanPolyData::New( );
+       cpd->SetInput( cntVTK->GetOutput( ) );
+       cpd->Update( );
+       cpd->UpdateInformation();
+
+       vtkPolyDataConnectivityFilter* conn = vtkPolyDataConnectivityFilter::New( );
+       conn->SetExtractionModeToLargestRegion( );
+       conn->SetInput( cpd->GetOutput( ) );
+       conn->Update( );
+       conn->UpdateInformation();
+
+       vtkStripper* vtkstripper = vtkStripper::New( );
+       vtkstripper->SetInput( conn->GetOutput() );
+       vtkstripper->Update();
+       vtkstripper->UpdateInformation();
+
+
+       vtkPolyData* polyDataResult =  cntVTK->GetOutput();
+       std::cout<<"Points "<<polyDataResult->GetNumberOfPoints()<<std::endl;
+       polyDataResult->Update( );
+       polyDataResult->UpdateInformation();
+
+//EED
+       /*
+ofstream myfile;
+myfile.open ("C:/Creatis/example.txt");
+myfile << "\n";
+polyDataResult->Print(myfile);
+myfile << "-------------------------------------\n";
+polyDataResult->GetLines()->Print(myfile);
+myfile.close();*/
+
+       cntVTK          -> Delete();
+       cpd                     -> Delete();
+       conn            -> Delete();
+
+
+//--Calculating control points
+
+       std::vector<double> vecX;
+       std::vector<double> vecY;
+       std::vector<double> vecXo;
+       std::vector<double> vecYo;
+       std::vector<double>::iterator vecXoi;
+       std::vector<double>::iterator vecYoi;
+       std::vector<double> vecZ;
+
+       std::vector<double> vecCtrlPointX;
+       std::vector<double> vecCtrlPointY;
+       std::vector<double> vecCtrlPointZ;
+
+
+       double *p;
+       double xAct=0;
+       double yAct=0;
+       int ii,size=polyDataResult->GetNumberOfPoints();
+ofstream myfile;
+myfile.open ("C:/Creatis/example2.txt");
+
+       size=polyDataResult->GetNumberOfPoints();
+       for (ii=0;ii<size;ii++)
+       {
+               if(ii==0)
+               {
+                       xAct=x;
+                       yAct=y;
+               }
+               p       = polyDataResult->GetPoint(ii);
+               double x=p[0];
+               double y=p[1];
+               /*if(fabs(yAct-y)>20)
+               {
+                       if((xAct-x)>1 || (xAct-x)<-1)
+                       {
+                       vecX.push_back( p[0] );
+                       vecY.push_back( p[1] );
+                       myfile <<p[0]<<","<<p[1]<<"\n";
+                       std::cout<<" x Anterior "<<xAct<<" x actual "<<x<<std::endl;
+               std::cout<<" y Anterior "<<yAct<<" y actual "<<y<<std::endl;
+               std::cout<<" x "<<p[0]<<" y "<<p[1]<<std::endl;
+                       vecZ.push_back( -900 );
+                       xAct=x;
+                       yAct=y;
+                       }
+                       else
+                       {
+                               vecXo.push_back(p[0]);
+                               vecYo.push_back(p[1]);
+                       }
+                       
+               }
+               else*/ if(fabs(xAct-x)>11)
+               {
+                       vecXo.push_back(p[0]);
+                       vecYo.push_back(p[1]);
+               }
+               else
+               {
+               vecX.push_back( p[0] );
+               myfile <<p[0]<<","<<p[1]<<"\n";
+               std::cout<<" x Anterior "<<xAct<<" x actual "<<x<<std::endl;
+               std::cout<<" y Anterior "<<yAct<<" y actual "<<y<<std::endl;
+               std::cout<<" x "<<p[0]<<" y "<<p[1]<<std::endl;
+               vecY.push_back( p[1] );
+               vecZ.push_back( -900 );
+               xAct=x;
+               yAct=y;
+               }
+               
+               
+       }
+
+       while(!vecXo.empty())
+       {
+               vecX.push_back(vecXo.back());
+               std::cout<<" x Siguiente "<<vecXo.back();
+               vecXo.pop_back();
+               vecZ.push_back( -900 );
+       }
+       while(!vecYo.empty())
+       {
+               vecY.push_back(vecYo.back());
+                       vecYo.pop_back();
+       }
+       myfile.close();
+
+       /*for(int l=0;l<vecX.size();l++)
+       {
+               if(l==0)
+               {
+            vecXo.push_back(p[0]);
+                       vecYo.push_back(p[1]);
+               }
+               else
+               {
+                       if(vecXoi[l-1]==)
+                       {
+                       }
+               }
+
+       }*/
+
+       ExtractControlPoints2D *extractcontrolpoints2d = new ExtractControlPoints2D();
+
+       extractcontrolpoints2d->SetContour( &vecX , &vecY , &vecZ );
+       
+int method=2;
+       if (method==0){
+               extractcontrolpoints2d->GetInitialControlPoints( &vecCtrlPointX , &vecCtrlPointY , &vecCtrlPointZ );
+       }
+       else if (method==1){
+               extractcontrolpoints2d->GetControlPoints(  &vecCtrlPointX , &vecCtrlPointY , &vecCtrlPointZ );
+       }
+       else if (method==2){
+               extractcontrolpoints2d->SetSamplingControlPoints( 15 );
+               extractcontrolpoints2d->GetSamplingControlPoints(  &vecCtrlPointX , &vecCtrlPointY , &vecCtrlPointZ );
+       }
+       //--Adding contour to the system
+
+       std::vector<int> actualInstantVector;
+       _instantPanel->getInstant( actualInstantVector );
+       actualInstantVector[1]=z;
+
+       int j,sizeCtrPt = vecCtrlPointX.size();
+       
+       manualContourModel *manModelContour =  kernelManager->factoryManualContourModel( typeofcontour );
+       manModelContour->SetNumberOfPointsSpline( ((sizeCtrPt/100)+1)*100 );
+       if (sizeCtrPt>=3){
+               for (j=0 ; j<sizeCtrPt ; j++)
+               {
+                       manModelContour->AddPoint( vecCtrlPointX[j] , vecCtrlPointY[j] , vecCtrlPointZ[j]  );
+               } // for
+               std::string theName;
+               //theName = _modelManager->createOutline( manModelContour, actualInstantVector );
+               theName = kernelManager->createOutline( manModelContour, actualInstantVector );
+               bool addedModel = theName.compare("") != 0;
+               if( addedModel )
+               {
+                       double spc[3];//Si no hay imagen pero hay contornos que spacing se pone por default
+                       _theViewPanel->getSpacing(spc); 
+                       //Adding the manualContourControler to interface objects structure
+                       //Adding the manualViewContour to interface objects structure           
+                       //_theViewPanel->getSceneManager()->setControlActiveStateOfALL( false );//This call is being done here because if the ROI is created underneath the previously created ROIS will still be active.
+                       _theViewPanel->configureViewControlTo(theName, manModelContour, spc, typeofcontour);
+                       //_theViewPanel->getSceneManager()->configureViewControlTo( theName, manModelContour,spc, typeofcontour ) ;
+               }       // if addedModel
+       } // if sizeCtrPt
+
+               
+                
+                WriterType::Pointer writer = WriterType::New();
+         CastFilterType3::Pointer caster = CastFilterType3::New();
+                
+                caster->SetInput( gradientMagnitude->GetOutput() );
+                writer->SetInput( caster->GetOutput() );
+                writer->SetFileName("Gradient Magnitude.png");
+                caster->SetOutputMinimum(   0 );
+                caster->SetOutputMaximum( 255 );
+                writer->Update();
+               
+                CastFilterType3::Pointer caster2 = CastFilterType3::New();
+                WriterType::Pointer writer2 = WriterType::New();
+
+                caster2->SetInput( sigmoid->GetOutput() );
+                writer2->SetInput( caster2->GetOutput() );
+                writer2->SetFileName("Sigmoid.png");
+                caster2->SetOutputMinimum(   0 );
+                caster2->SetOutputMaximum( 255 );
+                writer2->Update();
+
+                CastFilterType3::Pointer caster3 = CastFilterType3::New();
+                WriterType::Pointer writer3 = WriterType::New();
+
+                caster3->SetInput( fastMarching->GetOutput() );
+                writer3->SetInput( caster3->GetOutput() );
+                writer3->SetFileName("FastMarching.bmp");
+                caster3->SetOutputMinimum(   0 );
+                caster3->SetOutputMaximum( 255 );
+                writer3->Update();
+
+                CastFilterType3::Pointer caster4 = CastFilterType3::New();
+                WriterType::Pointer writer4 = WriterType::New();
+
+                caster4->SetInput( geodesicActiveContour->GetOutput() );
+                writer4->SetInput( caster4->GetOutput() );
+                writer4->SetFileName("GeodesicActiveContour.png");
+                caster4->SetOutputMinimum(   0 );
+                caster4->SetOutputMaximum( 255 );
+                writer4->Update();
+
+                CastFilterType3::Pointer caster5 = CastFilterType3::New();
+                WriterType::Pointer writer5 = WriterType::New();
+
+                caster5->SetInput( zeroCrossing->GetOutput() );
+                writer5->SetInput( caster5->GetOutput() );
+                writer5->SetFileName("ZeroCrossing.bmp");
+                caster5->SetOutputMinimum(   0 );
+                caster5->SetOutputMaximum( 255 );
+                writer5->Update();
+    }
+  catch( itk::ExceptionObject & excep )
+    {
+    std::cerr << "Exception caught !" << std::endl;
+    std::cerr << excep << std::endl;
+    }
+}
+
 void wxContourMainFrame::onSegmentationOneSlice(int isovalue,int sampling,int method){
        
        //JCP 20-10-08 Undo redo implementation