src/gdcmSerieHeader.cxx

   1 /*=========================================================================
   2
   3   Program:   gdcm
   4   Module:    $RCSfile: gdcmSerieHeader.cxx,v $
   5   Language:  C++
   6   Date:      $Date: 2005/01/31 12:19:34 $
   7   Version:   $Revision: 1.17 $
   8
   9   Copyright (c) CREATIS (Centre de Recherche et d'Applications en Traitement de
  10   l'Image). All rights reserved. See Doc/License.txt or
  11   http://www.creatis.insa-lyon.fr/Public/Gdcm/License.html for details.
  12
  13      This software is distributed WITHOUT ANY WARRANTY; without even
  14      the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  15      PURPOSE.  See the above copyright notices for more information.
  16
  17 =========================================================================*/
  18
  19 #include "gdcmSerieHeader.h"
  20 #include "gdcmDirList.h"
  21 #include "gdcmFile.h"
  22 #include "gdcmDebug.h"
  23
  24 #include <math.h>
  25 #include <algorithm>
  26 #include <vector>
  27
  28 namespace gdcm
  29 {
  30 typedef std::list<File* > GdcmFileList;
  31 typedef std::vector<File* > GdcmFileVector;
  32
  33 //-----------------------------------------------------------------------------
  34 // Constructor / Destructor
  35 /**
  36  * \brief   Constructor from a given SerieHeader
  37  */
  38 SerieHeader::SerieHeader()
  39 {
  40    //TODO
  41    //CoherentGdcmFileList.clear();
  42 }
  43
  44 /**
  45  * \brief   Canonical destructor.
  46  */
  47 SerieHeader::~SerieHeader()
  48 {
  49    // For all the Coherent File lists of the gdcm::Serie
  50    GdcmFileList *l = GetFirstCoherentFileList();
  51    while (l)
  52    {
  53       // For all the files of a Coherent File list
  54       for (GdcmFileList::iterator it = l->begin();
  55                                   it != l->end();
  56                                 ++it)
  57       {
  58          delete *it;
  59       }
  60       l->clear();
  61       l = GetNextCoherentFileList();
  62    }
  63 }
  64
  65 //-----------------------------------------------------------------------------
  66 // Print
  67
  68 /**
  69  * \brief   Canonical printer.
  70  */
  71 void SerieHeader::Print()
  72 {
  73    // For all the Coherent File lists of the gdcm::Serie
  74    CoherentFileListmap::iterator itl = CoherentGdcmFileListHT.begin();
  75    if ( itl == CoherentGdcmFileListHT.end() )
  76    {
  77       gdcmVerboseMacro( "No Coherent File list found" );
  78       return;
  79    }
  80    while (itl != CoherentGdcmFileListHT.end())
  81    {
  82       std::cout << "Serie UID :[" << itl->first << "]" << std::endl;
  83
  84       // For all the files of a Coherent File list
  85       for (GdcmFileList::iterator it =  (itl->second)->begin();
  86                                   it != (itl->second)->end();
  87                                 ++it)
  88       {
  89          std::cout << " --- " << (*it)->GetFileName() << std::endl;
  90       }
  91       ++itl;
  92    }
  93 }
  94
  95 //-----------------------------------------------------------------------------
  96 // Public
  97 /**
  98  * \brief add a gdcm::File to the list corresponding to its Serie UID
  99  * @param   filename Name of the file to deal with
 100  */
 101 void SerieHeader::AddFileName(std::string const &filename)
 102 {
 103    //directly use string and not const char*:
 104    File *header = new File( filename );
 105    if( header->IsReadable() )
 106    {
 107       // 0020 000e UI REL Series Instance UID
 108       std::string uid =  header->GetEntryValue (0x0020, 0x000e);
 109       // if uid == GDCM_UNFOUND then consistently we should find GDCM_UNFOUND
 110       // no need here to do anything special
 111
 112       if ( CoherentGdcmFileListHT.count(uid) == 0 )
 113       {
 114          gdcmVerboseMacro(" New Serie UID :[" << uid << "]");
 115          // create a std::list in 'uid' position
 116          CoherentGdcmFileListHT[uid] = new GdcmFileList;
 117       }
 118       // Current Serie UID and DICOM header seems to match add the file:
 119       CoherentGdcmFileListHT[uid]->push_back( header );
 120    }
 121    else
 122    {
 123       gdcmVerboseMacro("Could not read file: " << filename );
 124       delete header;
 125    }
 126 }
 127
 128 /**
 129  * \brief Sets the root Directory
 130  * @param   dir Name of the directory to deal with
 131  * @param recursive whether we want explore recursively the Directory
 132  */
 133 void SerieHeader::SetDirectory(std::string const &dir, bool recursive)
 134 {
 135    DirList dirList(dir, recursive); // OS specific
 136
 137    DirListType filenames_list = dirList.GetFilenames();
 138    for( DirListType::const_iterator it = filenames_list.begin();
 139         it != filenames_list.end(); ++it)
 140    {
 141       AddFileName( *it );
 142    }
 143 }
 144
 145 /**
 146  * \brief Sorts the given File List
 147  * \warning This could be implemented in a 'Strategy Pattern' approach
 148  *          But as I don't know how to do it, I leave it this way
 149  *          BTW, this is also a Strategy, I don't know this is the best approach :)
 150  */
 151 void SerieHeader::OrderGdcmFileList(GdcmFileList *CoherentGdcmFileList)
 152 {
 153    if( ImagePositionPatientOrdering( CoherentGdcmFileList ) )
 154    {
 155       return ;
 156    }
 157    else if( ImageNumberOrdering(CoherentGdcmFileList ) )
 158    {
 159       return ;
 160    }
 161    else
 162    {
 163       FileNameOrdering(CoherentGdcmFileList );
 164    }
 165 }
 166
 167 /**
 168  * \brief   Get the first List while visiting the CoherentFileListHT
 169  * @return  The first GdcmFileList if found, otherwhise NULL
 170  */
 171  std::list<File* > *SerieHeader::GetFirstCoherentFileList()
 172 // Why doesn't it compile ?!?
 173 //GdcmFileList *SerieHeader::GetFirstCoherentFileList()
 174 {
 175    ItListHt = CoherentGdcmFileListHT.begin();
 176    if( ItListHt != CoherentGdcmFileListHT.end() )
 177       return ItListHt->second;
 178    return NULL;
 179 }
 180
 181 /**
 182  * \brief   Get the next List while visiting the CoherentFileListHT
 183  * \note : meaningfull only if GetFirstCoherentFileList already called
 184  * @return  The next GdcmFileList if found, otherwhise NULL
 185  */
 186 std::list<File* > *SerieHeader::GetNextCoherentFileList()
 187 {
 188    gdcmAssertMacro (ItListHt != CoherentGdcmFileListHT.end());
 189
 190    ++ItListHt;
 191    if ( ItListHt != CoherentGdcmFileListHT.end() )
 192       return ItListHt->second;
 193    return NULL;
 194 }
 195
 196 /**
 197  * \brief   Get the Coherent Files list according to its Serie UID
 198  * @param SerieUID SerieUID
 199  * \return  pointer to the Coherent Filseslist if found, otherwhise NULL
 200  */
 201 GdcmFileList *SerieHeader::GetCoherentFileList(std::string SerieUID)
 202 {
 203    if ( CoherentGdcmFileListHT.count(SerieUID) == 0 )
 204       return 0;
 205    return CoherentGdcmFileListHT[SerieUID];
 206 }
 207
 208
 209 //-----------------------------------------------------------------------------
 210 // Protected
 211
 212 //-----------------------------------------------------------------------------
 213 // Private
 214 /**
 215  * \brief sorts the images, according to their Patient Position
 216  *  We may order, considering :
 217  *   -# Image Position Patient
 218  *   -# Image Number
 219  *   -# More to come :-)
 220  * @param CoherentGdcmFileList Coherent File list (same Serie UID) to sort
 221  * @return false only if the header is bugged !
 222  */
 223 bool SerieHeader::ImagePositionPatientOrdering(
 224                                        GdcmFileList *CoherentGdcmFileList )
 225 //based on Jolinda's algorithm
 226 {
 227    //iop is calculated based on the file file
 228    float cosines[6];
 229    float normal[3];
 230    float ipp[3];
 231    float dist;
 232    float min = 0, max = 0;
 233    bool first = true;
 234    int n=0;
 235    std::vector<float> distlist;
 236
 237    //!\todo rewrite this for loop.
 238    for ( GdcmFileList::const_iterator
 239          it = CoherentGdcmFileList->begin();
 240          it != CoherentGdcmFileList->end(); ++it )
 241    {
 242       if( first )
 243       {
 244          (*it)->GetImageOrientationPatient( cosines );
 245
 246          // You only have to do this once for all slices in the volume. Next,
 247          // for each slice, calculate the distance along the slice normal
 248          // using the IPP tag ("dist" is initialized to zero before reading
 249          // the first slice) :
 250          normal[0] = cosines[1]*cosines[5] - cosines[2]*cosines[4];
 251          normal[1] = cosines[2]*cosines[3] - cosines[0]*cosines[5];
 252          normal[2] = cosines[0]*cosines[4] - cosines[1]*cosines[3];
 253
 254          ipp[0] = (*it)->GetXOrigin();
 255          ipp[1] = (*it)->GetYOrigin();
 256          ipp[2] = (*it)->GetZOrigin();
 257
 258          dist = 0;
 259          for ( int i = 0; i < 3; ++i )
 260          {
 261             dist += normal[i]*ipp[i];
 262          }
 263
 264          if( dist == 0 )
 265          {
 266             return false;
 267          }
 268
 269          distlist.push_back( dist );
 270
 271          max = min = dist;
 272          first = false;
 273       }
 274       else
 275       {
 276          ipp[0] = (*it)->GetXOrigin();
 277          ipp[1] = (*it)->GetYOrigin();
 278          ipp[2] = (*it)->GetZOrigin();
 279
 280          dist = 0;
 281          for ( int i = 0; i < 3; ++i )
 282          {
 283             dist += normal[i]*ipp[i];
 284          }
 285
 286          if( dist == 0 )
 287          {
 288             return false;
 289          }
 290
 291          distlist.push_back( dist );
 292
 293          min = (min < dist) ? min : dist;
 294          max = (max > dist) ? max : dist;
 295       }
 296       ++n;
 297    }
 298
 299    // Then I order the slices according to the value "dist". Finally, once
 300    // I've read in all the slices, I calculate the z-spacing as the difference
 301    // between the "dist" values for the first two slices.
 302    GdcmFileVector CoherentGdcmFileVector(n);
 303    // CoherentGdcmFileVector.reserve( n );
 304    CoherentGdcmFileVector.resize( n );
 305    // gdcmAssertMacro( CoherentGdcmFileVector.capacity() >= n );
 306
 307    float step = (max - min)/(n - 1);
 308    int pos;
 309    n = 0;
 310
 311    //VC++ don't understand what scope is !! it -> it2
 312    for (GdcmFileList::const_iterator it2  = CoherentGdcmFileList->begin();
 313         it2 != CoherentGdcmFileList->end(); ++it2, ++n)
 314    {
 315       //2*n sort algo !!
 316       //Assumption: all files are present (no one missing)
 317       pos = (int)( fabs( (distlist[n]-min)/step) + .5 );
 318
 319       // a Dicom 'Serie' may contain scout views
 320       // and images may have differents directions
 321       // -> More than one may have the same 'pos'
 322       // Sorting has then NO meaning !
 323       if (CoherentGdcmFileVector[pos]==NULL)
 324          CoherentGdcmFileVector[pos] = *it2;
 325       else
 326          return false;
 327    }
 328
 329    CoherentGdcmFileList->clear();  // doesn't delete list elements, only node
 330
 331    //VC++ don't understand what scope is !! it -> it3
 332    for (GdcmFileVector::const_iterator it3  = CoherentGdcmFileVector.begin();
 333         it3 != CoherentGdcmFileVector.end(); ++it3)
 334    {
 335       CoherentGdcmFileList->push_back( *it3 );
 336    }
 337
 338    distlist.clear();
 339    CoherentGdcmFileVector.clear();
 340
 341    return true;
 342 }
 343
 344 /**
 345  * \brief sorts the images, according to their Image Number
 346  * @param CoherentGdcmFileList Coherent File list (same Serie UID) to sort
 347  * @return false only if the header is bugged !
 348  */
 349 bool SerieHeader::ImageNumberOrdering(GdcmFileList *CoherentGdcmFileList)
 350 {
 351    int min, max, pos;
 352    int n = 0;//CoherentGdcmFileList.size() is a O(N) operation
 353    unsigned char *partition;
 354
 355    GdcmFileList::const_iterator it = CoherentGdcmFileList->begin();
 356    min = max = (*it)->GetImageNumber();
 357
 358    for (; it != CoherentGdcmFileList->end(); ++it, ++n)
 359    {
 360       pos = (*it)->GetImageNumber();
 361
 362       //else
 363       min = (min < pos) ? min : pos;
 364       max = (max > pos) ? max : pos;
 365    }
 366
 367    // Find out if sorting worked:
 368    if( min == max || max == 0 || max > (n+min)) return false;
 369
 370    //bzeros(partition, n); //This function is deprecated, better use memset.
 371    partition = new unsigned char[n];
 372    memset(partition, 0, n);
 373
 374    GdcmFileVector CoherentGdcmFileVector(n);
 375
 376    //VC++ don't understand what scope is !! it -> it2
 377    for (GdcmFileList::const_iterator it2 = CoherentGdcmFileList->begin();
 378         it2 != CoherentGdcmFileList->end(); ++it2)
 379    {
 380       pos = (*it2)->GetImageNumber();
 381       CoherentGdcmFileVector[pos - min] = *it2;
 382       partition[pos - min]++;
 383    }
 384
 385    unsigned char mult = 1;
 386    for( int i=0; i<n ; i++ )
 387    {
 388       mult *= partition[i];
 389    }
 390
 391    //VC++ don't understand what scope is !! it -> it3
 392    CoherentGdcmFileList->clear();  // doesn't delete list elements, only node
 393    for ( GdcmFileVector::const_iterator it3 = CoherentGdcmFileVector.begin();
 394          it3 != CoherentGdcmFileVector.end(); ++it3 )
 395    {
 396       CoherentGdcmFileList->push_back( *it3 );
 397    }
 398    CoherentGdcmFileVector.clear();
 399
 400    delete[] partition;
 401
 402    return mult != 0;
 403 }
 404
 405 /**
 406  * \brief sorts the images, according to their File Name
 407  * @param CoherentGdcmFileList Coherent File list (same Serie UID) to sort
 408  * @return false only if the header is bugged !
 409  */
 410 bool SerieHeader::FileNameOrdering(GdcmFileList *)
 411 {
 412    //using the sort
 413    //sort(CoherentGdcmFileList.begin(), CoherentGdcmFileList.end());
 414    return true;
 415 }
 416
 417 } // end namespace gdcm
 418 //-----------------------------------------------------------------------------