src/gdcmHeaderHelper.cxx

   1 /*=========================================================================
   2
   3   Program:   gdcm
   4   Module:    $RCSfile: gdcmHeaderHelper.cxx,v $
   5   Language:  C++
   6   Date:      $Date: 2004/10/12 04:35:46 $
   7   Version:   $Revision: 1.43 $
   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 "gdcmHeaderHelper.h"
  20 #include "gdcmDirList.h"
  21 #include "gdcmDebug.h"
  22
  23 #include <math.h>
  24 #include <algorithm>
  25 #include <vector>
  26
  27 namespace gdcm
  28 {
  29
  30 typedef std::vector<Header* > GdcmHeaderVector;
  31 //-----------------------------------------------------------------------------
  32 // Constructor / Destructor
  33 SerieHeader::SerieHeader()
  34 {
  35    CoherentGdcmFileList.clear();
  36 }
  37
  38 SerieHeader::~SerieHeader()
  39 {
  40    /// \todo
  41    for ( GdcmHeaderList::const_iterator it = CoherentGdcmFileList.begin();
  42          it != CoherentGdcmFileList.end(); ++it)
  43    {
  44       delete *it;
  45    }
  46    CoherentGdcmFileList.clear();
  47 }
  48
  49 //-----------------------------------------------------------------------------
  50 // Print
  51
  52 //-----------------------------------------------------------------------------
  53 // Public
  54 /**
  55  * \brief add a File to the list based on file name
  56  * @param   filename Name of the file to deal with
  57  */
  58 void SerieHeader::AddFileName(std::string const & filename)
  59 {
  60    Header *header = new Header( filename );
  61    CoherentGdcmFileList.push_back( header );
  62 }
  63
  64 /**
  65  * \brief add a File to the list
  66  * @param   file Header to add
  67  */
  68 void SerieHeader::AddGdcmFile(Header *file)
  69 {
  70    CoherentGdcmFileList.push_back( file );
  71 }
  72
  73 /**
  74  * \brief Sets the Directory
  75  * @param   dir Name of the directory to deal with
  76  */
  77 void SerieHeader::SetDirectory(std::string const & dir)
  78 {
  79    DirList filenames_list(dir);  //OS specific
  80
  81    for( DirList::const_iterator it = filenames_list.begin();
  82         it != filenames_list.end(); ++it)
  83    {
  84       //use string and not const char*:
  85       Header *header = new Header( *it );
  86       CoherentGdcmFileList.push_back( header );
  87    }
  88 }
  89
  90 /**
  91  * \brief Sorts the File List
  92  * \warning This could be implemented in a 'Strategy Pattern' approach
  93  *          But as I don't know how to do it, I leave it this way
  94  *          BTW, this is also a Strategy, I don't know this is the best approach :)
  95  */
  96 void SerieHeader::OrderGdcmFileList()
  97 {
  98    if( ImagePositionPatientOrdering() )
  99    {
 100       return ;
 101    }
 102    else if( ImageNumberOrdering() )
 103    {
 104       return ;
 105    }
 106    else
 107    {
 108       FileNameOrdering();
 109    }
 110 }
 111
 112 //-----------------------------------------------------------------------------
 113 // Protected
 114
 115 //-----------------------------------------------------------------------------
 116 // Private
 117 /**
 118  * \ingroup Header
 119  * \brief sorts the images, according to their Patient Position
 120  *  We may order, considering :
 121  *   -# Image Number
 122  *   -# Image Position Patient
 123  *   -# More to come :)
 124  * @return false only if the header is bugged !
 125  */
 126 bool SerieHeader::ImagePositionPatientOrdering()
 127 //based on Jolinda's algorithm
 128 {
 129    //iop is calculated based on the file file
 130    float *cosines = new float[6];
 131    float normal[3];
 132    float ipp[3];
 133    float dist;
 134    float min = 0, max = 0;
 135    bool first = true;
 136    int n=0;
 137    std::vector<float> distlist;
 138
 139    //!\todo rewrite this for loop.
 140    for ( GdcmHeaderList::const_iterator
 141          it = CoherentGdcmFileList.begin();
 142          it != CoherentGdcmFileList.end(); ++it )
 143    {
 144       if( first )
 145       {
 146          (*it)->GetImageOrientationPatient( cosines );
 147
 148          //You only have to do this once for all slices in the volume. Next,
 149          // for each slice, calculate the distance along the slice normal
 150          // using the IPP tag ("dist" is initialized to zero before reading
 151          // the first slice) :
 152          normal[0] = cosines[1]*cosines[5] - cosines[2]*cosines[4];
 153          normal[1] = cosines[2]*cosines[3] - cosines[0]*cosines[5];
 154          normal[2] = cosines[0]*cosines[4] - cosines[1]*cosines[3];
 155
 156          ipp[0] = (*it)->GetXOrigin();
 157          ipp[1] = (*it)->GetYOrigin();
 158          ipp[2] = (*it)->GetZOrigin();
 159
 160          dist = 0;
 161          for ( int i = 0; i < 3; ++i )
 162          {
 163             dist += normal[i]*ipp[i];
 164          }
 165
 166          if( dist == 0 )
 167          {
 168             delete[] cosines;
 169             return false;
 170          }
 171
 172          distlist.push_back( dist );
 173
 174          max = min = dist;
 175          first = false;
 176       }
 177       else
 178       {
 179          ipp[0] = (*it)->GetXOrigin();
 180          ipp[1] = (*it)->GetYOrigin();
 181          ipp[2] = (*it)->GetZOrigin();
 182
 183          dist = 0;
 184          for ( int i = 0; i < 3; ++i )
 185          {
 186             dist += normal[i]*ipp[i];
 187          }
 188
 189          if( dist == 0 )
 190          {
 191             delete[] cosines;
 192             return false;
 193          }
 194
 195          distlist.push_back( dist );
 196
 197          min = (min < dist) ? min : dist;
 198          max = (max > dist) ? max : dist;
 199       }
 200       ++n;
 201    }
 202
 203    // Then I order the slices according to the value "dist". Finally, once
 204    // I've read in all the slices, I calculate the z-spacing as the difference
 205    // between the "dist" values for the first two slices.
 206    GdcmHeaderVector CoherentGdcmFileVector(n);
 207    // CoherentGdcmFileVector.reserve( n );
 208    CoherentGdcmFileVector.resize( n );
 209    // assert( CoherentGdcmFileVector.capacity() >= n );
 210
 211    float step = (max - min)/(n - 1);
 212    int pos;
 213    n = 0;
 214
 215    //VC++ don't understand what scope is !! it -> it2
 216    for (GdcmHeaderList::const_iterator it2  = CoherentGdcmFileList.begin();
 217         it2 != CoherentGdcmFileList.end(); ++it2, ++n)
 218    {
 219       //2*n sort algo !!
 220       //Assumption: all files are present (no one missing)
 221       pos = (int)( fabs( (distlist[n]-min)/step) + .5 );
 222
 223       CoherentGdcmFileVector[pos] = *it2;
 224    }
 225
 226    CoherentGdcmFileList.clear();  //this doesn't delete list's element, node only
 227
 228    //VC++ don't understand what scope is !! it -> it3
 229    for (GdcmHeaderVector::const_iterator it3  = CoherentGdcmFileVector.begin();
 230         it3 != CoherentGdcmFileVector.end(); ++it3)
 231    {
 232       CoherentGdcmFileList.push_back( *it3 );
 233    }
 234
 235    distlist.clear();
 236    CoherentGdcmFileVector.clear();
 237    delete[] cosines;
 238
 239    return true;
 240 }
 241
 242 /**
 243  * \ingroup Header
 244  * \brief sorts the images, according to their Image Number
 245  * @return false only if the header is bugged !
 246  */
 247
 248 bool SerieHeader::ImageNumberOrdering()
 249 {
 250    int min, max, pos;
 251    int n = 0;//CoherentGdcmFileList.size() is a O(N) operation !!
 252    unsigned char *partition;
 253
 254    GdcmHeaderList::const_iterator it = CoherentGdcmFileList.begin();
 255    min = max = (*it)->GetImageNumber();
 256
 257    for (; it != CoherentGdcmFileList.end(); ++it, ++n)
 258    {
 259       pos = (*it)->GetImageNumber();
 260
 261       //else
 262       min = (min < pos) ? min : pos;
 263    }
 264
 265    //bzeros(partition, n); //Cette fonction est déconseillée, utilisez plutôt memset.
 266    partition = new unsigned char[n];
 267    memset(partition, 0, n);
 268
 269    GdcmHeaderVector CoherentGdcmFileVector(n);
 270
 271    //VC++ don't understand what scope is !! it -> it2
 272    for (GdcmHeaderList::const_iterator it2 = CoherentGdcmFileList.begin();
 273         it2 != CoherentGdcmFileList.end(); ++it2)
 274    {
 275       pos = (*it2)->GetImageNumber();
 276       CoherentGdcmFileVector[pos - min] = *it2;
 277       partition[pos - min]++;
 278    }
 279
 280    unsigned char mult = 1;
 281    for( int i=0; i<n ; i++ )
 282    {
 283       mult *= partition[i];
 284    }
 285
 286    //VC++ don't understand what scope is !! it -> it3
 287    CoherentGdcmFileList.clear();  //this doesn't delete list's element, node only
 288    for ( GdcmHeaderVector::const_iterator it3 = CoherentGdcmFileVector.begin();
 289          it3 != CoherentGdcmFileVector.end(); ++it3 )
 290    {
 291       CoherentGdcmFileList.push_back( *it3 );
 292    }
 293    CoherentGdcmFileVector.clear();
 294
 295    delete[] partition;
 296
 297    return (mult != 0);
 298 }
 299
 300
 301 /**
 302  * \ingroup Header
 303  * \brief sorts the images, according to their File Name
 304  * @return false only if the header is bugged !
 305  */
 306 bool SerieHeader::FileNameOrdering()
 307 {
 308    //using the sort
 309    //sort(CoherentGdcmFileList.begin(), CoherentGdcmFileList.end());
 310    return true;
 311 }
 312
 313 } // end namespace gdcm
 314 //-----------------------------------------------------------------------------