src/gdcmHeaderHelper.cxx

   1 // gdcmHeaderHelper.cxx
   2 //-----------------------------------------------------------------------------
   3 #include "gdcmHeaderHelper.h"
   4 #include "gdcmDirList.h"
   5
   6 #include "gdcmDebug.h"
   7 #include <math.h>
   8 #include <algorithm>
   9 #include <vector>
  10
  11 //-----------------------------------------------------------------------------
  12 /**
  13  * \brief   constructor
  14  */
  15 gdcmHeaderHelper::gdcmHeaderHelper() : gdcmHeader( ) {
  16
  17 }
  18
  19 /**
  20  * \brief   constructor
  21  * @param   InFilename Name of the file to deal with
  22  * @param   exception_on_error
  23  * @param   enable_sequences = true to allow the header
  24  *          to be parsed *inside* the SeQuences,
  25  *          when they have an actual length
  26  * @param   ignore_shadow = true if user wants to skip shadow groups
  27  *           during parsing, to save memory space
  28  */
  29 gdcmHeaderHelper::gdcmHeaderHelper(const char *InFilename,
  30                                    bool exception_on_error,
  31                                    bool enable_sequences,
  32                                    bool ignore_shadow)
  33                  : gdcmHeader( InFilename,
  34                                exception_on_error,
  35                                enable_sequences,
  36                                ignore_shadow)
  37 {
  38 }
  39
  40 //-----------------------------------------------------------------------------
  41 // Print
  42
  43 //-----------------------------------------------------------------------------
  44 // Public
  45 /**
  46  * \brief   Returns the size (in bytes) of a single pixel of data.
  47  * @return  The size in bytes of a single pixel of data.
  48  *
  49  */
  50 int gdcmHeaderHelper::GetPixelSize() {
  51
  52      // 0028 0100 US IMG Bits Allocated
  53      // (in order no to be messed up by old RGB images)
  54    if (gdcmHeader::GetEntryByNumber(0x0028,0x0100) == "24")
  55       return 3;
  56
  57    std::string PixelType = GetPixelType();
  58    if (PixelType == "8U"  || PixelType == "8S")
  59       return 1;
  60    if (PixelType == "16U" || PixelType == "16S")
  61       return 2;
  62    if (PixelType == "32U" || PixelType == "32S")
  63       return 4;
  64    if (PixelType == "FD") // to help unfortunate users to manage DOUBLE
  65       return 8;
  66    dbg.Verbose(0, "gdcmHeader::GetPixelSize: Unknown pixel type");
  67    return 0;
  68 }
  69
  70 /**
  71  * \ingroup gdcmHeaderHelper
  72  * \brief   Build the Pixel Type of the image.
  73  *          Possible values are:
  74  *          - 8U  unsigned  8 bit,
  75  *          - 8S    signed  8 bit,
  76  *          - 16U unsigned 16 bit,
  77  *          - 16S   signed 16 bit,
  78  *          - 32U unsigned 32 bit,
  79  *          - 32S   signed 32 bit,
  80  *          - FD    Double,
  81  * \warning 12 bit images appear as 16 bit.
  82  *          24 bit images appear as 8 bit
  83  *          64 bit means 'DOUBLE' images
  84  *                 (no DOUBLE images in kosher DICOM,
  85  *                  but so usefull for people that miss them ;-)
  86  * @return
  87  */
  88 std::string gdcmHeaderHelper::GetPixelType() {
  89    std::string BitsAlloc;
  90    BitsAlloc = GetEntryByNumber(0x0028, 0x0100);
  91    if (BitsAlloc == GDCM_UNFOUND) { // Bits Allocated
  92       dbg.Verbose(0, "gdcmHeader::GetPixelType: unfound Bits Allocated");
  93       BitsAlloc = std::string("16");
  94    }
  95    if (BitsAlloc == "12")           // It will be unpacked
  96       BitsAlloc = std::string("16");
  97    else if (BitsAlloc == "24")      // (in order no to be messed up
  98       BitsAlloc = std::string("8"); // by old RGB images)
  99
 100    std::string Signed;
 101    Signed = GetEntryByNumber(0x0028, 0x0103);
 102    if (Signed == GDCM_UNFOUND) { // "Pixel Representation"
 103       dbg.Verbose(0, "gdcmHeader::GetPixelType: unfound Pixel Representation");
 104       BitsAlloc = std::string("0");
 105    }
 106    if (BitsAlloc == "64") // to help users that want to deal with DOUBLE
 107       return("FD");
 108
 109    if (Signed == "0")
 110       Signed = std::string("U");
 111    else
 112       Signed = std::string("S");
 113
 114    return( BitsAlloc + Signed);
 115 }
 116
 117 /**
 118   * \ingroup gdcmHeaderHelper
 119   * \brief gets the info from 0028,0030 : Pixel Spacing
 120   *             else 1.0
 121   * @return X dimension of a pixel
 122   */
 123 float gdcmHeaderHelper::GetXSpacing() {
 124     float xspacing, yspacing;
 125     std::string StrSpacing = GetEntryByNumber(0x0028,0x0030);
 126
 127    if (StrSpacing == GDCM_UNFOUND) {
 128       dbg.Verbose(0, "gdcmHeader::GetXSpacing: unfound Pixel Spacing (0028,0030)");
 129       return 1.;
 130     }
 131   int nbValues;
 132   if( (nbValues = sscanf( StrSpacing.c_str(), "%f\\%f", &yspacing, &xspacing)) != 2) {
 133     if (nbValues==1)  // if single value is found, xspacing is defaulted to yspacing
 134        return yspacing;
 135   }
 136   if (xspacing == 0.) {
 137     dbg.Verbose(0, "gdcmHeader::GetYSpacing: gdcmData/CT-MONO2-8-abdo.dcm problem");
 138     // seems to be a bug in the header ...
 139     sscanf( StrSpacing.c_str(), "%f\\0\\%f", &yspacing, &xspacing);
 140   }
 141   return xspacing;
 142 }
 143
 144 /**
 145   * \ingroup gdcmHeaderHelper
 146   * \brief gets the info from 0028,0030 : Pixel Spacing
 147   *             else 1.0
 148   * @return Y dimension of a pixel
 149   */
 150 float gdcmHeaderHelper::GetYSpacing() {
 151    float yspacing;
 152    std::string StrSpacing = GetEntryByNumber(0x0028,0x0030);
 153
 154    if (StrSpacing == GDCM_UNFOUND) {
 155       dbg.Verbose(0, "gdcmHeader::GetYSpacing: unfound Pixel Spacing (0028,0030)");
 156       return 1.;
 157     }
 158   sscanf( StrSpacing.c_str(), "%f", &yspacing);
 159   return yspacing;
 160 }
 161
 162 /**
 163   *\ingroup gdcmHeaderHelper
 164   *\brief gets the info from 0018,0088 : Space Between Slices
 165   *                else from 0018,0050 : Slice Thickness
 166    *                else 1.0
 167   * @return Z dimension of a voxel-to be
 168   */
 169 float gdcmHeaderHelper::GetZSpacing() {
 170    // Spacing Between Slices : distance entre le milieu de chaque coupe
 171    // Les coupes peuvent etre :
 172    //   jointives     (Spacing between Slices = Slice Thickness)
 173    //   chevauchantes (Spacing between Slices < Slice Thickness)
 174    //   disjointes    (Spacing between Slices > Slice Thickness)
 175    // Slice Thickness : epaisseur de tissus sur laquelle est acquis le signal
 176    //   ca interesse le physicien de l'IRM, pas le visualisateur de volumes ...
 177    //   Si le Spacing Between Slices est absent,
 178    //   on suppose que les coupes sont jointives
 179
 180    std::string StrSpacingBSlices = GetEntryByNumber(0x0018,0x0088);
 181
 182    if (StrSpacingBSlices == GDCM_UNFOUND) {
 183       dbg.Verbose(0, "gdcmHeader::GetZSpacing: unfound StrSpacingBSlices");
 184       std::string StrSliceThickness = GetEntryByNumber(0x0018,0x0050);
 185       if (StrSliceThickness == GDCM_UNFOUND)
 186          return 1.;
 187       else
 188          // if no 'Spacing Between Slices' is found,
 189          // we assume slices join together
 190          // (no overlapping, no interslice gap)
 191          // if they don't, we're fucked up
 192          return atof(StrSliceThickness.c_str());
 193    } else {
 194       return atof(StrSpacingBSlices.c_str());
 195    }
 196 }
 197
 198 /**
 199   *\ingroup gdcmHeaderHelper
 200   *\brief gets the info from 0028,1052 : Rescale Intercept
 201   * @return Rescale Intercept
 202  */
 203 float gdcmHeaderHelper::GetRescaleIntercept() {
 204   float resInter = 0.;
 205   std::string StrRescInter = GetEntryByNumber(0x0028,0x1052); //0028 1052 DS IMG Rescale Intercept
 206   if (StrRescInter != GDCM_UNFOUND) {
 207       if( sscanf( StrRescInter.c_str(), "%f", &resInter) != 1) {
 208          dbg.Verbose(0, "gdcmHeader::GetRescaleIntercept: Rescale Slope is empty");
 209            // bug in the element 0x0028,0x1052
 210       }
 211    }
 212   return resInter;
 213 }
 214
 215 /**
 216   *\ingroup gdcmHeaderHelper
 217   *\brief gets the info from 0028,1053 : Rescale Slope
 218   * @return Rescale Slope
 219  */
 220  float gdcmHeaderHelper::GetRescaleSlope() {
 221   float resSlope = 1.;
 222   std::string StrRescSlope = GetEntryByNumber(0x0028,0x1053); //0028 1053 DS IMG Rescale Slope
 223   if (StrRescSlope != GDCM_UNFOUND) {
 224       if( sscanf( StrRescSlope.c_str(), "%f", &resSlope) != 1) {
 225          dbg.Verbose(0, "gdcmHeader::GetRescaleSlope: Rescale Slope is empty");
 226            // bug in the element 0x0028,0x1053
 227       }
 228    }
 229    return resSlope;
 230 }
 231
 232 /**
 233   * \ingroup gdcmHeaderHelper
 234   * \brief This function is intended to user who doesn't want
 235   *   to have to manage a LUT and expects to get an RBG Pixel image
 236   *   (or a monochrome one ...)
 237   * \warning to be used with GetImagePixels()
 238   * @return 1 if Gray level, 3 if Color (RGB, YBR or PALETTE COLOR)
 239   */
 240 int gdcmHeaderHelper::GetNumberOfScalarComponents() {
 241    if (GetSamplesPerPixel() ==3)
 242       return 3;
 243
 244      // 0028 0100 US IMG Bits Allocated
 245      // (in order no to be messed up by old RGB images)
 246    if (gdcmHeader::GetEntryByNumber(0x0028,0x0100) == "24")
 247       return 3;
 248
 249    std::string PhotometricInterpretation =
 250                   gdcmHeader::GetEntryByNumber(0x0028,0x0004);
 251
 252    if ( ( PhotometricInterpretation == "PALETTE COLOR ") ) {
 253       if (HasLUT())   // PALETTE COLOR is NOT enough
 254          return 3;
 255       else
 256          return 1;
 257    }
 258
 259       //beware of trailing space at end of string
 260    if (PhotometricInterpretation.find(GDCM_UNFOUND) <
 261                            PhotometricInterpretation.length() ||
 262        PhotometricInterpretation.find("MONOCHROME1") <
 263                            PhotometricInterpretation.length() ||
 264        PhotometricInterpretation.find("MONOCHROME2") <
 265                            PhotometricInterpretation.length() )
 266        return 1;
 267     else
 268     // we assume that *all* kinds of YBR are dealt with
 269       return 3;
 270 }
 271
 272 /**
 273   * \ingroup gdcmHeaderHelper
 274   * \brief This function is intended to user that DOESN'T want
 275   *  to get RGB pixels image when it's stored as a PALETTE COLOR image
 276   *   - the (vtk) user is supposed to know how deal with LUTs -
 277   * \warning to be used with GetImagePixelsRaw()
 278   * @return 1 if Gray level, 3 if Color (RGB or YBR - NOT 'PALETTE COLOR' -)
 279   */
 280 int gdcmHeaderHelper::GetNumberOfScalarComponentsRaw() {
 281
 282      // 0028 0100 US IMG Bits Allocated
 283      // (in order no to be messed up by old RGB images)
 284    if (gdcmHeader::GetEntryByNumber(0x0028,0x0100) == "24")
 285       return 3;
 286
 287     // we assume that *all* kinds of YBR are dealt with
 288       return GetSamplesPerPixel();
 289 }
 290
 291 /**
 292   *\ingroup gdcmHeaderHelper
 293   *\brief gets the info from 0020,000d : Study Instance UID
 294   *\todo ? : return the ACR-NEMA element value if DICOM one is not found
 295   * @return Study Instance UID
 296  */
 297  std::string gdcmHeaderHelper::GetStudyUID(){
 298   return GetEntryByNumber(0x0020,0x000d); //0020 000d UI REL Study Instance UID
 299 }
 300
 301 /**
 302   *\ingroup gdcmHeaderHelper
 303   *\brief gets the info from 0020,000e : Series Instance UID
 304   *\todo ? : return the ACR-NEMA element value if DICOM one is not found
 305   * @return Series Instance UID
 306  */
 307  std::string gdcmHeaderHelper::GetSeriesUID(){
 308   return GetEntryByNumber(0x0020,0x000e); //0020 000e UI REL Series Instance UID
 309 }
 310
 311 /**
 312   *\ingroup gdcmHeaderHelper
 313   *\brief gets the info from 0008,0016 : SOP Class UID
 314   *\todo ? : return the ACR-NEMA element value if DICOM one is not found
 315   * @return SOP Class UID
 316  */
 317  std::string gdcmHeaderHelper::GetClassUID(){
 318   return GetEntryByNumber(0x0008,0x0016); //0008 0016 UI ID SOP Class UID
 319 }
 320
 321 /**
 322   *\brief gets the info from 0008,0018 : SOP Instance UID
 323   *\todo ? : return the ACR-NEMA element value if DICOM one is not found
 324   * @return SOP Instance UID
 325  */
 326  std::string gdcmHeaderHelper::GetInstanceUID(){
 327   return GetEntryByNumber(0x0008,0x0018); //0008 0018 UI ID SOP Instance UID
 328 }
 329 //
 330 // --------------  Remember ! ----------------------------------
 331 //
 332 // Image Position Patient                              (0020,0032):
 333 // If not found (ACR_NEMA) we try Image Position       (0020,0030)
 334 // If not found (ACR-NEMA), we consider Slice Location (0020,1041)
 335 //                                   or Location       (0020,0050)
 336 // as the Z coordinate,
 337 // 0. for all the coordinates if nothing is found
 338
 339 // \todo find a way to inform the caller nothing was found
 340 // \todo How to tell the caller a wrong number of values was found?
 341 //
 342 // ---------------------------------------------------------------
 343 //
 344
 345 /**
 346   * \brief gets the info from 0020,0032 : Image Position Patient
 347   *                 else from 0020,0030 : Image Position (RET)
 348   *                 else 0.
 349   * @return up-left image corner X position
 350   */
 351
 352 float gdcmHeaderHelper::GetXOrigin() {
 353     float xImPos, yImPos, zImPos;
 354     std::string StrImPos = GetEntryByNumber(0x0020,0x0032);
 355
 356     if (StrImPos == GDCM_UNFOUND) {
 357        dbg.Verbose(0, "gdcmHeader::GetXImagePosition: unfound Image Position Patient (0020,0032)");
 358        StrImPos = GetEntryByNumber(0x0020,0x0030); // For ACR-NEMA images
 359        if (StrImPos == GDCM_UNFOUND) {
 360           dbg.Verbose(0, "gdcmHeader::GetXImagePosition: unfound Image Position (RET) (0020,0030)");
 361           /// \todo How to tell the caller nothing was found ?
 362          return 0.;
 363        }
 364      }
 365    if( sscanf( StrImPos.c_str(), "%f\\%f\\%f", &xImPos, &yImPos, &zImPos) != 3)
 366      return 0.;
 367    return xImPos;
 368 }
 369
 370 /**
 371   * \brief gets the info from 0020,0032 : Image Position Patient
 372   *                 else from 0020,0030 : Image Position (RET)
 373   *                 else 0.
 374   * @return up-left image corner Y position
 375   */
 376 float gdcmHeaderHelper::GetYOrigin() {
 377     float xImPos, yImPos, zImPos;
 378     std::string StrImPos = GetEntryByNumber(0x0020,0x0032);
 379
 380     if (StrImPos == GDCM_UNFOUND) {
 381        dbg.Verbose(0, "gdcmHeader::GetYImagePosition: unfound Image Position Patient (0020,0032)");
 382        StrImPos = GetEntryByNumber(0x0020,0x0030); // For ACR-NEMA images
 383        if (StrImPos == GDCM_UNFOUND) {
 384           dbg.Verbose(0, "gdcmHeader::GetYImagePosition: unfound Image Position (RET) (0020,0030)");
 385           /// \todo How to tell the caller nothing was found ?
 386            return 0.;
 387        }
 388      }
 389    if( sscanf( StrImPos.c_str(), "%f\\%f\\%f", &xImPos, &yImPos, &zImPos) != 3)
 390      return 0.;
 391    return yImPos;
 392 }
 393
 394 /**
 395   * \brief gets the info from 0020,0032 : Image Position Patient
 396   * \               else from 0020,0030 : Image Position (RET)
 397   * \               else from 0020,1041 : Slice Location
 398   * \               else from 0020,0050 : Location
 399   * \               else 0.
 400   * @return up-left image corner Z position
 401   */
 402 float gdcmHeaderHelper::GetZOrigin() {
 403    float xImPos, yImPos, zImPos;
 404    std::string StrImPos = GetEntryByNumber(0x0020,0x0032);
 405    if (StrImPos != GDCM_UNFOUND) {
 406       if( sscanf( StrImPos.c_str(), "%f\\%f\\%f", &xImPos, &yImPos, &zImPos) != 3) {
 407          dbg.Verbose(0, "gdcmHeader::GetZImagePosition: wrong Image Position Patient (0020,0032)");
 408          return 0.;  // bug in the element 0x0020,0x0032
 409       } else {
 410          return zImPos;
 411       }
 412    }
 413    StrImPos = GetEntryByNumber(0x0020,0x0030); // For ACR-NEMA images
 414    if (StrImPos != GDCM_UNFOUND) {
 415       if( sscanf( StrImPos.c_str(), "%f\\%f\\%f", &xImPos, &yImPos, &zImPos) != 3) {
 416          dbg.Verbose(0, "gdcmHeader::GetZImagePosition: wrong Image Position (RET) (0020,0030)");
 417          return 0.;  // bug in the element 0x0020,0x0032
 418       } else {
 419          return zImPos;
 420       }
 421    }
 422    std::string StrSliceLocation = GetEntryByNumber(0x0020,0x1041);// for *very* old ACR-NEMA images
 423    if (StrSliceLocation != GDCM_UNFOUND) {
 424       if( sscanf( StrSliceLocation.c_str(), "%f", &zImPos) !=1) {
 425          dbg.Verbose(0, "gdcmHeader::GetZImagePosition: wrong Slice Location (0020,1041)");
 426          return 0.;  // bug in the element 0x0020,0x1041
 427       } else {
 428          return zImPos;
 429       }
 430    }
 431    dbg.Verbose(0, "gdcmHeader::GetZImagePosition: unfound Slice Location (0020,1041)");
 432    std::string StrLocation = GetEntryByNumber(0x0020,0x0050);
 433    if (StrLocation != GDCM_UNFOUND) {
 434       if( sscanf( StrLocation.c_str(), "%f", &zImPos) !=1) {
 435          dbg.Verbose(0, "gdcmHeader::GetZImagePosition: wrong Location (0020,0050)");
 436          return 0.;  // bug in the element 0x0020,0x0050
 437       } else {
 438          return zImPos;
 439       }
 440    }
 441    dbg.Verbose(0, "gdcmHeader::GetYImagePosition: unfound Location (0020,0050)");
 442    return 0.; // Hopeless
 443 }
 444
 445 /**
 446   * \brief gets the info from 0020,0013 : Image Number
 447   * \               else 0.
 448   * @return image number
 449   */
 450 int gdcmHeaderHelper::GetImageNumber() {
 451   // The function i atoi() takes the address of an area of memory as
 452   // parameter and converts the string stored at that location to an integer
 453   // using the external decimal to internal binary conversion rules. This may
 454   // be preferable to sscanf() since atoi() is a much smaller, simpler and
 455   // faster function. sscanf() can do all possible conversions whereas
 456   // atoi() can only do single decimal integer conversions.
 457   std::string StrImNumber = GetEntryByNumber(0x0020,0x0013); //0020 0013 IS REL Image Number
 458   if (StrImNumber != GDCM_UNFOUND) {
 459     return atoi( StrImNumber.c_str() );
 460   }
 461   return 0;   //Hopeless
 462 }
 463
 464 /**
 465   * \brief gets the info from 0008,0060 : Modality
 466   * @return Modality Type
 467   */
 468 ModalityType gdcmHeaderHelper::GetModality(void) {
 469   // 0008 0060 CS ID Modality
 470   std::string StrModality = GetEntryByNumber(0x0008,0x0060);
 471   if (StrModality != GDCM_UNFOUND) {
 472          if ( StrModality.find("AU") < StrModality.length()) return AU;
 473     else if ( StrModality.find("AS") < StrModality.length()) return AS;
 474     else if ( StrModality.find("BI") < StrModality.length()) return BI;
 475     else if ( StrModality.find("CF") < StrModality.length()) return CF;
 476     else if ( StrModality.find("CP") < StrModality.length()) return CP;
 477     else if ( StrModality.find("CR") < StrModality.length()) return CR;
 478     else if ( StrModality.find("CT") < StrModality.length()) return CT;
 479     else if ( StrModality.find("CS") < StrModality.length()) return CS;
 480     else if ( StrModality.find("DD") < StrModality.length()) return DD;
 481     else if ( StrModality.find("DF") < StrModality.length()) return DF;
 482     else if ( StrModality.find("DG") < StrModality.length()) return DG;
 483     else if ( StrModality.find("DM") < StrModality.length()) return DM;
 484     else if ( StrModality.find("DS") < StrModality.length()) return DS;
 485     else if ( StrModality.find("DX") < StrModality.length()) return DX;
 486     else if ( StrModality.find("ECG") < StrModality.length()) return ECG;
 487     else if ( StrModality.find("EPS") < StrModality.length()) return EPS;
 488     else if ( StrModality.find("FA") < StrModality.length()) return FA;
 489     else if ( StrModality.find("FS") < StrModality.length()) return FS;
 490     else if ( StrModality.find("HC") < StrModality.length()) return HC;
 491     else if ( StrModality.find("HD") < StrModality.length()) return HD;
 492     else if ( StrModality.find("LP") < StrModality.length()) return LP;
 493     else if ( StrModality.find("LS") < StrModality.length()) return LS;
 494     else if ( StrModality.find("MA") < StrModality.length()) return MA;
 495     else if ( StrModality.find("MR") < StrModality.length()) return MR;
 496     else if ( StrModality.find("NM") < StrModality.length()) return NM;
 497     else if ( StrModality.find("OT") < StrModality.length()) return OT;
 498     else if ( StrModality.find("PT") < StrModality.length()) return PT;
 499     else if ( StrModality.find("RF") < StrModality.length()) return RF;
 500     else if ( StrModality.find("RG") < StrModality.length()) return RG;
 501     else if ( StrModality.find("RTDOSE")  < StrModality.length()) return RTDOSE;
 502     else if ( StrModality.find("RTIMAGE") < StrModality.length()) return RTIMAGE;
 503     else if ( StrModality.find("RTPLAN")  < StrModality.length()) return RTPLAN;
 504     else if ( StrModality.find("RTSTRUCT")< StrModality.length()) return RTSTRUCT;
 505     else if ( StrModality.find("SM") < StrModality.length()) return SM;
 506     else if ( StrModality.find("ST") < StrModality.length()) return ST;
 507     else if ( StrModality.find("TG") < StrModality.length()) return TG;
 508     else if ( StrModality.find("US") < StrModality.length()) return US;
 509     else if ( StrModality.find("VF") < StrModality.length()) return VF;
 510     else if ( StrModality.find("XA") < StrModality.length()) return XA;
 511     else if ( StrModality.find("XC") < StrModality.length()) return XC;
 512
 513     else
 514     {
 515       /// \todo throw error return value ???
 516       /// specified <> unknow in our database
 517       return Unknow;
 518     }
 519   }
 520   return Unknow;
 521 }
 522
 523 /**
 524   * \brief gets the info from 0020,0037 : Image Orientation Patient
 525   * @param iop adress of the (6)float aray to receive values
 526   * @return cosines of image orientation patient
 527   */
 528 void gdcmHeaderHelper::GetImageOrientationPatient( float* iop ) {
 529
 530   //iop is supposed to be float[6]
 531   iop[0] = iop[1] = iop[2] = iop[3] = iop[4] = iop[5] = 0;
 532
 533   // 0020 0037 DS REL Image Orientation (Patient)
 534   std::string StrImOriPat = GetEntryByNumber(0x0020,0x0037);
 535   if (StrImOriPat != GDCM_UNFOUND) {
 536     if( sscanf( StrImOriPat.c_str(), "%f\\%f\\%f\\%f\\%f\\%f",
 537             &iop[0], &iop[1], &iop[2], &iop[3], &iop[4], &iop[5]) != 6) {
 538          dbg.Verbose(0, "gdcmHeader::GetImageOrientationPatient: wrong Image Orientation Patient (0020,0037)");
 539          return ;  // bug in the element 0x0020,0x0037
 540     }
 541     else
 542       return ;
 543   }
 544
 545   //For ACR-NEMA
 546   // 0020 0035 DS REL Image Orientation (RET)
 547   StrImOriPat = GetEntryByNumber(0x0020,0x0035);
 548   if (StrImOriPat != GDCM_UNFOUND) {
 549     if( sscanf( StrImOriPat.c_str(), "%f\\%f\\%f\\%f\\%f\\%f",
 550             &iop[0], &iop[1], &iop[2], &iop[3], &iop[4], &iop[5]) != 6) {
 551          dbg.Verbose(0, "gdcmHeader::GetImageOrientationPatient: wrong Image Orientation Patient (0020,0035)");
 552          return ;  // bug in the element 0x0020,0x0035
 553     }
 554     else
 555       return ;
 556   }
 557 }
 558
 559 //-----------------------------------------------------------------------------
 560 // Protected
 561
 562 //-----------------------------------------------------------------------------
 563 // Private
 564
 565 //-----------------------------------------------------------------------------
 566
 567
 568
 569 //-----------------------------------------------------------------------------
 570 // Constructor / Destructor
 571 gdcmSerieHeaderHelper::~gdcmSerieHeaderHelper(){
 572   /// \todo
 573   for (std::list<gdcmHeaderHelper*>::iterator it  = CoherentGdcmFileList.begin();
 574         it != CoherentGdcmFileList.end(); it++)
 575   {
 576     delete *it;
 577   }
 578   CoherentGdcmFileList.clear();
 579 }
 580
 581 //-----------------------------------------------------------------------------
 582 // Print
 583
 584 //-----------------------------------------------------------------------------
 585 // Public
 586 /**
 587  * \brief add a gdcmFile to the list based on file name
 588  * @param   filename Name of the file to deal with
 589  */
 590 void gdcmSerieHeaderHelper::AddFileName(std::string filename) {
 591   gdcmHeaderHelper *GdcmFile = new gdcmHeaderHelper( filename.c_str() );
 592   this->CoherentGdcmFileList.push_back( GdcmFile );
 593 }
 594
 595 /**
 596  * \brief add a gdcmFile to the list
 597  * @param   file gdcmHeaderHelper to add
 598  */
 599 void gdcmSerieHeaderHelper::AddGdcmFile(gdcmHeaderHelper *file){
 600   this->CoherentGdcmFileList.push_back( file );
 601 }
 602
 603 /**
 604  * \brief Sets the Directory
 605  * @param   dir Name of the directory to deal with
 606  */
 607 void gdcmSerieHeaderHelper::SetDirectory(std::string dir){
 608   gdcmDirList filenames_list(dir);  //OS specific
 609
 610   for(gdcmDirList::iterator it = filenames_list.begin();
 611       it !=filenames_list.end(); it++)
 612   {
 613     gdcmHeaderHelper *file = new gdcmHeaderHelper( it->c_str() );
 614     this->CoherentGdcmFileList.push_back( file );
 615   }
 616 }
 617
 618 /**
 619  * \brief Sorts the File List
 620  * \warning This could be implemented in a 'Strategy Pattern' approach
 621  *          But as I don't know how to do it, I leave it this way
 622  *          BTW, this is also a Strategy, I don't know this is the best approach :)
 623 */
 624 void gdcmSerieHeaderHelper::OrderGdcmFileList(){
 625   if( ImagePositionPatientOrdering() ) {
 626     return ;
 627   }
 628   else if( ImageNumberOrdering() ) {
 629     return ;
 630   } else  {
 631     FileNameOrdering();
 632   }
 633 }
 634
 635 /**
 636  * \brief Gets the *coherent* File List
 637  * @return the *coherent* File List
 638 */
 639 std::list<gdcmHeaderHelper*> &gdcmSerieHeaderHelper::GetGdcmFileList() {
 640   return CoherentGdcmFileList;
 641 }
 642
 643 //-----------------------------------------------------------------------------
 644 // Protected
 645
 646 //-----------------------------------------------------------------------------
 647 // Private
 648 /**
 649  * \ingroup gdcmHeaderHelper
 650  * \brief sorts the images, according to their Patient Position
 651  *  We may order, considering :
 652  *   -# Image Number
 653  *   -# Image Position Patient
 654  *   -# More to come :)
 655  * @return false only if the header is bugged !
 656  */
 657 bool gdcmSerieHeaderHelper::ImagePositionPatientOrdering()
 658 //based on Jolinda's algorithm
 659 {
 660   //iop is calculated based on the file file
 661   float *cosines = new float[6];
 662   float normal[3];
 663   float ipp[3];
 664   float dist;
 665   float min, max;
 666   bool first = true;
 667   int n=0;
 668   std::vector<float> distlist;
 669
 670   //!\todo rewrite this for loop.
 671   for (std::list<gdcmHeaderHelper*>::iterator it  = CoherentGdcmFileList.begin();
 672         it != CoherentGdcmFileList.end(); it++)
 673   {
 674     if(first) {
 675       (*it)->GetImageOrientationPatient(cosines);
 676
 677       //You only have to do this once for all slices in the volume. Next, for
 678       //each slice, calculate the distance along the slice normal using the IPP
 679       //tag ("dist" is initialized to zero before reading the first slice) :
 680       normal[0] = cosines[1]*cosines[5] - cosines[2]*cosines[4];
 681       normal[1] = cosines[2]*cosines[3] - cosines[0]*cosines[5];
 682       normal[2] = cosines[0]*cosines[4] - cosines[1]*cosines[3];
 683
 684       ipp[0] = (*it)->GetXOrigin();
 685       ipp[1] = (*it)->GetYOrigin();
 686       ipp[2] = (*it)->GetZOrigin();
 687
 688       dist = 0;
 689       for (int i = 0; i < 3; ++i)
 690           dist += normal[i]*ipp[i];
 691
 692       if( dist == 0 )
 693       {
 694         delete[] cosines;
 695         return false;
 696       }
 697
 698       distlist.push_back( dist );
 699
 700       max = min = dist;
 701       first = false;
 702     }
 703     else {
 704       ipp[0] = (*it)->GetXOrigin();
 705       ipp[1] = (*it)->GetYOrigin();
 706       ipp[2] = (*it)->GetZOrigin();
 707
 708       dist = 0;
 709       for (int i = 0; i < 3; ++i)
 710           dist += normal[i]*ipp[i];
 711
 712       if( dist == 0 )
 713       {
 714         delete[] cosines;
 715         return false;
 716       }
 717
 718       distlist.push_back( dist );
 719
 720       min = (min < dist) ? min : dist;
 721       max = (max > dist) ? max : dist;
 722     }
 723     n++;
 724   }
 725
 726     //Then I order the slices according to the value "dist". Finally, once
 727     //I've read in all the slices, I calculate the z-spacing as the difference
 728     //between the "dist" values for the first two slices.
 729     std::vector<gdcmHeaderHelper*> CoherentGdcmFileVector(n);
 730     //CoherentGdcmFileVector.reserve( n );
 731     CoherentGdcmFileVector.resize( n );
 732     //assert( CoherentGdcmFileVector.capacity() >= n );
 733
 734     float step = (max - min)/(n - 1);
 735     int pos;
 736     n = 0;
 737
 738     //VC++ don't understand what scope is !! it -> it2
 739     for (std::list<gdcmHeaderHelper*>::iterator it2  = CoherentGdcmFileList.begin();
 740         it2 != CoherentGdcmFileList.end(); it2++, n++)
 741     {
 742       //2*n sort algo !!
 743       //Assumption: all files are present (no one missing)
 744       pos = (int)( fabs( (distlist[n]-min)/step) + .5 );
 745
 746       CoherentGdcmFileVector[pos] = *it2;
 747     }
 748
 749   CoherentGdcmFileList.clear();  //this doesn't delete list's element, node only
 750
 751   //VC++ don't understand what scope is !! it -> it3
 752   for (std::vector<gdcmHeaderHelper*>::iterator it3  = CoherentGdcmFileVector.begin();
 753         it3 != CoherentGdcmFileVector.end(); it3++)
 754   {
 755     CoherentGdcmFileList.push_back( *it3 );
 756   }
 757
 758   distlist.clear();
 759   CoherentGdcmFileVector.clear();
 760   delete[] cosines;
 761
 762   return true;
 763 }
 764
 765 /**
 766  * \ingroup gdcmHeaderHelper
 767  * \brief sorts the images, according to their Image Number
 768  * @return false only if the header is bugged !
 769  */
 770
 771 bool gdcmSerieHeaderHelper::ImageNumberOrdering() {
 772   int min, max, pos;
 773   int n = 0;//CoherentGdcmFileList.size() is a O(N) operation !!
 774   unsigned char *partition;
 775
 776   std::list<gdcmHeaderHelper*>::iterator it  = CoherentGdcmFileList.begin();
 777   min = max = (*it)->GetImageNumber();
 778
 779   for (; it != CoherentGdcmFileList.end(); it++, n++)
 780   {
 781     pos = (*it)->GetImageNumber();
 782
 783     //else
 784     min = (min < pos) ? min : pos;
 785   }
 786
 787   //bzeros(partition, n); //Cette fonction est déconseillée, utilisez plutôt memset.
 788   partition = new unsigned char[n];
 789   memset(partition, 0, n);
 790
 791   std::vector<gdcmHeaderHelper*> CoherentGdcmFileVector(n);
 792
 793   //VC++ don't understand what scope is !! it -> it2
 794   for (std::list<gdcmHeaderHelper*>::iterator it2  = CoherentGdcmFileList.begin();
 795         it2 != CoherentGdcmFileList.end(); it2++)
 796   {
 797     pos = (*it2)->GetImageNumber();
 798     CoherentGdcmFileVector[pos - min] = *it2;
 799     partition[pos - min]++;
 800   }
 801
 802   unsigned char mult = 1;
 803   for(int i=0; i<n ; i++)
 804   {
 805     mult *= partition[i];
 806   }
 807
 808   //VC++ don't understand what scope is !! it -> it3
 809   CoherentGdcmFileList.clear();  //this doesn't delete list's element, node only
 810   for (std::vector<gdcmHeaderHelper*>::iterator it3  = CoherentGdcmFileVector.begin();
 811         it3 != CoherentGdcmFileVector.end(); it3++)
 812   {
 813     CoherentGdcmFileList.push_back( *it3 );
 814   }
 815   CoherentGdcmFileVector.clear();
 816
 817   delete[] partition;
 818   return (mult!=0);
 819 }
 820
 821
 822 /**
 823  * \ingroup gdcmHeaderHelper
 824  * \brief sorts the images, according to their File Name
 825  * @return false only if the header is bugged !
 826  */
 827  bool gdcmSerieHeaderHelper::FileNameOrdering() {
 828   //using the sort
 829   //sort(CoherentGdcmFileList.begin(), CoherentGdcmFileList.end());
 830   return true;
 831 }
 832
 833 //-----------------------------------------------------------------------------