common/clitkDicomRT_ROI.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 "clitkDicomRT_ROI.h"
  21 #include <vtkSmartPointer.h>
  22 #include <vtkAppendPolyData.h>
  23 #include <vtkImageClip.h>
  24 #include <vtkMarchingSquares.h>
  25 #include <vtkPolyDataWriter.h>
  26 #include <vtkVersion.h>
  27
  28 #if GDCM_MAJOR_VERSION >= 2
  29 #include "gdcmAttribute.h"
  30 #include "gdcmItem.h"
  31 #endif
  32
  33 //--------------------------------------------------------------------
  34 clitk::DicomRT_ROI::DicomRT_ROI()
  35 {
  36   mName = "NoName";
  37   mNumber = -1;
  38   mImage = NULL;
  39   mColor.resize(3);
  40   mColor[0] = mColor[1] = mColor[2] = 0;
  41   mMeshIsUpToDate = false;
  42   mBackgroundValue = 0;
  43   mForegroundValue = 1;
  44   SetDicomUptodateFlag(false);
  45   mFilename = "";
  46 }
  47 //--------------------------------------------------------------------
  48
  49
  50 //--------------------------------------------------------------------
  51 clitk::DicomRT_ROI::~DicomRT_ROI()
  52 {
  53 }
  54 //--------------------------------------------------------------------
  55
  56
  57 //--------------------------------------------------------------------
  58 void clitk::DicomRT_ROI::SetDisplayColor(double r, double v, double b)
  59 {
  60   mColor.resize(3);
  61   mColor[0] = r;
  62   mColor[1] = v;
  63   mColor[2] = b;
  64 }
  65 //--------------------------------------------------------------------
  66
  67
  68 //--------------------------------------------------------------------
  69 int clitk::DicomRT_ROI::GetROINumber() const
  70 {
  71   return mNumber;
  72 }
  73 //--------------------------------------------------------------------
  74
  75
  76 //--------------------------------------------------------------------
  77 void clitk::DicomRT_ROI::SetROINumber(int number)
  78 {
  79   mNumber = number;
  80 }
  81 //--------------------------------------------------------------------
  82
  83
  84 //--------------------------------------------------------------------
  85 const std::string & clitk::DicomRT_ROI::GetName() const
  86 {
  87   return mName;
  88 }
  89 //--------------------------------------------------------------------
  90
  91
  92 //--------------------------------------------------------------------
  93 const std::string & clitk::DicomRT_ROI::GetFilename() const
  94 {
  95   return mFilename;
  96 }
  97 //--------------------------------------------------------------------
  98
  99
 100 //--------------------------------------------------------------------
 101 const std::vector<double> & clitk::DicomRT_ROI::GetDisplayColor() const
 102 {
 103   return mColor;
 104 }
 105 //--------------------------------------------------------------------
 106
 107
 108 //--------------------------------------------------------------------
 109 void clitk::DicomRT_ROI::Print(std::ostream & os) const
 110 {
 111   os << "ROI " << mNumber << "\t" << mName
 112      << "\t(" << mColor[0] << " " << mColor[1] << " " << mColor[2] << ")"
 113      << "\t Contours = " << mListOfContours.size() << std::endl;
 114 }
 115 //--------------------------------------------------------------------
 116
 117
 118 //--------------------------------------------------------------------
 119 void clitk::DicomRT_ROI::SetBackgroundValueLabelImage(double bg)
 120 {
 121   mBackgroundValue = bg;
 122 }
 123 //--------------------------------------------------------------------
 124
 125
 126 //--------------------------------------------------------------------
 127 double clitk::DicomRT_ROI::GetBackgroundValueLabelImage() const
 128 {
 129   return mBackgroundValue;
 130 }
 131 //--------------------------------------------------------------------
 132
 133
 134 //--------------------------------------------------------------------
 135 void clitk::DicomRT_ROI::SetForegroundValueLabelImage(double bg)
 136 {
 137   mForegroundValue = bg;
 138 }
 139 //--------------------------------------------------------------------
 140
 141
 142 //--------------------------------------------------------------------
 143 double clitk::DicomRT_ROI::GetForegroundValueLabelImage() const
 144 {
 145   return mForegroundValue;
 146 }
 147 //--------------------------------------------------------------------
 148
 149
 150 //--------------------------------------------------------------------
 151 #if GDCM_MAJOR_VERSION >= 2
 152 bool clitk::DicomRT_ROI::Read(gdcm::Item * itemInfo, gdcm::Item * itemContour, double tol)
 153 {
 154   //FATAL("Error : compile vv with itk4 + external gdcm");
 155   // Keep dicom item
 156   mItemInfo = itemInfo;
 157   mItemContour = itemContour;
 158   // DD(mItemInfo);
 159
 160   // ROI number [Referenced ROI Number]
 161   const gdcm::DataSet & nesteddsInfo = mItemInfo->GetNestedDataSet();
 162   gdcm::Attribute<0x3006,0x0022> roinumber;
 163   roinumber.SetFromDataSet( nesteddsInfo );
 164   int nb1 = roinumber.GetValue();
 165
 166   // Check this is the same with the other item
 167   const gdcm::DataSet & nestedds = mItemContour->GetNestedDataSet();
 168   gdcm::Attribute<0x3006,0x0084> referencedroinumber;
 169   referencedroinumber.SetFromDataSet( nestedds );
 170   int nb2 = referencedroinumber.GetValue();
 171
 172   // Must never be different
 173   if (nb1 != nb2) {
 174     DD(nb2);
 175     DD(nb1);
 176     FATAL("nb1 must equal nb2" << std::endl);
 177   }
 178   mNumber = nb1;
 179
 180   // Retrieve ROI Name (in the info item)
 181   gdcm::Attribute<0x3006,0x26> roiname;
 182   roiname.SetFromDataSet( nesteddsInfo );
 183   mName = roiname.GetValue();
 184   // DD(mName);
 185
 186   // ROI Color [ROI Display Color]
 187   gdcm::Attribute<0x3006,0x002a> color = {};
 188   color.SetFromDataSet( nestedds );
 189   assert( color.GetNumberOfValues() == 3 );
 190   mColor[0] = color.GetValue(0);
 191   mColor[1] = color.GetValue(1);
 192   mColor[2] = color.GetValue(2);
 193
 194   // Read contours [Contour Sequence]
 195   gdcm::Tag tcsq(0x3006,0x0040);
 196   if( !nestedds.FindDataElement( tcsq ) )
 197     {
 198       std::cerr << "Warning. Could not read contour for structure <" << mName << ">, number" << mNumber << " ? I ignore it" << std::endl;
 199       SetDicomUptodateFlag(true);
 200       return false;
 201     }
 202   const gdcm::DataElement& csq = nestedds.GetDataElement( tcsq );
 203   mContoursSequenceOfItems = csq.GetValueAsSQ();
 204   gdcm::SmartPointer<gdcm::SequenceOfItems> & sqi2 = mContoursSequenceOfItems;
 205   if( !sqi2 || !sqi2->GetNumberOfItems() )
 206     {
 207     }
 208   unsigned int nitems = sqi2->GetNumberOfItems();
 209
 210   for(unsigned int i = 0; i < nitems; ++i)
 211     {
 212       gdcm::Item & j = sqi2->GetItem(i+1); // Item start at #1
 213       DicomRT_Contour::Pointer c = DicomRT_Contour::New();
 214       c->SetTolerance(tol);
 215       c->SetTransformMatrix(mTransformMatrix);
 216       bool b = c->Read(&j);
 217       if (b) {
 218         mListOfContours.push_back(c);
 219       }
 220     }
 221   SetDicomUptodateFlag(true);
 222   return true;
 223 }
 224 #else
 225 void clitk::DicomRT_ROI::Read(std::map<int, std::string> & rois, gdcm::SQItem * item, double tol)
 226 {
 227   // ROI number [Referenced ROI Number]
 228   mNumber = atoi(item->GetEntryValue(0x3006,0x0084).c_str());
 229
 230   // Retrieve ROI Name
 231   mName = rois[mNumber];
 232
 233   // ROI Color [ROI Display Color]
 234   mColor = clitk::parse_string<double>(item->GetEntryValue(0x3006,0x002a),'\\');
 235
 236   // Read contours [Contour Sequence]
 237   gdcm::SeqEntry * contours=item->GetSeqEntry(0x3006,0x0040);
 238   if (contours) {
 239     int i=0;
 240     for(gdcm::SQItem* j=contours->GetFirstSQItem(); j!=0; j=contours->GetNextSQItem()) {
 241       DicomRT_Contour::Pointer c = DicomRT_Contour::New();
 242       c->SetTolerance(tol);
 243       c->SetTransformMatrix(mTransformMatrix);
 244       bool b = c->Read(j);
 245       if (b) {
 246         mListOfContours.push_back(c);
 247       }
 248       ++i;
 249     }
 250   }
 251   else {
 252     std::cerr << "Warning. Could not read contour for structure <" << mName << ">, number" << mNumber << " ? I ignore it" << std::endl;
 253   }
 254   SetDicomUptodateFlag(true);
 255 }
 256 #endif
 257 //--------------------------------------------------------------------
 258
 259
 260 //--------------------------------------------------------------------
 261 void clitk::DicomRT_ROI::SetImage(vvImage::Pointer image)
 262 {
 263   mImage = image;
 264 }
 265 //--------------------------------------------------------------------
 266
 267
 268 //--------------------------------------------------------------------
 269 vtkPolyData * clitk::DicomRT_ROI::GetMesh()
 270 {
 271   if (!mMeshIsUpToDate) {
 272     ComputeMeshFromContour();
 273   }
 274   return mMesh;
 275 }
 276 //--------------------------------------------------------------------
 277
 278
 279 //--------------------------------------------------------------------
 280 clitk::DicomRT_Contour * clitk::DicomRT_ROI::GetContour(int n)
 281 {
 282   return mListOfContours[n];
 283 }
 284 //--------------------------------------------------------------------
 285
 286
 287 //--------------------------------------------------------------------
 288 void clitk::DicomRT_ROI::ComputeMeshFromContour()
 289 {
 290   vtkSmartPointer<vtkAppendPolyData> append = vtkSmartPointer<vtkAppendPolyData>::New();
 291   for(unsigned int i=0; i<mListOfContours.size(); i++) {
 292 #if VTK_MAJOR_VERSION <= 5
 293         append->AddInput(mListOfContours[i]->GetMesh());
 294 #else
 295         append->AddInputData(mListOfContours[i]->GetMesh());
 296 #endif
 297   }
 298   append->Update();
 299
 300   mMesh = vtkSmartPointer<vtkPolyData>::New();
 301   mMesh->DeepCopy(append->GetOutput());
 302   mMeshIsUpToDate = true;
 303 }
 304 //--------------------------------------------------------------------
 305
 306
 307 #if GDCM_MAJOR_VERSION >= 2
 308 //--------------------------------------------------------------------
 309 void clitk::DicomRT_ROI::UpdateDicomItem()
 310 {
 311   FATAL("Error : compile vv with itk4 + external gdcm");
 312
 313   if (GetDicomUptoDateFlag()) return;
 314   DD("ROI::UpdateDicomItem");
 315   DD(GetName());
 316
 317   // From now, only some item can be modified
 318
 319   // Set ROI Name 0x3006,0x26>
 320   gdcm::Attribute<0x3006,0x26> roiname;
 321   roiname.SetValue(GetName());
 322   gdcm::DataElement de = roiname.GetAsDataElement();
 323   gdcm::DataSet & ds = mItemInfo->GetNestedDataSet();
 324   ds.Replace(de);
 325
 326   // From MESH to CONTOURS
 327   ComputeContoursFromImage();
 328
 329   // Update contours
 330   DD(mListOfContours.size());
 331   for(uint i=0; i<mListOfContours.size(); i++) {
 332     DD(i);
 333     DicomRT_Contour::Pointer contour = mListOfContours[i];
 334     contour->UpdateDicomItem();//mItemContour);
 335   }
 336
 337   // Nb of contours
 338   unsigned int nitems = mContoursSequenceOfItems->GetNumberOfItems();
 339   DD(nitems);
 340
 341   // Write [Contour Sequence] = 0x3006,0x0040)
 342   gdcm::DataSet & dsc = mItemContour->GetNestedDataSet();
 343   gdcm::Tag tcsq(0x3006,0x0040);
 344   const gdcm::DataElement& csq = dsc.GetDataElement( tcsq );
 345   gdcm::DataElement dec(csq);
 346   dec.SetValue(*mContoursSequenceOfItems);
 347   dsc.Replace(dec);
 348
 349   gdcm::DataSet & a = mContoursSequenceOfItems->GetItem(1).GetNestedDataSet();
 350   gdcm::Attribute<0x3006,0x0050> at;
 351   gdcm::Tag tcontourdata(0x3006,0x0050);
 352   gdcm::DataElement contourdata = a.GetDataElement( tcontourdata );
 353   at.SetFromDataElement( contourdata );
 354   const double* points = at.GetValues();
 355   DD(points[0]);
 356
 357 }
 358 //--------------------------------------------------------------------
 359 #endif
 360
 361 //--------------------------------------------------------------------
 362 void clitk::DicomRT_ROI::SetFromBinaryImage(vvImage::Pointer image, int n,
 363                                             std::string name,
 364                                             std::vector<double> color,
 365                                             std::string filename)
 366 {
 367   // ROI number [Referenced ROI Number]
 368   mNumber = n;
 369
 370   // ROI Name
 371   mName = name;
 372   mFilename = filename;
 373
 374   // ROI Color [ROI Display Color]
 375   mColor = color;
 376
 377   // No contours [Contour Sequence]
 378   mListOfContours.clear();
 379
 380   // Set image
 381   mImage = image;
 382 }
 383 //--------------------------------------------------------------------
 384
 385
 386 //--------------------------------------------------------------------
 387 vvImage * clitk::DicomRT_ROI::GetImage() const
 388 {
 389   return mImage;
 390 }
 391 //--------------------------------------------------------------------
 392
 393
 394 //--------------------------------------------------------------------
 395 void clitk::DicomRT_ROI::SetTransformMatrix(vtkMatrix4x4* matrix)
 396 {
 397   mTransformMatrix = matrix;
 398 }
 399 //--------------------------------------------------------------------
 400
 401
 402 //--------------------------------------------------------------------
 403 void clitk::DicomRT_ROI::ComputeContoursFromImage()
 404 {
 405   FATAL("ComputeContoursFromImage should not be call. To be replace");
 406   DD("ComputeMeshFromImage");
 407
 408   // Check that an image is loaded
 409   if (!mImage) return;
 410
 411   // Only consider 3D here
 412   if (mImage->GetNumberOfDimensions() != 3) {
 413     FATAL("DicomRT_ROI::ComputeMeshFromImage only work with 3D images");
 414   }
 415
 416   // Get the VTK image
 417   vtkImageData * image = mImage->GetVTKImages()[0];
 418
 419   // Get initial extend for the clipping
 420   vtkSmartPointer<vtkImageClip> clipper = vtkSmartPointer<vtkImageClip>::New();
 421 #if VTK_MAJOR_VERSION <= 5
 422   clipper->SetInput(image);
 423 #else
 424   clipper->SetInputData(image);
 425 #endif
 426
 427   int* extent = image->GetExtent();
 428   DDV(extent, 6);
 429   //  std::vector<int> extend;
 430
 431
 432   // Loop on slice
 433   uint n = image->GetDimensions()[2];
 434   DD(n);
 435   DD(mListOfContours.size());
 436   mListOfContours.resize(n); /// ???FIXME
 437   DD(mListOfContours.size());
 438   std::vector<vtkSmartPointer<vtkPolyData> > contours;
 439   for(uint i=0; i<n; i++) {
 440     DD(i);
 441
 442     // FIXME     vtkDiscreteMarchingCubes INSTEAD
 443
 444
 445     vtkSmartPointer<vtkMarchingSquares> squares = vtkSmartPointer<vtkMarchingSquares>::New();
 446 #if VTK_MAJOR_VERSION <= 5
 447     squares->SetInput(image);
 448 #else
 449     squares->SetInputData(image);
 450 #endif
 451     squares->SetImageRange(extent[0], extent[1], extent[2], extent[3], i, i);
 452     squares->SetValue(1, 1.0);
 453     squares->Update();
 454     DD(squares->GetNumberOfContours());
 455
 456     //clitk::DicomRT_Contour * contour = new clitk::DicomRT_Contour();
 457     //mListOfContours[i]->SetMesh(squares->GetOutput());
 458
 459
 460     vtkSmartPointer<vtkPolyData> m = squares->GetOutput();
 461     contours.push_back(m);
 462
 463     /*
 464     // Clip to the current slice
 465     extent[4] = extent[5] = image->GetOrigin()[2]+i*image->GetSpacing()[2];
 466     DDV(extent, 6);
 467     // Prepare the marching squares
 468     vtkSmartPointer<vtkMarchingSquares> squares = vtkSmartPointer<vtkMarchingSquares>::New();
 469     clipper->SetOutputWholeExtent(extent[0],extent[1],extent[2],
 470                                 extent[3],extent[4],extent[5]);
 471
 472     squares->SetInput(clipper->GetOutput());
 473     squares->Update();
 474     DD(squares->GetNumberOfContours());
 475     mListOfContours[i]->SetMesh(squares->GetOutput());
 476     */
 477   }
 478   DD("done");
 479
 480   vtkSmartPointer<vtkAppendPolyData> append = vtkSmartPointer<vtkAppendPolyData>::New();
 481   for(unsigned int i=0; i<n; i++) {
 482 #if VTK_MAJOR_VERSION <= 5
 483     append->AddInput(contours[i]);
 484 #else
 485     append->AddInputData(contours[i]);
 486 #endif
 487   }
 488   append->Update();
 489
 490   mMesh = vtkSmartPointer<vtkPolyData>::New();
 491   mMesh->DeepCopy(append->GetOutput());
 492
 493   // Write vtk
 494   vtkPolyDataWriter * w = vtkPolyDataWriter::New();
 495 #if VTK_MAJOR_VERSION <= 5
 496   w->SetInput(mMesh);
 497 #else
 498   w->SetInputData(mMesh);
 499 #endif
 500   w->SetFileName("toto.vtk");
 501   w->Write();
 502
 503   DD("done");
 504 }
 505 //--------------------------------------------------------------------
 506
 507
 508 //--------------------------------------------------------------------
 509 #if CLITK_USE_SYSTEM_GDCM == 1
 510 void clitk::DicomRT_ROI::Read(vtkSmartPointer<vtkGDCMPolyDataReader> & reader, int roiindex, double tol)
 511 {
 512   vtkRTStructSetProperties * p = reader->GetRTStructSetProperties();
 513
 514   mName = p->GetStructureSetROIName(roiindex);
 515   mNumber = p->GetStructureSetROINumber(roiindex);
 516
 517   //mColor = //FIXME !!
 518
 519   // gdcm::Attribute<0x3006,0x002a> color = {};
 520
 521   // const gdcm::DataSet & nestedds = mItemContour->GetNestedDataSet();
 522   // color.SetFromDataSet( nestedds );
 523   // assert( color.GetNumberOfValues() == 3 );
 524   // mColor[0] = color.GetValue(0);
 525   // mColor[1] = color.GetValue(1);
 526   // mColor[2] = color.GetValue(2);
 527
 528
 529   SetDicomUptodateFlag(true);
 530   // Get the contour
 531   mMesh =  reader->GetOutput(roiindex);
 532   DicomRT_Contour::Pointer c = DicomRT_Contour::New();
 533   c->SetTolerance(tol);
 534   c->SetTransformMatrix(mTransformMatrix);
 535   c->SetMesh(mMesh); // FIXME no GetZ, not GetPoints
 536   mMeshIsUpToDate = true;
 537   mListOfContours.push_back(c);
 538 }
 539 #endif
 540 //--------------------------------------------------------------------
 541