src/gdcmSerieHelper.cxx

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