diff --git a/Modules/DICOMReader/mitkDICOMTagBasedSorter.cpp b/Modules/DICOMReader/mitkDICOMTagBasedSorter.cpp index 185814d2e6..5475550871 100644 --- a/Modules/DICOMReader/mitkDICOMTagBasedSorter.cpp +++ b/Modules/DICOMReader/mitkDICOMTagBasedSorter.cpp @@ -1,540 +1,540 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #include "mitkDICOMTagBasedSorter.h" #include #include mitk::DICOMTagBasedSorter::CutDecimalPlaces ::CutDecimalPlaces(unsigned int precision) :m_Precision(precision) { } mitk::DICOMTagBasedSorter::CutDecimalPlaces ::CutDecimalPlaces(const CutDecimalPlaces& other) :m_Precision(other.m_Precision) { } std::string mitk::DICOMTagBasedSorter::CutDecimalPlaces ::operator()(const std::string& input) const { // be a bit tolerant for tags such as image orientation orienatation, let only the first few digits matter (http://bugs.mitk.org/show_bug.cgi?id=12263) // iterate all fields, convert each to a number, cut this number as configured, then return a concatenated string with all cut-off numbers - static std::ostringstream resultString; + std::ostringstream resultString; resultString.str(std::string()); resultString.clear(); resultString.setf(std::ios::fixed, std::ios::floatfield); resultString.precision(m_Precision); - static std::stringstream ss(input); + std::stringstream ss(input); ss.str(input); ss.clear(); - static std::string item; + std::string item; + double number(0); + std::istringstream converter(item); while (std::getline(ss, item, '\\')) { - static std::istringstream converter(item); converter.str(item); converter.clear(); - static double number(0); if (converter >> number && converter.eof()) { // converted to double resultString << number; } else { // did not convert to double resultString << item; // just paste the unmodified string } if (!ss.eof()) { resultString << "\\"; } } return resultString.str(); } mitk::DICOMTagBasedSorter::TagValueProcessor* mitk::DICOMTagBasedSorter::CutDecimalPlaces ::Clone() const { return new CutDecimalPlaces(*this); } unsigned int mitk::DICOMTagBasedSorter::CutDecimalPlaces ::GetPrecision() const { return m_Precision; } mitk::DICOMTagBasedSorter ::DICOMTagBasedSorter() :DICOMDatasetSorter() ,m_StrictSorting(true) { } mitk::DICOMTagBasedSorter ::~DICOMTagBasedSorter() { for(TagValueProcessorMap::iterator ti = m_TagValueProcessor.begin(); ti != m_TagValueProcessor.end(); ++ti) { delete ti->second; } } mitk::DICOMTagBasedSorter ::DICOMTagBasedSorter(const DICOMTagBasedSorter& other ) :DICOMDatasetSorter(other) ,m_DistinguishingTags( other.m_DistinguishingTags ) ,m_SortCriterion( other.m_SortCriterion ) ,m_StrictSorting( other.m_StrictSorting ) { for(TagValueProcessorMap::const_iterator ti = other.m_TagValueProcessor.begin(); ti != other.m_TagValueProcessor.end(); ++ti) { m_TagValueProcessor[ti->first] = ti->second->Clone(); } } mitk::DICOMTagBasedSorter& mitk::DICOMTagBasedSorter ::operator=(const DICOMTagBasedSorter& other) { if (this != &other) { DICOMDatasetSorter::operator=(other); m_DistinguishingTags = other.m_DistinguishingTags; m_SortCriterion = other.m_SortCriterion; m_StrictSorting = other.m_StrictSorting; for(TagValueProcessorMap::const_iterator ti = other.m_TagValueProcessor.begin(); ti != other.m_TagValueProcessor.end(); ++ti) { m_TagValueProcessor[ti->first] = ti->second->Clone(); } } return *this; } bool mitk::DICOMTagBasedSorter ::operator==(const DICOMDatasetSorter& other) const { if (const DICOMTagBasedSorter* otherSelf = dynamic_cast(&other)) { if (this->m_StrictSorting != otherSelf->m_StrictSorting) return false; bool allTagsPresentAndEqual(true); if (this->m_DistinguishingTags.size() != otherSelf->m_DistinguishingTags.size()) return false; for (DICOMTagList::const_iterator myTag = this->m_DistinguishingTags.begin(); myTag != this->m_DistinguishingTags.end(); ++myTag) { allTagsPresentAndEqual &= (std::find( otherSelf->m_DistinguishingTags.begin(), otherSelf->m_DistinguishingTags.end(), *myTag ) != otherSelf->m_DistinguishingTags.end()); // other contains this tags // since size is equal, we don't need to check the inverse } if (!allTagsPresentAndEqual) return false; if (this->m_SortCriterion.IsNotNull() && otherSelf->m_SortCriterion.IsNotNull()) { return *(this->m_SortCriterion) == *(otherSelf->m_SortCriterion); } else { return this->m_SortCriterion.IsNull() && otherSelf->m_SortCriterion.IsNull(); } } else { return false; } } void mitk::DICOMTagBasedSorter ::PrintConfiguration(std::ostream& os, const std::string& indent) const { os << indent << "Tag based sorting:" << std::endl; for (DICOMTagList::const_iterator tagIter = m_DistinguishingTags.begin(); tagIter != m_DistinguishingTags.end(); ++tagIter) { os << indent << " Split on "; tagIter->Print(os); os << std::endl; } DICOMSortCriterion::ConstPointer crit = m_SortCriterion.GetPointer(); while (crit.IsNotNull()) { os << indent << " Sort by "; crit->Print(os); os << std::endl; crit = crit->GetSecondaryCriterion(); } } void mitk::DICOMTagBasedSorter ::SetStrictSorting(bool strict) { m_StrictSorting = strict; } bool mitk::DICOMTagBasedSorter ::GetStrictSorting() const { return m_StrictSorting; } mitk::DICOMTagList mitk::DICOMTagBasedSorter ::GetTagsOfInterest() { DICOMTagList allTags = m_DistinguishingTags; DICOMTagList sortingRelevantTags = m_SortCriterion->GetAllTagsOfInterest(); allTags.insert( allTags.end(), sortingRelevantTags.begin(), sortingRelevantTags.end() ); // append return allTags; } mitk::DICOMTagList mitk::DICOMTagBasedSorter ::GetDistinguishingTags() const { return m_DistinguishingTags; } const mitk::DICOMTagBasedSorter::TagValueProcessor* mitk::DICOMTagBasedSorter ::GetTagValueProcessorForDistinguishingTag(const DICOMTag& tag) const { TagValueProcessorMap::const_iterator loc = m_TagValueProcessor.find(tag); if (loc != m_TagValueProcessor.end()) { return loc->second; } else { return NULL; } } void mitk::DICOMTagBasedSorter ::AddDistinguishingTag( const DICOMTag& tag, TagValueProcessor* tagValueProcessor ) { m_DistinguishingTags.push_back(tag); m_TagValueProcessor[tag] = tagValueProcessor; } void mitk::DICOMTagBasedSorter ::SetSortCriterion( DICOMSortCriterion::ConstPointer criterion ) { m_SortCriterion = criterion; } mitk::DICOMSortCriterion::ConstPointer mitk::DICOMTagBasedSorter ::GetSortCriterion() const { return m_SortCriterion; } void mitk::DICOMTagBasedSorter ::Sort() { // 1. split // 2. sort each group GroupIDToListType groups = this->SplitInputGroups(); GroupIDToListType& sortedGroups = this->SortGroups( groups ); // 3. define output this->SetNumberOfOutputs(sortedGroups.size()); unsigned int outputIndex(0); for (GroupIDToListType::iterator groupIter = sortedGroups.begin(); groupIter != sortedGroups.end(); ++outputIndex, ++groupIter) { this->SetOutput(outputIndex, groupIter->second); } } std::string mitk::DICOMTagBasedSorter ::BuildGroupID( DICOMDatasetAccess* dataset ) { // just concatenate all tag values assert(dataset); std::stringstream groupID; groupID << "g"; for (DICOMTagList::iterator tagIter = m_DistinguishingTags.begin(); tagIter != m_DistinguishingTags.end(); ++tagIter) { groupID << tagIter->GetGroup() << tagIter->GetElement(); // make group/element part of the id to cover empty tags std::string rawTagValue = dataset->GetTagValueAsString(*tagIter); std::string processedTagValue; if ( m_TagValueProcessor[*tagIter] != NULL ) { processedTagValue = (*m_TagValueProcessor[*tagIter])(rawTagValue); } else { processedTagValue = rawTagValue; } groupID << processedTagValue; } // shorten ID? return groupID.str(); } mitk::DICOMTagBasedSorter::GroupIDToListType mitk::DICOMTagBasedSorter ::SplitInputGroups() { DICOMDatasetList input = GetInput(); // copy GroupIDToListType listForGroupID; for (DICOMDatasetList::iterator dsIter = input.begin(); dsIter != input.end(); ++dsIter) { DICOMDatasetAccess* dataset = *dsIter; assert(dataset); std::string groupID = this->BuildGroupID( dataset ); MITK_DEBUG << "Group ID for for " << dataset->GetFilenameIfAvailable() << ": " << groupID; listForGroupID[groupID].push_back(dataset); } MITK_DEBUG << "After tag based splitting: " << listForGroupID.size() << " groups"; return listForGroupID; } mitk::DICOMTagBasedSorter::GroupIDToListType& mitk::DICOMTagBasedSorter ::SortGroups(GroupIDToListType& groups) { if (m_SortCriterion.IsNotNull()) { /* Three steps here: 1. sort within each group - this may result in orders such as 1 2 3 4 6 7 8 10 12 13 14 2. create new groups by enforcing consecutive order within each group - resorts above example like 1 2 3 4 ; 6 7 8 ; 10 ; 12 13 14 3. sort all of the groups (not WITHIN each group) by their first frame - if earlier "distinguish" steps created groups like 6 7 8 ; 1 2 3 4 ; 10, then this step would sort them like 1 2 3 4 ; 6 7 8 ; 10 */ // Step 1: sort within the groups // for each output // sort by all configured tags, use secondary tags when equal or empty // make configurable: // - sorting order (ascending, descending) // - sort numerically // - ... ? unsigned int groupIndex(0); for (GroupIDToListType::iterator gIter = groups.begin(); gIter != groups.end(); ++groupIndex, ++gIter) { DICOMDatasetList& dsList = gIter->second; MITK_DEBUG << " --------------------------------------------------------------------------------"; MITK_DEBUG << " DICOMTagBasedSorter before sorting group : " << groupIndex; for (DICOMDatasetList::iterator oi = dsList.begin(); oi != dsList.end(); ++oi) { MITK_DEBUG << " INPUT : " << (*oi)->GetFilenameIfAvailable(); } std::sort( dsList.begin(), dsList.end(), ParameterizedDatasetSort( m_SortCriterion ) ); MITK_DEBUG << " --------------------------------------------------------------------------------"; MITK_DEBUG << " DICOMTagBasedSorter after sorting group : " << groupIndex; for (DICOMDatasetList::iterator oi = dsList.begin(); oi != dsList.end(); ++oi) { MITK_DEBUG << " OUTPUT : " << (*oi)->GetFilenameIfAvailable(); } MITK_DEBUG << " --------------------------------------------------------------------------------"; } GroupIDToListType consecutiveGroups; if (m_StrictSorting) { // Step 2: create new groups by enforcing consecutive order within each group unsigned int groupIndex(0); for (GroupIDToListType::iterator gIter = groups.begin(); gIter != groups.end(); ++gIter) { std::stringstream groupKey; groupKey << std::setfill('0') << std::setw(6) << groupIndex++; DICOMDatasetList& dsList = gIter->second; DICOMDatasetAccess* previousDS(NULL); unsigned int dsIndex(0); double constantDistance(0.0); bool constantDistanceInitialized(false); for (DICOMDatasetList::iterator dataset = dsList.begin(); dataset != dsList.end(); ++dsIndex, ++dataset) { if (dsIndex >0) // ignore the first dataset, we cannot check any distances yet.. { // for the second and every following dataset: // let the sorting criterion calculate a "distance" // if the distance is not 1, split off a new group! double currentDistance = m_SortCriterion->NumericDistance(previousDS, *dataset); if (constantDistanceInitialized) { if (fabs(currentDistance - constantDistance) < constantDistance * 0.01) // ok, deviation of up to 1% of distance is tolerated { // nothing to do, just ok } else if (currentDistance < mitk::eps) // close enough to 0 { // no numeric comparison possible? // rare case(?), just accept } else { // split! this is done by simply creating a new group (key) groupKey.str(std::string()); groupKey.clear(); groupKey << std::setfill('0') << std::setw(6) << groupIndex++; } } else { // second slice: learn about the expected distance! // heuristic: if distance is an integer, we check for a special case: // if the distance is integer and not 1/-1, then we assume // a missing slice right after the first slice // ==> split off slices // in all other cases: second dataset at this position, no need to split already, we are still learning about the images if ((currentDistance - (int)currentDistance == 0.0) && fabs(currentDistance) != 1.0) // exact comparison. An integer should not be expressed as 1.000000000000000000000000001! { groupKey.str(std::string()); groupKey.clear(); groupKey << std::setfill('0') << std::setw(6) << groupIndex++; } constantDistance = currentDistance; constantDistanceInitialized = true; } } consecutiveGroups[groupKey.str()].push_back(*dataset); previousDS = *dataset; } } } else { consecutiveGroups = groups; } // Step 3: sort all of the groups (not WITHIN each group) by their first frame /* build a list-1 of datasets with the first dataset one of each group sort this list-1 build a new result list-2: - iterate list-1, for each dataset - find the group that contains this dataset - add this group as the next element to list-2 return list-2 as the sorted output */ DICOMDatasetList firstSlices; for (GroupIDToListType::iterator gIter = consecutiveGroups.begin(); gIter != consecutiveGroups.end(); ++gIter) { assert(!gIter->second.empty()); firstSlices.push_back(gIter->second.front()); } std::sort( firstSlices.begin(), firstSlices.end(), ParameterizedDatasetSort( m_SortCriterion ) ); GroupIDToListType sortedResultBlocks; unsigned int groupKeyValue(0); for (DICOMDatasetList::iterator firstSlice = firstSlices.begin(); firstSlice != firstSlices.end(); ++firstSlice) { for (GroupIDToListType::iterator gIter = consecutiveGroups.begin(); gIter != consecutiveGroups.end(); ++groupKeyValue, ++gIter) { if (gIter->second.front() == *firstSlice) { std::stringstream groupKey; groupKey << std::setfill('0') << std::setw(6) << groupKeyValue; // try more than 999,999 groups and you are doomed (your application already is) sortedResultBlocks[groupKey.str()] = gIter->second; } } } groups = sortedResultBlocks; } return groups; } mitk::DICOMTagBasedSorter::ParameterizedDatasetSort ::ParameterizedDatasetSort(DICOMSortCriterion::ConstPointer criterion) :m_SortCriterion(criterion) { } bool mitk::DICOMTagBasedSorter::ParameterizedDatasetSort ::operator() (const mitk::DICOMDatasetAccess* left, const mitk::DICOMDatasetAccess* right) { assert(left); assert(right); assert(m_SortCriterion.IsNotNull()); return m_SortCriterion->IsLeftBeforeRight(left, right); } diff --git a/Modules/DICOMReader/mitkITKDICOMSeriesReaderHelper.cpp b/Modules/DICOMReader/mitkITKDICOMSeriesReaderHelper.cpp index a227eb355d..50cb1c8d22 100644 --- a/Modules/DICOMReader/mitkITKDICOMSeriesReaderHelper.cpp +++ b/Modules/DICOMReader/mitkITKDICOMSeriesReaderHelper.cpp @@ -1,219 +1,219 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ //#define MBILOG_ENABLE_DEBUG #include "mitkITKDICOMSeriesReaderHelper.h" #include "mitkITKDICOMSeriesReaderHelper.txx" #define switch3DCase(IOType, T) \ case IOType: return LoadDICOMByITK< T >(filenames, correctTilt, tiltInfo, io); bool mitk::ITKDICOMSeriesReaderHelper ::CanHandleFile(const std::string& filename) { MITK_DEBUG << "ITKDICOMSeriesReaderHelper::CanHandleFile " << filename; - static itk::GDCMImageIO::Pointer tester = itk::GDCMImageIO::New(); + itk::GDCMImageIO::Pointer tester = itk::GDCMImageIO::New(); if ( tester->CanReadFile(filename.c_str()) ) { tester->SetFileName( filename.c_str() ); tester->ReadImageInformation(); std::string numberOfFrames; if (tester->GetValueFromTag("0028|0008", numberOfFrames)) { std::istringstream converter(numberOfFrames); int i; if (converter >> i) { MITK_DEBUG << "Number of Frames for " << filename << ": " << numberOfFrames; return (i <= 1); // cannot handle multi-frame } else { return true; // we assume single-frame } } else { MITK_DEBUG << "No Number of Frames tag for " << filename; // friendly old single-frame file return true; } } MITK_DEBUG << "GDCMImageIO found: No DICOM in " << filename; // does not seem to be DICOM return false; } mitk::Image::Pointer mitk::ITKDICOMSeriesReaderHelper ::Load( const StringContainer& filenames, bool correctTilt, const GantryTiltInformation& tiltInfo ) { if( filenames.empty() ) { MITK_DEBUG << "Calling LoadDicomSeries with empty filename string container. Probably invalid application logic."; return NULL; // this is not actually an error but the result is very simple } typedef itk::GDCMImageIO DcmIoType; DcmIoType::Pointer io = DcmIoType::New(); try { if (io->CanReadFile(filenames.front().c_str())) { io->SetFileName(filenames.front().c_str()); io->ReadImageInformation(); if (io->GetPixelType() == itk::ImageIOBase::SCALAR) { switch (io->GetComponentType()) { switch3DCase(DcmIoType::UCHAR, unsigned char) switch3DCase(DcmIoType::CHAR, char) switch3DCase(DcmIoType::USHORT, unsigned short) switch3DCase(DcmIoType::SHORT, short) switch3DCase(DcmIoType::UINT, unsigned int) switch3DCase(DcmIoType::INT, int) switch3DCase(DcmIoType::ULONG, long unsigned int) switch3DCase(DcmIoType::LONG, long int) switch3DCase(DcmIoType::FLOAT, float) switch3DCase(DcmIoType::DOUBLE, double) default: MITK_ERROR << "Found unsupported DICOM scalar pixel type: (enum value) " << io->GetComponentType(); } } else if (io->GetPixelType() == itk::ImageIOBase::RGB) { switch (io->GetComponentType()) { switch3DCase(DcmIoType::UCHAR, itk::RGBPixel) switch3DCase(DcmIoType::CHAR, itk::RGBPixel) switch3DCase(DcmIoType::USHORT, itk::RGBPixel) switch3DCase(DcmIoType::SHORT, itk::RGBPixel) switch3DCase(DcmIoType::UINT, itk::RGBPixel) switch3DCase(DcmIoType::INT, itk::RGBPixel) switch3DCase(DcmIoType::ULONG, itk::RGBPixel) switch3DCase(DcmIoType::LONG, itk::RGBPixel) switch3DCase(DcmIoType::FLOAT, itk::RGBPixel) switch3DCase(DcmIoType::DOUBLE, itk::RGBPixel) default: MITK_ERROR << "Found unsupported DICOM scalar pixel type: (enum value) " << io->GetComponentType(); } } MITK_ERROR << "Unsupported DICOM pixel type"; return NULL; } } catch(itk::MemoryAllocationError& e) { MITK_ERROR << "Out of memory. Cannot load DICOM series: " << e.what(); } catch(std::exception& e) { MITK_ERROR << "Error encountered when loading DICOM series:" << e.what(); } catch(...) { MITK_ERROR << "Unspecified error encountered when loading DICOM series."; } return NULL; } #define switch3DnTCase(IOType, T) \ case IOType: return LoadDICOMByITK3DnT< T >(filenamesLists, correctTilt, tiltInfo, io); mitk::Image::Pointer mitk::ITKDICOMSeriesReaderHelper ::Load3DnT( const StringContainerList& filenamesLists, bool correctTilt, const GantryTiltInformation& tiltInfo ) { if( filenamesLists.empty() || filenamesLists.front().empty() ) { MITK_DEBUG << "Calling LoadDicomSeries with empty filename string container. Probably invalid application logic."; return NULL; // this is not actually an error but the result is very simple } typedef itk::GDCMImageIO DcmIoType; DcmIoType::Pointer io = DcmIoType::New(); try { if (io->CanReadFile(filenamesLists.front().front().c_str())) { io->SetFileName(filenamesLists.front().front().c_str()); io->ReadImageInformation(); if (io->GetPixelType() == itk::ImageIOBase::SCALAR) { switch (io->GetComponentType()) { switch3DnTCase(DcmIoType::UCHAR, unsigned char) switch3DnTCase(DcmIoType::CHAR, char) switch3DnTCase(DcmIoType::USHORT, unsigned short) switch3DnTCase(DcmIoType::SHORT, short) switch3DnTCase(DcmIoType::UINT, unsigned int) switch3DnTCase(DcmIoType::INT, int) switch3DnTCase(DcmIoType::ULONG, long unsigned int) switch3DnTCase(DcmIoType::LONG, long int) switch3DnTCase(DcmIoType::FLOAT, float) switch3DnTCase(DcmIoType::DOUBLE, double) default: MITK_ERROR << "Found unsupported DICOM scalar pixel type: (enum value) " << io->GetComponentType(); } } else if (io->GetPixelType() == itk::ImageIOBase::RGB) { switch (io->GetComponentType()) { switch3DnTCase(DcmIoType::UCHAR, itk::RGBPixel) switch3DnTCase(DcmIoType::CHAR, itk::RGBPixel) switch3DnTCase(DcmIoType::USHORT, itk::RGBPixel) switch3DnTCase(DcmIoType::SHORT, itk::RGBPixel) switch3DnTCase(DcmIoType::UINT, itk::RGBPixel) switch3DnTCase(DcmIoType::INT, itk::RGBPixel) switch3DnTCase(DcmIoType::ULONG, itk::RGBPixel) switch3DnTCase(DcmIoType::LONG, itk::RGBPixel) switch3DnTCase(DcmIoType::FLOAT, itk::RGBPixel) switch3DnTCase(DcmIoType::DOUBLE, itk::RGBPixel) default: MITK_ERROR << "Found unsupported DICOM scalar pixel type: (enum value) " << io->GetComponentType(); } } MITK_ERROR << "Unsupported DICOM pixel type"; return NULL; } } catch(itk::MemoryAllocationError& e) { MITK_ERROR << "Out of memory. Cannot load DICOM series: " << e.what(); } catch(std::exception& e) { MITK_ERROR << "Error encountered when loading DICOM series:" << e.what(); } catch(...) { MITK_ERROR << "Unspecified error encountered when loading DICOM series."; } return NULL; } diff --git a/Modules/DICOMReader/mitkNormalDirectionConsistencySorter.cpp b/Modules/DICOMReader/mitkNormalDirectionConsistencySorter.cpp index d152d03a6b..ccc16a9180 100644 --- a/Modules/DICOMReader/mitkNormalDirectionConsistencySorter.cpp +++ b/Modules/DICOMReader/mitkNormalDirectionConsistencySorter.cpp @@ -1,173 +1,171 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ //#define MBILOG_ENABLE_DEBUG #include "mitkNormalDirectionConsistencySorter.h" #include mitk::NormalDirectionConsistencySorter ::NormalDirectionConsistencySorter() :DICOMDatasetSorter() { } mitk::NormalDirectionConsistencySorter ::NormalDirectionConsistencySorter(const NormalDirectionConsistencySorter& other ) :DICOMDatasetSorter(other) { } mitk::NormalDirectionConsistencySorter ::~NormalDirectionConsistencySorter() { } void mitk::NormalDirectionConsistencySorter ::PrintConfiguration(std::ostream& os, const std::string& indent) const { os << indent << "NormalDirectionConsistencySorter" << std::endl; } mitk::NormalDirectionConsistencySorter& mitk::NormalDirectionConsistencySorter ::operator=(const NormalDirectionConsistencySorter& other) { if (this != &other) { DICOMDatasetSorter::operator=(other); } return *this; } bool mitk::NormalDirectionConsistencySorter ::operator==(const DICOMDatasetSorter& other) const { return dynamic_cast(&other) != NULL; } mitk::DICOMTagList mitk::NormalDirectionConsistencySorter ::GetTagsOfInterest() { DICOMTagList tags; tags.push_back( DICOMTag(0x0020, 0x0032) ); // ImagePositionPatient tags.push_back( DICOMTag(0x0020, 0x0037) ); // ImageOrientationPatient return tags; } void mitk::NormalDirectionConsistencySorter ::Sort() { DICOMDatasetList datasets = GetInput(); if (datasets.size() > 1) { // at some point in the code, there is the expectation that // the direction of the slice normals is the same as the direction between // first and last slice origin. We need to make this sure here, because // we want to feed the files into itk::ImageSeriesReader with the consistent // setting of ReverseOrderOff. static const DICOMTag tagImagePositionPatient = DICOMTag(0x0020,0x0032); // Image Position (Patient) static const DICOMTag tagImageOrientation = DICOMTag(0x0020, 0x0037); // Image Orientation DICOMDatasetAccess* firstDS = datasets.front(); DICOMDatasetAccess* lastDS = datasets.back(); // make sure here that the direction from slice to slice is the direction of // image normals... std::string imageOrientationString = firstDS->GetTagValueAsString( tagImageOrientation ); std::string imagePositionPatientFirst = firstDS->GetTagValueAsString( tagImagePositionPatient ); std::string imagePositionPatientLast = lastDS->GetTagValueAsString( tagImagePositionPatient ); - static Vector3D right; right.Fill(0.0); - static Vector3D up; up.Fill(0.0); - static bool hasOrientation(false); + Vector3D right; right.Fill(0.0); + Vector3D up; up.Fill(0.0); + bool hasOrientation(false); DICOMStringToOrientationVectors( imageOrientationString, right, up, hasOrientation ); - static Point3D firstOrigin; firstOrigin.Fill(0.0f); - static bool firstHasOrigin(false); + Point3D firstOrigin; firstOrigin.Fill(0.0f); + bool firstHasOrigin(false); firstOrigin = DICOMStringToPoint3D( imagePositionPatientFirst, firstHasOrigin ); - static Point3D lastOrigin; lastOrigin.Fill(0.0f); - static bool lastHasOrigin(false); + Point3D lastOrigin; lastOrigin.Fill(0.0f); + bool lastHasOrigin(false); lastOrigin = DICOMStringToPoint3D( imagePositionPatientLast, lastHasOrigin ); - static Vector3D normal; + Vector3D normal; normal[0] = right[1] * up[2] - right[2] * up[1]; normal[1] = right[2] * up[0] - right[0] * up[2]; normal[2] = right[0] * up[1] - right[1] * up[0]; normal.Normalize(); - static Vector3D directionOfSlices; + Vector3D directionOfSlices; directionOfSlices = lastOrigin - firstOrigin; directionOfSlices.Normalize(); - static double projection = 0.0; - projection = 0.0; - projection = normal * directionOfSlices; + double projection = normal * directionOfSlices; MITK_DEBUG << "Making sense of \norientation '" << imageOrientationString << "'\nfirst position '" << imagePositionPatientFirst << "'\nlast position '" << imagePositionPatientLast << "'"; MITK_DEBUG << "Normal: " << normal; MITK_DEBUG << "Direction of slices: " << directionOfSlices; MITK_DEBUG << "Projection of direction onto slice normal: " << projection; if ( projection < 0.0 ) { MITK_DEBUG << "Need to reverse filenames"; std::reverse( datasets.begin(), datasets.end() ); m_TiltInfo = GantryTiltInformation::MakeFromTagValues( imagePositionPatientLast, imagePositionPatientFirst, imageOrientationString, datasets.size() - 1 ); } else { m_TiltInfo = GantryTiltInformation::MakeFromTagValues( imagePositionPatientFirst, imagePositionPatientLast, imageOrientationString, datasets.size() - 1 ); } } else // just ONE dataset, do not forget to reset tilt information { m_TiltInfo = GantryTiltInformation(); // empty info } this->SetNumberOfOutputs(1); this->SetOutput(0, datasets); } mitk::GantryTiltInformation mitk::NormalDirectionConsistencySorter ::GetTiltInformation() const { return m_TiltInfo; } diff --git a/Modules/DICOMReader/mitkSortByImagePositionPatient.cpp b/Modules/DICOMReader/mitkSortByImagePositionPatient.cpp index a3b141811e..d228f7e86c 100644 --- a/Modules/DICOMReader/mitkSortByImagePositionPatient.cpp +++ b/Modules/DICOMReader/mitkSortByImagePositionPatient.cpp @@ -1,170 +1,168 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #include "mitkSortByImagePositionPatient.h" #include "mitkDICOMTag.h" mitk::SortByImagePositionPatient ::SortByImagePositionPatient(DICOMSortCriterion::Pointer secondaryCriterion) :DICOMSortCriterion(secondaryCriterion) { } mitk::SortByImagePositionPatient ::~SortByImagePositionPatient() { } mitk::SortByImagePositionPatient ::SortByImagePositionPatient(const SortByImagePositionPatient& other ) :DICOMSortCriterion(other) { } mitk::SortByImagePositionPatient& mitk::SortByImagePositionPatient ::operator=(const SortByImagePositionPatient& other) { if (this != &other) { DICOMSortCriterion::operator=(other); } return *this; } bool mitk::SortByImagePositionPatient ::operator==(const DICOMSortCriterion& other) const { return dynamic_cast(&other) != NULL; // same class } void mitk::SortByImagePositionPatient ::Print(std::ostream& os) const { os << "(0020,0032) Image Position (Patient) along normal of (0020,0037) Image Orientation (Patient)"; } mitk::DICOMTagList mitk::SortByImagePositionPatient ::GetTagsOfInterest() const { DICOMTagList tags; tags.push_back( DICOMTag(0x0020, 0x0032) ); // ImagePositionPatient tags.push_back( DICOMTag(0x0020, 0x0037) ); // ImageOrientationPatient return tags; } bool mitk::SortByImagePositionPatient ::IsLeftBeforeRight(const mitk::DICOMDatasetAccess* left, const mitk::DICOMDatasetAccess* right) const { bool possible(false); double distance = InternalNumericDistance(left, right, possible); // returns 0.0 if not possible if (possible) { return distance > 0.0; } else { return this->NextLevelIsLeftBeforeRight(left, right); } } double mitk::SortByImagePositionPatient ::InternalNumericDistance(const mitk::DICOMDatasetAccess* left, const mitk::DICOMDatasetAccess* right, bool& possible) const { // sort by distance to world origin, assuming (almost) equal orientation static const DICOMTag tagImagePositionPatient = DICOMTag(0x0020,0x0032); // Image Position (Patient) static const DICOMTag tagImageOrientation = DICOMTag(0x0020, 0x0037); // Image Orientation - static Vector3D leftRight; leftRight.Fill(0.0); - static Vector3D leftUp; leftUp.Fill(0.0); - static bool leftHasOrientation(false); + Vector3D leftRight; leftRight.Fill(0.0); + Vector3D leftUp; leftUp.Fill(0.0); + bool leftHasOrientation(false); DICOMStringToOrientationVectors( left->GetTagValueAsString( tagImageOrientation ), leftRight, leftUp, leftHasOrientation ); - static Vector3D rightRight; rightRight.Fill(0.0); - static Vector3D rightUp; rightUp.Fill(0.0); - static bool rightHasOrientation(false); + Vector3D rightRight; rightRight.Fill(0.0); + Vector3D rightUp; rightUp.Fill(0.0); + bool rightHasOrientation(false); DICOMStringToOrientationVectors( right->GetTagValueAsString( tagImageOrientation ), rightRight, rightUp, rightHasOrientation ); - static Point3D leftOrigin; leftOrigin.Fill(0.0f); - static bool leftHasOrigin(false); + Point3D leftOrigin; leftOrigin.Fill(0.0f); + bool leftHasOrigin(false); leftOrigin = DICOMStringToPoint3D( left->GetTagValueAsString( tagImagePositionPatient ), leftHasOrigin ); - static Point3D rightOrigin; rightOrigin.Fill(0.0f); - static bool rightHasOrigin(false); + Point3D rightOrigin; rightOrigin.Fill(0.0f); + bool rightHasOrigin(false); rightOrigin = DICOMStringToPoint3D( right->GetTagValueAsString( tagImagePositionPatient ), rightHasOrigin ); // we tolerate very small differences in image orientation, since we got to know about // acquisitions where these values change across a single series (7th decimal digit) // (http://bugs.mitk.org/show_bug.cgi?id=12263) // still, we want to check if our assumption of 'almost equal' orientations is valid for (unsigned int dim = 0; dim < 3; ++dim) { if ( fabs(leftRight[dim] - rightRight[dim]) > 0.0001 || fabs(leftUp[dim] - rightUp[dim]) > 0.0001) { MITK_ERROR << "Dicom images have different orientations."; throw std::logic_error("Dicom images have different orientations. Call GetSeries() first to separate images."); } } - static Vector3D normal; + Vector3D normal; normal[0] = leftRight[1] * leftUp[2] - leftRight[2] * leftUp[1]; normal[1] = leftRight[2] * leftUp[0] - leftRight[0] * leftUp[2]; normal[2] = leftRight[0] * leftUp[1] - leftRight[1] * leftUp[0]; - static double leftDistance = 0.0; - static double rightDistance = 0.0; - leftDistance = 0.0; - rightDistance = 0.0; + double leftDistance = 0.0; + double rightDistance = 0.0; // this computes the distance from world origin (0,0,0) ALONG THE NORMAL of the image planes for (unsigned int dim = 0; dim < 3; ++dim) { leftDistance += normal[dim] * leftOrigin[dim]; rightDistance += normal[dim] * rightOrigin[dim]; } // if we can sort by just comparing the distance, we do exactly that if ( fabs(leftDistance - rightDistance) >= mitk::eps) { possible = true; // default: compare position return rightDistance - leftDistance; // if (left < right> ==> diff > 0 } else { possible = false; return 0.0; } } double mitk::SortByImagePositionPatient ::NumericDistance(const mitk::DICOMDatasetAccess* from, const mitk::DICOMDatasetAccess* to) const { bool possible(false); double retVal = InternalNumericDistance(from, to, possible); // returns 0.0 if not possible return possible ? retVal : 0.0; }