]> Creatis software - creaMaracasVisu.git/blobdiff - lib/maracasVisuLib/src/interface/wxWindows/Contour/ContourExtractData.cxx
#3162 creaMaracasVisu Bug New Normal - Threshold layer
[creaMaracasVisu.git] / lib / maracasVisuLib / src / interface / wxWindows / Contour / ContourExtractData.cxx
index e64fa03a444370e831d5e8e3162e631fd35192d8..c08f5d977fa6307339d3bcfb6c9b189681d9a98d 100644 (file)
@@ -1,60 +1,90 @@
+/*# ---------------------------------------------------------------------
+#
+# Copyright (c) CREATIS (Centre de Recherche en Acquisition et Traitement de l'Image
+#                        pour la Sant�)
+# Authors : Eduardo Davila, Frederic Cervenansky, Claire Mouton
+# Previous Authors : Laurent Guigues, Jean-Pierre Roux
+# CreaTools website : www.creatis.insa-lyon.fr/site/fr/creatools_accueil
+#
+#  This software is governed by the CeCILL-B license under French law and
+#  abiding by the rules of distribution of free software. You can  use,
+#  modify and/ or redistribute the software under the terms of the CeCILL-B
+#  license as circulated by CEA, CNRS and INRIA at the following URL
+#  http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
+#  or in the file LICENSE.txt.
+#
+#  As a counterpart to the access to the source code and  rights to copy,
+#  modify and redistribute granted by the license, users are provided only
+#  with a limited warranty  and the software's author,  the holder of the
+#  economic rights,  and the successive licensors  have only  limited
+#  liability.
+#
+#  The fact that you are presently reading this means that you have had
+#  knowledge of the CeCILL-B license and that you accept its terms.
+# ------------------------------------------------------------------------ */
 
 #include "ContourExtractData.h"
 
 
 //----------------------------------------------------------------------
 ContourExtractData::ContourExtractData( bool okImagesResults)
-  {
-    this->imagedata                    = NULL;
+{
+       this->imagedata                 = NULL;
        imagedataValueResult    = NULL;
        imagedataMaskResult             = NULL;
        this->okImagesResults   = okImagesResults;
        _typeOperation                  = 0;
-  }
+}
 
-  // ------------------------------------------------------------------------
+// ------------------------------------------------------------------------
 
-  ContourExtractData::~ContourExtractData()
-  {
-  }
+ContourExtractData::~ContourExtractData()
+{
+}
 
 
 //----------------------------------------------------------------------
 void ContourExtractData::SetImage( vtkImageData* imagedata)
-  {
-    this->imagedata                    = imagedata;
-       if (this->okImagesResults==true){
-               InitVtkImagesResult();
-       }
-  }
+{
+       this->imagedata                 = imagedata;
+       this->imagedata->GetScalarRange(scalarRange);
+
+       // RaC 20-11-09 Changes in InitLstContoursLinesYPoints
+       int ext[6];
+       this->imagedata->GetWholeExtent(ext);
+       _sizeImageY = ext[3]-ext[2]+1;
+
+       // init vtk image result : valuesImage maskImage  
+       if (this->okImagesResults==true){  InitVtkImagesResult(); }
+}
 //----------------------------------------------------------------------
 void ContourExtractData::SetZtoBeAnalys( int z )
-  {
+{
        this->zImage                    = z;
-  }
+}
 
 //------------------------------------------------------------------------
-  void ContourExtractData::SetLstManualContourModel( std::vector<manualContourModel*> lstManConMod)
-  {
-         this->lstManConMod = lstManConMod;
-  }
+void ContourExtractData::SetLstManualContourModel( std::vector<manualBaseModel*> lstManConMod)
+{
+       this->lstManConMod = lstManConMod;
+}
 
 
 //------------------------------------------------------------------------
 void ContourExtractData::GetMinMaxPoint(int *minPoint, 
-                                                                                 int *maxPoint, 
-                                                                                 manualContourModel *manualcontourmodel
-                                                                                 )
+                                                                               int *maxPoint, 
+                                                                               manualBaseModel *manualcontourmodel
+                                                                               )
 {
        int i;
        //int   np              = manualcontourmodel->GetSizeLstPoints( );  // number of control points // JPRx
 
-// JSTG 26-02-08 ---------------------------------------------------------------------------------------
+       // JSTG 26-02-08 ---------------------------------------------------------------------------------------
        //int nps = manualviewbaseecontour->GetNumberOfPointsSpline(); // number of points in the spline
        int nps = manualcontourmodel->GetNumberOfPointsSpline(); // number of points in the spline
-//------------------------------------------------------------------------------------------------------
-       
-// JSTG 26-02-08 ---------------------------------------------------------------------------------------
+       //------------------------------------------------------------------------------------------------------
+
+       // JSTG 26-02-08 ---------------------------------------------------------------------------------------
        //double x,y,z,t;
        double x,y,z;
        //double delta=( double ) ( np  ) / ( double ) ( nps  );
@@ -69,72 +99,121 @@ void ContourExtractData::GetMinMaxPoint(int *minPoint,
                if (x>maxPoint[0]){ maxPoint[0]=(int)x; }
                if (y>maxPoint[1]){ maxPoint[1]=(int)y; }
        }
-//------------------------------------------------------------------------------------------------------
+       minPoint[0]--;
+       minPoint[1]--;
+       maxPoint[0]++;
+       maxPoint[1]++;
+       //------------------------------------------------------------------------------------------------------
 }
 
 //------------------------------------------------------------------------
 void ContourExtractData::GetMinMaxPoint_Of_LstManConMod(       int *minPoint, 
-                                                                                       int *maxPoint
-                                                                               )
+                                                                                                               int *maxPoint
+                                                                                                               )
 
 {
        int i,size = lstManConMod.size();
+
        for(i=0 ; i<size ; i++)
        {
                GetMinMaxPoint(minPoint,maxPoint,lstManConMod[i]);
        }
 }
-//------------------------------------------------------------------------
 
-int ContourExtractData::AnalisisContourInside(int x, 
-                                           int y, 
-                                           manualContourModel *manualcontourmodel
-                                           )
+
+//------------------------------------------------------------------------
+int ContourExtractData::AnalisisContourInsideV2(int x, int y, int iContour )
 {
+       bool inBorder=false;
        int result      = 0;
        int i;
-       //int   np              = manualcontourmodel->GetSizeLstPoints( );  // number of control points // JPRx
 
-// JSTG 26-02-08 ---------------------------------------------------------------------------------------
-       //int nps = manualviewbaseecontour->GetNumberOfPointsSpline(); // number of points in the spline
-       int nps = manualcontourmodel->GetNumberOfPointsSpline(); // number of points in the spline
-       //double x1,y1,z1,x2,y2,z2,t;
-       double x1,y1,z1,x2,y2,z2;
+       int nps=_lstlstlstVecX1[iContour][y].size();
+
+       double x1,y1,x2,y2;
+       double borderX, borderY;
        double xx1, yy1,xx2, yy2;
-       //double delta=( double ) ( np  ) / ( double ) ( nps  );
-       manualcontourmodel->UpdateSpline();
-//------------------------------------------------------------------------------------------------------
+       double xx=x, yy=y;
        double d;
-       bool ok;
-//     if (np>=2)
-//     {
-// JSTG 26-02-08 ---------------------------------------------------------------------------------------
-               nps--;
-               //manualcontourmodel->GetSplinePoint(0,x1,y1,z1);
-               manualcontourmodel->GetSpline_i_Point(0,&x1,&y1,&z1);
-               for (i=1; i<=nps; i++)
+
+       for (i=0; i<nps; i++)
+       {
+               x1=_lstlstlstVecX1[iContour][y][i];
+               y1=_lstlstlstVecY1[iContour][y][i];
+               x2=_lstlstlstVecX2[iContour][y][i];
+               y2=_lstlstlstVecY2[iContour][y][i];
+
+               borderX=x1;
+               borderY=y1;
+
+               if (y1<y2) 
+               { 
+                       xx1=x1; yy1=y1; xx2=x2; yy2=y2;
+               } else {
+                       xx1=x2; yy1=y2; xx2=x1; yy2=y1; 
+               } 
+
+               double difxx2xx1=fabs(xx2-xx1);
+               if (difxx2xx1==0) difxx2xx1=0.0000000001;
+
+               // Finding border looking in vertical direction  AND  verifing if pixel is at right of the line 
+               if  ( (yy>=yy1)&&(yy<=yy2) ) 
                {
-                       ok=false;
-                       //t= delta * (double)(i%nps);
-                       //manualcontourmodel->GetSplinePoint(t,x2,y2,z2);
-                       manualcontourmodel->GetSpline_i_Point(i,&x2,&y2,&z2);
-//------------------------------------------------------------------------------------------------------
                        //by triangle similarity
-                       if ( ((y1<y2)&&(y>=y1)&&(y<y2)) || ((y1>y2)&&(y<=y1)&&(y>y2)) )
+                       d = ( fabs(xx2-xx1)*(yy-yy1) ) / (yy2-yy1) ;
+                       if ( (xx1<=xx2)&&(x<(xx1+d)) )  {       result++;       }
+                       if ( (xx1>xx2)&&(x<(xx1-d)) )   {       result++;       } 
+
+                       if ( (yy2-yy1)/difxx2xx1 >= 1.0)
                        {
-                               if (y1<y2) { xx1=x1; yy1=y1; xx2=x2; yy2=y2;} else { xx1=x2; yy1=y2; xx2=x1; yy2=y1; } 
-                               d = ( fabs(xx2-xx1)*(y-yy1) ) / (yy2-yy1) ;
-                               if (  ((xx1<xx2)&&(x<(xx1+d)))  ||  ((xx1>xx2)&&(x<(xx1-d)))  ) { result++; }
-                       } // if
-                       x1=x2; y1=y2; z1=z2;
-               } // for i
-//     } //if
+                               if (xx1<=xx2)   
+                               { 
+                                       borderX = xx1+d;
+                                       borderY = y;
+                               } else { 
+                                       borderX = xx1-d;
+                                       borderY = y;
+                               } 
+                       } 
+               } // if point inside y 
+
+
+               // Finding border looking in vertical direction
+               if ( ((xx1<=xx2)&&(xx>=xx1)&&(xx<xx2)) || ((xx1>xx2)&&(xx>=xx2)&&(xx<xx1)) )
+               {
+                       if ( (yy2-yy1)/difxx2xx1 <= 1.0)
+                       {
+                               //by triangle similarity
+                               d = ( fabs(xx1-xx)*(yy2-yy1) ) / difxx2xx1;
+                               if (yy1+d<=yy2) {
+                                       borderX=x;
+                                       borderY=yy1+d;
+                               }
+                       }                       
+               } // if point inside x                  
+
+
+               //Border verification
+               if (   (x==(int)borderX) && (y==(int)borderY)   )  { inBorder=true; }// if point in border
+
+               // Verification : border in horizontal line 
+               if ( ((int)y1==(int)y2)  &&  ((int)y1==y) && (x1<x2) && (x>=x1) && (x<=x2))             { inBorder=true;        }
+               if ( ((int)y1==(int)y2)  &&  ((int)y1==y) && (x2<x1) && (x>=x2) && (x<=x1))             { inBorder=true;        }
+
+               if (inBorder==true){ i=nps; }           
+       } // for i
+
+       if (inBorder==true)     { result=1;     }
+
        return result;
 }
 
-//------------------------------------------------------------------------
 
 
+//------------------------------------------------------------------------
+// typeOperation=0 AND
+// typeOperation=1 OR
+// typeOperation=2 XOR
 bool ContourExtractData::isInside(int x, int y, int typeOperation)
 {
        bool result                             = false;
@@ -142,44 +221,90 @@ bool ContourExtractData::isInside(int x, int y, int typeOperation)
        int i,size                              = this->lstManConMod.size();
        int numberInside        = 0;
 
+       /* RaC 20-11-09 (C1) Changes to use the method without the image.
        int ext[6];
        imagedata->GetExtent(ext);
 
        if ((x>=0) && (x<=ext[1]) && (y>=0) && (y<=ext[3]))
        {
+       */
 
-               if (typeOperation==0)  // AND  Intersection
+       if (typeOperation==0)  // AND  Intersection
+       {
+               for (i=0;i<size;i++)
                {
-                       for (i=0;i<size;i++)
+                       // To process the statistics of the Points contour the procedure is different
+                       // RaC 19-09-09
+                       manualBaseModel *mbm = lstManConMod[i];
+                       if(mbm->GetTypeModel()==7)
                        {
-                               numberLeft =  AnalisisContourInside(x,y, lstManConMod[i] );
-                               if ( (numberLeft % 2) ==1){         numberInside++;  }
-                       }
-                       if ( numberInside == (size) ){ result=true; }
-               } // AND  Intersection
+                               if(mbm->IsPoint(x,y)==true)
+                               {
+                                       numberLeft=1;
+                               }
+                       }//if
+                       else
+                       {
+                               numberLeft =  AnalisisContourInsideV2(x,y, i );
+                       }//else
+
+                       if ( (numberLeft % 2) ==1){         numberInside++;  }
+               }
+               if ( numberInside == (size) ){ result=true; }
+       } // AND  Intersection
 
+       numberLeft=0;
 
-               if (typeOperation==1)  // OR  All
+       if (typeOperation==1)  // OR  All
+       {
+               for (i=0;i<size;i++)
                {
-                       for (i=0;i<size;i++)
+                       // To process the statistics of the Points contour the procedure is different
+                       // RaC 19-09-09
+                       manualBaseModel *mbm = lstManConMod[i];
+                       if(mbm->GetTypeModel()==7)
                        {
-                               numberLeft =  AnalisisContourInside(x,y, lstManConMod[i] );
-                               if ( (numberLeft % 2) ==1){ result=true;  }
-                       }
-               } // OR  All
+                               if(mbm->IsPoint(x,y)==true)
+                               {
+                                       numberLeft=1;
+                               }
+                       }//if
+                       else
+                       {
+                               numberLeft =  AnalisisContourInsideV2(x,y, i );
+                       }//else
+                       if ( (numberLeft % 2) ==1){ result=true;  }
+               }
+       } // OR  All
+
+       numberLeft=0;
 
-               if (typeOperation==2)  // XOR  crown
+       if (typeOperation==2)  // XOR  crown
+       {
+               for (i=0;i<size;i++)
                {
-                       for (i=0;i<size;i++)
+                       // To process the statistics of the Points contour the procedure is different
+                       // RaC 19-09-09
+                       manualBaseModel *mbm = lstManConMod[i];
+                       if(mbm->GetTypeModel()==7)
+                       {
+                               if(mbm->IsPoint(x,y)==true)
+                               {
+                                       numberLeft=1;
+                               }
+                       }//if
+                       else
                        {
-                               numberLeft = numberLeft + AnalisisContourInside(x,y, lstManConMod[i] );
-                       }
-                       if ( numberLeft % 2 ==1){       result = true;  } 
-               }// XOR  crown
+                               numberLeft =  numberLeft + AnalisisContourInsideV2(x,y, i );
+                       }//else
 
+               }
+               if ( numberLeft % 2 ==1){       result = true;  } 
+       }// XOR  crown
 
 
-       }
+
+       //RaC }
 
        return result;
 }
@@ -188,31 +313,71 @@ bool ContourExtractData::isInside(int x, int y, int typeOperation)
 
 double ContourExtractData::GetDataValue(int x, int y, int z)
 {
-//     wxVtk2DBaseView *wxvtk2dbaseview = (wxVtk2DBaseView*)wxvtkbaseview;
-//     int z = (int)wxvtk2dbaseview->GetVtkBaseData()->GetZ();
+       //      wxVtk2DBaseView *wxvtk2dbaseview = (wxVtk2DBaseView*)wxvtkbaseview;
+       //      int z = (int)wxvtk2dbaseview->GetVtkBaseData()->GetZ();
        //JSTG 13-03-08-----
-//EED OJO avec JS      _zz = z;
+       //EED OJO avec JS       _zz = z;
        //------------------
        double result;
        void *p;
        p = imagedata->GetScalarPointer(x,y,z);
 
-       if (imagedata->GetScalarType()==VTK_UNSIGNED_CHAR)
-       {       unsigned char *pp = (unsigned char*)p;
+       if (imagedata->GetScalarType()==VTK_CHAR)
+       {
+               char *pp = (char*)p;
                result = (double)(*pp);
        }
-       if (imagedata->GetScalarType()==VTK_FLOAT)
-       {       float *pp = (float*)p;
+       else if (imagedata->GetScalarType()==VTK_SIGNED_CHAR)
+       {
+               signed char *pp = (signed char*)p;
+               result = (double)(*pp);
+       }
+       else if (imagedata->GetScalarType()==VTK_UNSIGNED_CHAR)
+       {
+               unsigned char *pp = (unsigned char*)p;
                result = (double)(*pp);
        }
-       if (imagedata->GetScalarType()==VTK_SHORT)
-       {       short *pp = (short*)p;
+       else if (imagedata->GetScalarType()==VTK_SHORT)
+       {
+               short *pp = (short*)p;
                result = (double)(*pp);
        }
-       if (imagedata->GetScalarType()==VTK_UNSIGNED_SHORT)
-       {       unsigned short *pp = (unsigned short*)p;
+       else if (imagedata->GetScalarType()==VTK_UNSIGNED_SHORT)
+       {
+               unsigned short *pp = (unsigned short*)p;
                result = (double)(*pp);
        }
+       else if (imagedata->GetScalarType()==VTK_INT)
+       {
+               int *pp = (int*)p;
+               result = (double)(*pp);
+       }
+       else if (imagedata->GetScalarType()==VTK_UNSIGNED_INT)
+       {
+               unsigned int *pp = (unsigned int*)p;
+               result = (double)(*pp);
+       }
+       else if (imagedata->GetScalarType()==VTK_LONG)
+       {
+               long *pp = (long*)p;
+               result = (double)(*pp);
+       }
+       else if (imagedata->GetScalarType()==VTK_UNSIGNED_LONG)
+       {
+               unsigned long *pp = (unsigned long*)p;
+               result = (double)(*pp);
+       }
+       else if (imagedata->GetScalarType()==VTK_FLOAT)
+       {       
+               float *pp = (float*)p;
+               result = (double)(*pp);
+       }
+       else if (imagedata->GetScalarType()==VTK_DOUBLE)
+       {
+               double *pp = (double*)p;
+               result = (double)(*pp);
+       }
+
        return result;
 }
 
@@ -222,10 +387,13 @@ void ContourExtractData::PutVtkImageDataResultValue( int x, int y, int z, double
 {
        unsigned short *pValue;
        unsigned short *pMask;
-       pValue  = (unsigned short *)imagedataValueResult->GetScalarPointer(x,y,z);
-       pMask   = (unsigned short *)imagedataMaskResult->GetScalarPointer(x,y,z);
-       *pMask  = 255;
-       *pValue = (unsigned short)value;
+//EED 2017-12-18
+       imagedataValueResult->SetScalarComponentFromDouble(x,y,z,0,value);
+       imagedataMaskResult->SetScalarComponentFromDouble(x,y,z,0,255);
+//     pValue  = (unsigned short *)imagedataValueResult->GetScalarPointer(x,y,z);
+//     *pValue = (unsigned short)value;
+//     pMask   = (unsigned char *)imagedataMaskResult->GetScalarPointer(x,y,z);
+//     *pMask  = 255;
 }
 
 //------------------------------------------------------------------------
@@ -237,18 +405,13 @@ void ContourExtractData::ResetImageResult(int z)
                unsigned short *pMask;
                pValue  = (unsigned short *)imagedataValueResult->GetScalarPointer(0,0,z);
                pMask   = (unsigned short *)imagedataMaskResult->GetScalarPointer(0,0,z);
-               
+
                int ext[6];
                imagedataValueResult->GetExtent(ext);
 
-               int i,size = (ext[1]-ext[0]+1) * (ext[3]-ext[2]+1); 
-               for(i=0; i<size; i++)
-               {
-                       *pMask  = 0;
-                       *pValue = 0;
-                       pMask++;
-                       pValue++;
-               }// for
+               int size = (ext[1]-ext[0]+1) * (ext[3]-ext[2]+1); 
+               memset(pValue,0, size*imagedataValueResult->GetScalarSize() );
+               memset(pMask,0, size*imagedataMaskResult->GetScalarSize() );
        } // if
 }
 
@@ -271,27 +434,44 @@ void ContourExtractData::CalculateImageResult()
                maxPoint[1] = -999999;
 
                GetMinMaxPoint_Of_LstManConMod(minPoint,maxPoint);
-               for (j=minPoint[1]; j<maxPoint[1]; j++)
+               InitLstContoursLinesYPoints();
+
+               for (j=minPoint[1]; j<=maxPoint[1]; j++)
                {
-                       for (i=minPoint[0]; i<maxPoint[0]; i++)
+                       for (i=minPoint[0]; i<=maxPoint[0]; i++)
                        {
-                               if (isInside(i,j,_typeOperation)==true)
+
+                               //RaC 20-11-09 Changes specified in isInside (C1)
+                               int ext[6];
+                               imagedata->GetExtent(ext);
+
+                               if ((i>=0) && (i<=ext[1]) && (j>=0) && (j<=ext[3]))
                                {
-                                       value = GetDataValue(i,j,zImage);
-                                       PutVtkImageDataResultValue(i,j,zImage,  value );
-                               } // if
+                  if (isInside(i,j,_typeOperation)==true)
+                                       {
+                                         value = GetDataValue(i,j,zImage);
+                                         if ( (value>=scalarRange[0]) && (value<=scalarRange[1]) )
+                                         {
+                                               PutVtkImageDataResultValue(i,j,zImage,  value );
+                                         } // scalarRange
+                                       } // if isInside
+
+                               } // if ext
                        } // for i
                } // for j
 
 
                imagedataValueResult->Modified();
                imagedataMaskResult->Modified();
-       }
+               imagedataValueResult->Update();
+               imagedataMaskResult->Update();
+       } // if
 
 }
 
 //------------------------------------------------------------------------
-void ContourExtractData::GetValuesInsideCrown(std::vector<double> *pLstValue,
+void ContourExtractData::GetValuesInsideCrown( int *numberOfPixels, 
+                                                                                               std::vector<double> *pLstValue,
                                                                                                std::vector<double> *pLstValuePosX,
                                                                                                std::vector<double> *pLstValuePosY,
                                                                                                std::vector<double> *pLstValuePosZ)
@@ -301,16 +481,16 @@ void ContourExtractData::GetValuesInsideCrown(std::vector<double> *pLstValue,
        pLstValuePosY->clear();
        pLstValuePosZ->clear();
 
-//     if (okImagesResults==true)
-//     {
-//             ResetImageResult(zImage);
-//     }
+       //      if (okImagesResults==true)
+       //      {
+       //              ResetImageResult(zImage);
+       //      }
 
        int minPoint[2];
        int maxPoint[2];
        int i,j;
        double value;
-
+       int acum=0;
 
        minPoint[0] = 999999;
        minPoint[1] = 999999;
@@ -318,37 +498,36 @@ void ContourExtractData::GetValuesInsideCrown(std::vector<double> *pLstValue,
        maxPoint[1] = -999999;
 
        GetMinMaxPoint_Of_LstManConMod(minPoint,maxPoint);
+       InitLstContoursLinesYPoints();
 
-
-       for (j=minPoint[1]; j<maxPoint[1]; j++)
+       for (j=minPoint[1]; j<=maxPoint[1]; j++)
        {
-               for (i=minPoint[0]; i<maxPoint[0]; i++)
+               for (i=minPoint[0]; i<=maxPoint[0]; i++)
                {
-                       if (isInside(i,j,_typeOperation)==true)
+                       //RaC 20-11-09 Changes specified in isInside (C1)
+                       int ext[6];
+                       imagedata->GetExtent(ext);
+
+                       if ((i>=0) && (i<=ext[1]) && (j>=0) && (j<=ext[3]))
                        {
-                               value = GetDataValue(i,j,zImage);
-
-// Borrame
-//                             if (okImagesResults==true){
-//                                     PutVtkImageDataResultValue(i,j,zImage,  value );
-//                             }
-
-                               pLstValue               -> push_back( value );
-                               pLstValuePosX   -> push_back( i );
-                               pLstValuePosY   -> push_back( j );
-                               pLstValuePosZ   -> push_back( -1 );
-                       } // if
+                   if (isInside(i,j,_typeOperation)==true)
+                               {
+
+                                       acum++;
+                                       value = GetDataValue(i,j,zImage);
+                                       if ( (value>=scalarRange[0]) && (value<=scalarRange[1]) )
+                                       {
+                                               pLstValue               -> push_back( value );
+                                               pLstValuePosX   -> push_back( i );
+                                               pLstValuePosY   -> push_back( j );
+                                               pLstValuePosZ   -> push_back( -1 );
+                                       } // scalarRange
+                               } // if isInside
+                       } // if ext
                } // for
        } // for
 
-
-// Borrame
-//     if (this->okImagesResults==true){
-//             imagedataValueResult->Modified();
-//             imagedataMaskResult->Modified();
-//     }
-
-
+       *numberOfPixels = acum;
 }
 
 //------------------------------------------------------------------------
@@ -365,41 +544,93 @@ vtkImageData *ContourExtractData::GetVtkImageMaskResult()
 // ------------------------------------------------------------------------
 void ContourExtractData::InitVtkImagesResult()
 {
-         int ext[6];
-         int newDim[3];
-         double spc[3];
-         int scalartype;
-
-         imagedata->GetSpacing(spc);
-         imagedata->GetExtent(ext);
-         newDim[0]=ext[1]-ext[0]+1;
-         newDim[1]=ext[3]-ext[2]+1;
-         newDim[2]=ext[5]-ext[4]+1;
-         scalartype = imagedata->GetScalarType();
-
-         if (imagedataValueResult!=NULL)
-         {
-                 imagedataValueResult->Delete();
-         }
-         imagedataValueResult = vtkImageData::New();
-//       imagedataValueResult->SetScalarType(scalartype);
-         imagedataValueResult->SetScalarTypeToUnsignedShort();
-         imagedataValueResult->SetSpacing(spc);
-         imagedataValueResult->SetDimensions( newDim );
-         imagedataValueResult->AllocateScalars();
-
-         if (imagedataMaskResult!=NULL)
-         {
-                 imagedataMaskResult->Delete();
-         }
-         imagedataMaskResult  = vtkImageData::New();
-//       imagedataMaskResult->SetScalarType(scalartype);
-         imagedataMaskResult->SetScalarTypeToUnsignedShort();
-         imagedataMaskResult->SetSpacing(spc);
-         imagedataMaskResult->SetDimensions( newDim );
-         imagedataMaskResult->AllocateScalars();
+       int ext[6];
+       int newDim[3];
+       double spc[3];
+       int scalartype;
+
+       imagedata->GetSpacing(spc);
+       imagedata->GetExtent(ext);
+       newDim[0]=ext[1]-ext[0]+1;
+       newDim[1]=ext[3]-ext[2]+1;
+       newDim[2]=ext[5]-ext[4]+1;
+       scalartype = imagedata->GetScalarType();
+
+       if (imagedataValueResult!=NULL)
+       {
+               imagedataValueResult->Delete();
+       }
+       imagedataValueResult = vtkImageData::New();
+
+//EED 2017-12-18
+       imagedataValueResult->SetScalarType(scalartype);
+       //imagedataValueResult->SetScalarTypeToUnsignedShort();
+
+       imagedataValueResult->SetSpacing(spc);
+       imagedataValueResult->SetDimensions( newDim );
+       imagedataValueResult->AllocateScalars();
+
+       if (imagedataMaskResult!=NULL)
+       {
+               imagedataMaskResult->Delete();
+       }
+       imagedataMaskResult  = vtkImageData::New();
+
+//EED 2017-12-18
+//     imagedataMaskResult->SetScalarType(scalartype);
+//     imagedataMaskResult->SetScalarTypeToUnsignedShort();
+       imagedataMaskResult->SetScalarTypeToUnsignedChar();
+
+       imagedataMaskResult->SetSpacing(spc);
+       imagedataMaskResult->SetDimensions( newDim );
+       imagedataMaskResult->AllocateScalars();
+}
+
+
+//------------------------------------------------------------------------
+void ContourExtractData::InitVolumeStatistics()
+{
+       vol_rCountRange                         = 0;
+       vol_rsize                                       = 0;
+       vol_minValue                            = 9999999;
+       vol_maxValue                            =-9999999;
+       vol_acum_average                        = 0;
+       vol_acum_standardeviation       = 0;
 }
 
+//------------------------------------------------------------------------
+void ContourExtractData::SetVolumeStatistics(int rCountRange, 
+                                                                                        int rsize,
+                                                                                        double minValue,
+                                                                                        double maxValue,
+                                                                                        double acum_average,
+                                                                                        double acum_standardeviation)
+{
+       vol_rCountRange                         = vol_rCountRange + rCountRange; 
+       vol_rsize                                       = vol_rsize + rsize; 
+
+       if (minValue<vol_minValue){ vol_minValue = minValue;  }
+       if (maxValue>vol_maxValue){ vol_maxValue = maxValue;  }
+
+       vol_acum_average                        = vol_acum_average + acum_average; 
+       vol_acum_standardeviation       = vol_acum_standardeviation + acum_standardeviation; 
+}
+
+//------------------------------------------------------------------------
+void ContourExtractData::GetVolumeStatistics(int *vol_rCountRange, 
+                                                                                        int *vol_rsize,
+                                                                                        double *vol_minValue,
+                                                                                        double *vol_maxValue,
+                                                                                        double *vol_average,
+                                                                                        double *vol_standardeviation)
+{
+       *vol_rCountRange                = this->vol_rCountRange; 
+       *vol_rsize                              = this->vol_rsize; 
+       *vol_minValue                   = this->vol_minValue; 
+       *vol_maxValue                   = this->vol_maxValue; 
+       *vol_average                    = this->vol_acum_average / this->vol_rsize;
+       *vol_standardeviation   = sqrt(this->vol_acum_standardeviation / this->vol_rsize);
+}
 
 
 //------------------------------------------------------------------------
@@ -411,55 +642,59 @@ void ContourExtractData::Statistics( std::vector<double> *inputLstValue,
                                                                        double  *rmin, 
                                                                        double  *rmax,
                                                                        double  *raverage,
-                                                                       double  *rstandardeviation)
+                                                                       double  *rstandardeviation
+                                                                       )
 {
-         double min                            = 0;
-         double max                            = 0;
-         double average                        = 0;
-         double standardeviation       = 0;
-         double acum                           = 0;
-         int    size                           = 0;
-         int    countRange                     = 0;
-         double ng;
-
-         if (inputLstValue!=NULL)
-         {
+       double min                                              = 0;
+       double max                                              = 0;
+       double average                                  = 0;
+       double standardeviation                 = 0;
+       double acum_average                             = 0;
+       double acum_standardeviation            = 0;
+       int      size                                           = 0;
+       int      countRange                                     = 0;
+       double ng;
+
+       if (inputLstValue!=NULL)
+       {
                size=inputLstValue->size();
                if (size>0){
-              max=(*inputLstValue)[0];
-                  min=(*inputLstValue)[0];
-     // Average , countRange
+                       max=(*inputLstValue)[0];
+                       min=(*inputLstValue)[0];
+                       // Average , countRange
                        int i;
                        for ( i=0; i<size; i++ )
                        {
                                ng=(*inputLstValue)[i];
-                               acum = acum + ng;
+                               acum_average = acum_average + ng;
                                if (max<ng) max=ng;             // Max
                                if (min>ng) min=ng;     // Min
                                if ((ng>=grayRangeMin) && (ng<=grayRangeMax)) countRange++;  // countRange
-                       }
-                       average = acum / size;
+                       } // for average
+                       average = acum_average / size;
 
-         // Standar Deviation
-                       acum=0;
+                       // Standar Deviation
+                       acum_standardeviation=0;
                        double tmp;
-                       for ( i=0; i<size; i++ )
+                       for ( i=0; i<size; i++ )
                        {
-                tmp = (*inputLstValue)[i] - average;
-                               acum = acum + tmp*tmp;
-                       }
-                       standardeviation = sqrt(acum/size);
-
-               }
-         }
-
-         // OUTPUT
-               *rsize                          = size; 
-               *rCountRange            = countRange;
-               *rmin                           = min; 
-               *rmax                           = max;
-               *raverage                       = average;
-               *rstandardeviation      = standardeviation;
+                               tmp = (*inputLstValue)[i] - average;
+                               acum_standardeviation = acum_standardeviation + tmp*tmp;
+                       } // for standar deviation
+                       standardeviation = sqrt(acum_standardeviation/size);
+                       SetVolumeStatistics(countRange, (*rsize),
+                                                               min,max,
+                                                           acum_average,acum_standardeviation);
+               } // if size
+       } // if NULL
+
+       // OUTPUT
+       *rsize                          = size; 
+       *rCountRange            = countRange;
+       *rmin                           = min; 
+       *rmax                           = max;
+       *raverage                       = average;
+       *rstandardeviation      = standardeviation;
 }
 
 //------------------------------------------------------------------------
@@ -467,3 +702,80 @@ void ContourExtractData::SetTypeOperation(int type)
 {
        _typeOperation=type;
 }
+
+//------------------------------------------------------------------------
+void ContourExtractData::Fill_lstlstlstVecXY(int iContour, int sizeY)
+{
+       int i,y;
+       double x1,y1,z1,x2,y2,z2;
+       manualBaseModel *manualcontourmodel= lstManConMod[iContour];    
+       int nps = manualcontourmodel->GetNumberOfPointsSpline(); // number of points in the spline
+       manualcontourmodel->UpdateSpline();
+       //------------------------------------------------------------------------------------------------------
+
+       for (y=0;y<sizeY;y++)
+       {
+               manualcontourmodel->GetSpline_i_Point(0,&x1,&y1,&z1);
+               x1=x1+0.5; y1=y1+0.5;
+               for (i=1; i<nps; i++)
+               {
+                       manualcontourmodel->GetSpline_i_Point(i,&x2,&y2,&z2);
+                       x2=x2+0.5; y2=y2+0.5;
+                       if ( ((y1<y2)&&(y>=y1)&&(y<=y2)) || ((y1>y2)&&(y>=y2)&&(y<=y1)) || ((int)y1==y) || ((int)y2==y)  )
+                       {
+                               _lstlstlstVecX1[iContour][y].push_back(x1);
+                               _lstlstlstVecY1[iContour][y].push_back(y1);
+                               _lstlstlstVecX2[iContour][y].push_back(x2);
+                               _lstlstlstVecY2[iContour][y].push_back(y2);
+                       } 
+                       x1=x2; y1=y2; z1=z2;
+               } // for i  Points in spline
+       } // y
+
+}
+
+void ContourExtractData::InitLstContoursLinesYPoints()
+{
+       // init InInside Optimisation  
+       int i;
+
+       _lstlstlstVecX1.clear();
+       _lstlstlstVecY1.clear();
+       _lstlstlstVecX2.clear();
+       _lstlstlstVecY2.clear();
+
+       /* RaC 20-11-09
+       int ext[6];
+       this->imagedata->GetWholeExtent(ext);
+       int sizeY = ext[3]-ext[2]+1;
+       */
+       std::vector<double> vecDouble;
+       std::vector< std::vector<double> > vecVecDouble;
+       for ( i=0 ; i<_sizeImageY ; i++ )
+       {
+               vecVecDouble.push_back( vecDouble );
+       }
+
+       //Fill structure with points
+       int sizeContours = lstManConMod.size();
+       for( i=0 ; i<sizeContours ; i++ )
+       {
+               _lstlstlstVecX1.push_back( vecVecDouble );
+               _lstlstlstVecY1.push_back( vecVecDouble );
+               _lstlstlstVecX2.push_back( vecVecDouble );
+               _lstlstlstVecY2.push_back( vecVecDouble );
+               Fill_lstlstlstVecXY(i,_sizeImageY);
+       }
+
+}
+
+void ContourExtractData::SetScalarRange(double min, double max)
+{
+       scalarRange[0]=min;
+       scalarRange[1]=max;
+}
+
+void ContourExtractData::SetSizeImageY(int pSizeImageY)
+{
+       _sizeImageY=pSizeImageY;
+}