lib/Substraction.cxx

   1 /*
   2 # ---------------------------------------------------------------------
   3 #
   4 # Copyright (c) CREATIS (Centre de Recherche en Acquisition et Traitement de l'Image
   5 #                        pour la Santé)
   6 # Authors : Eduardo Davila, Frederic Cervenansky, Claire Mouton
   7 # Previous Authors : Laurent Guigues, Jean-Pierre Roux
   8 # CreaTools website : www.creatis.insa-lyon.fr/site/fr/creatools_accueil
   9 #
  10 #  This software is governed by the CeCILL-B license under French law and
  11 #  abiding by the rules of distribution of free software. You can  use,
  12 #  modify and/ or redistribute the software under the terms of the CeCILL-B
  13 #  license as circulated by CEA, CNRS and INRIA at the following URL
  14 #  http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
  15 #  or in the file LICENSE.txt.
  16 #
  17 #  As a counterpart to the access to the source code and  rights to copy,
  18 #  modify and redistribute granted by the license, users are provided only
  19 #  with a limited warranty  and the software's author,  the holder of the
  20 #  economic rights,  and the successive licensors  have only  limited
  21 #  liability.
  22 #
  23 #  The fact that you are presently reading this means that you have had
  24 #  knowledge of the CeCILL-B license and that you accept its terms.
  25 # ------------------------------------------------------------------------      */
  26
  27 #include "Substraction.h"
  28
  29 Substraction::Substraction(vtkImageData* imageData1,vtkImageData* imageData2, int uZLevel,int lZLevel, std::vector<double> uColor, std::vector<double> lColor, std::vector<double> mColor)
  30 {
  31         imageResult= vtkImageData::New();
  32         sizeImage=0;
  33         uZeroLevel=uZLevel;
  34         lZeroLevel=lZLevel;
  35         if(uColor.size() != 0)
  36         {
  37                 upperColor[0] = uColor[0];
  38                 upperColor[1] = uColor[1];
  39                 upperColor[2] = uColor[2];
  40         }
  41         else
  42         {
  43                 upperColor[0] = 255;
  44                 upperColor[1] = 255;
  45                 upperColor[2] = 255;
  46         }
  47         if(mColor.size() != 0)
  48         {
  49                 mediumColor[0] = mColor[0];
  50                 mediumColor[1] = mColor[1];
  51                 mediumColor[2] = mColor[2];
  52         }
  53         else
  54         {
  55                 mediumColor[0] = 125;
  56                 mediumColor[1] = 125;
  57                 mediumColor[2] = 125;
  58         }
  59         if(lColor.size() != 0)
  60         {
  61                 lowerColor[0] = lColor[0];
  62                 lowerColor[1] = lColor[1];
  63                 lowerColor[2] = lColor[2];
  64         }
  65         else
  66         {
  67                 lowerColor[0] = 0;
  68                 lowerColor[1] = 0;
  69                 lowerColor[2] = 0;
  70         }
  71
  72         //Original image type this case is an unsigned char (3)
  73         //int t=imageData1->GetScalarType(); // JPR : unused
  74         //substracting the image
  75         substractImage(imageData1, imageData2);
  76 }
  77
  78 Substraction::~Substraction()
  79 {
  80         if(imageResult!=NULL)imageResult->Delete();
  81 }
  82
  83 //----------------------------------------------------------------------------
  84 // Methods
  85 //----------------------------------------------------------------------------
  86
  87 /*
  88 Calculate the new image and save it in the attribute imageResult
  89 it is used if the user had given the imageData
  90 */
  91 void Substraction::substractImage(vtkImageData* imageData1, vtkImageData* imageData2)
  92 {
  93         //dimensions of the image (extent)
  94     int ext[6];
  95         //setting the dimensionality (1d or 2d or 3d )
  96     int newDim[3];
  97         //image spacing
  98     double spc[3];
  99
 100         //getting the information from the original image
 101     imageData1->GetSpacing(spc);
 102     imageData1->GetExtent(ext);
 103
 104         //this a 2d image
 105     newDim[0]=ext[1]-ext[0]+1;
 106     newDim[1]=ext[3]-ext[2]+1;
 107     newDim[2]=1;// in general it is ext[5]-ext[4]+1
 108
 109
 110     imageType = imageData1->GetScalarType();
 111         //initializing the image that represents the substracted image
 112     initialize(newDim,spc);
 113         //Time to substract
 114     substract(imageData1, imageData2);
 115 }
 116
 117 /*
 118   getting ready the imageResult
 119 */
 120 void Substraction::initialize(int dimensions[], double spacing[])
 121 {
 122         //setting image data of the imageResult
 123         imageResult->SetSpacing(spacing);
 124         imageResult->SetDimensions(dimensions);
 125 //EED 2017-01-01 Migration VTK7
 126 #if VTK_MAJOR_VERSION <= 5
 127         imageResult->SetScalarType(imageType);
 128         imageResult->AllocateScalars();
 129         imageResult->Update();
 130 #else
 131         imageResult->AllocateScalars(imageType,1);
 132 #endif
 133
 134 }
 135
 136 /*
 137          Setting the values for the
 138 */
 139 void Substraction::substract(vtkImageData* imageData1, vtkImageData* imageData2)
 140 {
 141         if(imageType == VTK_CHAR)
 142         {
 143                 /*
 144                 images pointers
 145                 */
 146                 // pointers to get into the image
 147                 char* dataImagePointer1=NULL;
 148                 char* dataImagePointer2=NULL;
 149                 char* dataImageResultPointer=NULL;
 150
 151                 substractByType(dataImagePointer1, dataImagePointer2, dataImageResultPointer, imageData1, imageData2);
 152         }
 153         else if(imageType == VTK_UNSIGNED_CHAR)
 154         {
 155                 /*
 156                 images pointers
 157                 */
 158                 // pointers to get into the image
 159                 unsigned char* dataImagePointer1=NULL;
 160                 unsigned char* dataImagePointer2=NULL;
 161                 unsigned char* dataImageResultPointer=NULL;
 162
 163                 substractByType(dataImagePointer1, dataImagePointer2, dataImageResultPointer, imageData1, imageData2);
 164         }
 165         else if(imageType == VTK_SIGNED_CHAR)
 166         {
 167                 /*
 168                 images pointers
 169                 */
 170                 // pointers to get into the image
 171                 signed char* dataImagePointer1=NULL;
 172                 signed char* dataImagePointer2=NULL;
 173                 signed char* dataImageResultPointer=NULL;
 174
 175                 substractByType(dataImagePointer1, dataImagePointer2, dataImageResultPointer, imageData1, imageData2);
 176         }
 177         else if(imageType == VTK_SHORT)
 178         {
 179                 /*
 180                 images pointers
 181                 */
 182                 // pointers to get into the image
 183                 short* dataImagePointer1=NULL;
 184                 short* dataImagePointer2=NULL;
 185                 short* dataImageResultPointer=NULL;
 186
 187                 substractByType(dataImagePointer1, dataImagePointer2, dataImageResultPointer, imageData1, imageData2);
 188         }
 189         else if(imageType == VTK_UNSIGNED_SHORT)
 190         {
 191                 /*
 192                 images pointers
 193                 */
 194                 // pointers to get into the image
 195                 unsigned short* dataImagePointer1=NULL;
 196                 unsigned short* dataImagePointer2=NULL;
 197                 unsigned short* dataImageResultPointer=NULL;
 198
 199                 substractByType(dataImagePointer1, dataImagePointer2, dataImageResultPointer, imageData1, imageData2);
 200         }
 201         else if(imageType == VTK_INT)
 202         {
 203                 /*
 204                 images pointers
 205                 */
 206                 // pointers to get into the image
 207                 int* dataImagePointer1=NULL;
 208                 int* dataImagePointer2=NULL;
 209                 int* dataImageResultPointer=NULL;
 210
 211                 substractByType(dataImagePointer1, dataImagePointer2, dataImageResultPointer, imageData1, imageData2);
 212         }
 213         else if(imageType == VTK_UNSIGNED_INT)
 214         {
 215                 /*
 216                 images pointers
 217                 */
 218                 // pointers to get into the image
 219                 unsigned int* dataImagePointer1=NULL;
 220                 unsigned int* dataImagePointer2=NULL;
 221                 unsigned int* dataImageResultPointer=NULL;
 222
 223                 substractByType(dataImagePointer1, dataImagePointer2, dataImageResultPointer, imageData1, imageData2);
 224         }
 225         else if(imageType == VTK_LONG)
 226         {
 227                 /*
 228                 images pointers
 229                 */
 230                 // pointers to get into the image
 231                 long* dataImagePointer1=NULL;
 232                 long* dataImagePointer2=NULL;
 233                 long* dataImageResultPointer=NULL;
 234
 235                 substractByType(dataImagePointer1, dataImagePointer2, dataImageResultPointer, imageData1, imageData2);
 236         }
 237         else if(imageType == VTK_UNSIGNED_LONG)
 238         {
 239                 /*
 240                 images pointers
 241                 */
 242                 // pointers to get into the image
 243                 unsigned long* dataImagePointer1=NULL;
 244                 unsigned long* dataImagePointer2=NULL;
 245                 unsigned long* dataImageResultPointer=NULL;
 246
 247                 substractByType(dataImagePointer1, dataImagePointer2, dataImageResultPointer, imageData1, imageData2);
 248         }
 249         else if(imageType == VTK_FLOAT)
 250         {
 251                 /*
 252                 images pointers
 253                 */
 254                 // pointers to get into the image
 255                 float* dataImagePointer1=NULL;
 256                 float* dataImagePointer2=NULL;
 257                 float* dataImageResultPointer=NULL;
 258
 259                 substractByType(dataImagePointer1, dataImagePointer2, dataImageResultPointer, imageData1, imageData2);
 260         }
 261         else if(imageType == VTK_DOUBLE)
 262         {
 263                 /*
 264                 images pointers
 265                 */
 266                 // pointers to get into the image
 267                 double* dataImagePointer1=NULL;
 268                 double* dataImagePointer2=NULL;
 269                 double* dataImageResultPointer=NULL;
 270
 271                 substractByType(dataImagePointer1, dataImagePointer2, dataImageResultPointer, imageData1, imageData2);
 272         }
 273 }
 274
 275 template <class T>
 276 void Substraction::substractByType(T* dataImagePointer1, T* dataImagePointer2, T* dataImageResultPointer, vtkImageData *imageData1, vtkImageData *imageData2)
 277 {
 278         // we start where the  image starts
 279         dataImagePointer1=(T*)imageData1->GetScalarPointer(0,0,0);
 280         dataImagePointer2=(T*)imageData2->GetScalarPointer(0,0,0);
 281         dataImageResultPointer=(T*)imageResult->GetScalarPointer(0,0,0);
 282         /*
 283         Image Size
 284         */
 285         int ext[6];
 286         imageData1->GetExtent(ext);
 287         int sx,sy,sz;
 288         sx=ext[1]-ext[0]+1;
 289         sy=ext[3]-ext[2]+1;
 290         sz=ext[5]-ext[4]+1;
 291
 292         sizeImage=sx*sy*sz;
 293         //-----------------
 294         //A3
 295         //-----------------
 296         //walking in the image
 297         int i=0,j=0,k=0,counter=0,nU=0,nL=0,nZ=0;
 298         double sum1=0,sum2=0;
 299         for(i=0;i<sx;i++)
 300         {
 301                 for(j=0;j<sy;j++)
 302                 {
 303                         for(k=0;k<sz;k++)
 304                         {
 305
 306                                 // this is for getting just the grey level in that position
 307                                 //originalValue=(short)imageData->GetScalarComponentAsFloat(i,j,k,0);
 308
 309                                 // we get the pointer to the position (i,j,k)y that way we can get the
 310                                 //grey level and we can change it
 311
 312                                 dataImagePointer1=(T*)imageData1->GetScalarPointer(i,j,k);
 313                                 dataImagePointer2=(T*)imageData2->GetScalarPointer(i,j,k);
 314                                 dataImageResultPointer=(T*)imageResult->GetScalarPointer(i,j,k);
 315
 316                                 sum1=(int)(dataImagePointer1[0]) + (int)(dataImagePointer1[1]) + (int)(dataImagePointer1[2]);
 317                                 sum1=sum1/3;
 318                                 sum2=(int)(dataImagePointer2[0]) + (int)(dataImagePointer2[1]) + (int)(dataImagePointer2[2]);
 319                                 sum2=sum2/3;
 320                                 if((sum1 - sum2) < lZeroLevel)
 321                                 {
 322                                         dataImageResultPointer[0] =(T) lowerColor[0];
 323                                         dataImageResultPointer[1] =(T) lowerColor[1];
 324                                         dataImageResultPointer[2] =(T) lowerColor[2];
 325                                         nL++;
 326                                 }
 327                                 else if((sum1 - sum2) > uZeroLevel)
 328                                 {
 329                                         dataImageResultPointer[0] =(T) upperColor[0];
 330                                         dataImageResultPointer[1] =(T) upperColor[1];
 331                                         dataImageResultPointer[2] =(T) upperColor[2];
 332                                         nU++;
 333                                 }
 334                                 else
 335                                 {
 336                                         dataImageResultPointer[0] =(T) mediumColor[0];
 337                                         dataImageResultPointer[1] =(T) mediumColor[1];
 338                                         dataImageResultPointer[2] =(T) mediumColor[2];
 339                                         nZ++;
 340                                 }
 341                                 counter++;
 342                         }
 343                 }
 344         }
 345 }
 346
 347 /*
 348 Returns the filtered image
 349 */
 350 vtkImageData* Substraction::getSubstractedImage()
 351 {
 352         return imageResult;
 353 }
 354
 355 /*
 356         Get Image Size
 357 */
 358 int Substraction::getImageSize()
 359 {
 360         return sizeImage;
 361 }