src/gdcmHeaderHelper.cxx

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