src/gdcmSerieHelper.cxx

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