common/clitkDicomRT_Contour.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_Contour.h"
  21 #include <vtkCellArray.h>
  22
  23 #if GDCM_MAJOR_VERSION == 2
  24 #include "gdcmAttribute.h"
  25 #include "gdcmItem.h"
  26 #endif
  27
  28 //--------------------------------------------------------------------
  29 clitk::DicomRT_Contour::DicomRT_Contour()
  30 {
  31   mMeshIsUpToDate = false;
  32   mNbOfPoints = 0;
  33   mZ = -1;
  34 }
  35 //--------------------------------------------------------------------
  36
  37
  38 //--------------------------------------------------------------------
  39 clitk::DicomRT_Contour::~DicomRT_Contour()
  40 {
  41
  42 }
  43 //--------------------------------------------------------------------
  44
  45
  46
  47 //--------------------------------------------------------------------
  48 void clitk::DicomRT_Contour::Print(std::ostream & os) const
  49 {
  50   DD("TODO : print Contours");
  51 }
  52 //--------------------------------------------------------------------
  53
  54
  55 //--------------------------------------------------------------------
  56 void clitk::DicomRT_Contour::UpdateDicomItem()
  57 {
  58   DD("DicomRT_Contour::UpdateDicomItem");
  59
  60   gdcm::DataSet & nestedds = mItem->GetNestedDataSet();
  61
  62   //NON CONSIDER CONTOUR ITEM NOT SEQ ITEM ?
  63
  64   // Contour type [Contour Geometric Type]
  65   gdcm::Attribute<0x3006,0x0042> contgeotype;
  66   contgeotype.SetFromDataSet( nestedds );
  67
  68   // Number of points [Number of Contour Points]
  69   gdcm::Attribute<0x3006,0x0046> numcontpoints;
  70   numcontpoints.SetFromDataSet( nestedds );
  71   DD(mNbOfPoints);
  72   mNbOfPoints = numcontpoints.GetValue();
  73   DD(mNbOfPoints);
  74
  75   // Contour dicom values from DataPoints
  76   //ComputeDataPointsFromMesh();
  77   uint nb = mData->GetNumberOfPoints();
  78   std::vector<double> points;
  79   points.resize(mNbOfPoints*3);
  80   for(unsigned int i=0; i<nb; i++) {
  81     double * p = mData->GetPoint(i);
  82     points[i*3] = p[0];
  83     points[i*3+1] = p[1];
  84     points[i*3+2] = p[2];
  85   }
  86
  87   // Get attribute
  88   gdcm::Attribute<0x3006,0x0050> at;
  89   gdcm::Tag tcontourdata(0x3006,0x0050);
  90   gdcm::DataElement contourdata = nestedds.GetDataElement( tcontourdata );
  91   at.SetFromDataElement( contourdata );
  92
  93   // Set attribute
  94   at.SetValues(&points[0], points.size(), false);
  95   DD(at.GetValues()[0]);
  96
  97   DD("replace");
  98   nestedds.Replace(at.GetAsDataElement());
  99   DD("done");
 100
 101   // Change Number of points [Number of Contour Points]
 102   numcontpoints.SetValue(nb);
 103   nestedds.Replace(numcontpoints.GetAsDataElement());
 104
 105   // Test
 106   gdcm::DataElement aa = nestedds.GetDataElement( tcontourdata );
 107   at.SetFromDataElement( aa );
 108   const double* bb = at.GetValues();
 109   DD(bb[0]);
 110
 111 }
 112 //--------------------------------------------------------------------
 113
 114
 115 //--------------------------------------------------------------------
 116 #if GDCM_MAJOR_VERSION == 2
 117 bool clitk::DicomRT_Contour::Read(gdcm::Item * item)
 118 {
 119   mItem = item;
 120
 121   const gdcm::DataSet& nestedds2 = item->GetNestedDataSet();
 122
 123   // Contour type [Contour Geometric Type]
 124   gdcm::Attribute<0x3006,0x0042> contgeotype;
 125   contgeotype.SetFromDataSet( nestedds2 );
 126
 127   if (contgeotype.GetValue() != "CLOSED_PLANAR " && contgeotype.GetValue() != "POINT ") { ///WARNING to the space after the name ...
 128     //std::cerr << "Skip this contour : type=" << mType << std::endl;
 129     return false;
 130   }
 131   if (contgeotype.GetValue() == "POINT ") {
 132     std::cerr << "Warning: POINT type not fully supported. (don't use GetMesh() with this!)"
 133       << std::endl;
 134   }
 135
 136   gdcm::Attribute<0x3006,0x0046> numcontpoints;
 137   numcontpoints.SetFromDataSet( nestedds2 );
 138   // Number of points [Number of Contour Points]
 139   mNbOfPoints = numcontpoints.GetValue();
 140   // DD(mNbOfPoints);
 141
 142   gdcm::Attribute<0x3006,0x0050> at;
 143   gdcm::Tag tcontourdata(0x3006,0x0050);
 144   const gdcm::DataElement & contourdata = nestedds2.GetDataElement( tcontourdata );
 145   at.SetFromDataElement( contourdata );
 146   const double* points = at.GetValues();
 147 <<<<<<< Updated upstream
 148 <<<<<<< Updated upstream
 149 <<<<<<< variant A
 150 >>>>>>> variant B
 151   //  unsigned int npts = at.GetNumberOfValues() / 3;
 152 ======= end
 153 =======
 154   //  unsigned int npts = at.GetNumberOfValues() / 3;
 155 >>>>>>> Stashed changes
 156 =======
 157   //  unsigned int npts = at.GetNumberOfValues() / 3;
 158 >>>>>>> Stashed changes
 159
 160   assert(at.GetNumberOfValues() == static_cast<unsigned int>(mNbOfPoints)*3);
 161
 162   // Organize values
 163   mData = vtkSmartPointer<vtkPoints>::New();
 164   mData->SetDataTypeToDouble();
 165   mData->SetNumberOfPoints(mNbOfPoints);
 166   for(unsigned int i=0; i<mNbOfPoints; i++) {
 167     double p[3];
 168     p[0] = points[i*3];
 169     p[1] = points[i*3+1];
 170     p[2] = points[i*3+2];
 171     mData->SetPoint(i, p);
 172     if (mZ == -1) mZ = p[2];
 173     if (p[2] != mZ) {
 174       DD(i);
 175       DD(p[2]);
 176       DD(mZ);
 177       std::cout << "ERROR ! contour not in the same slice" << std::endl;
 178       assert(p[2] == mZ);
 179     }
 180   }
 181
 182   return true;
 183
 184 }
 185 #else
 186 bool clitk::DicomRT_Contour::Read(gdcm::SQItem * item)
 187 {
 188
 189   // Contour type [Contour Geometric Type]
 190   mType = item->GetEntryValue(0x3006,0x0042);
 191   // DD(mType);
 192   if (mType != "CLOSED_PLANAR " && mType != "POINT ") { ///WARNING to the space after the name ...
 193     //std::cerr << "Skip this contour : type=" << mType << std::endl;
 194     return false;
 195   }
 196   if (mType == "POINT ") {
 197     std::cerr << "Warning: POINT type not fully supported. (don't use GetMesh() with this!)"
 198       << std::endl;
 199   }
 200
 201   // Number of points [Number of Contour Points]
 202   mNbOfPoints = parse_value<int>(item->GetEntryValue(0x3006,0x0046));
 203   // DD(mNbOfPoints);
 204
 205   // Read values [Contour Data]
 206   std::vector<float> points = parse_string<float>(item->GetEntryValue(0x3006,0x0050),'\\');
 207   assert(points.size() == static_cast<unsigned int>(mNbOfPoints)*3);
 208
 209   // Organize values
 210   mData = vtkSmartPointer<vtkPoints>::New();
 211   mData->SetDataTypeToDouble();
 212   mData->SetNumberOfPoints(mNbOfPoints);
 213   for(unsigned int i=0; i<mNbOfPoints; i++) {
 214     double p[3];
 215     p[0] = points[i*3];
 216     p[1] = points[i*3+1];
 217     p[2] = points[i*3+2];
 218     mData->SetPoint(i, p);
 219     if (mZ == -1) mZ = p[2];
 220     if (p[2] != mZ) {
 221       DD(i);
 222       DD(p[2]);
 223       DD(mZ);
 224       std::cout << "ERROR ! contour not in the same slice" << std::endl;
 225       assert(p[2] == mZ);
 226     }
 227   }
 228
 229   return true;
 230 }
 231 #endif
 232 //--------------------------------------------------------------------
 233
 234
 235 //--------------------------------------------------------------------
 236 vtkPolyData * clitk::DicomRT_Contour::GetMesh()
 237 {
 238   if (!mMeshIsUpToDate) {
 239     ComputeMeshFromDataPoints();
 240   }
 241   return mMesh;
 242 }
 243 //--------------------------------------------------------------------
 244
 245
 246 //--------------------------------------------------------------------
 247 void clitk::DicomRT_Contour::SetMesh(vtkPolyData * mesh)
 248 {
 249   mMesh = mesh;
 250 }
 251 //--------------------------------------------------------------------
 252
 253
 254 //--------------------------------------------------------------------
 255 void clitk::DicomRT_Contour::ComputeMeshFromDataPoints()
 256 {
 257 //  DD("ComputeMesh Contour");
 258   mMesh = vtkSmartPointer<vtkPolyData>::New();
 259   mMesh->Allocate(); //for cell structures
 260   mPoints = vtkSmartPointer<vtkPoints>::New();
 261   mMesh->SetPoints(mPoints);
 262   vtkIdType ids[2];
 263   for (unsigned int idx=0 ; idx<mNbOfPoints ; idx++) {
 264     mMesh->GetPoints()->InsertNextPoint(mData->GetPoint(idx)[0],
 265                                         mData->GetPoint(idx)[1],
 266                                         mData->GetPoint(idx)[2]);
 267     ids[0]=idx;
 268     ids[1]=(ids[0]+1) % mNbOfPoints; //0-1,1-2,...,n-1-0
 269     mMesh->GetLines()->InsertNextCell(2,ids);
 270   }
 271   mMeshIsUpToDate = true;
 272 }
 273 //--------------------------------------------------------------------
 274
 275
 276 //--------------------------------------------------------------------
 277 void clitk::DicomRT_Contour::ComputeDataPointsFromMesh()
 278 {
 279   DD("ComputeDataPointsFromMesh");
 280
 281
 282   /*todo
 283
 284   mMesh = vtkSmartPointer<vtkPolyData>::New();
 285   mMesh->Allocate(); //for cell structures
 286   mPoints = vtkSmartPointer<vtkPoints>::New();
 287   mMesh->SetPoints(mPoints);
 288   vtkIdType ids[2];
 289   for (unsigned int idx=0 ; idx<mNbOfPoints ; idx++) {
 290     mMesh->GetPoints()->InsertNextPoint(mData->GetPoint(idx)[0],
 291                                         mData->GetPoint(idx)[1],
 292                                         mData->GetPoint(idx)[2]);
 293     ids[0]=idx;
 294     ids[1]=(ids[0]+1) % mNbOfPoints; //0-1,1-2,...,n-1-0
 295     mMesh->GetLines()->InsertNextCell(2,ids);
 296   }
 297   mMeshIsUpToDate = true;
 298 */
 299 }
 300 //--------------------------------------------------------------------