X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=blobdiff_plain;f=vtk%2FvtkGdcm4DSplitter.cxx;h=e3e27d813c4eb72ee04cbe72971da146689d38be;hb=refs%2Fheads%2Fmaster;hp=b039a8b45173ae33c5bb0ce05e26a5c6415011d9;hpb=0046646fb08bd47ba81a0892d2f253eaadaefc9e;p=gdcm.git diff --git a/vtk/vtkGdcm4DSplitter.cxx b/vtk/vtkGdcm4DSplitter.cxx index b039a8b4..e3e27d81 100644 --- a/vtk/vtkGdcm4DSplitter.cxx +++ b/vtk/vtkGdcm4DSplitter.cxx @@ -3,8 +3,8 @@ Program: gdcm Module: $RCSfile: vtkGdcm4DSplitter.cxx,v $ Language: C++ - Date: $Date: 2011/03/29 15:45:38 $ - Version: $Revision: 1.4 $ + Date: $Date: 2011/10/25 20:36:21 $ + Version: $Revision: 1.16 $ Copyright (c) CREATIS (Centre de Recherche et d'Applications en Traitement de l'Image). All rights reserved. See Doc/License.txt or @@ -38,14 +38,22 @@ vtkGdcm4DSplitter A 'Dicom Serie' doesn't mean always the same thing : a given Slice along the time a given Volume at a given time -Sometimes, an image within a serie is so artefacted than user decides to replace -it by an other image. +Sometimes, an image within a serie is so artefacted than user decides to replace it + by an other image. User needs to be aware, *only him* knows want he wants to do. -vtkGdcm4DSplitter class does the job for hom +vtkGdcm4DSplitter class does the job for him (despite its name, it works on 3D or 2D+T images too) -User will have to specify some points +==> To (try to) understand how the file store is organised, user is deeply encourage to use something like: + +PrintFile dirin=Root_Directory_holding_the_images rec > 1.txt + +open it with any text editor, and/or (Linux users) +grep a_supposed_to_be_string_of_interest 1.txt + + +Aware user will have to specify some points : . Choose input data ------------------- @@ -56,6 +64,8 @@ User will have to specify some points bool setVectDirName(std::vector &vectDirName); - a list of files bool setVectFileName(std::vector &vectFileName); +- a list of gdcm::File* + bool setVectGdcmFile(std::vector &vectGdcmFile); - Recursive directory exploration void setRecursive(bool recursive); @@ -70,81 +80,122 @@ User will have to specify some points - User choosen tag void setSplitOnTag(unsigned short splitGroup, unsigned short splitElem); void setSplitConvertToFloat(bool conv); - + - UserDefined Function + void setSortOnUserFunction (FoncComp f); + + for 'true' 3D image sets : + - if you want to get a single 3D vtkImageData, use SplitOnOrientation -i.e. no split- + - if you want to get a vector of 2D vtkImageData, use SplitOnPosition -i.e. one slice in each 'XCoherent fileset'- + + for 'true' 4D multi-orientation image sets (i.e. a stack of n axial + m sagital + p coronal images, at different instants ...) + --> this is 5D, right? + (almost) nothing done, yet : + . use setSplitOnly() + . Use a first time vtkGdcm4DSplitter with setSplitOnOrientation(); + . Get the VectGdcmFileLists (a std::vector of 'XCoherent fileset') + . use vtkGdcm4DSplitter, with as many setVectGdcmFile(std::vector &vectGdcmFile) you need + one per element of std::vector + think on 'spliting' and 'sorting' it, according to your needs. + . Choose 'sort' criterion : -------------------------- - ImagePositionPatient - void setSortOnPosition(); - - User choosen tag - void setSortOnTag(unsigned short sortGroup, unsigned short sortElem); - void setSortConvertToFloat(bool conv) - -. Execute + void setSortOnPosition(); + - ImageOrientationPatient + ==> Only in your dreams! + ==> or, please, write a IOP sorter ... + - UserDefined Function + void setSortOnUserFunction (FoncComp f); + - File name + void setSortOnFileName() + +- Choose 'flip' option +---------------------- +VTK reader 'flips upside down' the images to display them correctly ?!? +When a 'math guy' gets the image, the first thing he does is flipping it again. +Avoid loosing CPU time using : + void SetFlipY(bool); // Wrong = no flip + +. Execute : ----------- bool Go(); . Get the result ---------------- - -a single vtkImageData: - vtkImageData *GetImageData(); +// -a single vtkImageData: +// vtkImageData *GetImageData(); - a vector of vtkImageData std::vector *GetImageDataVector(); ===================================================================== */ -#include "gdcmSerieHelper.h" - #include "vtkGdcmReader.h" #include "vtkGdcm4DSplitter.h" #include #include "gdcmSerieHelper.h" // for ImagePositionPatientOrdering() #include // for atof +// Constructor / Destructor +/** + * \brief Constructor from a given vtkGdcm4DSplitter + */ vtkGdcm4DSplitter::vtkGdcm4DSplitter() : SplitOnPosition(false), SplitOnOrientation(false), SplitOnTag(false), - SplitGroup(0), SplitElem(0), + SplitGroup(0), SplitElem(0), SplitConvertToFloat(false), SplitOnly(false), - SortOnPosition(false), SortOnOrientation(false), SortOnTag(false), - SortGroup(0), SortElem(0), SortConvertToFloat(false), + SortOnPosition(false), SortOnOrientation(false), /* SortOnTag(false),*/ SortOnFileName(false), SortOnUserFunction(false), + /*SortGroup(0), SortElem(0), SortConvertToFloat(false),*/ - Recursive(false), TypeDir(0), - verbose(true) + Recursive(false),FlipY(true), TypeDir(0), + verbose(false) { } - std::vector * vtkGdcm4DSplitter::GetImageDataVector() - { - if (verbose) std::cout << "TypeDir " << TypeDir << std::endl; - if (TypeResult == 2) - return ImageDataVector; - else - if (TypeResult == 1) - { - std::vector *t = new std::vector; - t->push_back( ImageData ); - return t; - } - else - return (std::vector*) NULL; - } - - vtkImageData *vtkGdcm4DSplitter::GetImageData() +/** + * \brief Canonical destructor. + */ + vtkGdcm4DSplitter::~vtkGdcm4DSplitter() { - if (TypeResult == 1) - return ImageData; - else - if (TypeResult == 1) - { - return (*ImageDataVector)[0]; - } - else - return (vtkImageData*) NULL; + /// \TODO : delete everything that must be! + + GDCM_NAME_SPACE::XCoherentFileSetmap::iterator it; + std::vector::iterator it2; + + for ( it = xcm.begin(); // for each std::vector + it != xcm.end(); + ++it) + { + for ( it2 = (*it).second->begin(); // for each gdcm::File* + it2 != (*it).second->end(); + ++it2) + { + (*it2)->Delete(); // delete gdcm::File + } + delete (*it).second; // delete the now empty std::vector + } + + // VectGdcmFileLists is a vector of pointers on elements of xcm : nothing to do! + } - - + + // Locate Data to process + // ====================== +/** + * \brief sets the directories exploration mode + * @param recursive whether we want explore recursively the root Directory + */ +void vtkGdcm4DSplitter::setRecursive(bool recursive) +{ + Recursive=recursive; +} + +/** + * \brief Sets the root Directory to get the images from + * @param dirName name of the directory to deal with + */ bool vtkGdcm4DSplitter::setDirName(std::string &dirName) { if ( ! GDCM_NAME_SPACE::DirList::IsDirectory(dirName) ) @@ -156,7 +207,11 @@ User will have to specify some points TypeDir=1; return true; } - + +/** + * \brief Sets a list of Directories to get the images from + * @param vectDirName vector of directory names to deal with + */ bool vtkGdcm4DSplitter::setVectDirName(std::vector &vectDirName) { int nbDir = vectDirName.size(); @@ -173,7 +228,11 @@ User will have to specify some points TypeDir=2; return true; } - + +/** + * \brief Sets a list of files names to read + * @param vectFileName vector of file names to deal with + */ bool vtkGdcm4DSplitter::setVectFileName(std::vector &vectFileName) { if ( vectFileName.size() == 0) @@ -186,29 +245,187 @@ User will have to specify some points return true; } - bool vtkGdcm4DSplitter::CompareOnSortTagConvertToFloat(GDCM_NAME_SPACE::File *file1, GDCM_NAME_SPACE::File *file2) - { - if (verbose) printf ("%04x %04x\n", this->SortGroup,this->SortElem); - if (verbose) std :: cout << file1->GetEntryString(SortGroup,SortElem).c_str() << " : " - << atof(file1->GetEntryString(SortGroup,SortElem).c_str()) - << std::endl; - return atof(file1->GetEntryString(SortGroup,SortElem).c_str()) < atof(file2->GetEntryString(SortGroup,SortElem).c_str()); +/** + * \brief Sets an already filled std::vector of gdcm::File * + * @param vectGdcmFileName std::vector of gdcm::File * + */ + + bool vtkGdcm4DSplitter::setVectGdcmFile(GDCM_NAME_SPACE::FileList *vectGdcmFile) + { + if ( vectGdcmFile->size() == 0) + { + std::cout << "[ vectGdcmFile ] : empty list" << std::endl; + return false; + } + TypeDir=4; + VectGdcmFile = vectGdcmFile; } - bool vtkGdcm4DSplitter::CompareOnSortTag(GDCM_NAME_SPACE::File *file1, GDCM_NAME_SPACE::File *file2) + + // Split criterion + // =============== + // + +/** + * \brief asks for splitting Filesets according to the Position + */ + void vtkGdcm4DSplitter::setSplitOnPosition() { - return file1->GetEntryString(SortGroup,SortElem) < file2->GetEntryString(SortGroup,SortElem); + SplitOnPosition=true; + SplitOnOrientation=false; + SplitOnTag=false; } - - bool vtkGdcm4DSplitter::Go() + +/** + * \brief asks for splitting Filesets according to the Orientation + */ + void vtkGdcm4DSplitter::setSplitOnOrientation() { - if (!SplitOnPosition && !SplitOnOrientation && !SplitOnTag ) - { - ///\TODO (?) Throw an exception "Choose Splitting mode before!" - std::cout << "Choose Splitting mode before!" << std::endl; - return false; - } + SplitOnPosition=false; + SplitOnOrientation=true; + SplitOnTag=false; + } + +/** + * \brief asks for splitting Filesets according to the value of a given Tag + * @param group group number of the target Element + * @param element element number of the target Element + */ + void vtkGdcm4DSplitter::setSplitOnTag(unsigned short int splitGroup, unsigned short int splitElem) + { + SplitOnPosition=false; + SplitOnOrientation=false; + SplitOnTag=true; + SplitGroup=splitGroup; + SplitElem=splitElem; + } + +/** + * \brief asks for converting to 'float' the tag values used as a splitting criteria (lexicographic order may not be suitable) + */ + void vtkGdcm4DSplitter::setSplitConvertToFloat(bool conv) {SplitConvertToFloat=conv;} + +/** + * \brief asks for splitting Filesets according to what was asked for (no sorting, no reading data) + */ + void vtkGdcm4DSplitter::setSplitOnly(bool s) + { + SplitOnly = s; + } + + // Sort criterion + // ============== + // + +/** + * \brief asks for IPP sorting each XCoherent gdcm::FILE set + */ + void vtkGdcm4DSplitter::setSortOnPosition() + { + SortOnPosition=true; + SortOnOrientation=false; + //SortOnTag=false; + SortOnFileName=false; + SortOnUserFunction=false; + SortOnPosition=true; + } + + // use setSortOnUserFunction, instead! + // void setSortOnTag(unsigned short int sortGroup, unsigned short int sortElem) + // { + // SortOnPosition=false; + // SortOnOrientation=false; + // SortOnTag=true; + // SortOnFileName=false; + // SortOnUserFunction=false; + // SortGroup=sortGroup; SortElem=sortElem; + // } + + +/** + * \brief sets a user supplied function (comparison) + * @param f comparison function + */ + void vtkGdcm4DSplitter::setSortOnUserFunction (FoncComp f) + { + UserCompareFunction=f; + SortOnPosition=false; + SortOnOrientation=false; + //SortOnTag=false; + SortOnFileName=false; + SortOnUserFunction=true; + } + + // void setSortConvertToFloat(bool conv) + // { + // SortConvertToFloat=conv; + // } +/** + * \brief asks for sorting each XCoherent gdcm::FILE set, according to the File names + */ + void vtkGdcm4DSplitter::setSortOnFileName() + { + SortOnPosition=false; + SortOnOrientation=false; + //SortOnTag=false; + SortOnFileName=true; + SortOnUserFunction=false; + } + +/** + * \brief returns a std::vector of gdcm::FileList* (gdcm::FileList : actually, a std::vector of gdcm::File*) + */ + std::vector *vtkGdcm4DSplitter::GetVectGdcmFileLists() + { + if (SplitOnly) + return &VectGdcmFileLists; + + GDCM_NAME_SPACE::XCoherentFileSetmap::iterator it; + for ( it = xcm.begin(); + it != xcm.end(); + ++it) + { + VectGdcmFileLists.push_back((*it).second); + } + return &VectGdcmFileLists; + } + +/** + * \brief returns a std::vector of [2D/3D, depending on what was passed] vtkImageData* + */ + std::vector * vtkGdcm4DSplitter::GetImageDataVector() + { + if (SplitOnly) + return NULL; + return ImageDataVector; + } + +/** + * \brief when user _knows_ only _one_ vtkImageData* is returned he may be interested in not getting a vector... + */ + vtkImageData *vtkGdcm4DSplitter::GetImageData() + { + if (SplitOnly) + return NULL; + return (*ImageDataVector)[0]; + } + + +// bool vtkGdcm4DSplitter::CompareOnSortTagConvertToFloat(GDCM_NAME_SPACE::File *file1, GDCM_NAME_SPACE::File *file2) +// { +// return atof(file1->GetEntryString(vtkGdcm4DSplitter::SortGroup,vtkGdcm4DSplitter::SortElem).c_str()) < atof(file2->GetEntryString(vtkGdcm4DSplitter::SortGroup,vtkGdcm4DSplitter::SortElem).c_str()); +// return atof(file1->GetEntryString(SortGroup,SortElem).c_str()) < atof(file2->GetEntryString(SortGroup,SortElem).c_str()); +// } + +// bool vtkGdcm4DSplitter::CompareOnSortTag(GDCM_NAME_SPACE::File *file1, GDCM_NAME_SPACE::File *file2) +// { +// return file1->GetEntryString(vtkGdcm4DSplitter::SortGroup,vtkGdcm4DSplitter::SortElem) < file2->GetEntryString(vtkGdcm4DSplitter::SortGroup,vtkGdcm4DSplitter::SortElem); +// } + + + bool vtkGdcm4DSplitter::Go() + { /// How To : /* entree nom de directory / Vecteur de noms? @@ -220,72 +437,45 @@ User will have to specify some points trier chaque Coherent file set passer chacun a un vtkGcdmReader retourner le (vecteur de) vtkImageData - */ - - GDCM_NAME_SPACE::SerieHelper *s; - s = GDCM_NAME_SPACE::SerieHelper::New(); + */ - GDCM_NAME_SPACE::File *f; - GDCM_NAME_SPACE::DirListType fileNames; - - if (TypeDir == 0 ) - { - ///\TODO (?) Throw an exception "Set input Directory name(s) / file names before!" - std::cout << "Set input Directory name(s) / file names before!" << std::endl; - return false; - } - else if (TypeDir == 1) + if (!SplitOnPosition && !SplitOnOrientation && !SplitOnTag) { - GDCM_NAME_SPACE::DirList dirlist(DirName, Recursive); // NO recursive exploration - fileNames = dirlist.GetFilenames(); // all the file names - } - - else if (TypeDir == 2) - { - int nbDir = VectDirName.size(); - GDCM_NAME_SPACE::DirListType tmpFileNames; - for (int iDir=0; iDirSetMaxSizeLoadEntry(maxSize); - f->SetFileName( *it ); - if (f->Load()) - l->push_back(f); - else - std::cout << " Fail to load [" << *it << "]" << std::endl; - } + GDCM_NAME_SPACE::FileList *l; - GDCM_NAME_SPACE::XCoherentFileSetmap xcm; + GDCM_NAME_SPACE::SerieHelper *s; + s = GDCM_NAME_SPACE::SerieHelper::New(); +// +// Load the gdcm::File* set, according to user's requierements +// ------------------------------------------------------------ +// + l = getGdcmFileList(); +std::cout << l->size() << " gdcm::File read" << std::endl; +// +// Split the gdcm::File* set, according to user's requierements +// ------------------------------------------------------------ +// if (SplitOnOrientation) { - s->SetDropDuplicatePositions(false); - xcm = s->SplitOnOrientation(l); + s->SetDropDuplicatePositions(false); + xcm = s->SplitOnOrientation(l); } else if (SplitOnPosition) { - s->SetDropDuplicatePositions(true); - xcm = s->SplitOnPosition(l); + s->SetDropDuplicatePositions(true); + xcm = s->SplitOnPosition(l); + + // the key of xcm follows lexicographical order + // (that may be different than the 'distance' order) + // we have to reorganize it! + + reorgXCoherentFileSetmap(xcm); } else if (SplitOnTag) { @@ -297,41 +487,49 @@ User will have to specify some points xcm = s->SplitOnTagValue(l, SplitGroup, SplitElem); else { - xcm = s->SplitOnTagValueConvertToFloat(l, SplitGroup, SplitElem); + xcm = s->SplitOnTagValueConvertToFloat(l, SplitGroup, SplitElem); } - } - + } + if (xcm.size() == 0) { - if(verbose) std::cout << "Empty XCoherent File Set after 'split' ?!?" << std::endl; + if(verbose) + std::cout << "Empty XCoherent File Set after 'split' ?!?" << std::endl; return false; } - else if (xcm.size() == 1) - TypeResult=1; else - TypeResult=2; - - ImageDataVector = new std::vector; - // vtkGdcmReader *reader = vtkGdcmReader::New(); // move inside the loop, or be clever using vtk! - - for (GDCM_NAME_SPACE::XCoherentFileSetmap::iterator i = xcm.begin(); - i != xcm.end(); - ++i) { - if (verbose) - std::cout << "--- xCoherentName = [" << (*i).first << "]" << std::endl; + std::cout << xcm.size() << " XCoherent entries found" << std::endl; } - // XCoherentFileSetmap map < critère de split, FileList (= std::vector de gdcm::File*) > + + // put here, to avoid segfault when unaware user sets SplitOnly to true, and the asks for ImageDataVector + ImageDataVector = new std::vector; + + std::cout <<"SplitOnly " << SplitOnly << std::endl; + if(SplitOnly) + return true; +// +// +// ------------------------------------------------------------ +// +// ImageDataVector = new std::vector; + + /// \TODO move inside the loop, or be clever using vtk! + vtkGdcmReader *reader = vtkGdcmReader::New(); // move inside the loop, or be clever using vtk! + + // XCoherentFileSetmap map < critère de split, FileList (= std::vector de gdcm::File*) > for (GDCM_NAME_SPACE::XCoherentFileSetmap::iterator i = xcm.begin(); - i != xcm.end(); - ++i) + i != xcm.end(); + ++i) { - - vtkGdcmReader *reader = vtkGdcmReader::New(); /// \FIXME : unable to delete! - + //vtkGdcmReader *reader = vtkGdcmReader::New(); /// \TODO FIXME : unable to delete! + + reader->SetFlipY(FlipY); + // better user SetFileLowerLeft() + /// \TODO : modify vtkGdcmReader ! if (verbose) - std::cout << "==========================================xCoherentName = [" << (*i).first << "]" << std::endl; + std::cout << " --- xCoherentName = [" << (*i).first << "]" << std::endl; if (SortOnPosition) { @@ -344,15 +542,7 @@ User will have to specify some points else if (SortOnOrientation) { if (verbose) std::cout << "SortOnOrientation" << std::endl; - /// \TODO SortOnOrientation() - // Within a 'just to see' program, - // OrderFileList() causes trouble, since some files - // (eg:MIP views) don't have 'Position', now considered as mandatory - // --> Activated on user demand. - - // Information is in : - // 0020,0032 : Image Position Patient - // 0020,0030 : Image Position (RET) + /// \TODO (?) SortOnOrientation() // we still miss an algo to sort an Orientation, given by 6 cosines! // Anything like this, in GDCM2? @@ -372,38 +562,445 @@ User will have to specify some points if (verbose) std::cout << "Out of SortOnFileName" << std::endl; } - else if (SortOnTag) - { - if (verbose) std::cout << "SortOnTag" << std::endl; - printf ("--> %04x %04x\n", SortGroup,SortElem); - if ( SortConvertToFloat ) - s->SetUserLessThanFunction( reinterpret_cast - ( &vtkGdcm4DSplitter::CompareOnSortTagConvertToFloat)); - else - s->SetUserLessThanFunction( reinterpret_cast - ( &vtkGdcm4DSplitter::CompareOnSortTag)); +// else if (SortOnTag) +// { +// if (verbose) std::cout << "SortOnTag" << std::endl; +// printf ("--> %04x %04x\n", SortGroup,SortElem); +// std::cout << "Sorry, troubles not solved yet; use SortOnUserFunction, right now!" << std::endl; + + /* ==> WARNING : This one has troubles; do NOT use it, right now! + // a pointer to fonction cannot be casted as a pointer to member function! + // Use SortOnUserFunction, instead! + + // if ( SortConvertToFloat ) + // s->SetUserLessThanFunction( reinterpret_cast + // ( &vtkGdcm4DSplitter::CompareOnSortTagConvertToFloat)); + // else + // s->SetUserLessThanFunction( reinterpret_cast + ( &vtkGdcm4DSplitter::CompareOnSortTag)); // Anything like this, in GDCM2? + // s->UserOrdering((*i).second); + */ + +// //if (verbose) std::cout << "Out of SortOnTag" << std::endl; +// std::cout << "NO ordering performed :-( " << std::endl; +// } + + else if (SortOnUserFunction) + { + if (verbose) std::cout << "SortOnUserFunction" << std::endl; + s->SetUserLessThanFunction( UserCompareFunction ); + // Anything like this, in GDCM2? s->UserOrdering((*i).second); - if (verbose) std::cout << "Out of SortOnTag" << std::endl; + if (verbose) std::cout << "Out of SortOnUserFunction" << std::endl; } reader->SetCoherentFileList((*i).second); reader->Update(); /// \TODO : remove the following - //if (verbose) reader->GetOutput()->PrintSelf(std::cout, vtkIndent(2)); + if (verbose) { + // std::cout << "reader->GetOutput() :" << std::endl; + // reader->GetOutput()->PrintSelf(std::cout, vtkIndent(2)); + } ImageDataVector->push_back(reader->GetOutput() ); - - std::cout << std::endl; + + std::vector::iterator it; + //if (verbose) // JPR + + for(it=ImageDataVector->begin(); it!=ImageDataVector->end(); ++it) { + // std::cout << "-in vtkGdcm4DSplitter --------------------------" << std::endl; + // (*it)->PrintSelf(std::cout, vtkIndent(2)); + } + //std::cout << std::endl; } - //reader->Delete(); // \TODO : fix + reader->Delete(); // \TODO : fix s->Delete(); - f->Delete(); + // f->Delete(); delete l; return true; } + +/** + * \brief Load the gdcm::File* set, according to user's requierements + * returns a std::vector of gdcm::File* (gdcm::FileList : actually, a std::vector of gdcm::File*) + */ + + +GDCM_NAME_SPACE::FileList *vtkGdcm4DSplitter::getGdcmFileList() +{ + + GDCM_NAME_SPACE::File *f; + GDCM_NAME_SPACE::DirListType fileNames; + + // + // Fill fileNames with the user supplied file names (in any) + // ------------------------------------------------ + // + if (TypeDir == 0 ) // Nothing was set as input... + { + ///\TODO (?) Throw an exception "Set input Directory name(s) / file names before!" + std::cout << "Set input Directory name(s) / file names before!" << std::endl; + return false; + } + else if (TypeDir == 1) // A root directory name was set as input + { + GDCM_NAME_SPACE::DirList dirlist(DirName, Recursive); // NO recursive exploration + fileNames = dirlist.GetFilenames(); // all the file names + } + + else if (TypeDir == 2) // a std::vector of directory names was set as input + { + int nbDir = VectDirName.size(); + GDCM_NAME_SPACE::DirListType tmpFileNames; + for (int iDir=0; iDirSetMaxSizeLoadEntry(maxSize); + f->SetFileName( *it ); + if (f->Load()) + l->push_back(f); + else + std::cout << " Fail to load [" << *it << "]" << std::endl; + } + } + return l; +} + + + + void vtkGdcm4DSplitter::reorgXCoherentFileSetmap(GDCM_NAME_SPACE::XCoherentFileSetmap &xcm) + { + /* + the key of the 'XCoherentFileSetmap', is a std::string, used as if it was found in the Dicom header + Normaly(?), it's suitable for almost everything ... + ... but the 'Image Position Patient'. + We need to order the 'XCoherentFileSetmap' (NOT the content of each XCoherentFileSet!) according to the IPP, + using Jolinda Smith's algorithm. + (we use a subset of the one defined in gdcm::SerieHelper) + +*/ + + ELEM e; + std::vector vectElem; + +/* just to remember : + typedef struct + { + std::string strIPP; + double dist; + GDCM_NAME_SPACE::File *file; + } ELEM; +*/ + + bool Debug=true; + + for (GDCM_NAME_SPACE::XCoherentFileSetmap::iterator i = xcm.begin(); + i != xcm.end(); + ++i) + { + if (verbose) + std::cout << "--- xCoherentName = [" << (*i).first << "]" << std::endl; + + e.strIPP = (*i).first; + e.file = *(((*i).second)->begin()); // all the gdcm::File of a given xcm item *have* the same IPP; first one is enough + e.dist=0.0; + vectElem.push_back(e); + } + sortVectElem(&vectElem); + + // now, create the final std::map ! + // final_xcm + // xcm = final_xcm; + // check what we need to free ! + + int dist; + char c_dist[100]; + std::string str_dist; + int lgr=vectElem.size(); + GDCM_NAME_SPACE::XCoherentFileSetmap final_xcm; + for (int i2=0; i2 *fileList) +{ +//based on Jolinda Smith's algorithm +// +// NOTE : if you need to use Jolinda Smith's algorithm, get the one inside gdcm::SerieHelper +// this one is a light version. + +//Tags always use the same coordinate system, where "x" is left +//to right, "y" is posterior to anterior, and "z" is foot to head (RAH). + + //iop is calculated based on the file file + float cosines[6]; + double normal[3]; + double ipp[3]; + double dist; + double min = 0, max = 0; + bool first = true; + + //double ZSpacing; // useless here! // JPR + bool DirectOrder = true; // remove it! + + // ZSpacing = -1.0; // will be updated if process doesn't fail + + //std::multimap distmultimap; // JPR + std::multimap distmultimap; // JPR + + // Use a multimap to sort the distances from 0,0,0 + //for ( FileList::const_iterator // JPR + for ( std::vector::iterator // JPR + it = fileList->begin(); + it != fileList->end(); ++it ) + { + //gdcmDebugMacro("deal with " << (*it)->file->GetFileName() ); + if ( first ) + { + (*it).file->GetImageOrientationPatient( cosines ); + + // The "Image Orientation Patient" tag gives the direction cosines + // for the rows and columns for the three axes defined above. + // Typical axial slices will have a value 1/0/0/0/1/0: + // rows increase from left to right, + // columns increase from posterior to anterior. This is your everyday + // "looking up from the bottom of the head with the eyeballs up" image. + + // The "Image Position Patient" tag gives the coordinates of the first + // voxel in the image in the "RAH" coordinate system, relative to some + // origin. + + // First, calculate the slice normal from IOP : + + // You only have to do this once for all slices in the volume. Next, + // for each slice, calculate the distance along the slice normal + // using the IPP ("Image Position Patient") tag. + // ("dist" is initialized to zero before reading the first slice) : + + normal[0] = cosines[1]*cosines[5] - cosines[2]*cosines[4]; + normal[1] = cosines[2]*cosines[3] - cosines[0]*cosines[5]; + normal[2] = cosines[0]*cosines[4] - cosines[1]*cosines[3]; + + // For each slice (here : the first), calculate the distance along + // the slice normal using the IPP tag + + ipp[0] = (*it).file->GetXOrigin(); + ipp[1] = (*it).file->GetYOrigin(); + ipp[2] = (*it).file->GetZOrigin(); + + dist = 0; + for ( int i = 0; i < 3; ++i ) + { + dist += normal[i]*ipp[i]; + } + + //gdcmDebugMacro("dist : " << dist); + distmultimap.insert(std::pair(dist, *it)); + + max = min = dist; + first = false; + } + else + { + // Next, for each slice, calculate the distance along the slice normal + // using the IPP tag + ipp[0] = (*it).file->GetXOrigin(); + ipp[1] = (*it).file->GetYOrigin(); + ipp[2] = (*it).file->GetZOrigin(); + + dist = 0; + for ( int i = 0; i < 3; ++i ) + { + dist += normal[i]*ipp[i]; + } + + (*it).dist = dist; // JPR + + distmultimap.insert(std::pair(dist, *it)); + //gdcmDebugMacro("dist : " << dist); + min = (min < dist) ? min : dist; + max = (max > dist) ? max : dist; + } + } + + // gdcmDebugMacro("After parsing vector, nb of elements : " << fileList->size() ); + +/*Useless here. // JPR + + // Find out if min/max are coherent + if ( min == max ) + { + gdcmWarningMacro("Looks like all images have the exact same image position. " + << "No PositionPatientOrdering sort performed. " + << "No 'ZSpacing' calculated! "); + return false; + } +*/ + +/* Useless here, 'split' already done. // JPR + + // Check to see if image shares a common position + bool ok = true; + for (std::multimap::iterator it2 = distmultimap.begin(); + it2 != distmultimap.end(); + ++it2) + { + + gdcmDebugMacro("Check if image shares a common position : " << ((*it2).second).file->GetFileName() ); + + if (distmultimap.count((*it2).first) != 1) + { + gdcmWarningMacro("File: [" + << ((*it2).second->GetFileName()) + << "] : more than ONE file at distance: '" + << (*it2).first + << " (position is not unique!) " + << "No PositionPatientOrdering sort performed. " + << "No 'ZSpacing' calculated! "); + + ok = false; + } + } + if (!ok) + { + if (! DropDuplicatePositions) + return false; + } + +*/ + +// Now, we can calculate Z Spacing as the difference +// between the "dist" values for the first two slices. + +// The following (un)-commented out code is let here +// to be re-used by whomsoever is interested... + + //std::multimap::iterator it5 = distmultimap.begin(); + //double d1 = (*it5).first; + //it5++; + //double d2 = (*it5).first; + //ZSpacing = d1-d2; + //if (ZSpacing < 0.0) + // ZSpacing = - ZSpacing; + + fileList->clear(); // doesn't delete list elements, only nodes + +// Acording to user requierement, we sort direct order or reverse order. + if (DirectOrder) + { + for (std::multimap::iterator it3 = distmultimap.begin(); + it3 != distmultimap.end(); + ++it3) + { + fileList->push_back( (*it3).second ); + +/* useless here! // JPR + + if (DropDuplicatePositions) + { + // ImagePositionPatientOrdering wrong duplicates are found ??? + // --> fixed. See comment + + it3 = distmultimap.upper_bound((*it3).first); // skip all duplicates + // the upper_bound function increments the iterator to the next non-duplicate entry + // The for loop iteration also increments the iterator, which causes the code to skip every other image + // --> decrement the iterator after the upper_bound function call + it3--; + if (it3 == distmultimap.end() ) // if last image, stop iterate + break; + } +*/ + } + } + else // user asked for reverse order + { + std::multimap::const_iterator it4; + it4 = distmultimap.end(); + do + { + it4--; + fileList->push_back( (*it4).second ); + +/* useless here // JPR + + if (DropDuplicatePositions) // skip all duplicates + { + // lower_bound finds the next element that is + // less than or *equal to* the current value! + //it4 = distmultimap.lower_bound((*it4).first); + + // David Feng's fix + std::multimap::const_iterator itPrev = it4; + while (itPrev->first == it4->first) + --itPrev; + it4 = itPrev; + + if (it4 == distmultimap.begin() ) // if first image, stop iterate + break; + } +*/ + } while (it4 != distmultimap.begin() ); + } + + distmultimap.clear(); + + return true; + +}