1 /*# ---------------------------------------------------------------------
3 # Copyright (c) CREATIS (Centre de Recherche en Acquisition et Traitement de l'Image
5 # Authors : Eduardo Davila, Frederic Cervenansky, Claire Mouton
6 # Previous Authors : Laurent Guigues, Jean-Pierre Roux
7 # CreaTools website : www.creatis.insa-lyon.fr/site/fr/creatools_accueil
9 # This software is governed by the CeCILL-B license under French law and
10 # abiding by the rules of distribution of free software. You can use,
11 # modify and/ or redistribute the software under the terms of the CeCILL-B
12 # license as circulated by CEA, CNRS and INRIA at the following URL
13 # http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
14 # or in the file LICENSE.txt.
16 # As a counterpart to the access to the source code and rights to copy,
17 # modify and redistribute granted by the license, users are provided only
18 # with a limited warranty and the software's author, the holder of the
19 # economic rights, and the successive licensors have only limited
22 # The fact that you are presently reading this means that you have had
23 # knowledge of the CeCILL-B license and that you accept its terms.
24 # ------------------------------------------------------------------------ */
26 #include "ContourExtractData.h"
29 //----------------------------------------------------------------------
30 ContourExtractData::ContourExtractData( bool okImagesResults)
32 this->imagedata = NULL;
33 imagedataValueResult = NULL;
34 imagedataMaskResult = NULL;
35 this->okImagesResults = okImagesResults;
39 // ------------------------------------------------------------------------
41 ContourExtractData::~ContourExtractData()
46 //----------------------------------------------------------------------
47 void ContourExtractData::SetImage( vtkImageData* imagedata)
49 this->imagedata = imagedata;
50 this->imagedata->GetScalarRange(scalarRange);
52 // RaC 20-11-09 Changes in InitLstContoursLinesYPoints
55 //EED 2017-01-01 Migration VTK7
56 #if VTK_MAJOR_VERSION <= 5
57 this->imagedata->GetWholeExtent(ext);
59 this->imagedata->GetExtent(ext);
62 _sizeImageY = ext[3]-ext[2]+1;
64 // init vtk image result : valuesImage maskImage
65 if (this->okImagesResults==true){ InitVtkImagesResult(); }
67 //----------------------------------------------------------------------
68 void ContourExtractData::SetZtoBeAnalys( int z )
73 //------------------------------------------------------------------------
74 void ContourExtractData::SetLstManualContourModel( std::vector<manualBaseModel*> lstManConMod)
76 this->lstManConMod = lstManConMod;
80 //------------------------------------------------------------------------
81 void ContourExtractData::GetMinMaxPoint(int *minPoint,
83 manualBaseModel *manualcontourmodel
87 //int np = manualcontourmodel->GetSizeLstPoints( ); // number of control points // JPRx
89 // JSTG 26-02-08 ---------------------------------------------------------------------------------------
90 //int nps = manualviewbaseecontour->GetNumberOfPointsSpline(); // number of points in the spline
91 int nps = manualcontourmodel->GetNumberOfPointsSpline(); // number of points in the spline
92 //------------------------------------------------------------------------------------------------------
94 // JSTG 26-02-08 ---------------------------------------------------------------------------------------
97 //double delta=( double ) ( np ) / ( double ) ( nps );
98 manualcontourmodel->UpdateSpline();
101 //t= delta * (double)i;
102 //manualcontourmodel->GetSplinePoint(t,x,y,z);
103 manualcontourmodel->GetSpline_i_Point(i,&x,&y,&z);
104 if (x<minPoint[0]){ minPoint[0]=(int)x; }
105 if (y<minPoint[1]){ minPoint[1]=(int)y; }
106 if (x>maxPoint[0]){ maxPoint[0]=(int)x; }
107 if (y>maxPoint[1]){ maxPoint[1]=(int)y; }
113 //------------------------------------------------------------------------------------------------------
116 //------------------------------------------------------------------------
117 void ContourExtractData::GetMinMaxPoint_Of_LstManConMod( int *minPoint,
122 int i,size = lstManConMod.size();
124 for(i=0 ; i<size ; i++)
126 GetMinMaxPoint(minPoint,maxPoint,lstManConMod[i]);
131 //------------------------------------------------------------------------
132 int ContourExtractData::AnalisisContourInsideV2(int x, int y, int iContour )
138 int nps=_lstlstlstVecX1[iContour][y].size();
141 double borderX, borderY;
142 double xx1, yy1,xx2, yy2;
146 for (i=0; i<nps; i++)
148 x1=_lstlstlstVecX1[iContour][y][i];
149 y1=_lstlstlstVecY1[iContour][y][i];
150 x2=_lstlstlstVecX2[iContour][y][i];
151 y2=_lstlstlstVecY2[iContour][y][i];
158 xx1=x1; yy1=y1; xx2=x2; yy2=y2;
160 xx1=x2; yy1=y2; xx2=x1; yy2=y1;
163 double difxx2xx1=fabs(xx2-xx1);
164 if (difxx2xx1==0) difxx2xx1=0.0000000001;
166 // Finding border looking in vertical direction AND verifing if pixel is at right of the line
167 if ( (yy>=yy1)&&(yy<=yy2) )
169 //by triangle similarity
170 d = ( fabs(xx2-xx1)*(yy-yy1) ) / (yy2-yy1) ;
171 if ( (xx1<=xx2)&&(x<(xx1+d)) ) { result++; }
172 if ( (xx1>xx2)&&(x<(xx1-d)) ) { result++; }
174 if ( (yy2-yy1)/difxx2xx1 >= 1.0)
185 } // if point inside y
188 // Finding border looking in vertical direction
189 if ( ((xx1<=xx2)&&(xx>=xx1)&&(xx<xx2)) || ((xx1>xx2)&&(xx>=xx2)&&(xx<xx1)) )
191 if ( (yy2-yy1)/difxx2xx1 <= 1.0)
193 //by triangle similarity
194 d = ( fabs(xx1-xx)*(yy2-yy1) ) / difxx2xx1;
200 } // if point inside x
203 //Border verification
204 if ( (x==(int)borderX) && (y==(int)borderY) ) { inBorder=true; }// if point in border
206 // Verification : border in horizontal line
207 if ( ((int)y1==(int)y2) && ((int)y1==y) && (x1<x2) && (x>=x1) && (x<=x2)) { inBorder=true; }
208 if ( ((int)y1==(int)y2) && ((int)y1==y) && (x2<x1) && (x>=x2) && (x<=x1)) { inBorder=true; }
210 if (inBorder==true){ i=nps; }
213 if (inBorder==true) { result=1; }
220 //------------------------------------------------------------------------
221 // typeOperation=0 AND
222 // typeOperation=1 OR
223 // typeOperation=2 XOR
224 bool ContourExtractData::isInside(int x, int y, int typeOperation)
228 int i,size = this->lstManConMod.size();
229 int numberInside = 0;
231 /* RaC 20-11-09 (C1) Changes to use the method without the image.
233 imagedata->GetExtent(ext);
235 if ((x>=0) && (x<=ext[1]) && (y>=0) && (y<=ext[3]))
239 if (typeOperation==0) // AND Intersection
243 // To process the statistics of the Points contour the procedure is different
245 manualBaseModel *mbm = lstManConMod[i];
246 if(mbm->GetTypeModel()==7)
248 if(mbm->IsPoint(x,y)==true)
255 numberLeft = AnalisisContourInsideV2(x,y, i );
258 if ( (numberLeft % 2) ==1){ numberInside++; }
260 if ( numberInside == (size) ){ result=true; }
261 } // AND Intersection
265 if (typeOperation==1) // OR All
269 // To process the statistics of the Points contour the procedure is different
271 manualBaseModel *mbm = lstManConMod[i];
272 if(mbm->GetTypeModel()==7)
274 if(mbm->IsPoint(x,y)==true)
281 numberLeft = AnalisisContourInsideV2(x,y, i );
283 if ( (numberLeft % 2) ==1){ result=true; }
289 if (typeOperation==2) // XOR crown
293 // To process the statistics of the Points contour the procedure is different
295 manualBaseModel *mbm = lstManConMod[i];
296 if(mbm->GetTypeModel()==7)
298 if(mbm->IsPoint(x,y)==true)
305 numberLeft = numberLeft + AnalisisContourInsideV2(x,y, i );
309 if ( numberLeft % 2 ==1){ result = true; }
319 //------------------------------------------------------------------------
321 double ContourExtractData::GetDataValue(int x, int y, int z)
323 // wxVtk2DBaseView *wxvtk2dbaseview = (wxVtk2DBaseView*)wxvtkbaseview;
324 // int z = (int)wxvtk2dbaseview->GetVtkBaseData()->GetZ();
326 //EED OJO avec JS _zz = z;
330 p = imagedata->GetScalarPointer(x,y,z);
332 if (imagedata->GetScalarType()==VTK_CHAR)
335 result = (double)(*pp);
337 else if (imagedata->GetScalarType()==VTK_SIGNED_CHAR)
339 signed char *pp = (signed char*)p;
340 result = (double)(*pp);
342 else if (imagedata->GetScalarType()==VTK_UNSIGNED_CHAR)
344 unsigned char *pp = (unsigned char*)p;
345 result = (double)(*pp);
347 else if (imagedata->GetScalarType()==VTK_SHORT)
349 short *pp = (short*)p;
350 result = (double)(*pp);
352 else if (imagedata->GetScalarType()==VTK_UNSIGNED_SHORT)
354 unsigned short *pp = (unsigned short*)p;
355 result = (double)(*pp);
357 else if (imagedata->GetScalarType()==VTK_INT)
360 result = (double)(*pp);
362 else if (imagedata->GetScalarType()==VTK_UNSIGNED_INT)
364 unsigned int *pp = (unsigned int*)p;
365 result = (double)(*pp);
367 else if (imagedata->GetScalarType()==VTK_LONG)
370 result = (double)(*pp);
372 else if (imagedata->GetScalarType()==VTK_UNSIGNED_LONG)
374 unsigned long *pp = (unsigned long*)p;
375 result = (double)(*pp);
377 else if (imagedata->GetScalarType()==VTK_FLOAT)
379 float *pp = (float*)p;
380 result = (double)(*pp);
382 else if (imagedata->GetScalarType()==VTK_DOUBLE)
384 double *pp = (double*)p;
385 result = (double)(*pp);
391 //------------------------------------------------------------------------
393 void ContourExtractData::PutVtkImageDataResultValue( int x, int y, int z, double value )
395 unsigned short *pValue;
396 unsigned short *pMask;
398 imagedataValueResult->SetScalarComponentFromDouble(x,y,z,0,value);
399 imagedataMaskResult->SetScalarComponentFromDouble(x,y,z,0,255);
400 // pValue = (unsigned short *)imagedataValueResult->GetScalarPointer(x,y,z);
401 // *pValue = (unsigned short)value;
402 // pMask = (unsigned char *)imagedataMaskResult->GetScalarPointer(x,y,z);
406 //------------------------------------------------------------------------
407 void ContourExtractData::ResetImageResult(int z)
409 if (okImagesResults==true)
411 unsigned short *pValue;
412 unsigned short *pMask;
413 pValue = (unsigned short *)imagedataValueResult->GetScalarPointer(0,0,z);
414 pMask = (unsigned short *)imagedataMaskResult->GetScalarPointer(0,0,z);
417 imagedataValueResult->GetExtent(ext);
419 int size = (ext[1]-ext[0]+1) * (ext[3]-ext[2]+1);
420 memset(pValue,0, size*imagedataValueResult->GetScalarSize() );
421 memset(pMask,0, size*imagedataMaskResult->GetScalarSize() );
426 //------------------------------------------------------------------------
427 void ContourExtractData::CalculateImageResult()
429 if (okImagesResults==true)
431 ResetImageResult(zImage);
438 minPoint[0] = 999999;
439 minPoint[1] = 999999;
440 maxPoint[0] = -999999;
441 maxPoint[1] = -999999;
443 GetMinMaxPoint_Of_LstManConMod(minPoint,maxPoint);
444 InitLstContoursLinesYPoints();
446 for (j=minPoint[1]; j<=maxPoint[1]; j++)
448 for (i=minPoint[0]; i<=maxPoint[0]; i++)
451 //RaC 20-11-09 Changes specified in isInside (C1)
453 imagedata->GetExtent(ext);
455 if ((i>=0) && (i<=ext[1]) && (j>=0) && (j<=ext[3]))
457 if (isInside(i,j,_typeOperation)==true)
459 value = GetDataValue(i,j,zImage);
460 if ( (value>=scalarRange[0]) && (value<=scalarRange[1]) )
462 PutVtkImageDataResultValue(i,j,zImage, value );
471 imagedataValueResult->Modified();
472 imagedataMaskResult->Modified();
474 //EED 2017-01-01 Migration VTK7
475 #if VTK_MAJOR_VERSION <= 5
476 imagedataValueResult->Update();
477 imagedataMaskResult->Update();
486 //------------------------------------------------------------------------
487 void ContourExtractData::GetValuesInsideCrown( int *numberOfPixels,
488 std::vector<double> *pLstValue,
489 std::vector<double> *pLstValuePosX,
490 std::vector<double> *pLstValuePosY,
491 std::vector<double> *pLstValuePosZ)
494 pLstValuePosX->clear();
495 pLstValuePosY->clear();
496 pLstValuePosZ->clear();
498 // if (okImagesResults==true)
500 // ResetImageResult(zImage);
509 minPoint[0] = 999999;
510 minPoint[1] = 999999;
511 maxPoint[0] = -999999;
512 maxPoint[1] = -999999;
514 GetMinMaxPoint_Of_LstManConMod(minPoint,maxPoint);
515 InitLstContoursLinesYPoints();
517 for (j=minPoint[1]; j<=maxPoint[1]; j++)
519 for (i=minPoint[0]; i<=maxPoint[0]; i++)
521 //RaC 20-11-09 Changes specified in isInside (C1)
523 imagedata->GetExtent(ext);
525 if ((i>=0) && (i<=ext[1]) && (j>=0) && (j<=ext[3]))
527 if (isInside(i,j,_typeOperation)==true)
531 value = GetDataValue(i,j,zImage);
532 if ( (value>=scalarRange[0]) && (value<=scalarRange[1]) )
534 pLstValue -> push_back( value );
535 pLstValuePosX -> push_back( i );
536 pLstValuePosY -> push_back( j );
537 pLstValuePosZ -> push_back( -1 );
544 *numberOfPixels = acum;
547 //------------------------------------------------------------------------
549 vtkImageData *ContourExtractData::GetVtkImageValueResult()
551 return imagedataValueResult;
553 //------------------------------------------------------------------------
554 vtkImageData *ContourExtractData::GetVtkImageMaskResult()
556 return imagedataMaskResult;
558 // ------------------------------------------------------------------------
559 void ContourExtractData::InitVtkImagesResult()
566 imagedata->GetSpacing(spc);
567 imagedata->GetExtent(ext);
568 newDim[0]=ext[1]-ext[0]+1;
569 newDim[1]=ext[3]-ext[2]+1;
570 newDim[2]=ext[5]-ext[4]+1;
571 scalartype = imagedata->GetScalarType();
573 if (imagedataValueResult!=NULL)
575 imagedataValueResult->Delete();
577 imagedataValueResult = vtkImageData::New();
578 imagedataValueResult->SetSpacing(spc);
579 imagedataValueResult->SetDimensions( newDim );
581 //EED 2017-01-01 Migration VTK7
582 #if VTK_MAJOR_VERSION <= 5
585 imagedataValueResult->SetScalarType(scalartype);
586 //imagedataValueResult->SetScalarTypeToUnsignedShort();
587 imagedataValueResult->AllocateScalars();
589 imagedataValueResult->AllocateScalars(scalartype,1);
593 if (imagedataMaskResult!=NULL)
595 imagedataMaskResult->Delete();
597 imagedataMaskResult = vtkImageData::New();
600 imagedataMaskResult->SetSpacing(spc);
601 imagedataMaskResult->SetDimensions( newDim );
603 //EED 2017-01-01 Migration VTK7
604 #if VTK_MAJOR_VERSION <= 5
606 //imagedataMaskResult->SetScalarTypeToUnsignedShort();
607 imagedataMaskResult->SetScalarTypeToUnsignedChar();
608 imagedataMaskResult->AllocateScalars();
610 imagedataMaskResult->AllocateScalars(VTK_UNSIGNED_CHAR,1);
616 //------------------------------------------------------------------------
617 void ContourExtractData::InitVolumeStatistics()
621 vol_minValue = 9999999;
622 vol_maxValue =-9999999;
623 vol_acum_average = 0;
624 vol_acum_standardeviation = 0;
627 //------------------------------------------------------------------------
628 void ContourExtractData::SetVolumeStatistics(int rCountRange,
633 double acum_standardeviation)
635 vol_rCountRange = vol_rCountRange + rCountRange;
636 vol_rsize = vol_rsize + rsize;
638 if (minValue<vol_minValue){ vol_minValue = minValue; }
639 if (maxValue>vol_maxValue){ vol_maxValue = maxValue; }
641 vol_acum_average = vol_acum_average + acum_average;
642 vol_acum_standardeviation = vol_acum_standardeviation + acum_standardeviation;
645 //------------------------------------------------------------------------
646 void ContourExtractData::GetVolumeStatistics(int *vol_rCountRange,
648 double *vol_minValue,
649 double *vol_maxValue,
651 double *vol_standardeviation)
653 *vol_rCountRange = this->vol_rCountRange;
654 *vol_rsize = this->vol_rsize;
655 *vol_minValue = this->vol_minValue;
656 *vol_maxValue = this->vol_maxValue;
657 *vol_average = this->vol_acum_average / this->vol_rsize;
658 *vol_standardeviation = sqrt(this->vol_acum_standardeviation / this->vol_rsize);
662 //------------------------------------------------------------------------
663 void ContourExtractData::Statistics( std::vector<double> *inputLstValue,
671 double *rstandardeviation
677 double standardeviation = 0;
678 double acum_average = 0;
679 double acum_standardeviation = 0;
684 if (inputLstValue!=NULL)
686 size=inputLstValue->size();
688 max=(*inputLstValue)[0];
689 min=(*inputLstValue)[0];
690 // Average , countRange
692 for ( i=0; i<size; i++ )
694 ng=(*inputLstValue)[i];
695 acum_average = acum_average + ng;
696 if (max<ng) max=ng; // Max
697 if (min>ng) min=ng; // Min
698 if ((ng>=grayRangeMin) && (ng<=grayRangeMax)) countRange++; // countRange
700 average = acum_average / size;
703 acum_standardeviation=0;
705 for ( i=0; i<size; i++ )
707 tmp = (*inputLstValue)[i] - average;
708 acum_standardeviation = acum_standardeviation + tmp*tmp;
709 } // for standar deviation
710 standardeviation = sqrt(acum_standardeviation/size);
711 SetVolumeStatistics(countRange, (*rsize),
713 acum_average,acum_standardeviation);
719 *rCountRange = countRange;
723 *rstandardeviation = standardeviation;
726 //------------------------------------------------------------------------
727 void ContourExtractData::SetTypeOperation(int type)
732 //------------------------------------------------------------------------
733 void ContourExtractData::Fill_lstlstlstVecXY(int iContour, int sizeY)
736 double x1,y1,z1,x2,y2,z2;
737 manualBaseModel *manualcontourmodel= lstManConMod[iContour];
738 int nps = manualcontourmodel->GetNumberOfPointsSpline(); // number of points in the spline
739 manualcontourmodel->UpdateSpline();
740 //------------------------------------------------------------------------------------------------------
742 for (y=0;y<sizeY;y++)
744 manualcontourmodel->GetSpline_i_Point(0,&x1,&y1,&z1);
745 x1=x1+0.5; y1=y1+0.5;
746 for (i=1; i<nps; i++)
748 manualcontourmodel->GetSpline_i_Point(i,&x2,&y2,&z2);
749 x2=x2+0.5; y2=y2+0.5;
750 if ( ((y1<y2)&&(y>=y1)&&(y<=y2)) || ((y1>y2)&&(y>=y2)&&(y<=y1)) || ((int)y1==y) || ((int)y2==y) )
752 _lstlstlstVecX1[iContour][y].push_back(x1);
753 _lstlstlstVecY1[iContour][y].push_back(y1);
754 _lstlstlstVecX2[iContour][y].push_back(x2);
755 _lstlstlstVecY2[iContour][y].push_back(y2);
758 } // for i Points in spline
763 void ContourExtractData::InitLstContoursLinesYPoints()
765 // init InInside Optimisation
768 _lstlstlstVecX1.clear();
769 _lstlstlstVecY1.clear();
770 _lstlstlstVecX2.clear();
771 _lstlstlstVecY2.clear();
775 this->imagedata->GetWholeExtent(ext);
776 int sizeY = ext[3]-ext[2]+1;
778 std::vector<double> vecDouble;
779 std::vector< std::vector<double> > vecVecDouble;
780 for ( i=0 ; i<_sizeImageY ; i++ )
782 vecVecDouble.push_back( vecDouble );
785 //Fill structure with points
786 int sizeContours = lstManConMod.size();
787 for( i=0 ; i<sizeContours ; i++ )
789 _lstlstlstVecX1.push_back( vecVecDouble );
790 _lstlstlstVecY1.push_back( vecVecDouble );
791 _lstlstlstVecX2.push_back( vecVecDouble );
792 _lstlstlstVecY2.push_back( vecVecDouble );
793 Fill_lstlstlstVecXY(i,_sizeImageY);
798 void ContourExtractData::SetScalarRange(double min, double max)
804 void ContourExtractData::SetSizeImageY(int pSizeImageY)
806 _sizeImageY=pSizeImageY;