common/clitkDicomRTStruct2ImageFilter.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
  23 // clitk
  24 #include "clitkDicomRTStruct2ImageFilter.h"
  25 #include "clitkImageCommon.h"
  26
  27 // vtk
  28 #include <vtkPolyDataToImageStencil.h>
  29 #include <vtkSmartPointer.h>
  30 #include <vtkImageStencil.h>
  31 #include <vtkLinearExtrusionFilter.h>
  32 #include <vtkMetaImageWriter.h>
  33
  34
  35 //--------------------------------------------------------------------
  36 clitk::DicomRTStruct2ImageFilter::DicomRTStruct2ImageFilter()
  37 {
  38   mROI = NULL;
  39   mWriteOutput = false;
  40   mCropMask = true;
  41 }
  42 //--------------------------------------------------------------------
  43
  44
  45 //--------------------------------------------------------------------
  46 clitk::DicomRTStruct2ImageFilter::~DicomRTStruct2ImageFilter()
  47 {
  48
  49 }
  50 //--------------------------------------------------------------------
  51
  52
  53 //--------------------------------------------------------------------
  54 bool clitk::DicomRTStruct2ImageFilter::ImageInfoIsSet() const
  55 {
  56   return mSize.size() && mSpacing.size() && mOrigin.size();
  57 }
  58 //--------------------------------------------------------------------
  59
  60
  61 //--------------------------------------------------------------------
  62 void clitk::DicomRTStruct2ImageFilter::SetWriteOutputFlag(bool b)
  63 {
  64   mWriteOutput = b;
  65 }
  66 //--------------------------------------------------------------------
  67
  68
  69 //--------------------------------------------------------------------
  70 void clitk::DicomRTStruct2ImageFilter::SetROI(clitk::DicomRT_ROI * roi)
  71 {
  72   mROI = roi;
  73 }
  74 //--------------------------------------------------------------------
  75
  76
  77 //--------------------------------------------------------------------
  78 void clitk::DicomRTStruct2ImageFilter::SetCropMaskEnabled(bool b)
  79 {
  80   mCropMask = b;
  81 }
  82 //--------------------------------------------------------------------
  83
  84
  85 //--------------------------------------------------------------------
  86 void clitk::DicomRTStruct2ImageFilter::SetOutputImageFilename(std::string s)
  87 {
  88   mOutputFilename = s;
  89   mWriteOutput = true;
  90 }
  91 //--------------------------------------------------------------------
  92
  93
  94 //--------------------------------------------------------------------
  95 void clitk::DicomRTStruct2ImageFilter::SetImage(vvImage::Pointer image)
  96 {
  97   if (image->GetNumberOfDimensions() != 3) {
  98     std::cerr << "Error. Please provide a 3D image." << std::endl;
  99     exit(0);
 100   }
 101   mSpacing.resize(3);
 102   mOrigin.resize(3);
 103   mSize.resize(3);
 104   for(unsigned int i=0; i<3; i++) {
 105     mSpacing[i] = image->GetSpacing()[i];
 106     mOrigin[i] = image->GetOrigin()[i];
 107     mSize[i] = image->GetSize()[i];
 108   }
 109 }
 110 //--------------------------------------------------------------------
 111
 112
 113 //--------------------------------------------------------------------
 114 void clitk::DicomRTStruct2ImageFilter::SetImageFilename(std::string f)
 115 {
 116   itk::ImageIOBase::Pointer header = clitk::readImageHeader(f);
 117   if (header->GetNumberOfDimensions() < 3) {
 118     std::cerr << "Error. Please provide a 3D image instead of " << f << std::endl;
 119     exit(0);
 120   }
 121   if (header->GetNumberOfDimensions() > 3) {
 122     std::cerr << "Warning dimension > 3 are ignored" << std::endl;
 123   }
 124   mSpacing.resize(3);
 125   mOrigin.resize(3);
 126   mSize.resize(3);
 127   for(unsigned int i=0; i<3; i++) {
 128     mSpacing[i] = header->GetSpacing(i);
 129     mOrigin[i] = header->GetOrigin(i);
 130     mSize[i] = header->GetDimensions(i);
 131   }
 132 }
 133 //--------------------------------------------------------------------
 134
 135
 136 //--------------------------------------------------------------------
 137 void clitk::DicomRTStruct2ImageFilter::SetOutputOrigin(const double* origin)
 138 {
 139   std::copy(origin,origin+3,std::back_inserter(mOrigin));
 140 }
 141 //--------------------------------------------------------------------
 142
 143
 144 //--------------------------------------------------------------------
 145 void clitk::DicomRTStruct2ImageFilter::SetOutputSpacing(const double* spacing)
 146 {
 147   std::copy(spacing,spacing+3,std::back_inserter(mSpacing));
 148 }
 149 //--------------------------------------------------------------------
 150
 151
 152 //--------------------------------------------------------------------
 153 void clitk::DicomRTStruct2ImageFilter::SetOutputSize(const unsigned long* size)
 154 {
 155   std::copy(size,size+3,std::back_inserter(mSize));
 156 }
 157 //--------------------------------------------------------------------
 158
 159
 160 //--------------------------------------------------------------------
 161 void clitk::DicomRTStruct2ImageFilter::Update()
 162 {
 163   if (!mROI) {
 164     std::cerr << "Error. No ROI set, please use SetROI." << std::endl;
 165     exit(0);
 166   }
 167   if (!ImageInfoIsSet()) {
 168     std::cerr << "Error. Please provide image info (spacing/origin) with SetImageFilename" << std::endl;
 169     exit(0);
 170   }
 171
 172   // Get Mesh
 173   vtkPolyData * mesh = mROI->GetMesh();
 174
 175   // Get bounds
 176   double *bounds=mesh->GetBounds();
 177   // for(int i=0; i<6; i++){
 178 //     DD(bounds[i]);
 179 //   }
 180
 181   // Compute origin
 182   std::vector<double> origin;
 183   origin.resize(3);
 184   origin[0] = floor((bounds[0]-mOrigin[0])/mSpacing[0]-2)*mSpacing[0]+mOrigin[0];
 185   origin[1] = floor((bounds[2]-mOrigin[1])/mSpacing[1]-2)*mSpacing[1]+mOrigin[1];
 186   origin[2] = floor((bounds[4]-mOrigin[2])/mSpacing[2]-2)*mSpacing[2]+mOrigin[2];
 187
 188   // Compute extend
 189   std::vector<double> extend;
 190   extend.resize(3);
 191   extend[0] = ceil((bounds[1]-origin[0])/mSpacing[0]+4);
 192   extend[1] = ceil((bounds[3]-origin[1])/mSpacing[1]+4);
 193   extend[2] = ceil((bounds[5]-origin[2])/mSpacing[2]+4);
 194
 195   // If no crop, set initial image size/origin
 196   if (!mCropMask) {
 197     for(int i=0; i<3; i++) {
 198       origin[i] = mOrigin[i];
 199       extend[i] = mSize[i]-1;
 200     }
 201   }
 202
 203   // Create new output image
 204   mBinaryImage = vtkSmartPointer<vtkImageData>::New();
 205   mBinaryImage->SetScalarTypeToUnsignedChar();
 206   mBinaryImage->SetOrigin(&origin[0]);
 207   mBinaryImage->SetSpacing(&mSpacing[0]);
 208   mBinaryImage->SetExtent(0, extend[0],
 209                           0, extend[1],
 210                           0, extend[2]);
 211   mBinaryImage->AllocateScalars();
 212
 213   // for(int i=0; i<3; i++){
 214   //     DD(origin[i]);
 215   //     DD(extend[i]);
 216   //     DD(mBinaryImage->GetDimensions()[i]);
 217   //   }
 218   memset(mBinaryImage->GetScalarPointer(), 0,
 219          mBinaryImage->GetDimensions()[0]*mBinaryImage->GetDimensions()[1]*mBinaryImage->GetDimensions()[2]*sizeof(unsigned char));
 220
 221   // Extrude
 222   vtkSmartPointer<vtkLinearExtrusionFilter> extrude=vtkSmartPointer<vtkLinearExtrusionFilter>::New();
 223   extrude->SetInput(mesh);
 224   ///We extrude in the -slice_spacing direction to respect the FOCAL convention (NEEDED !)
 225   extrude->SetVector(0, 0, -mSpacing[2]);
 226
 227   // Binarization
 228   vtkSmartPointer<vtkPolyDataToImageStencil> sts=vtkSmartPointer<vtkPolyDataToImageStencil>::New();
 229   //The following line is extremely important
 230   //http://www.nabble.com/Bug-in-vtkPolyDataToImageStencil--td23368312.html#a23370933
 231   sts->SetTolerance(0);
 232   sts->SetInformationInput(mBinaryImage);
 233   sts->SetInput(extrude->GetOutput());
 234   //sts->SetInput(mesh);
 235
 236   vtkSmartPointer<vtkImageStencil> stencil=vtkSmartPointer<vtkImageStencil>::New();
 237   stencil->SetStencil(sts->GetOutput());
 238   stencil->SetInput(mBinaryImage);
 239   stencil->ReverseStencilOn();
 240   stencil->Update();
 241
 242   /*
 243   vtkSmartPointer<vtkMetaImageWriter> w = vtkSmartPointer<vtkMetaImageWriter>::New();
 244   w->SetInput(stencil->GetOutput());
 245   w->SetFileName("binary2.mhd");
 246   w->Write();
 247   */
 248
 249   mBinaryImage->ShallowCopy(stencil->GetOutput());
 250
 251   if (mWriteOutput) {
 252     typedef itk::Image<unsigned char, 3> ImageType;
 253     typedef itk::VTKImageToImageFilter<ImageType> ConnectorType;
 254     ConnectorType::Pointer connector = ConnectorType::New();
 255     connector->SetInput(GetOutput());
 256     connector->Update();
 257     clitk::writeImage<ImageType>(connector->GetOutput(), mOutputFilename);
 258   }
 259 }
 260 //--------------------------------------------------------------------
 261
 262
 263
 264 //--------------------------------------------------------------------
 265 vtkImageData * clitk::DicomRTStruct2ImageFilter::GetOutput()
 266 {
 267   assert(mBinaryImage);
 268   return mBinaryImage;
 269 }
 270 //--------------------------------------------------------------------
 271