common/clitkDicomRT_ROI_ConvertToImageFilter.cxx

   1 /*=========================================================================
   2   Program:         vv http://www.creatis.insa-lyon.fr/rio/vv
   3   Main authors :   XX XX XX
   4
   5   Authors belongs to:
   6   - University of LYON           http://www.universite-lyon.fr/
   7   - Léon Bérard cancer center    http://www.centreleonberard.fr
   8   - CREATIS CNRS laboratory      http://www.creatis.insa-lyon.fr
   9
  10   This software is distributed WITHOUT ANY WARRANTY; without even
  11   the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  12   PURPOSE.  See the copyright notices for more information.
  13
  14   It is distributed under dual licence
  15   - BSD       http://www.opensource.org/licenses/bsd-license.php
  16   - CeCILL-B  http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
  17
  18   =========================================================================*/
  19
  20 #include <iterator>
  21 #include <algorithm>
  22 #include "clitkDicomRT_ROI_ConvertToImageFilter.h"
  23 #include <vtkPolyDataToImageStencil.h>
  24 #include <vtkSmartPointer.h>
  25 #include <vtkImageStencil.h>
  26 #include <vtkLinearExtrusionFilter.h>
  27 #include "clitkImageCommon.h"
  28
  29 //--------------------------------------------------------------------
  30 clitk::DicomRT_ROI_ConvertToImageFilter::DicomRT_ROI_ConvertToImageFilter()
  31 {
  32   mROI = NULL;
  33   mWriteOutput = false;
  34   mCropMask = true;
  35 }
  36 //--------------------------------------------------------------------
  37
  38
  39 //--------------------------------------------------------------------
  40 clitk::DicomRT_ROI_ConvertToImageFilter::~DicomRT_ROI_ConvertToImageFilter()
  41 {
  42
  43 }
  44 //--------------------------------------------------------------------
  45
  46 bool clitk::DicomRT_ROI_ConvertToImageFilter::ImageInfoIsSet() const
  47 {
  48   return mSize.size() && mSpacing.size() && mOrigin.size();
  49 }
  50
  51 //--------------------------------------------------------------------
  52 void clitk::DicomRT_ROI_ConvertToImageFilter::SetROI(clitk::DicomRT_ROI * roi)
  53 {
  54   mROI = roi;
  55 }
  56 //--------------------------------------------------------------------
  57
  58
  59 //--------------------------------------------------------------------
  60 void clitk::DicomRT_ROI_ConvertToImageFilter::SetCropMaskEnabled(bool b)
  61 {
  62   mCropMask = b;
  63 }
  64 //--------------------------------------------------------------------
  65
  66
  67 //--------------------------------------------------------------------
  68 void clitk::DicomRT_ROI_ConvertToImageFilter::SetOutputImageFilename(std::string s)
  69 {
  70   mOutputFilename = s;
  71   mWriteOutput = true;
  72 }
  73 //--------------------------------------------------------------------
  74
  75
  76 //--------------------------------------------------------------------
  77 void clitk::DicomRT_ROI_ConvertToImageFilter::SetImageFilename(std::string f)
  78 {
  79   itk::ImageIOBase::Pointer header = clitk::readImageHeader(f);
  80   if (header->GetNumberOfDimensions() < 3) {
  81     std::cerr << "Error. Please provide a 3D image instead of " << f << std::endl;
  82     exit(0);
  83   }
  84   if (header->GetNumberOfDimensions() > 3) {
  85     std::cerr << "Warning dimension > 3 are ignored" << std::endl;
  86   }
  87   mSpacing.resize(3);
  88   mOrigin.resize(3);
  89   mSize.resize(3);
  90   for(unsigned int i=0; i<3; i++) {
  91     mSpacing[i] = header->GetSpacing(i);
  92     mOrigin[i] = header->GetOrigin(i);
  93     mSize[i] = header->GetDimensions(i);
  94   }
  95 }
  96 //--------------------------------------------------------------------
  97
  98 void clitk::DicomRT_ROI_ConvertToImageFilter::SetOutputOrigin(const double* origin)
  99 {
 100   std::copy(origin,origin+3,std::back_inserter(mOrigin));
 101 }
 102 //--------------------------------------------------------------------
 103 void clitk::DicomRT_ROI_ConvertToImageFilter::SetOutputSpacing(const double* spacing)
 104 {
 105   std::copy(spacing,spacing+3,std::back_inserter(mSpacing));
 106 }
 107 //--------------------------------------------------------------------
 108 void clitk::DicomRT_ROI_ConvertToImageFilter::SetOutputSize(const unsigned long* size)
 109 {
 110   std::copy(size,size+3,std::back_inserter(mSize));
 111 }
 112
 113 //--------------------------------------------------------------------
 114 void clitk::DicomRT_ROI_ConvertToImageFilter::Update()
 115 {
 116   if (!mROI) {
 117     std::cerr << "Error. No ROI set, please use SetROI." << std::endl;
 118     exit(0);
 119   }
 120   if (!ImageInfoIsSet()) {
 121     std::cerr << "Error. Please provide image info (spacing/origin) with SetImageFilename" << std::endl;
 122     exit(0);
 123   }
 124   // DD("Update");
 125
 126   // Get Mesh
 127   vtkPolyData * mesh = mROI->GetMesh();
 128   DD(mesh->GetNumberOfCells());
 129
 130   // Get bounds
 131   double *bounds=mesh->GetBounds();
 132   // for(int i=0; i<6; i++){
 133 //     DD(bounds[i]);
 134 //   }
 135
 136   // Compute origin
 137   std::vector<double> origin;
 138   origin.resize(3);
 139   origin[0] = floor((bounds[0]-mOrigin[0])/mSpacing[0]-2)*mSpacing[0]+mOrigin[0];
 140   origin[1] = floor((bounds[2]-mOrigin[1])/mSpacing[1]-2)*mSpacing[1]+mOrigin[1];
 141   origin[2] = floor((bounds[4]-mOrigin[2])/mSpacing[2]-2)*mSpacing[2]+mOrigin[2];
 142
 143   // Compute extend
 144   std::vector<double> extend;
 145   extend.resize(3);
 146   extend[0] = ceil((bounds[1]-origin[0])/mSpacing[0]+4);
 147   extend[1] = ceil((bounds[3]-origin[1])/mSpacing[1]+4);
 148   extend[2] = ceil((bounds[5]-origin[2])/mSpacing[2]+4);
 149
 150   // If no crop, set initial image size/origin
 151   if (!mCropMask) {
 152     for(int i=0; i<3; i++) {
 153       origin[i] = mOrigin[i];
 154       extend[i] = mSize[i]-1;
 155     }
 156   }
 157
 158   // Create new output image
 159   mBinaryImage = vtkSmartPointer<vtkImageData>::New();
 160   mBinaryImage->SetScalarTypeToUnsignedChar();
 161   mBinaryImage->SetOrigin(&origin[0]);
 162   mBinaryImage->SetSpacing(&mSpacing[0]);
 163   mBinaryImage->SetExtent(0, extend[0],
 164                           0, extend[1],
 165                           0, extend[2]);
 166   mBinaryImage->AllocateScalars();
 167
 168   // for(int i=0; i<3; i++){
 169   //     DD(origin[i]);
 170   //     DD(extend[i]);
 171   //     DD(mBinaryImage->GetDimensions()[i]);
 172   //   }
 173   memset(mBinaryImage->GetScalarPointer(), 0,
 174          mBinaryImage->GetDimensions()[0]*mBinaryImage->GetDimensions()[1]*mBinaryImage->GetDimensions()[2]*sizeof(unsigned char));
 175
 176   // Extrude
 177   vtkSmartPointer<vtkLinearExtrusionFilter> extrude=vtkSmartPointer<vtkLinearExtrusionFilter>::New();
 178   extrude->SetInput(mesh);
 179   extrude->SetVector(0, 0, mROI->GetContourSpacing());
 180
 181   // Binarization
 182   vtkSmartPointer<vtkPolyDataToImageStencil> sts=vtkSmartPointer<vtkPolyDataToImageStencil>::New();
 183   //The following line is extremely important
 184   //http://www.nabble.com/Bug-in-vtkPolyDataToImageStencil--td23368312.html#a23370933
 185   sts->SetTolerance(0);
 186   sts->SetInformationInput(mBinaryImage);
 187   sts->SetInput(extrude->GetOutput());
 188   //sts->SetInput(mesh);
 189
 190   vtkSmartPointer<vtkImageStencil> stencil=vtkSmartPointer<vtkImageStencil>::New();
 191   stencil->SetStencil(sts->GetOutput());
 192   stencil->SetInput(mBinaryImage);
 193   stencil->ReverseStencilOn();
 194   stencil->Update();
 195   mBinaryImage->ShallowCopy(stencil->GetOutput());
 196
 197   if (mWriteOutput) {
 198     typedef itk::Image<unsigned char, 3> ImageType;
 199     typedef itk::VTKImageToImageFilter<ImageType> ConnectorType;
 200     ConnectorType::Pointer connector = ConnectorType::New();
 201     connector->SetInput(GetOutput());
 202     connector->Update();
 203     clitk::writeImage<ImageType>(connector->GetOutput(), mOutputFilename);
 204   }
 205 }
 206 //--------------------------------------------------------------------
 207
 208
 209
 210 //--------------------------------------------------------------------
 211 vtkImageData * clitk::DicomRT_ROI_ConvertToImageFilter::GetOutput()
 212 {
 213   assert(mBinaryImage);
 214   return mBinaryImage;
 215 }
 216 //--------------------------------------------------------------------
 217