]> Creatis software - creaRigidRegistration.git/commitdiff
Added Substraction and Convolution Libraries
authortrillos <trillos>
Thu, 15 Oct 2009 08:32:41 +0000 (08:32 +0000)
committertrillos <trillos>
Thu, 15 Oct 2009 08:32:41 +0000 (08:32 +0000)
lib/Convolution.cxx [new file with mode: 0644]
lib/Convolution.h [new file with mode: 0644]
lib/Substraction.cxx [new file with mode: 0644]
lib/Substraction.h [new file with mode: 0644]

diff --git a/lib/Convolution.cxx b/lib/Convolution.cxx
new file mode 100644 (file)
index 0000000..08b32cf
--- /dev/null
@@ -0,0 +1,47 @@
+#include "Convolution.h"
+
+/*
+* Constructor
+*/
+//------------------------------------------------------------
+Convolution::Convolution()
+{
+       _image=NULL;
+       _convolve = vtkImageConvolve::New();
+       _cast = vtkImageCast::New();
+}
+
+/*
+* Destructor
+*/
+//------------------------------------------------------------
+Convolution::~Convolution()
+{
+       if (_convolve != NULL ) { _convolve->Delete(); }
+       if (_cast != NULL ) { _cast->Delete(); }
+}
+
+vtkImageData *Convolution::getImage()
+{
+       return _cast->GetOutput();
+}
+
+void Convolution::setImage(vtkImageData *image)
+{
+       _image = image;
+       _convolve->SetInput(_image);
+}
+
+void Convolution::setFactor(double factor)
+{
+       _factor = (factor/100.0)*5.0;   
+}
+
+void Convolution::Run()
+{
+       double kernel[] = {0.0,1.0,0.0,1.0,-_factor,1.0,0.0,1.0,0.0};
+       _convolve->SetKernel3x3(kernel);
+       _convolve->Update();
+       _cast->SetInput(_convolve->GetOutput());
+       _cast->SetOutputScalarTypeToDouble();
+}
\ No newline at end of file
diff --git a/lib/Convolution.h b/lib/Convolution.h
new file mode 100644 (file)
index 0000000..b9e752c
--- /dev/null
@@ -0,0 +1,37 @@
+#include "vtkImageData.h"
+#include "vtkImageConvolve.h"
+#include "vtkImageCast.h"
+#include <iostream>
+#include <fstream>
+#include <string>
+#include <vector>
+
+class Convolution
+{
+       public: 
+               Convolution();
+               ~Convolution();
+
+               //Gets the result
+               vtkImageData* getImage();
+               
+               void setImage(vtkImageData *image);
+               void setFactor(double factor);
+
+               void Run();
+               
+       // --- Atributes --- //
+       private: 
+
+               //Original Image
+               vtkImageData *_image;
+
+               //Convolution factor
+               double _factor;
+
+               //Convolution Filter
+               vtkImageConvolve *_convolve;
+
+               //Casting Filter
+               vtkImageCast *_cast;
+};
\ No newline at end of file
diff --git a/lib/Substraction.cxx b/lib/Substraction.cxx
new file mode 100644 (file)
index 0000000..ee25e70
--- /dev/null
@@ -0,0 +1,983 @@
+#include "Substraction.h"
+
+Substraction::Substraction(vtkImageData* imageData1,vtkImageData* imageData2, int uZLevel,int lZLevel, std::vector<double> uColor, std::vector<double> lColor, std::vector<double> mColor)
+{
+       imageResult= vtkImageData::New();
+       sizeImage=0;
+       uZeroLevel=uZLevel;
+       lZeroLevel=lZLevel;
+       if(uColor.size() != NULL)
+       {
+               upperColor[0] = uColor[0];
+               upperColor[1] = uColor[1];
+               upperColor[2] = uColor[2];
+       }
+       else
+       {
+               upperColor[0] = 255;
+               upperColor[1] = 255;
+               upperColor[2] = 255;
+       }       
+       if(mColor.size() != NULL)
+       {
+               mediumColor[0] = mColor[0];
+               mediumColor[1] = mColor[1];
+               mediumColor[2] = mColor[2];
+       }
+       else
+       {
+               mediumColor[0] = 125;
+               mediumColor[1] = 125;
+               mediumColor[2] = 125;
+       }       
+       if(lColor.size() != NULL)
+       {
+               lowerColor[0] = lColor[0];
+               lowerColor[1] = lColor[1];
+               lowerColor[2] = lColor[2];
+       }
+       else
+       {
+               lowerColor[0] = 0;
+               lowerColor[1] = 0;
+               lowerColor[2] = 0;
+       }
+       
+       //Original image type this case is an unsigned char (3)
+       int t=imageData1->GetScalarType();
+       //substracting the image
+       substractImage(imageData1, imageData2);
+}
+
+Substraction::~Substraction()
+{
+       if(imageResult!=NULL)imageResult->Delete();
+}
+
+//----------------------------------------------------------------------------
+// Methods
+//----------------------------------------------------------------------------
+
+
+/*
+Calculate the new image and save it in the attribute imageResult
+it is used if the user had given the imageData
+*/
+void Substraction::substractImage(vtkImageData* imageData1, vtkImageData* imageData2)
+{
+       //dimensions of the image (extent)
+       int ext[6];
+       //setting the dimensionality (1d or 2d or 3d )
+    int newDim[3];
+       //image spacing
+    double spc[3];
+  
+       //getting the information from the original image
+       imageData1->GetSpacing(spc);
+       imageData1->GetExtent(ext);
+       
+       //this a 2d image
+       newDim[0]=ext[1]-ext[0]+1;
+    newDim[1]=ext[3]-ext[2]+1;
+    newDim[2]=1;// in general it is ext[5]-ext[4]+1
+
+
+       imageType = imageData1->GetScalarType();
+       //initializing the image that represents the substracted image
+       initialize(newDim,spc);
+       //Time to substract
+       substract(imageData1, imageData2);      
+}
+
+/*
+  getting ready the imageResult
+*/
+void Substraction::initialize(int dimensions[], double spacing[])
+{
+       //setting image data of the imageResult
+       imageResult->SetScalarType(imageType);
+       imageResult->SetSpacing(spacing);
+       imageResult->SetDimensions(dimensions);
+       imageResult->AllocateScalars();
+       imageResult->Update();
+}
+
+/*
+        Setting the values for the
+*/
+void Substraction::substract(vtkImageData* imageData1, vtkImageData* imageData2)
+{
+       if(imageType == VTK_CHAR)
+       {
+               /*
+               images pointers
+               */
+               // pointers to get into the image
+               char* dataImagePointer1=NULL;
+               char* dataImagePointer2=NULL;
+               char* dataImageResultPointer=NULL;
+               // we start where the  image starts
+               dataImagePointer1=(char*)imageData1->GetScalarPointer(0,0,0);
+               dataImagePointer2=(char*)imageData2->GetScalarPointer(0,0,0);
+               dataImageResultPointer=(char*)imageResult->GetScalarPointer(0,0,0);
+       
+               /*
+               Image Size
+               */
+               int ext[6];
+               imageData1->GetExtent(ext);
+               int sx,sy,sz;
+               sx=ext[1]-ext[0]+1;
+               sy=ext[3]-ext[2]+1;
+               sz=ext[5]-ext[4]+1;
+
+               sizeImage=sx*sy*sz;
+
+               //-----------------
+               //A3    
+               //-----------------
+               //walking in the image
+               int i=0,j=0,k=0,counter=0,nU=0,nL=0,nZ=0;
+               double sum1=0,sum2=0;
+               for(i=0;i<sx;i++)
+               {
+                       for(j=0;j<sy;j++)
+                       {
+                               for(k=0;k<sz;k++)
+                               {
+                                       
+                                       // this is for getting just the grey level in that position
+                                       //originalValue=(short)imageData->GetScalarComponentAsFloat(i,j,k,0);
+                               
+                                       // we get the pointer to the position (i,j,k)y that way we can get the 
+                                       //grey level and we can change it
+                                       dataImagePointer1=(char*)imageData1->GetScalarPointer(i,j,k);
+                                       dataImagePointer2=(char*)imageData2->GetScalarPointer(i,j,k);
+                                       dataImageResultPointer=(char*)imageResult->GetScalarPointer(i,j,k);
+
+                                       sum1=(int)(dataImagePointer1[0]) + (int)(dataImagePointer1[1]) + (int)(dataImagePointer1[2]);
+                                       sum1=sum1/3;
+                                       sum2=(int)(dataImagePointer2[0]) + (int)(dataImagePointer2[1]) + (int)(dataImagePointer2[2]);
+                                       sum2=sum2/3;                            
+                                       if((sum1 - sum2) < lZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(char) lowerColor[0];
+                                               dataImageResultPointer[1] =(char) lowerColor[1];
+                                               dataImageResultPointer[2] =(char) lowerColor[2];
+                                               nL++;
+                                       }
+                                       else if((sum1 - sum2) > uZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(char) upperColor[0];
+                                               dataImageResultPointer[1] =(char) upperColor[1];
+                                               dataImageResultPointer[2] =(char) upperColor[2];
+                                               nU++;
+                                       }
+                                       else
+                                       {
+                                               dataImageResultPointer[0] =(char) mediumColor[0];
+                                               dataImageResultPointer[1] =(char) mediumColor[1];
+                                               dataImageResultPointer[2] =(char) mediumColor[2];
+                                               nZ++;
+                                       }                               
+                                       counter++;
+                               }
+                       }
+               }
+       }
+       else if(imageType == VTK_UNSIGNED_CHAR)
+       {
+               /*
+               images pointers
+               */
+               // pointers to get into the image
+               unsigned char* dataImagePointer1=NULL;
+               unsigned char* dataImagePointer2=NULL;
+               unsigned char* dataImageResultPointer=NULL;
+               // we start where the  image starts
+               dataImagePointer1=(unsigned char*)imageData1->GetScalarPointer(0,0,0);
+               dataImagePointer2=(unsigned char*)imageData2->GetScalarPointer(0,0,0);
+               dataImageResultPointer=(unsigned char*)imageResult->GetScalarPointer(0,0,0);
+       
+               /*
+               Image Size
+               */
+               int ext[6];
+               imageData1->GetExtent(ext);
+               int sx,sy,sz;
+               sx=ext[1]-ext[0]+1;
+               sy=ext[3]-ext[2]+1;
+               sz=ext[5]-ext[4]+1;
+
+               sizeImage=sx*sy*sz;
+
+               //-----------------
+               //A3    
+               //-----------------
+               //walking in the image
+               int i=0,j=0,k=0,counter=0,nU=0,nL=0,nZ=0;
+               double sum1=0,sum2=0;
+               for(i=0;i<sx;i++)
+               {
+                       for(j=0;j<sy;j++)
+                       {
+                               for(k=0;k<sz;k++)
+                               {
+                                       
+                                       // this is for getting just the grey level in that position
+                                       //originalValue=(short)imageData->GetScalarComponentAsFloat(i,j,k,0);
+                               
+                                       // we get the pointer to the position (i,j,k)y that way we can get the 
+                                       //grey level and we can change it
+                                       dataImagePointer1=(unsigned char*)imageData1->GetScalarPointer(i,j,k);
+                                       dataImagePointer2=(unsigned char*)imageData2->GetScalarPointer(i,j,k);
+                                       dataImageResultPointer=(unsigned char*)imageResult->GetScalarPointer(i,j,k);
+
+                                       sum1=(int)(dataImagePointer1[0]) + (int)(dataImagePointer1[1]) + (int)(dataImagePointer1[2]);
+                                       sum1=sum1/3;
+                                       sum2=(int)(dataImagePointer2[0]) + (int)(dataImagePointer2[1]) + (int)(dataImagePointer2[2]);
+                                       sum2=sum2/3;                            
+                                       if((sum1 - sum2) < lZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(unsigned char) lowerColor[0];
+                                               dataImageResultPointer[1] =(unsigned char) lowerColor[1];
+                                               dataImageResultPointer[2] =(unsigned char) lowerColor[2];
+                                               nL++;
+                                       }
+                                       else if((sum1 - sum2) > uZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(unsigned char) upperColor[0];
+                                               dataImageResultPointer[1] =(unsigned char) upperColor[1];
+                                               dataImageResultPointer[2] =(unsigned char) upperColor[2];
+                                               nU++;
+                                       }
+                                       else
+                                       {
+                                               dataImageResultPointer[0] =(unsigned char) mediumColor[0];
+                                               dataImageResultPointer[1] =(unsigned char) mediumColor[1];
+                                               dataImageResultPointer[2] =(unsigned char) mediumColor[2];
+                                               nZ++;
+                                       }                               
+                                       counter++;
+                               }
+                       }
+               }               
+       }
+       else if(imageType == VTK_SIGNED_CHAR)
+       {
+               /*
+               images pointers
+               */
+               // pointers to get into the image
+               signed char* dataImagePointer1=NULL;
+               signed char* dataImagePointer2=NULL;
+               signed char* dataImageResultPointer=NULL;
+               // we start where the  image starts
+               dataImagePointer1=(signed char*)imageData1->GetScalarPointer(0,0,0);
+               dataImagePointer2=(signed char*)imageData2->GetScalarPointer(0,0,0);
+               dataImageResultPointer=(signed char*)imageResult->GetScalarPointer(0,0,0);
+       
+               /*
+               Image Size
+               */
+               int ext[6];
+               imageData1->GetExtent(ext);
+               int sx,sy,sz;
+               sx=ext[1]-ext[0]+1;
+               sy=ext[3]-ext[2]+1;
+               sz=ext[5]-ext[4]+1;
+
+               sizeImage=sx*sy*sz;
+
+               //-----------------
+               //A3    
+               //-----------------
+               //walking in the image
+               int i=0,j=0,k=0,counter=0,nU=0,nL=0,nZ=0;
+               double sum1=0,sum2=0;
+               for(i=0;i<sx;i++)
+               {
+                       for(j=0;j<sy;j++)
+                       {
+                               for(k=0;k<sz;k++)
+                               {
+                                       
+                                       // this is for getting just the grey level in that position
+                                       //originalValue=(short)imageData->GetScalarComponentAsFloat(i,j,k,0);
+                               
+                                       // we get the pointer to the position (i,j,k)y that way we can get the 
+                                       //grey level and we can change it
+                                       dataImagePointer1=(signed char*)imageData1->GetScalarPointer(i,j,k);
+                                       dataImagePointer2=(signed char*)imageData2->GetScalarPointer(i,j,k);
+                                       dataImageResultPointer=(signed char*)imageResult->GetScalarPointer(i,j,k);
+
+                                       sum1=(int)(dataImagePointer1[0]) + (int)(dataImagePointer1[1]) + (int)(dataImagePointer1[2]);
+                                       sum1=sum1/3;
+                                       sum2=(int)(dataImagePointer2[0]) + (int)(dataImagePointer2[1]) + (int)(dataImagePointer2[2]);
+                                       sum2=sum2/3;                            
+                                       if((sum1 - sum2) < lZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(signed char) lowerColor[0];
+                                               dataImageResultPointer[1] =(signed char) lowerColor[1];
+                                               dataImageResultPointer[2] =(signed char) lowerColor[2];
+                                               nL++;
+                                       }
+                                       else if((sum1 - sum2) > uZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(signed char) upperColor[0];
+                                               dataImageResultPointer[1] =(signed char) upperColor[1];
+                                               dataImageResultPointer[2] =(signed char) upperColor[2];
+                                               nU++;
+                                       }
+                                       else
+                                       {
+                                               dataImageResultPointer[0] =(signed char) mediumColor[0];
+                                               dataImageResultPointer[1] =(signed char) mediumColor[1];
+                                               dataImageResultPointer[2] =(signed char) mediumColor[2];
+                                               nZ++;
+                                       }                               
+                                       counter++;
+                               }
+                       }
+               }
+       }
+       else if(imageType == VTK_SHORT)
+       {
+               /*
+               images pointers
+               */
+               // pointers to get into the image
+               short* dataImagePointer1=NULL;
+               short* dataImagePointer2=NULL;
+               short* dataImageResultPointer=NULL;
+               // we start where the  image starts
+               dataImagePointer1=(short*)imageData1->GetScalarPointer(0,0,0);
+               dataImagePointer2=(short*)imageData2->GetScalarPointer(0,0,0);
+               dataImageResultPointer=(short*)imageResult->GetScalarPointer(0,0,0);
+       
+               /*
+               Image Size
+               */
+               int ext[6];
+               imageData1->GetExtent(ext);
+               int sx,sy,sz;
+               sx=ext[1]-ext[0]+1;
+               sy=ext[3]-ext[2]+1;
+               sz=ext[5]-ext[4]+1;
+
+               sizeImage=sx*sy*sz;
+
+               //-----------------
+               //A3    
+               //-----------------
+               //walking in the image
+               int i=0,j=0,k=0,counter=0,nU=0,nL=0,nZ=0;
+               double sum1=0,sum2=0;
+               for(i=0;i<sx;i++)
+               {
+                       for(j=0;j<sy;j++)
+                       {
+                               for(k=0;k<sz;k++)
+                               {
+                                       
+                                       // this is for getting just the grey level in that position
+                                       //originalValue=(short)imageData->GetScalarComponentAsFloat(i,j,k,0);
+                               
+                                       // we get the pointer to the position (i,j,k)y that way we can get the 
+                                       //grey level and we can change it
+                                       dataImagePointer1=(short*)imageData1->GetScalarPointer(i,j,k);
+                                       dataImagePointer2=(short*)imageData2->GetScalarPointer(i,j,k);
+                                       dataImageResultPointer=(short*)imageResult->GetScalarPointer(i,j,k);
+
+                                       sum1=(int)(dataImagePointer1[0]) + (int)(dataImagePointer1[1]) + (int)(dataImagePointer1[2]);
+                                       sum1=sum1/3;
+                                       sum2=(int)(dataImagePointer2[0]) + (int)(dataImagePointer2[1]) + (int)(dataImagePointer2[2]);
+                                       sum2=sum2/3;                            
+                                       if((sum1 - sum2) < lZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(short) lowerColor[0];
+                                               dataImageResultPointer[1] =(short) lowerColor[1];
+                                               dataImageResultPointer[2] =(short) lowerColor[2];
+                                               nL++;
+                                       }
+                                       else if((sum1 - sum2) > uZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(short) upperColor[0];
+                                               dataImageResultPointer[1] =(short) upperColor[1];
+                                               dataImageResultPointer[2] =(short) upperColor[2];
+                                               nU++;
+                                       }
+                                       else
+                                       {
+                                               dataImageResultPointer[0] =(short) mediumColor[0];
+                                               dataImageResultPointer[1] =(short) mediumColor[1];
+                                               dataImageResultPointer[2] =(short) mediumColor[2];
+                                               nZ++;
+                                       }                               
+                                       counter++;
+                               }
+                       }
+               }
+       }
+       else if(imageType == VTK_UNSIGNED_SHORT)
+       {
+               /*
+               images pointers
+               */
+               // pointers to get into the image
+               unsigned short* dataImagePointer1=NULL;
+               unsigned short* dataImagePointer2=NULL;
+               unsigned short* dataImageResultPointer=NULL;
+               // we start where the  image starts
+               dataImagePointer1=(unsigned short*)imageData1->GetScalarPointer(0,0,0);
+               dataImagePointer2=(unsigned short*)imageData2->GetScalarPointer(0,0,0);
+               dataImageResultPointer=(unsigned short*)imageResult->GetScalarPointer(0,0,0);
+       
+               /*
+               Image Size
+               */
+               int ext[6];
+               imageData1->GetExtent(ext);
+               int sx,sy,sz;
+               sx=ext[1]-ext[0]+1;
+               sy=ext[3]-ext[2]+1;
+               sz=ext[5]-ext[4]+1;
+
+               sizeImage=sx*sy*sz;
+
+               //-----------------
+               //A3    
+               //-----------------
+               //walking in the image
+               int i=0,j=0,k=0,counter=0,nU=0,nL=0,nZ=0;
+               double sum1=0,sum2=0;
+               for(i=0;i<sx;i++)
+               {
+                       for(j=0;j<sy;j++)
+                       {
+                               for(k=0;k<sz;k++)
+                               {
+                                       
+                                       // this is for getting just the grey level in that position
+                                       //originalValue=(short)imageData->GetScalarComponentAsFloat(i,j,k,0);
+                               
+                                       // we get the pointer to the position (i,j,k)y that way we can get the 
+                                       //grey level and we can change it
+                                       dataImagePointer1=(unsigned short*)imageData1->GetScalarPointer(i,j,k);
+                                       dataImagePointer2=(unsigned short*)imageData2->GetScalarPointer(i,j,k);
+                                       dataImageResultPointer=(unsigned short*)imageResult->GetScalarPointer(i,j,k);
+
+                                       sum1=(int)(dataImagePointer1[0]) + (int)(dataImagePointer1[1]) + (int)(dataImagePointer1[2]);
+                                       sum1=sum1/3;
+                                       sum2=(int)(dataImagePointer2[0]) + (int)(dataImagePointer2[1]) + (int)(dataImagePointer2[2]);
+                                       sum2=sum2/3;                            
+                                       if((sum1 - sum2) < lZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(unsigned short) lowerColor[0];
+                                               dataImageResultPointer[1] =(unsigned short) lowerColor[1];
+                                               dataImageResultPointer[2] =(unsigned short) lowerColor[2];
+                                               nL++;
+                                       }
+                                       else if((sum1 - sum2) > uZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(unsigned short) upperColor[0];
+                                               dataImageResultPointer[1] =(unsigned short) upperColor[1];
+                                               dataImageResultPointer[2] =(unsigned short) upperColor[2];
+                                               nU++;
+                                       }
+                                       else
+                                       {
+                                               dataImageResultPointer[0] =(unsigned short) mediumColor[0];
+                                               dataImageResultPointer[1] =(unsigned short) mediumColor[1];
+                                               dataImageResultPointer[2] =(unsigned short) mediumColor[2];
+                                               nZ++;
+                                       }                               
+                                       counter++;
+                               }
+                       }
+               }               
+       }
+       else if(imageType == VTK_INT)
+       {
+               /*
+               images pointers
+               */
+               // pointers to get into the image
+               int* dataImagePointer1=NULL;
+               int* dataImagePointer2=NULL;
+               int* dataImageResultPointer=NULL;
+               // we start where the  image starts
+               dataImagePointer1=(int*)imageData1->GetScalarPointer(0,0,0);
+               dataImagePointer2=(int*)imageData2->GetScalarPointer(0,0,0);
+               dataImageResultPointer=(int*)imageResult->GetScalarPointer(0,0,0);
+       
+               /*
+               Image Size
+               */
+               int ext[6];
+               imageData1->GetExtent(ext);
+               int sx,sy,sz;
+               sx=ext[1]-ext[0]+1;
+               sy=ext[3]-ext[2]+1;
+               sz=ext[5]-ext[4]+1;
+
+               sizeImage=sx*sy*sz;
+
+               //-----------------
+               //A3    
+               //-----------------
+               //walking in the image
+               int i=0,j=0,k=0,counter=0,nU=0,nL=0,nZ=0;
+               double sum1=0,sum2=0;
+               for(i=0;i<sx;i++)
+               {
+                       for(j=0;j<sy;j++)
+                       {
+                               for(k=0;k<sz;k++)
+                               {
+                                       
+                                       // this is for getting just the grey level in that position
+                                       //originalValue=(short)imageData->GetScalarComponentAsFloat(i,j,k,0);
+                               
+                                       // we get the pointer to the position (i,j,k)y that way we can get the 
+                                       //grey level and we can change it
+                                       dataImagePointer1=(int*)imageData1->GetScalarPointer(i,j,k);
+                                       dataImagePointer2=(int*)imageData2->GetScalarPointer(i,j,k);
+                                       dataImageResultPointer=(int*)imageResult->GetScalarPointer(i,j,k);
+
+                                       sum1=(int)(dataImagePointer1[0]) + (int)(dataImagePointer1[1]) + (int)(dataImagePointer1[2]);
+                                       sum1=sum1/3;
+                                       sum2=(int)(dataImagePointer2[0]) + (int)(dataImagePointer2[1]) + (int)(dataImagePointer2[2]);
+                                       sum2=sum2/3;                            
+                                       if((sum1 - sum2) < lZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(int) lowerColor[0];
+                                               dataImageResultPointer[1] =(int) lowerColor[1];
+                                               dataImageResultPointer[2] =(int) lowerColor[2];
+                                               nL++;
+                                       }
+                                       else if((sum1 - sum2) > uZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(int) upperColor[0];
+                                               dataImageResultPointer[1] =(int) upperColor[1];
+                                               dataImageResultPointer[2] =(int) upperColor[2];
+                                               nU++;
+                                       }
+                                       else
+                                       {
+                                               dataImageResultPointer[0] =(int) mediumColor[0];
+                                               dataImageResultPointer[1] =(int) mediumColor[1];
+                                               dataImageResultPointer[2] =(int) mediumColor[2];
+                                               nZ++;
+                                       }                               
+                                       counter++;
+                               }
+                       }
+               }
+       }
+       else if(imageType == VTK_UNSIGNED_INT)
+       {
+               /*
+               images pointers
+               */
+               // pointers to get into the image
+               unsigned int* dataImagePointer1=NULL;
+               unsigned int* dataImagePointer2=NULL;
+               unsigned int* dataImageResultPointer=NULL;
+               // we start where the  image starts
+               dataImagePointer1=(unsigned int*)imageData1->GetScalarPointer(0,0,0);
+               dataImagePointer2=(unsigned int*)imageData2->GetScalarPointer(0,0,0);
+               dataImageResultPointer=(unsigned int*)imageResult->GetScalarPointer(0,0,0);
+       
+               /*
+               Image Size
+               */
+               int ext[6];
+               imageData1->GetExtent(ext);
+               int sx,sy,sz;
+               sx=ext[1]-ext[0]+1;
+               sy=ext[3]-ext[2]+1;
+               sz=ext[5]-ext[4]+1;
+
+               sizeImage=sx*sy*sz;
+
+               //-----------------
+               //A3    
+               //-----------------
+               //walking in the image
+               int i=0,j=0,k=0,counter=0,nU=0,nL=0,nZ=0;
+               double sum1=0,sum2=0;
+               for(i=0;i<sx;i++)
+               {
+                       for(j=0;j<sy;j++)
+                       {
+                               for(k=0;k<sz;k++)
+                               {
+                                       
+                                       // this is for getting just the grey level in that position
+                                       //originalValue=(short)imageData->GetScalarComponentAsFloat(i,j,k,0);
+                               
+                                       // we get the pointer to the position (i,j,k)y that way we can get the 
+                                       //grey level and we can change it
+                                       dataImagePointer1=(unsigned int*)imageData1->GetScalarPointer(i,j,k);
+                                       dataImagePointer2=(unsigned int*)imageData2->GetScalarPointer(i,j,k);
+                                       dataImageResultPointer=(unsigned int*)imageResult->GetScalarPointer(i,j,k);
+
+                                       sum1=(int)(dataImagePointer1[0]) + (int)(dataImagePointer1[1]) + (int)(dataImagePointer1[2]);
+                                       sum1=sum1/3;
+                                       sum2=(int)(dataImagePointer2[0]) + (int)(dataImagePointer2[1]) + (int)(dataImagePointer2[2]);
+                                       sum2=sum2/3;                            
+                                       if((sum1 - sum2) < lZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(unsigned int) lowerColor[0];
+                                               dataImageResultPointer[1] =(unsigned int) lowerColor[1];
+                                               dataImageResultPointer[2] =(unsigned int) lowerColor[2];
+                                               nL++;
+                                       }
+                                       else if((sum1 - sum2) > uZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(unsigned int) upperColor[0];
+                                               dataImageResultPointer[1] =(unsigned int) upperColor[1];
+                                               dataImageResultPointer[2] =(unsigned int) upperColor[2];
+                                               nU++;
+                                       }
+                                       else
+                                       {
+                                               dataImageResultPointer[0] =(unsigned int) mediumColor[0];
+                                               dataImageResultPointer[1] =(unsigned int) mediumColor[1];
+                                               dataImageResultPointer[2] =(unsigned int) mediumColor[2];
+                                               nZ++;
+                                       }                               
+                                       counter++;
+                               }
+                       }
+               }               
+       }
+       else if(imageType == VTK_LONG)
+       {
+               /*
+               images pointers
+               */
+               // pointers to get into the image
+               long* dataImagePointer1=NULL;
+               long* dataImagePointer2=NULL;
+               long* dataImageResultPointer=NULL;
+               // we start where the  image starts
+               dataImagePointer1=(long*)imageData1->GetScalarPointer(0,0,0);
+               dataImagePointer2=(long*)imageData2->GetScalarPointer(0,0,0);
+               dataImageResultPointer=(long*)imageResult->GetScalarPointer(0,0,0);
+       
+               /*
+               Image Size
+               */
+               int ext[6];
+               imageData1->GetExtent(ext);
+               int sx,sy,sz;
+               sx=ext[1]-ext[0]+1;
+               sy=ext[3]-ext[2]+1;
+               sz=ext[5]-ext[4]+1;
+
+               sizeImage=sx*sy*sz;
+
+               //-----------------
+               //A3    
+               //-----------------
+               //walking in the image
+               int i=0,j=0,k=0,counter=0,nU=0,nL=0,nZ=0;
+               double sum1=0,sum2=0;
+               for(i=0;i<sx;i++)
+               {
+                       for(j=0;j<sy;j++)
+                       {
+                               for(k=0;k<sz;k++)
+                               {
+                                       
+                                       // this is for getting just the grey level in that position
+                                       //originalValue=(short)imageData->GetScalarComponentAsFloat(i,j,k,0);
+                               
+                                       // we get the pointer to the position (i,j,k)y that way we can get the 
+                                       //grey level and we can change it
+                                       dataImagePointer1=(long*)imageData1->GetScalarPointer(i,j,k);
+                                       dataImagePointer2=(long*)imageData2->GetScalarPointer(i,j,k);
+                                       dataImageResultPointer=(long*)imageResult->GetScalarPointer(i,j,k);
+
+                                       sum1=(int)(dataImagePointer1[0]) + (int)(dataImagePointer1[1]) + (int)(dataImagePointer1[2]);
+                                       sum1=sum1/3;
+                                       sum2=(int)(dataImagePointer2[0]) + (int)(dataImagePointer2[1]) + (int)(dataImagePointer2[2]);
+                                       sum2=sum2/3;                            
+                                       if((sum1 - sum2) < lZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(long) lowerColor[0];
+                                               dataImageResultPointer[1] =(long) lowerColor[1];
+                                               dataImageResultPointer[2] =(long) lowerColor[2];
+                                               nL++;
+                                       }
+                                       else if((sum1 - sum2) > uZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(long) upperColor[0];
+                                               dataImageResultPointer[1] =(long) upperColor[1];
+                                               dataImageResultPointer[2] =(long) upperColor[2];
+                                               nU++;
+                                       }
+                                       else
+                                       {
+                                               dataImageResultPointer[0] =(long) mediumColor[0];
+                                               dataImageResultPointer[1] =(long) mediumColor[1];
+                                               dataImageResultPointer[2] =(long) mediumColor[2];
+                                               nZ++;
+                                       }                               
+                                       counter++;
+                               }
+                       }
+               }
+       }
+       else if(imageType == VTK_UNSIGNED_LONG)
+       {
+               /*
+               images pointers
+               */
+               // pointers to get into the image
+               unsigned long* dataImagePointer1=NULL;
+               unsigned long* dataImagePointer2=NULL;
+               unsigned long* dataImageResultPointer=NULL;
+               // we start where the  image starts
+               dataImagePointer1=(unsigned long*)imageData1->GetScalarPointer(0,0,0);
+               dataImagePointer2=(unsigned long*)imageData2->GetScalarPointer(0,0,0);
+               dataImageResultPointer=(unsigned long*)imageResult->GetScalarPointer(0,0,0);
+       
+               /*
+               Image Size
+               */
+               int ext[6];
+               imageData1->GetExtent(ext);
+               int sx,sy,sz;
+               sx=ext[1]-ext[0]+1;
+               sy=ext[3]-ext[2]+1;
+               sz=ext[5]-ext[4]+1;
+
+               sizeImage=sx*sy*sz;
+
+               //-----------------
+               //A3    
+               //-----------------
+               //walking in the image
+               int i=0,j=0,k=0,counter=0,nU=0,nL=0,nZ=0;
+               double sum1=0,sum2=0;
+               for(i=0;i<sx;i++)
+               {
+                       for(j=0;j<sy;j++)
+                       {
+                               for(k=0;k<sz;k++)
+                               {
+                                       
+                                       // this is for getting just the grey level in that position
+                                       //originalValue=(short)imageData->GetScalarComponentAsFloat(i,j,k,0);
+                               
+                                       // we get the pointer to the position (i,j,k)y that way we can get the 
+                                       //grey level and we can change it
+                                       dataImagePointer1=(unsigned long*)imageData1->GetScalarPointer(i,j,k);
+                                       dataImagePointer2=(unsigned long*)imageData2->GetScalarPointer(i,j,k);
+                                       dataImageResultPointer=(unsigned long*)imageResult->GetScalarPointer(i,j,k);
+
+                                       sum1=(int)(dataImagePointer1[0]) + (int)(dataImagePointer1[1]) + (int)(dataImagePointer1[2]);
+                                       sum1=sum1/3;
+                                       sum2=(int)(dataImagePointer2[0]) + (int)(dataImagePointer2[1]) + (int)(dataImagePointer2[2]);
+                                       sum2=sum2/3;                            
+                                       if((sum1 - sum2) < lZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(unsigned long) lowerColor[0];
+                                               dataImageResultPointer[1] =(unsigned long) lowerColor[1];
+                                               dataImageResultPointer[2] =(unsigned long) lowerColor[2];
+                                               nL++;
+                                       }
+                                       else if((sum1 - sum2) > uZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(unsigned long) upperColor[0];
+                                               dataImageResultPointer[1] =(unsigned long) upperColor[1];
+                                               dataImageResultPointer[2] =(unsigned long) upperColor[2];
+                                               nU++;
+                                       }
+                                       else
+                                       {
+                                               dataImageResultPointer[0] =(unsigned long) mediumColor[0];
+                                               dataImageResultPointer[1] =(unsigned long) mediumColor[1];
+                                               dataImageResultPointer[2] =(unsigned long) mediumColor[2];
+                                               nZ++;
+                                       }                               
+                                       counter++;
+                               }
+                       }
+               }               
+       }
+       else if(imageType == VTK_FLOAT)
+       {
+               /*
+               images pointers
+               */
+               // pointers to get into the image
+               float* dataImagePointer1=NULL;
+               float* dataImagePointer2=NULL;
+               float* dataImageResultPointer=NULL;
+               // we start where the  image starts
+               dataImagePointer1=(float*)imageData1->GetScalarPointer(0,0,0);
+               dataImagePointer2=(float*)imageData2->GetScalarPointer(0,0,0);
+               dataImageResultPointer=(float*)imageResult->GetScalarPointer(0,0,0);
+       
+               /*
+               Image Size
+               */
+               int ext[6];
+               imageData1->GetExtent(ext);
+               int sx,sy,sz;
+               sx=ext[1]-ext[0]+1;
+               sy=ext[3]-ext[2]+1;
+               sz=ext[5]-ext[4]+1;
+
+               sizeImage=sx*sy*sz;
+
+               //-----------------
+               //A3    
+               //-----------------
+               //walking in the image
+               int i=0,j=0,k=0,counter=0,nU=0,nL=0,nZ=0;
+               double sum1=0,sum2=0;
+               for(i=0;i<sx;i++)
+               {
+                       for(j=0;j<sy;j++)
+                       {
+                               for(k=0;k<sz;k++)
+                               {
+                                       
+                                       // this is for getting just the grey level in that position
+                                       //originalValue=(short)imageData->GetScalarComponentAsFloat(i,j,k,0);
+                               
+                                       // we get the pointer to the position (i,j,k)y that way we can get the 
+                                       //grey level and we can change it
+                                       dataImagePointer1=(float*)imageData1->GetScalarPointer(i,j,k);
+                                       dataImagePointer2=(float*)imageData2->GetScalarPointer(i,j,k);
+                                       dataImageResultPointer=(float*)imageResult->GetScalarPointer(i,j,k);
+
+                                       sum1=(int)(dataImagePointer1[0]) + (int)(dataImagePointer1[1]) + (int)(dataImagePointer1[2]);
+                                       sum1=sum1/3;
+                                       sum2=(int)(dataImagePointer2[0]) + (int)(dataImagePointer2[1]) + (int)(dataImagePointer2[2]);
+                                       sum2=sum2/3;                            
+                                       if((sum1 - sum2) < lZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(float) lowerColor[0];
+                                               dataImageResultPointer[1] =(float) lowerColor[1];
+                                               dataImageResultPointer[2] =(float) lowerColor[2];
+                                               nL++;
+                                       }
+                                       else if((sum1 - sum2) > uZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(float) upperColor[0];
+                                               dataImageResultPointer[1] =(float) upperColor[1];
+                                               dataImageResultPointer[2] =(float) upperColor[2];
+                                               nU++;
+                                       }
+                                       else
+                                       {
+                                               dataImageResultPointer[0] =(float) mediumColor[0];
+                                               dataImageResultPointer[1] =(float) mediumColor[1];
+                                               dataImageResultPointer[2] =(float) mediumColor[2];
+                                               nZ++;
+                                       }                               
+                                       counter++;
+                               }
+                       }
+               }
+       }
+       else if(imageType == VTK_DOUBLE)
+       {
+               /*
+               images pointers
+               */
+               // pointers to get into the image
+               double* dataImagePointer1=NULL;
+               double* dataImagePointer2=NULL;
+               double* dataImageResultPointer=NULL;
+               // we start where the  image starts
+               dataImagePointer1=(double*)imageData1->GetScalarPointer(0,0,0);
+               dataImagePointer2=(double*)imageData2->GetScalarPointer(0,0,0);
+               dataImageResultPointer=(double*)imageResult->GetScalarPointer(0,0,0);
+       
+               /*
+               Image Size
+               */
+               int ext[6];
+               imageData1->GetExtent(ext);
+               int sx,sy,sz;
+               sx=ext[1]-ext[0]+1;
+               sy=ext[3]-ext[2]+1;
+               sz=ext[5]-ext[4]+1;
+
+               sizeImage=sx*sy*sz;
+
+               //-----------------
+               //A3    
+               //-----------------
+               //walking in the image
+               int i=0,j=0,k=0,counter=0,nU=0,nL=0,nZ=0;
+               double sum1=0,sum2=0;
+               for(i=0;i<sx;i++)
+               {
+                       for(j=0;j<sy;j++)
+                       {
+                               for(k=0;k<sz;k++)
+                               {
+                                       
+                                       // this is for getting just the grey level in that position
+                                       //originalValue=(short)imageData->GetScalarComponentAsFloat(i,j,k,0);
+                               
+                                       // we get the pointer to the position (i,j,k)y that way we can get the 
+                                       //grey level and we can change it
+                                       dataImagePointer1=(double*)imageData1->GetScalarPointer(i,j,k);
+                                       dataImagePointer2=(double*)imageData2->GetScalarPointer(i,j,k);
+                                       dataImageResultPointer=(double*)imageResult->GetScalarPointer(i,j,k);
+
+                                       sum1=(int)(dataImagePointer1[0]) + (int)(dataImagePointer1[1]) + (int)(dataImagePointer1[2]);
+                                       sum1=sum1/3;
+                                       sum2=(int)(dataImagePointer2[0]) + (int)(dataImagePointer2[1]) + (int)(dataImagePointer2[2]);
+                                       sum2=sum2/3;                            
+                                       if((sum1 - sum2) < lZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(double) lowerColor[0];
+                                               dataImageResultPointer[1] =(double) lowerColor[1];
+                                               dataImageResultPointer[2] =(double) lowerColor[2];
+                                               nL++;
+                                       }
+                                       else if((sum1 - sum2) > uZeroLevel)
+                                       {
+                                               dataImageResultPointer[0] =(double) upperColor[0];
+                                               dataImageResultPointer[1] =(double) upperColor[1];
+                                               dataImageResultPointer[2] =(double) upperColor[2];
+                                               nU++;
+                                       }
+                                       else
+                                       {
+                                               dataImageResultPointer[0] =(double) mediumColor[0];
+                                               dataImageResultPointer[1] =(double) mediumColor[1];
+                                               dataImageResultPointer[2] =(double) mediumColor[2];
+                                               nZ++;
+                                       }                               
+                                       counter++;
+                               }
+                       }
+               }
+       }
+}
+/*
+Returns the filtered image
+*/
+vtkImageData* Substraction::getSubstractedImage()
+{
+       return imageResult;
+}
+
+/*
+       Get Image Size
+*/
+int Substraction::getImageSize()
+{ 
+       return sizeImage;
+}
\ No newline at end of file
diff --git a/lib/Substraction.h b/lib/Substraction.h
new file mode 100644 (file)
index 0000000..b8f5cfb
--- /dev/null
@@ -0,0 +1,84 @@
+#include "vtkImageData.h"
+
+#include <string>
+#include <vector>
+
+class Substraction 
+{
+
+       //----------------------------------------------------------------------------------------
+       // Methods definition
+       //----------------------------------------------------------------------------------------
+       public:
+       //--------------------------
+       //Constructor & Destructor
+       //--------------------------            
+               Substraction(vtkImageData* imageData1, vtkImageData* imageData2, int uZLevel,int lZLevel, std::vector<double> uColor, std::vector<double> lColor, std::vector<double> mColor);
+               ~Substraction();
+       //--------------------------
+       //Methods
+       //--------------------------
+               /*
+               getting ready the points
+               */
+               void initialize(int dimensions[],double spacing[]);
+               /*
+               Calculate the new image and save it in the attribute imageResult
+               it is used if the user had given the imageData
+               */
+               void substractImage(vtkImageData* imageData1, vtkImageData* imageData2);
+               
+               /*
+               Returns the ImageResult
+               */
+               vtkImageData* getSubstractedImage();
+               
+               /*
+                Get Image Size
+               */
+               int getImageSize();
+               
+               /*
+               constructing image substract
+               */
+               void substract(vtkImageData* imageData1, vtkImageData* imageData2);
+               
+
+
+       //----------------------------------------------------------------------------------------
+       // Attributes declaration
+       //----------------------------------------------------------------------------------------
+       private: 
+               /*
+                Substracted Image
+               */
+               vtkImageData* imageResult;
+               /*
+                image size dimx*dimy*dimz
+               */
+               int sizeImage;
+               /*
+               upper zero level for doing the Substraction
+               */
+               int uZeroLevel;
+               /*
+               lower zero level for doing the Substraction
+               */
+               int lZeroLevel;
+               /*
+               Color for the upper threshold
+               */
+               int upperColor[3];
+               /*
+               Color for the lower threshold
+               */
+               int lowerColor[3];
+               /*
+               Color for the medium threshold
+               */
+               int mediumColor[3];
+               /*
+               Image type
+               */
+               int imageType;
+};
\ No newline at end of file