src/gdcmHeaderHelper.cxx

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