bbtk/src/bbmaracasvisuAxeVolume.cxx

   1 /*# ---------------------------------------------------------------------
   2 #
   3 # Copyright (c) CREATIS (Centre de Recherche en Acquisition et Traitement de l'Image
   4 #                        pour la Sant�)
   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
   8 #
   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.
  15 #
  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
  20 #  liability.
  21 #
  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 # ------------------------------------------------------------------------ */
  25
  26 #include "bbmaracasvisuAxeVolume.h"
  27 #include "bbcreaMaracasVisuPackage.h"
  28
  29 namespace bbcreaMaracasVisu
  30 {
  31
  32
  33 BBTK_ADD_BLACK_BOX_TO_PACKAGE(creaMaracasVisu,AxeVolume)
  34 BBTK_BLACK_BOX_IMPLEMENTATION(AxeVolume,bbtk::AtomicBlackBox);
  35
  36 int AxeVolume::GetTypeFormat( std::string formatStr , vtkImageData* image )
  37 {
  38         int outputformat = VTK_UNSIGNED_CHAR;
  39         if (formatStr=="SAME")
  40         {
  41                 if (image!=NULL) outputformat = image->GetScalarType();
  42         }
  43         else if (formatStr=="VTK_BIT")                          outputformat = VTK_BIT;                         // 1
  44         else if (formatStr=="VTK_CHAR")                         outputformat = VTK_CHAR;                        // 2
  45         else if (formatStr=="VTK_SIGNED_CHAR")          outputformat = VTK_SIGNED_CHAR;         // 15
  46         else if (formatStr=="VTK_UNSIGNED_CHAR")        outputformat = VTK_UNSIGNED_CHAR;       // 3
  47         else if (formatStr=="VTK_SHORT")                        outputformat = VTK_SHORT;                       // 4
  48         else if (formatStr=="VTK_UNSIGNED_SHORT")       outputformat = VTK_UNSIGNED_SHORT;      // 5
  49         else if (formatStr=="VTK_INT")                          outputformat = VTK_INT;                 // 6
  50         else if (formatStr=="VTK_UNSIGNED_INT")         outputformat = VTK_UNSIGNED_INT;        // 7
  51         else if (formatStr=="VTK_LONG")                         outputformat = VTK_LONG;                // 8
  52         else if (formatStr=="VTK_UNSIGNED_LONG")        outputformat = VTK_UNSIGNED_LONG;       // 9
  53         else if (formatStr=="VTK_FLOAT")                        outputformat = VTK_FLOAT;               // 10
  54         else if (formatStr=="VTK_DOUBLE")                       outputformat = VTK_DOUBLE;              // 11
  55         else if (formatStr=="MET_CHAR")                         outputformat = VTK_CHAR;                        // 2
  56         else if (formatStr=="MET_UCHAR")                        outputformat = VTK_UNSIGNED_CHAR;       // 3
  57         else if (formatStr=="MET_SHORT")                        outputformat = VTK_SHORT;                       // 4
  58         else if (formatStr=="MET_USHORT")                       outputformat = VTK_UNSIGNED_SHORT;      // 5
  59         else if (formatStr=="MET_SHORT")                        outputformat = VTK_SHORT;               // 4
  60         else if (formatStr=="MET_FLOAT")                        outputformat = VTK_FLOAT;               // 10
  61         else if (formatStr=="MET_DOUBLE")                       outputformat = VTK_DOUBLE;              // 11
  62
  63     return outputformat;
  64 }
  65
  66
  67 void AxeVolume::Process()
  68 {
  69         if ( mimage!=NULL )
  70         {
  71                 mimage->Delete();
  72         }
  73
  74         int ext[6];
  75         bbGetInputIn()->GetExtent(ext);
  76         int sizeX=ext[1]-ext[0]+1;
  77         int sizeY=ext[3]-ext[2]+1;
  78         int sizeZ=ext[5]-ext[4]+1;
  79
  80     double spc[3];
  81         bbGetInputIn()->GetSpacing(spc);
  82     double invSpc[3];
  83     invSpc[0] = 1/spc[0];
  84     invSpc[1] = 1/spc[1];
  85     invSpc[2] = 1/spc[2];
  86
  87         int outputformat = GetTypeFormat( bbGetInputOutputFormat() , bbGetInputIn() );
  88
  89         mimage = vtkImageData::New();
  90         mimage->SetSpacing(bbGetInputIn()->GetSpacing());
  91         mimage->SetDimensions(bbGetInputIn()->GetDimensions());
  92         mimage->SetExtent(bbGetInputIn()->GetExtent());
  93         mimage->SetOrigin(bbGetInputIn()->GetOrigin());
  94         mimage->SetScalarType( outputformat );
  95         mimage->AllocateScalars();
  96
  97     int sizeLstPointR   = bbGetInputlstPointR().size();
  98         int iAxe,sizeAxe        = bbGetInputlstPointX().size();
  99         int ii;
 100         int sizeImage           = sizeX*sizeY*sizeZ;
 101         unsigned short *p;
 102
 103         // Clean image
 104         p = (unsigned short*)mimage->GetScalarPointer (0, 0, 0);
 105         for ( ii=0 ; ii<sizeImage ; ii++)
 106         {
 107                 *p = 0;
 108                 p++;
 109         }
 110
 111
 112 #pragma omp parallel for
 113         for (iAxe=0 ; iAxe<sizeAxe; iAxe++)
 114          {
 115                 int i,j,k;
 116                 double rx,ry,rz;
 117                 double r,rr;
 118                 double px,py,pz;
 119                 double px1,py1,pz1;
 120                 double px2,py2,pz2;
 121                 printf("AxeVolume %d/%d\n ",iAxe,sizeAxe);
 122                 if (iAxe<sizeLstPointR)
 123                 {
 124                         r = bbGetInputlstPointR()[ iAxe ]* invSpc[0];
 125                 } else {
 126                         if (bbGetInputlstPointR().size()>=1)
 127                         {
 128                                 r = bbGetInputlstPointR()[ bbGetInputlstPointR().size()-1 ]    * invSpc[0];
 129                         } else {
 130                                 r = 1;
 131                         }
 132                 }
 133                 px = bbGetInputlstPointX()[iAxe] * invSpc[0];
 134                 py = bbGetInputlstPointY()[iAxe] * invSpc[1];
 135                 pz = bbGetInputlstPointZ()[iAxe] * invSpc[2];
 136                 px1 = px - r;
 137                 py1 = py - r;
 138                 pz1 = pz - r;
 139                 px2 = px + r;
 140                 py2 = py + r;
 141                 pz2 = pz + r;
 142                 rr=r*r;
 143
 144                 for ( i=px1 ; i<=px2 ; i++ )
 145                 {
 146                         rx      =       i - px;
 147                         rx      =       rx*rx;
 148                         for ( j=py1 ; j<py2 ; j++ )
 149                         {
 150                                 ry      =       j - py;
 151                                 ry      =       ry*ry;
 152                                 for ( k=pz1 ; k<pz2 ; k++ )
 153                                 {
 154                                         if ( (i>=0) && (i<sizeX) && (j>=0) && (j<sizeY) &&(k>=0) && (k<sizeZ) )
 155                                         {
 156 //                                              p = (unsigned short*)mimage->GetScalarPointer (i, j, k);
 157 //                                              if (*p==0)
 158                                                 if ( mimage->GetScalarComponentAsDouble(i,j,k,0)==0 )
 159                                                 {
 160                                                         rz      =       k - pz;
 161                                                         rz      =       rz*rz;
 162                                                         if ( rx + ry + rz <= rr )
 163                                                         {
 164 //                                                              *p=255;
 165                                                                 mimage->SetScalarComponentFromDouble (i,j,k,0, bbGetInputValue() );
 166                                                         }
 167                                                 } // *p==0
 168                                         } // if inside point
 169                                 } //for k
 170                         } //for j
 171                 } //for i
 172         } // for iAxe
 173     bbSetOutputOut( mimage );
 174 }
 175
 176
 177         //-----------------------------------------------------------------
 178         void AxeVolume::bbUserSetDefaultValues()
 179         {
 180                 mimage=NULL;
 181                 bbSetInputOutputFormat("SAME");
 182                 bbSetInputValue(255);
 183         }
 184
 185         //-----------------------------------------------------------------
 186         void AxeVolume::bbUserInitializeProcessing()
 187         {
 188         }
 189
 190         //-----------------------------------------------------------------
 191         void AxeVolume::bbUserFinalizeProcessing()
 192         {
 193         }
 194
 195         //-----------------------------------------------------------------
 196
 197 }
 198 // EO namespace bbcreaMaracasVisu
 199
 200