tools/clitkDicomRTPlan2Gate.cxx

   1 /*=========================================================================
   2   Program:   vv                     http://www.creatis.insa-lyon.fr/rio/vv
   3
   4   Authors belong to:
   5   - University of LYON              http://www.universite-lyon.fr/
   6   - Léon Bérard cancer center       http://www.centreleonberard.fr
   7   - CREATIS CNRS laboratory         http://www.creatis.insa-lyon.fr
   8
   9   This software is distributed WITHOUT ANY WARRANTY; without even
  10   the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  11   PURPOSE.  See the copyright notices for more information.
  12
  13   It is distributed under dual licence
  14
  15   - BSD        See included LICENSE.txt file
  16   - CeCILL-B   http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
  17 ===========================================================================**/
  18 #ifndef CLITKDICOMRTPLAN2GATE_CXX
  19 #define CLITKDICOMRTPLAN2GATE_CXX
  20
  21 #include "clitkDicomRTPlan2Gate_ggo.h"
  22
  23 #include "clitkIO.h"
  24 #include "clitkCommon.h"
  25
  26 #include "gdcmFile.h"
  27 #include "gdcmReader.h"
  28 #include "gdcmAttribute.h"
  29
  30 int main(int argc, char * argv[])
  31 {
  32
  33   // init command line
  34   GGO(clitkDicomRTPlan2Gate, args_info);
  35   CLITK_INIT;
  36   char* input_file = args_info.input_arg;
  37   char* output_file = args_info.output_arg;
  38
  39   gdcm::Reader reader;
  40   reader.SetFileName(input_file);
  41   if(!reader.Read()) {
  42     std::cerr << "Failed to read: " << input_file << std::endl;
  43     return 1;
  44   }
  45
  46   const gdcm::DataSet & ds = reader.GetFile().GetDataSet();
  47
  48   //Check file type
  49   gdcm::MediaStorage ms;
  50   ms.SetFromFile(reader.GetFile());
  51   if(ms != gdcm::MediaStorage::RTIonPlanStorage) {
  52     std::cerr << "The file <" << input_file << "> is not a RT Ion Plan." << std::endl;
  53     exit(0);
  54   }
  55
  56   std::ofstream output;
  57   output.open(output_file);
  58   output << "#TREATMENT-PLAN-DESCRIPTION" << "\n";
  59
  60   gdcm::Attribute<0x300a,0x2> label;
  61   label.SetFromDataSet(ds);
  62   gdcm::Attribute<0x300a,0x3> planName;
  63   planName.SetFromDataSet(ds);
  64
  65   output << "#PlanName\n"<< planName.GetValue() << "\n";
  66
  67   gdcm::Tag tFractionGroupSQ(0x300a,0x70);
  68   gdcm::Attribute<0x300a,0x71> fractionID;
  69   gdcm::Attribute<0x300a,0x78> numberOfFractions;
  70   gdcm::Attribute<0x300a,0x80> numberOfBeams;
  71
  72   gdcm::Tag tReferencedBeamSQ(0x300c,0x4);
  73   gdcm::Attribute<0x300a,0x86> meterSet;
  74   gdcm::Attribute<0x300c,0x6> beamID;
  75
  76   const gdcm::DataElement &FractionGroupSQ = ds.GetDataElement( tFractionGroupSQ );
  77   gdcm::SmartPointer<gdcm::SequenceOfItems> sqf = FractionGroupSQ.GetValueAsSQ();
  78   gdcm::Item &ReferencedBeamItem = sqf->GetItem(1);
  79   gdcm::DataSet &ReferencedBeamData = ReferencedBeamItem.GetNestedDataSet();
  80   numberOfFractions.SetFromDataSet(ReferencedBeamData);
  81
  82   output << "#NumberOfFractions\n"<< numberOfFractions.GetValue() << "\n";
  83
  84   for(unsigned int fidx = 1; fidx <= sqf->GetNumberOfItems(); ++fidx) {
  85     gdcm::Item &ReferencedBeamItem = sqf->GetItem(fidx);
  86     gdcm::DataSet &ReferencedBeamData = ReferencedBeamItem.GetNestedDataSet();
  87
  88     fractionID.SetFromDataSet(ReferencedBeamData);
  89     numberOfFractions.SetFromDataSet(ReferencedBeamData);
  90     numberOfBeams.SetFromDataSet(ReferencedBeamData);
  91     output << "##FractionID\n"<<  fractionID.GetValue() << "\n";
  92     output << "##NumberOfFields\n"<<  numberOfBeams.GetValue() << "\n";
  93
  94     const gdcm::DataElement &ReferencedBeamSQ = ReferencedBeamData.GetDataElement(tReferencedBeamSQ);
  95     gdcm::SmartPointer<gdcm::SequenceOfItems> sqb = ReferencedBeamSQ.GetValueAsSQ();
  96     for(unsigned int bidx = 1; bidx <= sqb->GetNumberOfItems(); ++bidx) {
  97       gdcm::Item &ReferencedBeamGroupItem = sqb->GetItem(bidx);
  98       gdcm::DataSet &ReferencedBeamGroupData = ReferencedBeamGroupItem.GetNestedDataSet();
  99
 100       meterSet.SetFromDataSet(ReferencedBeamGroupData);
 101       beamID.SetFromDataSet(ReferencedBeamGroupData);
 102       output << "###FieldsID\n"<<  beamID.GetValue() << "\n";
 103     }
 104   }
 105
 106   gdcm::Tag tPatientSetupSQ(0x300a,0x180);
 107   gdcm::Attribute<0x0018,0x5100> position;
 108   gdcm::Attribute<0x300a,0x1b0> setupTechnique;
 109   gdcm::Attribute<0x300a,0x182> setupID;
 110
 111   const gdcm::DataElement &PatientSetupSQ = ds.GetDataElement( tPatientSetupSQ );
 112   gdcm::SmartPointer<gdcm::SequenceOfItems> sqs = PatientSetupSQ.GetValueAsSQ();
 113   for(unsigned int sidx = 1; sidx <= sqs->GetNumberOfItems(); ++sidx) {
 114     gdcm::Item &PatientSetupItem = sqs->GetItem(sidx);
 115     gdcm::DataSet &PatientSetupData = PatientSetupItem.GetNestedDataSet();
 116
 117     position.SetFromDataSet(PatientSetupData);
 118     setupTechnique.SetFromDataSet(PatientSetupData);
 119     setupID.SetFromDataSet(PatientSetupData);
 120   }
 121
 122   gdcm::Tag tIonBeamSQ(0x300a,0x3a2);
 123   gdcm::Attribute<0x300a,0xb2> machine;
 124   gdcm::Attribute<0x300a,0xc6> radiationType;
 125   gdcm::Attribute<0x300a,0xce> treatmentType;
 126   gdcm::Attribute<0x300a,0xc0> ionBeamID;
 127   gdcm::Attribute<0x300a,0x10e> finalCumulativeMeterSetWeight;
 128   gdcm::Attribute<0x300a,0x110> numberOfControlPoints;
 129   gdcm::Attribute<0x300c,0x6a> patientSetupID;
 130   gdcm::Attribute<0x300a,0xd0> numberOfWedges;
 131   gdcm::Attribute<0x300a,0xe0> numberOfCompensators;
 132   gdcm::Attribute<0x300a,0xed> numberOfBoli;
 133   gdcm::Attribute<0x300a,0xf0> numberOfBlocks;
 134   gdcm::Attribute<0x300a,0x312> numberOfRangeShifters;
 135   gdcm::Attribute<0x300a,0x330> numberOfLateralSpreadingDevices;
 136   gdcm::Attribute<0x300a,0x340> numberOfRangeModulators;
 137
 138   //gdcm::Tag tSnoutSQ(0x300a,0x30c);
 139   //gdcm::Attribute<0x300a,0x30f> snoutID;
 140
 141   gdcm::Tag tControlPointSQ(0x300a,0x3a8);
 142   gdcm::Attribute<0x300a,0x15> energyUnit;
 143   gdcm::Attribute<0x300a,0x114> energyValue;
 144   gdcm::Attribute<0x300a,0x112> controlPointIndex;
 145   gdcm::Attribute<0x300a,0x134> cumulativeMetersetWeight;
 146   gdcm::Attribute<0x300a,0x392> numberOfScannedSpots;
 147   gdcm::Attribute<0x300a,0x390> spotTunnedID;
 148
 149   gdcm::Tag tSpotPositionsData(0x300a,0x394);
 150   gdcm::Attribute<0x300a,0x394> spotPositions;
 151   gdcm::Tag tSpotWeightsData(0x300a,0x396);
 152   gdcm::Attribute<0x300a,0x396> spotWeights;
 153
 154   gdcm::Attribute<0x300a,0x11e> gantryAngle;
 155   gdcm::Attribute<0x300a,0x122> patientSupportAngle;
 156   gdcm::Attribute<0x300a,0x12c> isocenterPosition;
 157
 158   long totalCumMeterSet = 0;
 159
 160   const gdcm::DataElement &IonBeamSQ = ds.GetDataElement( tIonBeamSQ );
 161   gdcm::SmartPointer<gdcm::SequenceOfItems> sqi = IonBeamSQ.GetValueAsSQ();
 162
 163   for(unsigned int iidx = 1; iidx <= sqi->GetNumberOfItems(); ++iidx) {
 164     gdcm::Item &IonBeamItem = sqi->GetItem(iidx);
 165     gdcm::DataSet &IonBeamData = IonBeamItem.GetNestedDataSet();
 166     finalCumulativeMeterSetWeight.SetFromDataSet(IonBeamData);
 167     totalCumMeterSet += finalCumulativeMeterSetWeight.GetValue();
 168   }
 169
 170   output << "#TotalMetersetWeightOfAllFields\n" <<  totalCumMeterSet << "\n";
 171
 172   for(unsigned int iidx = 1; iidx <= sqi->GetNumberOfItems(); ++iidx) {
 173     gdcm::Item &IonBeamItem = sqi->GetItem(iidx);
 174     gdcm::DataSet &IonBeamData = IonBeamItem.GetNestedDataSet();
 175
 176     machine.SetFromDataSet(IonBeamData);
 177     radiationType.SetFromDataSet(IonBeamData);
 178     treatmentType.SetFromDataSet(IonBeamData);
 179     ionBeamID.SetFromDataSet(IonBeamData);
 180     finalCumulativeMeterSetWeight.SetFromDataSet(IonBeamData);
 181     numberOfControlPoints.SetFromDataSet(IonBeamData);
 182     patientSetupID.SetFromDataSet(IonBeamData);
 183     numberOfWedges.SetFromDataSet(IonBeamData);
 184     numberOfCompensators.SetFromDataSet(IonBeamData);
 185     numberOfBoli.SetFromDataSet(IonBeamData);
 186     numberOfBlocks.SetFromDataSet(IonBeamData);
 187     numberOfRangeShifters.SetFromDataSet(IonBeamData);
 188     numberOfLateralSpreadingDevices.SetFromDataSet(IonBeamData);
 189     numberOfRangeModulators.SetFromDataSet(IonBeamData);
 190
 191     /*
 192     const gdcm::DataElement &SnoutSQ = IonBeamData.GetDataElement(tSnoutSQ);
 193     gdcm::SmartPointer<gdcm::SequenceOfItems> sqsn = SnoutSQ.GetValueAsSQ();
 194     for(unsigned int snidx = 1; snidx <= sqsn->GetNumberOfItems(); ++snidx) {
 195       gdcm::Item &SnoutItem = sqsn->GetItem(snidx);
 196       gdcm::DataSet &SnoutData = SnoutItem.GetNestedDataSet();
 197
 198       snoutID.SetFromDataSet(SnoutData);
 199     }
 200     */
 201
 202     const gdcm::DataElement &ControlPointSQ = IonBeamData.GetDataElement(tControlPointSQ);
 203     gdcm::SmartPointer<gdcm::SequenceOfItems> sqcp = ControlPointSQ.GetValueAsSQ();
 204     gdcm::Item &ControlPointItem = sqcp->GetItem(1);
 205     gdcm::DataSet &ControlPointData = ControlPointItem.GetNestedDataSet();
 206
 207     gantryAngle.SetFromDataSet(ControlPointData);
 208     patientSupportAngle.SetFromDataSet(ControlPointData);
 209     isocenterPosition.SetFromDataSet(ControlPointData);
 210
 211     output << "\n#FIELD-DESCRIPTION" << "\n";
 212     output << "###FieldID\n"<< ionBeamID.GetValue() << "\n";
 213     output << "###FinalCumulativeMeterSetWeight\n"<<  finalCumulativeMeterSetWeight.GetValue() << "\n";
 214     output << "###GantryAngle\n" << gantryAngle.GetValue() << "\n";
 215     output << "###PatientSupportAngle\n" << patientSupportAngle.GetValue() << "\n";
 216     output << "###IsocenterPosition\n" << isocenterPosition.GetValue(0);
 217     output << " " << isocenterPosition.GetValue(1);
 218     output << " " << isocenterPosition.GetValue(2) << "\n";
 219     output << "###NumberOfControlPoints\n" << numberOfControlPoints.GetValue() << "\n\n";
 220     output << "#SPOTS-DESCRIPTION" << "\n";
 221
 222     for(unsigned int cpidx = 1; cpidx <= sqcp->GetNumberOfItems(); ++cpidx) {
 223       gdcm::Item &ControlPointItem = sqcp->GetItem(cpidx);
 224       gdcm::DataSet &ControlPointData = ControlPointItem.GetNestedDataSet();
 225
 226       if(cpidx == 1) {
 227         gantryAngle.SetFromDataSet(ControlPointData);
 228         patientSupportAngle.SetFromDataSet(ControlPointData);
 229         isocenterPosition.SetFromDataSet(ControlPointData);
 230       }
 231
 232       if(cpidx < sqcp->GetNumberOfItems()) {
 233         energyUnit.SetFromDataSet(ControlPointData);
 234         energyValue.SetFromDataSet(ControlPointData);
 235         controlPointIndex.SetFromDataSet(ControlPointData);
 236         cumulativeMetersetWeight.SetFromDataSet(ControlPointData);
 237         numberOfScannedSpots.SetFromDataSet(ControlPointData);
 238         spotTunnedID.SetFromDataSet(ControlPointData);
 239
 240         const gdcm::DataElement &spotPositionsData = ControlPointData.GetDataElement(tSpotPositionsData);
 241         spotPositions.SetFromDataElement(spotPositionsData);
 242         const gdcm::DataElement &spotWeightsData = ControlPointData.GetDataElement(tSpotWeightsData);
 243         spotWeights.SetFromDataElement(spotWeightsData);
 244       }
 245
 246       if(cpidx == sqcp->GetNumberOfItems()) {
 247         energyUnit.SetValue( "0" );
 248         energyValue.SetValue( 0 );
 249         controlPointIndex.SetFromDataSet(ControlPointData);
 250         cumulativeMetersetWeight.SetFromDataSet(ControlPointData);
 251         numberOfScannedSpots.SetValue( 0 );
 252         spotTunnedID.SetValue( "0" );
 253       }
 254
 255       output << "####ControlPointIndex\n" << controlPointIndex.GetValue() << "\n";
 256       output << "####SpotTunnedID\n" << spotTunnedID.GetValue() << "\n";
 257       output << "####CumulativeMetersetWeight\n" << cumulativeMetersetWeight.GetValue() << "\n";
 258       output << "####Energy (MeV)\n"<< energyValue.GetValue() << "\n";
 259       output << "####NbOfScannedSpots\n" << numberOfScannedSpots.GetValue() << "\n";
 260       output << "####X Y Weight\n";
 261       for( unsigned int sp=0; sp < numberOfScannedSpots.GetValue(); ++sp)
 262         output << spotPositions.GetValue(2*sp) << " " << spotPositions.GetValue(2*sp+1) << " " << spotWeights.GetValue(sp) << "\n";
 263     }
 264   }
 265
 266   return 0;
 267 }
 268
 269 #endif