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