From c112f1284e728cbaf463e04857c489b3e51577d1 Mon Sep 17 00:00:00 2001 From: Andrey Fedorov Date: Fri, 16 Sep 2016 14:49:23 -0400 Subject: [PATCH 1/2] BUG: fixes to improve handling of single frame SEG --- include/ConverterBase.h | 68 +++++++++++++++++++++++------------- libsrc/ImageSEGConverter.cpp | 4 +++ 2 files changed, 47 insertions(+), 25 deletions(-) diff --git a/include/ConverterBase.h b/include/ConverterBase.h index e19f73f6..2a52ef40 100644 --- a/include/ConverterBase.h +++ b/include/ConverterBase.h @@ -76,12 +76,17 @@ namespace dcmqi { vnl_vector sliceDirection = vnl_cross_3d(rowDirection, colDirection); sliceDirection.normalize(); + cout << "Row direction: " << rowDirection << endl; + cout << "Col direction: " << colDirection << endl; + for(int i=0;i<3;i++){ dir[i][0] = rowDirection[i]; dir[i][1] = colDirection[i]; dir[i][2] = sliceDirection[i]; } + cout << "Z direction: " << sliceDirection << endl; + return 0; } @@ -100,6 +105,8 @@ namespace dcmqi { map frame2overlap; double minDistance; + sliceSpacing = 0; + unsigned numFrames = fgInterface.getNumberOfFrames(); /* Framesorter is to be moved to DCMTK at some point @@ -190,35 +197,45 @@ namespace dcmqi { } } - // sort all unique distances, this will be used to check consistency of - // slice spacing, and also to locate the slice position from ImagePositionPatient - // later when we read the segments - sort(originDistances.begin(), originDistances.end()); - - sliceSpacing = fabs(originDistances[0]-originDistances[1]); + // it IS possible to have a segmentation object containing just one frame! + if(numFrames>1){ + // WARNING: this should be improved further. Spacing should be calculated for + // consecutive frames of the individual segment. Right now, all frames are considered + // indiscriminately. Question is whether it should be computed at all, considering we do + // not have any information about whether the 2 frames are adjacent or not, so perhaps we should + // always rely on the declared spacing, and not even try to compute it? + // TODO: discuss this with the QIICR team! + + // sort all unique distances, this will be used to check consistency of + // slice spacing, and also to locate the slice position from ImagePositionPatient + // later when we read the segments + sort(originDistances.begin(), originDistances.end()); + + sliceSpacing = fabs(originDistances[0]-originDistances[1]); + + for(int i=1;i 0.001){ + cerr << "WARNING: Inter-slice distance " << originDistances[i] << + " difference exceeded threshold: " << delta << endl; + } + } - for(int i=1;i 0.001){ - cerr << "WARNING: Inter-slice distance " << originDistances[i] << - " difference exceeded threshold: " << delta << endl; + sliceExtent = fabs(originDistances[0]-originDistances[originDistances.size()-1]); + unsigned overlappingFramesCnt = 0; + for(map::const_iterator it=frame2overlap.begin(); + it!=frame2overlap.end();++it){ + if(it->second>1) + overlappingFramesCnt++; } - } - sliceExtent = fabs(originDistances[0]-originDistances[originDistances.size()-1]); - unsigned overlappingFramesCnt = 0; - for(map::const_iterator it=frame2overlap.begin(); - it!=frame2overlap.end();++it){ - if(it->second>1) - overlappingFramesCnt++; + cout << "Total frames: " << numFrames << endl; + cout << "Total frames with unique IPP: " << originDistances.size() << endl; + cout << "Total overlapping frames: " << overlappingFramesCnt << endl; + cout << "Origin: " << imageOrigin << endl; } - cout << "Total frames: " << numFrames << endl; - cout << "Total frames with unique IPP: " << originDistances.size() << endl; - cout << "Total overlapping frames: " << overlappingFramesCnt << endl; - cout << "Origin: " << imageOrigin << endl; - return 0; } @@ -241,7 +258,8 @@ namespace dcmqi { } else if(pixelMeasures->getSliceThickness(spacingFloat,0).good() && spacingFloat != 0){ // SliceThickness can be carried forward from the source images, and may not be what we need // As an example, this ePAD example has 1.25 carried from CT, but true computed thickness is 1! - cerr << "WARNING: SliceThickness is present and is " << spacingFloat << ". NOT using it!" << endl; + cerr << "WARNING: SliceThickness is present and is " << spacingFloat << ". using it!" << endl; + spacing[2] = spacingFloat; } return 0; } diff --git a/libsrc/ImageSEGConverter.cpp b/libsrc/ImageSEGConverter.cpp index 5cba429c..b867b61c 100644 --- a/libsrc/ImageSEGConverter.cpp +++ b/libsrc/ImageSEGConverter.cpp @@ -459,6 +459,10 @@ namespace dcmqi { } const double tolerance = 1e-5; + if(!imageSpacing[2] && !computedSliceSpacing){ + cerr << "FATAL ERROR: No sufficient information to derive slice spacing! Unable to interpret the data." << endl; + throw -1; + } if(!imageSpacing[2]){ imageSpacing[2] = computedSliceSpacing; } From 26e806b6fceda031e1345b1e24f0474026a6d2af Mon Sep 17 00:00:00 2001 From: Andrey Fedorov Date: Sat, 15 Oct 2016 11:16:32 -0400 Subject: [PATCH 2/2] BUG: do not consider computed slice thickness for SEG As discussed in #89, in the general case, we cannot rely on the slice thickness computed from the positions of the individual frames, since frames may not be contiguous. When source image references are available, and the source images are also available, additional checks could be possible, but not in dcmqi. --- libsrc/ImageSEGConverter.cpp | 10 +--------- 1 file changed, 1 insertion(+), 9 deletions(-) diff --git a/libsrc/ImageSEGConverter.cpp b/libsrc/ImageSEGConverter.cpp index b867b61c..d3d284a9 100644 --- a/libsrc/ImageSEGConverter.cpp +++ b/libsrc/ImageSEGConverter.cpp @@ -458,18 +458,10 @@ namespace dcmqi { throw -1; } - const double tolerance = 1e-5; - if(!imageSpacing[2] && !computedSliceSpacing){ + if(!imageSpacing[2]){ cerr << "FATAL ERROR: No sufficient information to derive slice spacing! Unable to interpret the data." << endl; throw -1; } - if(!imageSpacing[2]){ - imageSpacing[2] = computedSliceSpacing; - } - else if(fabs(imageSpacing[2]-computedSliceSpacing)>tolerance){ - cerr << "WARNING: Declared slice spacing is significantly different from the one declared in DICOM!" << - " Declared = " << imageSpacing[2] << " Computed = " << computedSliceSpacing << endl; - } // Region size ImageType::SizeType imageSize;