common/clitkDicomRT_StructureSet.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_StructureSet.h"
  21 #include <vtksys/SystemTools.hxx>
  22 #include "gdcmFile.h"
  23
  24 //--------------------------------------------------------------------
  25 clitk::DicomRT_StructureSet::DicomRT_StructureSet()
  26 {
  27   mStudyID = "NoStudyID";
  28   mStudyTime = "NoStudyTime";
  29   mStudyDate = "NoStudyDate";
  30   mLabel = "NoLabel";
  31   mName = "NoName";
  32   mDate = "NoDate";
  33   mTime = "NoTime";
  34   mFile = NULL;
  35 }
  36 //--------------------------------------------------------------------
  37
  38
  39 //--------------------------------------------------------------------
  40 clitk::DicomRT_StructureSet::~DicomRT_StructureSet()
  41 {
  42 }
  43 //--------------------------------------------------------------------
  44
  45
  46 //--------------------------------------------------------------------
  47 const std::string & clitk::DicomRT_StructureSet::GetStudyID() const
  48 {
  49   return mStudyID;
  50 }
  51 //--------------------------------------------------------------------
  52
  53
  54 //--------------------------------------------------------------------
  55 const std::string & clitk::DicomRT_StructureSet::GetStudyTime() const
  56 {
  57   return mStudyTime;
  58 }
  59 //--------------------------------------------------------------------
  60
  61
  62 //--------------------------------------------------------------------
  63 const std::string & clitk::DicomRT_StructureSet::GetStudyDate() const
  64 {
  65   return mStudyDate;
  66 }
  67 //--------------------------------------------------------------------
  68
  69
  70 //--------------------------------------------------------------------
  71 const std::string & clitk::DicomRT_StructureSet::GetLabel() const
  72 {
  73   return mLabel;
  74 }
  75 //--------------------------------------------------------------------
  76
  77
  78 //--------------------------------------------------------------------
  79 const std::string & clitk::DicomRT_StructureSet::GetName() const
  80 {
  81   return mName;
  82 }
  83 //--------------------------------------------------------------------
  84
  85
  86 //--------------------------------------------------------------------
  87 const std::string & clitk::DicomRT_StructureSet::GetDate() const
  88 {
  89   return mDate;
  90 }
  91 //--------------------------------------------------------------------
  92
  93
  94 //--------------------------------------------------------------------
  95 const std::string & clitk::DicomRT_StructureSet::GetTime() const
  96 {
  97   return mTime;
  98 }
  99 //--------------------------------------------------------------------
 100
 101
 102 //--------------------------------------------------------------------
 103 // const std::vector<clitk::DicomRT_ROI::Pointer> & clitk::DicomRT_StructureSet::GetListOfROI() const
 104 // {
 105 //   return mListOfROI;
 106 // }
 107 //--------------------------------------------------------------------
 108
 109
 110 //--------------------------------------------------------------------
 111 clitk::DicomRT_ROI* clitk::DicomRT_StructureSet::GetROIFromROINumber(int n)
 112 {
 113   if (mROIs.find(n) == mROIs.end()) {
 114     std::cerr << "No ROI number " << n << std::endl;
 115     return NULL;
 116   }
 117   return mROIs[n];
 118 }
 119 //--------------------------------------------------------------------
 120
 121
 122 //--------------------------------------------------------------------
 123 void clitk::DicomRT_StructureSet::Print(std::ostream & os) const
 124 {
 125   os << "Study ID      = " << mStudyID << std::endl
 126      << "Study Date    = " << mStudyDate << std::endl
 127      << "Study Time    = " << mStudyTime << std::endl
 128      << "Struct Label  = " << mLabel << std::endl
 129      << "Struct Name   = " << mName << std::endl
 130      << "Struct Time   = " << mTime << std::endl
 131      << "Number of ROI = " << mROIs.size() << std::endl;
 132   for(ROIConstIteratorType iter = mROIs.begin(); iter != mROIs.end(); iter++) {
 133     iter->second->Print(os);
 134   }
 135 }
 136 //--------------------------------------------------------------------
 137
 138
 139 #if GDCM_MAJOR_VERSION == 2
 140 //--------------------------------------------------------------------
 141 int clitk::DicomRT_StructureSet::ReadROINumber(const gdcm::Item & item)
 142 {
 143   // 0x3006,0x0022 = [ROI Number]
 144   const gdcm::DataSet & nestedds = item.GetNestedDataSet();
 145   gdcm::Attribute<0x3006,0x0022> roinumber;
 146   roinumber.SetFromDataSet( nestedds );
 147   return roinumber.GetValue();
 148 }
 149 //--------------------------------------------------------------------
 150 #endif
 151
 152 //--------------------------------------------------------------------
 153 void clitk::DicomRT_StructureSet::Write(const std::string & filename)
 154 {
 155 #if GDCM_MAJOR_VERSION == 2
 156   DD("DCM RT Writer");
 157
 158   // Assert that the gdcm file is still open (we can write only if it was readed)
 159   if (mFile == NULL) {
 160     //assert(mFile != NULL);
 161     FATAL("Sorry, I can write DICOM only if it was read first from a file with 'Read' function");
 162   }
 163
 164   // Loop and update each ROI
 165   for(ROIIteratorType iter = mROIs.begin(); iter != mROIs.end(); iter++) {
 166     iter->second->UpdateDicomItem();
 167   }
 168
 169   // Write [ Structure Set ROI Sequence ] = 0x3006,0x0020
 170   gdcm::DataSet & ds = mFile->GetDataSet();
 171   gdcm::Tag tssroisq(0x3006,0x0020);
 172   const gdcm::DataElement &ssroisq = ds.GetDataElement( tssroisq );
 173   gdcm::DataElement de(ssroisq);
 174   de.SetValue(*mROIInfoSequenceOfItems);
 175   ds.Replace(de);
 176
 177   // Write [ ROI Contour Sequence ] = 0x3006,0x0039
 178   DD("ici");
 179   gdcm::Tag troicsq(0x3006,0x0039);
 180   const gdcm::DataElement &roicsq = ds.GetDataElement( troicsq );
 181   gdcm::DataElement de2(roicsq);
 182   de2.SetValue(*mROIContoursSequenceOfItems);
 183   ds.Replace(de2);
 184
 185   //DEBUG
 186   gdcm::DataSet & a = mROIContoursSequenceOfItems->GetItem(1).GetNestedDataSet();
 187   gdcm::Tag tcsq(0x3006,0x0040);
 188   const gdcm::DataElement& csq = a.GetDataElement( tcsq );
 189   gdcm::SmartPointer<gdcm::SequenceOfItems> sqi2 = csq.GetValueAsSQ();
 190   gdcm::Item & j = sqi2->GetItem(1);
 191   gdcm::DataSet & b = j.GetNestedDataSet();
 192   gdcm::Attribute<0x3006,0x0050> at;
 193   gdcm::Tag tcontourdata(0x3006,0x0050);
 194   gdcm::DataElement contourdata = b.GetDataElement( tcontourdata );
 195   at.SetFromDataElement( contourdata );
 196   const double* points = at.GetValues();
 197   DD(points[0]);
 198
 199
 200   // Write dicom
 201   gdcm::Writer writer;
 202   //writer.CheckFileMetaInformationOff();
 203   writer.SetFileName(filename.c_str());
 204   writer.SetFile(*mFile);
 205   DD("before write");
 206   writer.Write();
 207   DD("End write");
 208 #else
 209   FATAL("Sorry not compatible with GDCM1, use GDCM2");
 210 #endif
 211 }
 212 //--------------------------------------------------------------------
 213
 214
 215 //--------------------------------------------------------------------
 216 void clitk::DicomRT_StructureSet::Read(const std::string & filename)
 217 {
 218   // Open DICOM
 219 #if GDCM_MAJOR_VERSION == 2
 220   // Read gdcm file
 221   mReader = new gdcm::Reader;
 222   mReader->SetFileName(filename.c_str());
 223   mReader->Read();
 224   mFile = &(mReader->GetFile());
 225   const gdcm::DataSet & ds = mFile->GetDataSet();
 226
 227   // Check file type
 228   //Verify if the file is a RT-Structure-Set dicom file
 229   gdcm::MediaStorage ms;
 230   ms.SetFromFile(*mFile);
 231   if( ms != gdcm::MediaStorage::RTStructureSetStorage )
 232     {
 233     std::cerr << "Error. the file " << filename
 234               << " is not a Dicom Struct ? (must have a SOP Class UID [0008|0016] = 1.2.840.10008.5.1.4.1.1.481.3 ==> [RT Structure Set Storage])"
 235               << std::endl;
 236     exit(0);
 237     }
 238
 239   gdcm::Attribute<0x8,0x60> modality;
 240   modality.SetFromDataSet( ds );
 241   if( modality.GetValue() != "RTSTRUCT" )
 242     {
 243     std::cerr << "Error. the file " << filename
 244               << " is not a Dicom Struct ? (must have 0x0008,0x0060 = RTSTRUCT [RT Structure Set Storage])"
 245               << std::endl;
 246     exit(0);
 247     }
 248
 249   // Read global info
 250   gdcm::Attribute<0x20,0x10> studyid;
 251   studyid.SetFromDataSet( ds );
 252   gdcm::Attribute<0x8,0x20> studytime;
 253   studytime.SetFromDataSet( ds );
 254   gdcm::Attribute<0x8,0x30> studydate;
 255   studydate.SetFromDataSet( ds );
 256   gdcm::Attribute<0x3006,0x02> label;
 257   label.SetFromDataSet( ds );
 258   gdcm::Attribute<0x3006,0x04> atname;
 259   atname.SetFromDataSet( ds );
 260   gdcm::Attribute<0x3006,0x09> time;
 261   time.SetFromDataSet( ds );
 262
 263   mStudyID   = studyid.GetValue();
 264   mStudyTime = studytime.GetValue();
 265   mStudyDate = studydate.GetValue();
 266   mLabel     = label.GetValue();
 267   mName      = atname.GetValue();
 268   mTime      = time.GetValue();
 269
 270   // Temporary store the list of items
 271   std::map<int, gdcm::Item*> mMapOfROIInfo;
 272   std::map<int, gdcm::Item*> mMapOfROIContours;
 273
 274   //----------------------------------
 275   // Read all ROI Names and number
 276   // 0x3006,0x0020 = [ Structure Set ROI Sequence ]
 277   gdcm::Tag tssroisq(0x3006,0x0020);
 278   const gdcm::DataElement &ssroisq = ds.GetDataElement( tssroisq );
 279   mROIInfoSequenceOfItems = ssroisq.GetValueAsSQ();
 280   gdcm::SmartPointer<gdcm::SequenceOfItems> & roi_seq = mROIInfoSequenceOfItems;
 281   assert(roi_seq); // TODO error message
 282   for(unsigned int ridx = 0; ridx < roi_seq->GetNumberOfItems(); ++ridx)
 283     {
 284     gdcm::Item & item = roi_seq->GetItem( ridx + 1); // Item starts at 1
 285     const gdcm::DataSet& nestedds = item.GetNestedDataSet();
 286
 287     gdcm::Attribute<0x3006,0x26> roiname;
 288     roiname.SetFromDataSet( nestedds );
 289     std::string name = roiname.GetValue(); // 0x3006,0x0026 = [ROI Name]
 290
 291     // 0x3006,0x0022 = [ROI Number]
 292     int nb = ReadROINumber(item);
 293
 294     // Store the item
 295     mMapOfROIInfo[nb] = &item;
 296
 297     // Check if such a number already exist
 298     if (mMapOfROIName.find(nb) != mMapOfROIName.end()) {
 299       std::cerr << "WARNING. A Roi already exist with the number "
 300         << nb << ". I replace." << std::endl;
 301     }
 302     // Add in map
 303     mMapOfROIName[nb] = name;
 304     }
 305
 306   //----------------------------------
 307   // Read all ROI item
 308   // 0x3006,0x0039 = [ ROI Contour Sequence ]
 309   gdcm::Tag troicsq(0x3006,0x0039);
 310   const gdcm::DataElement &roicsq = ds.GetDataElement( troicsq );
 311   gdcm::SmartPointer<gdcm::SequenceOfItems> roi_contour_seq = roicsq.GetValueAsSQ();
 312   mROIContoursSequenceOfItems = roi_contour_seq;
 313   assert(roi_contour_seq); // TODO error message
 314   for(unsigned int ridx = 0; ridx < roi_contour_seq->GetNumberOfItems(); ++ridx) {
 315     gdcm::Item & item = roi_contour_seq->GetItem( ridx + 1); // Item starts at 1
 316     // ROI number [Referenced ROI Number]
 317     const gdcm::DataSet& nestedds = item.GetNestedDataSet();
 318     gdcm::Attribute<0x3006,0x0084> referencedroinumber;
 319     referencedroinumber.SetFromDataSet( nestedds );
 320     int nb = referencedroinumber.GetValue();
 321     // Store the item
 322     mMapOfROIContours[nb] = &item;
 323   }
 324
 325   //----------------------------------
 326   // Create the ROIs
 327   for(std::map<int, gdcm::Item*>::iterator i = mMapOfROIInfo.begin(); i != mMapOfROIInfo.end(); i++) {
 328     int nb = i->first;//ReadROINumber(i);//mROIIndex[i];
 329     // Create the roi
 330     DicomRT_ROI::Pointer roi = DicomRT_ROI::New();
 331     roi->Read(mMapOfROIInfo[nb], mMapOfROIContours[nb]);
 332     //    mListOfROI.push_back(roi);
 333     //    mMapOfROIIndex[nb] = i;
 334     mROIs[nb] = roi;
 335   }
 336
 337   //----------------------------------------------------------------------------------------
 338   //----------------------------------------------------------------------------------------
 339   //----------------------------------------------------------------------------------------
 340 #else
 341   mFile = new gdcm::File;
 342   mFile->SetFileName(filename.c_str());
 343   mFile->SetMaxSizeLoadEntry(16384); // Needed ...
 344   mFile->SetLoadMode(gdcm::LD_NOSHADOW); // don't load shadow tags (in order to save memory)
 345   mFile->Load();
 346
 347   // Check file type
 348   //Verify if the file is a RT-Structure-Set dicom file
 349   if (!gdcm::Util::DicomStringEqual(mFile->GetEntryValue(0x0008,0x0016),"1.2.840.10008.5.1.4.1.1.481.3")) {  //SOP clas UID
 350     std::cerr << "Error. the file " << filename
 351               << " is not a Dicom Struct ? (must have a SOP Class UID [0008|0016] = 1.2.840.10008.5.1.4.1.1.481.3 ==> [RT Structure Set Storage])"
 352               << std::endl;
 353     exit(0);
 354   }
 355   if (!gdcm::Util::DicomStringEqual(mFile->GetEntryValue(0x0008,0x0060),"RTSTRUCT")) {  //SOP clas UID
 356     std::cerr << "Error. the file " << filename
 357               << " is not a Dicom Struct ? (must have 0x0008,0x0060 = RTSTRUCT [RT Structure Set Storage])"
 358               << std::endl;
 359     exit(0);
 360   }
 361
 362   // Read global info
 363   mStudyID   = mFile->GetValEntry(0x0020,0x0010)->GetValue();
 364   mStudyTime = mFile->GetValEntry(0x008,0x0020)->GetValue();
 365   mStudyDate = mFile->GetValEntry(0x008,0x0030)->GetValue();
 366   mLabel     = mFile->GetValEntry(0x3006,0x002)->GetValue();
 367   if (!mFile->GetValEntry(0x3006,0x004)) {
 368     mName = "Anonymous";
 369   }
 370   else {
 371     mName = mFile->GetValEntry(0x3006,0x004)->GetValue();
 372   }
 373   mTime      = mFile->GetValEntry(0x3006,0x009)->GetValue();
 374
 375   //----------------------------------
 376   // Read all ROI Names and number
 377   // 0x3006,0x0020 = [ Structure Set ROI Sequence ]
 378   gdcm::SeqEntry * roi_seq=mFile->GetSeqEntry(0x3006,0x0020);
 379   assert(roi_seq); // TODO error message
 380   for (gdcm::SQItem* r=roi_seq->GetFirstSQItem(); r!=0; r=roi_seq->GetNextSQItem()) {
 381     std::string name = r->GetEntryValue(0x3006,0x0026);      // 0x3006,0x0026 = [ROI Name]
 382     int nb = atoi(r->GetEntryValue(0x3006,0x0022).c_str());  // 0x3006,0x0022 = [ROI Number]
 383     // Check if such a number already exist
 384     if (mMapOfROIName.find(nb) != mMapOfROIName.end()) {
 385       std::cerr << "WARNING. A Roi already exist with the number "
 386                 << nb << ". I replace." << std::endl;
 387     }
 388     // Add in map
 389     mMapOfROIName[nb] = name;
 390   }
 391
 392   //----------------------------------
 393   // Read all ROI
 394   // 0x3006,0x0039 = [ ROI Contour Sequence ]
 395   gdcm::SeqEntry * roi_contour_seq=mFile->GetSeqEntry(0x3006,0x0039);
 396   assert(roi_contour_seq); // TODO error message
 397   int n=0;
 398   for (gdcm::SQItem* r=roi_contour_seq->GetFirstSQItem(); r!=0; r=roi_contour_seq->GetNextSQItem()) {
 399     DicomRT_ROI::Pointer roi = DicomRT_ROI::New();
 400     roi->Read(mMapOfROIName, r);
 401     mROIs[roi->GetROINumber()] = roi;
 402     n++;
 403   }
 404
 405 #endif
 406 }
 407 //--------------------------------------------------------------------
 408
 409
 410 //--------------------------------------------------------------------
 411 bool clitk::DicomRT_StructureSet::IsDicomRTStruct(const std::string & filename)
 412 {
 413   // Open DICOM
 414 #if GDCM_MAJOR_VERSION == 2
 415   // Read gdcm file
 416   mReader = new gdcm::Reader;
 417   mReader->SetFileName(filename.c_str());
 418   mReader->Read();
 419   mFile = &(mReader->GetFile());
 420   const gdcm::DataSet & ds = mFile->GetDataSet();
 421
 422   // Check file type
 423   //Verify if the file is a RT-Structure-Set dicom file
 424   gdcm::MediaStorage ms;
 425   ms.SetFromFile(*mFile);
 426   if( ms != gdcm::MediaStorage::RTStructureSetStorage ) return false;
 427
 428   gdcm::Attribute<0x8,0x60> modality;
 429   modality.SetFromDataSet( ds );
 430   if( modality.GetValue() != "RTSTRUCT" ) return false;
 431
 432   return true;
 433
 434   //----------------------------------------------------------------------------------------
 435 #else
 436   mFile = new gdcm::File;
 437   mFile->SetFileName(filename.c_str());
 438   mFile->SetMaxSizeLoadEntry(16384); // Needed ...
 439   mFile->SetLoadMode(gdcm::LD_NOSHADOW); // don't load shadow tags (in order to save memory)
 440   mFile->Load();
 441
 442   // Check file type
 443   //Verify if the file is a RT-Structure-Set dicom file
 444   if (!gdcm::Util::DicomStringEqual(mFile->GetEntryValue(0x0008,0x0016),"1.2.840.10008.5.1.4.1.1.481.3"))
 445     return false;
 446   if (!gdcm::Util::DicomStringEqual(mFile->GetEntryValue(0x0008,0x0060),"RTSTRUCT"))
 447     return false;
 448
 449   return true;
 450
 451 #endif
 452 }
 453 //--------------------------------------------------------------------
 454
 455
 456 //--------------------------------------------------------------------
 457 int clitk::DicomRT_StructureSet::AddBinaryImageAsNewROI(vvImage * im, std::string n)
 458 {
 459   // Search max ROI number
 460   int max = -1;
 461   for(ROIConstIteratorType iter = mROIs.begin(); iter != mROIs.end(); iter++) {
 462     //  for(unsigned int i=0; i<mListOfROI.size(); i++) {
 463     clitk::DicomRT_ROI::Pointer roi = iter->second;
 464     if (roi->GetROINumber() > max)
 465       max = roi->GetROINumber();
 466   }
 467   ++max;
 468
 469   // Compute name
 470   std::ostringstream oss;
 471   oss << vtksys::SystemTools::GetFilenameName(vtksys::SystemTools::GetFilenameWithoutLastExtension(n));
 472   mMapOfROIName[max] = oss.str();
 473
 474   // Set color
 475   std::vector<double> color;
 476   color.push_back(1);
 477   color.push_back(0);
 478   color.push_back(0);
 479
 480   // Create ROI
 481   DicomRT_ROI::Pointer roi = DicomRT_ROI::New();
 482   roi->SetFromBinaryImage(im, max, oss.str(), color, n);
 483   mROIs[max] = roi;
 484   return max;
 485 }
 486 //--------------------------------------------------------------------
 487
 488