2 #include "ContourExtractData.h"
5 //----------------------------------------------------------------------
6 ContourExtractData::ContourExtractData( bool okImagesResults)
8 this->imagedata = NULL;
9 imagedataValueResult = NULL;
10 imagedataMaskResult = NULL;
11 this->okImagesResults = okImagesResults;
15 // ------------------------------------------------------------------------
17 ContourExtractData::~ContourExtractData()
22 //----------------------------------------------------------------------
23 void ContourExtractData::SetImage( vtkImageData* imagedata)
25 this->imagedata = imagedata;
27 // init vtk image result : valuesImage maskImage
28 if (this->okImagesResults==true){ InitVtkImagesResult(); }
30 //----------------------------------------------------------------------
31 void ContourExtractData::SetZtoBeAnalys( int z )
36 //------------------------------------------------------------------------
37 void ContourExtractData::SetLstManualContourModel( std::vector<manualBaseModel*> lstManConMod)
39 this->lstManConMod = lstManConMod;
43 //------------------------------------------------------------------------
44 void ContourExtractData::GetMinMaxPoint(int *minPoint,
46 manualBaseModel *manualcontourmodel
50 //int np = manualcontourmodel->GetSizeLstPoints( ); // number of control points // JPRx
52 // JSTG 26-02-08 ---------------------------------------------------------------------------------------
53 //int nps = manualviewbaseecontour->GetNumberOfPointsSpline(); // number of points in the spline
54 int nps = manualcontourmodel->GetNumberOfPointsSpline(); // number of points in the spline
55 //------------------------------------------------------------------------------------------------------
57 // JSTG 26-02-08 ---------------------------------------------------------------------------------------
60 //double delta=( double ) ( np ) / ( double ) ( nps );
61 manualcontourmodel->UpdateSpline();
64 //t= delta * (double)i;
65 //manualcontourmodel->GetSplinePoint(t,x,y,z);
66 manualcontourmodel->GetSpline_i_Point(i,&x,&y,&z);
67 if (x<minPoint[0]){ minPoint[0]=(int)x; }
68 if (y<minPoint[1]){ minPoint[1]=(int)y; }
69 if (x>maxPoint[0]){ maxPoint[0]=(int)x; }
70 if (y>maxPoint[1]){ maxPoint[1]=(int)y; }
76 //------------------------------------------------------------------------------------------------------
79 //------------------------------------------------------------------------
80 void ContourExtractData::GetMinMaxPoint_Of_LstManConMod( int *minPoint,
85 int i,size = lstManConMod.size();
87 for(i=0 ; i<size ; i++)
89 GetMinMaxPoint(minPoint,maxPoint,lstManConMod[i]);
94 //------------------------------------------------------------------------
95 int ContourExtractData::AnalisisContourInsideV2(int x, int y, int iContour )
101 int nps=_lstlstlstVecX1[iContour][y].size();
104 double borderX, borderY;
105 double xx1, yy1,xx2, yy2;
109 for (i=0; i<nps; i++)
111 x1=_lstlstlstVecX1[iContour][y][i];
112 y1=_lstlstlstVecY1[iContour][y][i];
113 x2=_lstlstlstVecX2[iContour][y][i];
114 y2=_lstlstlstVecY2[iContour][y][i];
121 xx1=x1; yy1=y1; xx2=x2; yy2=y2;
123 xx1=x2; yy1=y2; xx2=x1; yy2=y1;
126 double difxx2xx1=fabs(xx2-xx1);
127 if (difxx2xx1==0) difxx2xx1=0.0000000001;
129 // Finding border looking in vertical direction AND verifing if pixel is at right of the line
130 if ( (yy>=yy1)&&(yy<=yy2) )
132 //by triangle similarity
133 d = ( fabs(xx2-xx1)*(yy-yy1) ) / (yy2-yy1) ;
134 if ( (xx1<=xx2)&&(x<(xx1+d)) ) { result++; }
135 if ( (xx1>xx2)&&(x<(xx1-d)) ) { result++; }
137 if ( (yy2-yy1)/difxx2xx1 >= 1.0)
148 } // if point inside y
151 // Finding border looking in vertical direction
152 if ( ((xx1<=xx2)&&(xx>=xx1)&&(xx<xx2)) || ((xx1>xx2)&&(xx>=xx2)&&(xx<xx1)) )
154 if ( (yy2-yy1)/difxx2xx1 <= 1.0)
156 //by triangle similarity
157 d = ( fabs(xx1-xx)*(yy2-yy1) ) / difxx2xx1;
163 } // if point inside x
166 //Border verification
167 if ( (x==(int)borderX) && (y==(int)borderY) ) { inBorder=true; }// if point in border
169 // Verification : border in horizontal line
170 if ( ((int)y1==(int)y2) && ((int)y1==y) && (x1<x2) && (x>=x1) && (x<=x2)) { inBorder=true; }
171 if ( ((int)y1==(int)y2) && ((int)y1==y) && (x2<x1) && (x>=x2) && (x<=x1)) { inBorder=true; }
173 if (inBorder==true){ i=nps; }
176 if (inBorder==true) { result=1; }
183 //------------------------------------------------------------------------
184 bool ContourExtractData::isInside(int x, int y, int typeOperation)
188 int i,size = this->lstManConMod.size();
189 int numberInside = 0;
192 imagedata->GetExtent(ext);
194 if ((x>=0) && (x<=ext[1]) && (y>=0) && (y<=ext[3]))
197 if (typeOperation==0) // AND Intersection
201 // To process the statistics of the Points contour the procedure is different
203 manualBaseModel *mbm = lstManConMod[i];
204 if(mbm->GetTypeModel()==7)
206 if(mbm->IsPoint(x,y)==true)
213 numberLeft = AnalisisContourInsideV2(x,y, i );
216 if ( (numberLeft % 2) ==1){ numberInside++; }
218 if ( numberInside == (size) ){ result=true; }
219 } // AND Intersection
223 if (typeOperation==1) // OR All
227 // To process the statistics of the Points contour the procedure is different
229 manualBaseModel *mbm = lstManConMod[i];
230 if(mbm->GetTypeModel()==7)
232 if(mbm->IsPoint(x,y)==true)
239 numberLeft = AnalisisContourInsideV2(x,y, i );
241 if ( (numberLeft % 2) ==1){ result=true; }
247 if (typeOperation==2) // XOR crown
251 // To process the statistics of the Points contour the procedure is different
253 manualBaseModel *mbm = lstManConMod[i];
254 if(mbm->GetTypeModel()==7)
256 if(mbm->IsPoint(x,y)==true)
263 numberLeft = numberLeft + AnalisisContourInsideV2(x,y, i );
267 if ( numberLeft % 2 ==1){ result = true; }
277 //------------------------------------------------------------------------
279 double ContourExtractData::GetDataValue(int x, int y, int z)
281 // wxVtk2DBaseView *wxvtk2dbaseview = (wxVtk2DBaseView*)wxvtkbaseview;
282 // int z = (int)wxvtk2dbaseview->GetVtkBaseData()->GetZ();
284 //EED OJO avec JS _zz = z;
288 p = imagedata->GetScalarPointer(x,y,z);
290 if (imagedata->GetScalarType()==VTK_CHAR)
293 result = (double)(*pp);
295 else if (imagedata->GetScalarType()==VTK_SIGNED_CHAR)
297 signed char *pp = (signed char*)p;
298 result = (double)(*pp);
300 else if (imagedata->GetScalarType()==VTK_UNSIGNED_CHAR)
302 unsigned char *pp = (unsigned char*)p;
303 result = (double)(*pp);
305 else if (imagedata->GetScalarType()==VTK_SHORT)
307 short *pp = (short*)p;
308 result = (double)(*pp);
310 else if (imagedata->GetScalarType()==VTK_UNSIGNED_SHORT)
312 unsigned short *pp = (unsigned short*)p;
313 result = (double)(*pp);
315 else if (imagedata->GetScalarType()==VTK_INT)
318 result = (double)(*pp);
320 else if (imagedata->GetScalarType()==VTK_UNSIGNED_INT)
322 unsigned int *pp = (unsigned int*)p;
323 result = (double)(*pp);
325 else if (imagedata->GetScalarType()==VTK_LONG)
328 result = (double)(*pp);
330 else if (imagedata->GetScalarType()==VTK_UNSIGNED_LONG)
332 unsigned long *pp = (unsigned long*)p;
333 result = (double)(*pp);
335 else if (imagedata->GetScalarType()==VTK_FLOAT)
337 float *pp = (float*)p;
338 result = (double)(*pp);
340 else if (imagedata->GetScalarType()==VTK_DOUBLE)
342 double *pp = (double*)p;
343 result = (double)(*pp);
349 //------------------------------------------------------------------------
351 void ContourExtractData::PutVtkImageDataResultValue( int x, int y, int z, double value )
353 unsigned short *pValue;
354 unsigned short *pMask;
355 pValue = (unsigned short *)imagedataValueResult->GetScalarPointer(x,y,z);
356 pMask = (unsigned short *)imagedataMaskResult->GetScalarPointer(x,y,z);
358 *pValue = (unsigned short)value;
361 //------------------------------------------------------------------------
362 void ContourExtractData::ResetImageResult(int z)
364 if (okImagesResults==true)
366 unsigned short *pValue;
367 unsigned short *pMask;
368 pValue = (unsigned short *)imagedataValueResult->GetScalarPointer(0,0,z);
369 pMask = (unsigned short *)imagedataMaskResult->GetScalarPointer(0,0,z);
372 imagedataValueResult->GetExtent(ext);
374 int size = (ext[1]-ext[0]+1) * (ext[3]-ext[2]+1);
375 memset(pValue,0,size*2);
376 memset(pMask,0,size*2);
381 //------------------------------------------------------------------------
382 void ContourExtractData::CalculateImageResult()
384 if (okImagesResults==true)
386 ResetImageResult(zImage);
393 minPoint[0] = 999999;
394 minPoint[1] = 999999;
395 maxPoint[0] = -999999;
396 maxPoint[1] = -999999;
398 GetMinMaxPoint_Of_LstManConMod(minPoint,maxPoint);
399 InitLstContoursLinesYPoints();
401 for (j=minPoint[1]; j<=maxPoint[1]; j++)
403 for (i=minPoint[0]; i<=maxPoint[0]; i++)
405 if (isInside(i,j,_typeOperation)==true)
407 value = GetDataValue(i,j,zImage);
408 PutVtkImageDataResultValue(i,j,zImage, value );
414 imagedataValueResult->Modified();
415 imagedataMaskResult->Modified();
416 imagedataValueResult->Update();
417 imagedataMaskResult->Update();
422 //------------------------------------------------------------------------
423 void ContourExtractData::GetValuesInsideCrown(std::vector<double> *pLstValue,
424 std::vector<double> *pLstValuePosX,
425 std::vector<double> *pLstValuePosY,
426 std::vector<double> *pLstValuePosZ)
429 pLstValuePosX->clear();
430 pLstValuePosY->clear();
431 pLstValuePosZ->clear();
433 // if (okImagesResults==true)
435 // ResetImageResult(zImage);
444 minPoint[0] = 999999;
445 minPoint[1] = 999999;
446 maxPoint[0] = -999999;
447 maxPoint[1] = -999999;
449 GetMinMaxPoint_Of_LstManConMod(minPoint,maxPoint);
450 InitLstContoursLinesYPoints();
452 for (j=minPoint[1]; j<=maxPoint[1]; j++)
454 for (i=minPoint[0]; i<=maxPoint[0]; i++)
456 if (isInside(i,j,_typeOperation)==true)
458 value = GetDataValue(i,j,zImage);
461 // if (okImagesResults==true){
462 // PutVtkImageDataResultValue(i,j,zImage, value );
465 pLstValue -> push_back( value );
466 pLstValuePosX -> push_back( i );
467 pLstValuePosY -> push_back( j );
468 pLstValuePosZ -> push_back( -1 );
474 // if (this->okImagesResults==true){
475 // imagedataValueResult->Modified();
476 // imagedataMaskResult->Modified();
482 //------------------------------------------------------------------------
484 vtkImageData *ContourExtractData::GetVtkImageValueResult()
486 return imagedataValueResult;
488 //------------------------------------------------------------------------
489 vtkImageData *ContourExtractData::GetVtkImageMaskResult()
491 return imagedataMaskResult;
493 // ------------------------------------------------------------------------
494 void ContourExtractData::InitVtkImagesResult()
501 imagedata->GetSpacing(spc);
502 imagedata->GetExtent(ext);
503 newDim[0]=ext[1]-ext[0]+1;
504 newDim[1]=ext[3]-ext[2]+1;
505 newDim[2]=ext[5]-ext[4]+1;
506 scalartype = imagedata->GetScalarType();
508 if (imagedataValueResult!=NULL)
510 imagedataValueResult->Delete();
512 imagedataValueResult = vtkImageData::New();
513 // imagedataValueResult->SetScalarType(scalartype);
514 imagedataValueResult->SetScalarTypeToUnsignedShort();
515 imagedataValueResult->SetSpacing(spc);
516 imagedataValueResult->SetDimensions( newDim );
517 imagedataValueResult->AllocateScalars();
519 if (imagedataMaskResult!=NULL)
521 imagedataMaskResult->Delete();
523 imagedataMaskResult = vtkImageData::New();
524 // imagedataMaskResult->SetScalarType(scalartype);
525 imagedataMaskResult->SetScalarTypeToUnsignedShort();
526 imagedataMaskResult->SetSpacing(spc);
527 imagedataMaskResult->SetDimensions( newDim );
528 imagedataMaskResult->AllocateScalars();
532 //------------------------------------------------------------------------
533 void ContourExtractData::InitVolumeStatistics()
537 vol_minValue = 9999999;
538 vol_maxValue =-9999999;
539 vol_acum_average = 0;
540 vol_acum_standardeviation = 0;
543 //------------------------------------------------------------------------
544 void ContourExtractData::SetVolumeStatistics(int rCountRange,
549 double acum_standardeviation)
551 vol_rCountRange = vol_rCountRange + rCountRange;
552 vol_rsize = vol_rsize + rsize;
554 if (minValue<vol_minValue){ vol_minValue = minValue; }
555 if (maxValue>vol_maxValue){ vol_maxValue = maxValue; }
557 vol_acum_average = vol_acum_average + acum_average;
558 vol_acum_standardeviation = vol_acum_standardeviation + acum_standardeviation;
561 //------------------------------------------------------------------------
562 void ContourExtractData::GetVolumeStatistics(int *vol_rCountRange,
564 double *vol_minValue,
565 double *vol_maxValue,
567 double *vol_standardeviation)
569 *vol_rCountRange = this->vol_rCountRange;
570 *vol_rsize = this->vol_rsize;
571 *vol_minValue = this->vol_minValue;
572 *vol_maxValue = this->vol_maxValue;
573 *vol_average = this->vol_acum_average / this->vol_rsize;
574 *vol_standardeviation = sqrt(this->vol_acum_standardeviation / this->vol_rsize);
578 //------------------------------------------------------------------------
579 void ContourExtractData::Statistics( std::vector<double> *inputLstValue,
587 double *rstandardeviation
593 double standardeviation = 0;
594 double acum_average = 0;
595 double acum_standardeviation = 0;
600 if (inputLstValue!=NULL)
602 size=inputLstValue->size();
604 max=(*inputLstValue)[0];
605 min=(*inputLstValue)[0];
606 // Average , countRange
608 for ( i=0; i<size; i++ )
610 ng=(*inputLstValue)[i];
611 acum_average = acum_average + ng;
612 if (max<ng) max=ng; // Max
613 if (min>ng) min=ng; // Min
614 if ((ng>=grayRangeMin) && (ng<=grayRangeMax)) countRange++; // countRange
616 average = acum_average / size;
619 acum_standardeviation=0;
621 for ( i=0; i<size; i++ )
623 tmp = (*inputLstValue)[i] - average;
624 acum_standardeviation = acum_standardeviation + tmp*tmp;
625 } // for standar deviation
626 standardeviation = sqrt(acum_standardeviation/size);
627 SetVolumeStatistics(countRange, size,
629 acum_average,acum_standardeviation);
635 *rCountRange = countRange;
639 *rstandardeviation = standardeviation;
642 //------------------------------------------------------------------------
643 void ContourExtractData::SetTypeOperation(int type)
648 //------------------------------------------------------------------------
649 void ContourExtractData::Fill_lstlstlstVecXY(int iContour, int sizeY)
652 double x1,y1,z1,x2,y2,z2;
653 manualBaseModel *manualcontourmodel= lstManConMod[iContour];
654 int nps = manualcontourmodel->GetNumberOfPointsSpline(); // number of points in the spline
655 manualcontourmodel->UpdateSpline();
656 //------------------------------------------------------------------------------------------------------
658 for (y=0;y<sizeY;y++)
660 manualcontourmodel->GetSpline_i_Point(0,&x1,&y1,&z1);
661 x1=x1+0.5; y1=y1+0.5;
662 for (i=1; i<nps; i++)
664 manualcontourmodel->GetSpline_i_Point(i,&x2,&y2,&z2);
665 x2=x2+0.5; y2=y2+0.5;
666 if ( ((y1<y2)&&(y>=y1)&&(y<=y2)) || ((y1>y2)&&(y>=y2)&&(y<=y1)) || ((int)y1==y) || ((int)y2==y) )
668 _lstlstlstVecX1[iContour][y].push_back(x1);
669 _lstlstlstVecY1[iContour][y].push_back(y1);
670 _lstlstlstVecX2[iContour][y].push_back(x2);
671 _lstlstlstVecY2[iContour][y].push_back(y2);
674 } // for i Points in spline
679 void ContourExtractData::InitLstContoursLinesYPoints()
681 // init InInside Optimisation
684 _lstlstlstVecX1.clear();
685 _lstlstlstVecY1.clear();
686 _lstlstlstVecX2.clear();
687 _lstlstlstVecY2.clear();
690 this->imagedata->GetWholeExtent(ext);
691 int sizeY = ext[3]-ext[2]+1;
692 std::vector<double> vecDouble;
693 std::vector< std::vector<double> > vecVecDouble;
694 for ( i=0 ; i<sizeY ; i++ )
696 vecVecDouble.push_back( vecDouble );
699 //Fill structure with points
700 int sizeContours = lstManConMod.size();
701 for( i=0 ; i<sizeContours ; i++ )
703 _lstlstlstVecX1.push_back( vecVecDouble );
704 _lstlstlstVecY1.push_back( vecVecDouble );
705 _lstlstlstVecX2.push_back( vecVecDouble );
706 _lstlstlstVecY2.push_back( vecVecDouble );
707 Fill_lstlstlstVecXY(i,sizeY);