Skip to content

Commit

Permalink
be tolerant of new segment numbers
Browse files Browse the repository at this point in the history
  • Loading branch information
@spalte
spalte committed Feb 1, 2024
1 parent d30215c commit 4143036
Show file tree
Hide file tree
Showing 3 changed files with 71 additions and 65 deletions.
19 changes: 10 additions & 9 deletions include/dcmqi/OverlapUtil.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -27,6 +27,7 @@
#include "dcmtk/ofstd/oftypes.h"
#include "dcmtk/ofstd/ofvector.h"
#include <set>
#include <map>

class DcmSegmentation;

Expand Down Expand Up @@ -91,7 +92,7 @@ class OverlapUtil

/// Lists frames for each segment where segment with index i is represented by the vector at index i,
/// and index 0 is unused. I.e. index i is segment number, value is vector of physical frame numbers.
typedef OFVector<OFVector<Uint32>> FramesForSegment;
typedef std::map<Uint32, std::set<Uint32>> FramesForSegment;

// Set of segments present on each frame.
typedef OFVector<std::set<Uint32>> SegmentsForFrame;
Expand Down Expand Up @@ -119,7 +120,7 @@ class OverlapUtil

/// Matrix of N x N segment numbers, where N is the number of segments.
/// Value is 1 at x,y if x and y overlap, 0 if they don't overlap, and -1 if not initialized.
typedef OFVector<OFVector<Sint8>> OverlapMatrix;
typedef std::map<Uint32, std::map<Uint32, Sint8>> OverlapMatrix;

/// Group of non-overlapping segments (each represented by its segment number)
typedef OFVector<OFVector<Uint32>> SegmentGroups;
Expand All @@ -137,7 +138,7 @@ class OverlapUtil
, m_frameNumber(f)
{
}
/// Segment number as used in DICOM segmentation object (1-n)
/// Segment number as used in DICOM segmentation object
Uint16 m_segmentNumber;
/// Logical frame number (number of frame in DistinctFramePositions vector)
Uint16 m_frameNumber;
Expand Down Expand Up @@ -190,15 +191,15 @@ class OverlapUtil
OFCondition getSegmentsByPosition(SegmentsByPosition& result);

/** Get phyiscal frames for a specific segment by its segment number
* @param segmentNumber Segment number to get frames for (1..n)
* @param frames Resulting vector of physical frame numbers (first frame is frame 0)
* @param segmentNumber Segment number for which to get frames
* @param frames Resulting vector of physical frame numbers
* @return EC_Normal if successful, error otherwise
*/
OFCondition getFramesForSegment(const Uint32 segmentNumber, OFVector<Uint32>& frames);
OFCondition getFramesForSegment(const Uint32 segmentNumber, std::set<Uint32>& frames);

/** Get the all the segments present on a specified frame
* @param frameNumber The frame number for which to get segments
* @param segments Resulting set of segment numbers (1..n)
* @param segments Resulting set of segment numbers
* @return EC_Normal if successful, error otherwise
*/
OFCondition getSegmentsForFrame(const Uint32 frameNumber, std::set<Uint32>& segments);
Expand Down Expand Up @@ -247,12 +248,12 @@ class OverlapUtil
/** Get the list of segments within a label map segmentation frame.
* Does not cache the result. To be used exclusively on a label map frame.
* @param frameNumber The frame number for which to get labels
* @param segments The resulting set of segments on the frame (1..n)
* @param segments The resulting set of segments on the frame
* @return EC_Normal if successful, error otherwise
*/
OFCondition getSegmentsForLabelMapFrame(const Uint32 frameNumber, std::set<Uint32>& segments);

/** Get the segment number (1..n) for a binary or fractional segmentation frame.
/** Get the segment number for a binary or fractional segmentation frame.
* Does not cache the result. To be used exclusively on a binary or fractional
* segmentation frame.
* @param frameNumber The frame number for which to get labels
Expand Down
15 changes: 8 additions & 7 deletions libsrc/Dicom2ItkConverter.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -124,32 +124,34 @@ itk::SmartPointer<ShortImageType> Dicom2ItkConverter::nextResult()
{
// Iterate over frames for this segment, and copy the data into the ITK image.
// Afterwards, the ITK image will have the complete data belonging to that segment
OverlapUtil::FramesForSegment::value_type framesForSegment;
OverlapUtil::FramesForSegment::mapped_type framesForSegment;
m_overlapUtil.getFramesForSegment(*segNum, framesForSegment);
for (size_t frameIndex = 0; frameIndex < framesForSegment.size(); frameIndex++)
for (std::set<Uint32>::iterator it = framesForSegment.begin();
it != framesForSegment.end();
it++)
{
// Copy the data from the frame into the ITK image
ShortImageType::PointType frameOriginPoint;
ShortImageType::IndexType frameOriginIndex;
result = getITKImageOrigin(framesForSegment[frameIndex], frameOriginPoint);
result = getITKImageOrigin(*it, frameOriginPoint);
if (result.bad())
{
cerr << "ERROR: Failed to get origin for frame " << framesForSegment[frameIndex] << " of segment "
cerr << "ERROR: Failed to get origin for frame " << *it << " of segment "
<< *segNum << endl;
m_groupIterator = m_segmentGroups.end();
return nullptr;
}
if (!itkImage->TransformPhysicalPointToIndex(frameOriginPoint, frameOriginIndex))
{
cerr << "ERROR: Frame " << framesForSegment[frameIndex] << " origin " << frameOriginPoint
cerr << "ERROR: Frame " << *it << " origin " << frameOriginPoint
<< " is outside image geometry!" << frameOriginIndex << endl;
cerr << "Image size: " << itkImage->GetBufferedRegion().GetSize() << endl;
m_groupIterator = m_segmentGroups.end();
return nullptr;
}
// Handling differs depending on whether the segmentation is binary or fractional
// (we have to unpack binary frames before copying them into the ITK image)
const DcmIODTypes::Frame* rawFrame = m_segDoc->getFrame(framesForSegment[frameIndex]);
const DcmIODTypes::Frame* rawFrame = m_segDoc->getFrame(*it);
const DcmIODTypes::Frame* unpackedFrame = NULL;
if (m_segDoc->getSegmentationType() == DcmSegTypes::ST_BINARY)
{
Expand Down Expand Up @@ -372,7 +374,6 @@ OFCondition Dicom2ItkConverter::getNonOverlappingSegmentGroups(const bool mergeS
cout << "Will not merge segments: Splitting segments into " << segmentGroups.size() << " groups" << endl;
}
return result;
;
}

// -------------------------------------------------------------------------------------
Expand Down
102 changes: 53 additions & 49 deletions libsrc/OverlapUtil.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -52,7 +52,7 @@ OverlapUtil::OverlapUtil()
, m_segmentsForFrame()
, m_logicalFramePositions()
, m_segmentsByPosition()
, m_segmentOverlapMatrix(0)
, m_segmentOverlapMatrix()
, m_nonOverlappingSegments()
, m_seg()
{
Expand Down Expand Up @@ -99,18 +99,13 @@ OFCondition OverlapUtil::getFramesByPosition(DistinctFramePositions& result)
return cond;
}

OFCondition OverlapUtil::getFramesForSegment(const Uint32 segmentNumber, OFVector<Uint32>& frames)
OFCondition OverlapUtil::getFramesForSegment(const Uint32 segmentNumber, std::set<Uint32>& frames)
{
if (!m_seg)
{
DCMSEG_ERROR("getFramesForSegment(): No segmentation object set");
return EC_IllegalCall;
}
if ((segmentNumber == 0) || (segmentNumber > m_seg->getNumberOfSegments()))
{
DCMSEG_ERROR("getFramesForSegment(): Segment number " << segmentNumber << " is out of range");
return EC_IllegalParameter;
}
if (m_framesForSegment.empty())
{
FGInterface& fg = m_seg->getFunctionalGroups();
Expand All @@ -121,9 +116,6 @@ OFCondition OverlapUtil::getFramesForSegment(const Uint32 segmentNumber, OFVecto
return EC_IllegalParameter;
}
Uint32 numFrames = static_cast<Uint32>(m_seg->getNumberOfFrames());
// use sets for performance and deduplication when adding frames to each segment index
OFVector<std::set<Uint32>> frameSetsForSegment;
frameSetsForSegment.resize(m_seg->getNumberOfSegments());
// Use the helper function getSegmentsForFrame() to pivot the data
// from segments on each frame, to the list of frames with contain each segment
for (Uint32 f = 0; f < numFrames; f++)
Expand All @@ -137,7 +129,7 @@ OFCondition OverlapUtil::getFramesForSegment(const Uint32 segmentNumber, OFVecto
itr != segments.end();
itr++)
{
frameSetsForSegment[(*itr) - 1].insert(f);
m_framesForSegment[*itr].insert(f);
}
}
else
Expand All @@ -146,19 +138,8 @@ OFCondition OverlapUtil::getFramesForSegment(const Uint32 segmentNumber, OFVecto
return result;
}
}
// copy the sets back into the OFVectors in m_framesForSegment
m_framesForSegment.resize(frameSetsForSegment.size());
for (Uint32 s = 0; s < frameSetsForSegment.size(); s++) {
m_framesForSegment[s].reserve(frameSetsForSegment[s].size());
for (std::set<Uint32>::iterator itr = frameSetsForSegment[s].begin();
itr != frameSetsForSegment[s].end();
itr++)
{
m_framesForSegment[s].push_back(*itr);
}
}
}
frames = m_framesForSegment[segmentNumber - 1];
frames = m_framesForSegment[segmentNumber];
return EC_Normal;
}

Expand Down Expand Up @@ -337,7 +318,6 @@ OFCondition OverlapUtil::getSegmentsByPosition(SegmentsByPosition& result)
{
return cond;
}
size_t numSegments = m_seg->getNumberOfSegments();
if (m_logicalFramePositions.empty())
{
cond = getFramesByPosition(m_logicalFramePositions);
Expand Down Expand Up @@ -438,19 +418,21 @@ OFCondition OverlapUtil::getNonOverlappingSegments(SegmentGroups& segmentGroups)
// grouped segments. If not, add them to the same group. If yes, create a new group
// and add them there.
m_nonOverlappingSegments.push_back(OFVector<Uint32>());
for (size_t i = 0; i < m_segmentOverlapMatrix.size(); ++i)
for (OverlapMatrix::iterator matrixItr = m_segmentOverlapMatrix.begin();
matrixItr != m_segmentOverlapMatrix.end();
matrixItr++)
{
// Loop over all groups and check whether the current segment overlaps with any of them
OFBool overlaps = OFFalse;
for (size_t j = 0; j < m_nonOverlappingSegments.size(); ++j)
{
// Loop over all segments in the current group
for (OFVector<Uint32>::iterator it = m_nonOverlappingSegments[j].begin();
for (SegmentGroups::value_type::iterator it = m_nonOverlappingSegments[j].begin();
it != m_nonOverlappingSegments[j].end();
++it)
{
// Check if the current segment overlaps with the segment in the current group
if (m_segmentOverlapMatrix[i][(*it) - 1] == 1)
if (m_segmentOverlapMatrix[matrixItr->first][*it] == 1)
{
overlaps = OFTrue;
break;
Expand All @@ -459,15 +441,15 @@ OFCondition OverlapUtil::getNonOverlappingSegments(SegmentGroups& segmentGroups)
if (!overlaps)
{
// Add segment to current group
m_nonOverlappingSegments[j].push_back(i + 1);
m_nonOverlappingSegments[j].push_back(matrixItr->first);
break;
}
}
if (overlaps)
{
// Create new group and add segment to it
m_nonOverlappingSegments.push_back(OFVector<Uint32>());
m_nonOverlappingSegments.back().push_back(i + 1);
m_nonOverlappingSegments.back().push_back(matrixItr->first);
}
}
}
Expand Down Expand Up @@ -495,7 +477,7 @@ void OverlapUtil::printSegmentsByPosition(OFStringStream& ss)
for (size_t i = 0; i < m_segmentsByPosition.size(); ++i)
{
OFStringStream tempSS;
for (std::set<SegNumAndFrameNum>::iterator it = m_segmentsByPosition[i].begin();
for (SegmentsByPosition::value_type::iterator it = m_segmentsByPosition[i].begin();
it != m_segmentsByPosition[i].end();
++it)
{
Expand Down Expand Up @@ -530,7 +512,7 @@ void OverlapUtil::printNonOverlappingSegments(OFStringStream& ss)
for (size_t i = 0; i < m_nonOverlappingSegments.size(); ++i)
{
ss << "Group #" << i << ": ";
for (OFVector<Uint32>::iterator it = m_nonOverlappingSegments[i].begin();
for (SegmentGroups::value_type::iterator it = m_nonOverlappingSegments[i].begin();
it != m_nonOverlappingSegments[i].end();
++it)
{
Expand Down Expand Up @@ -560,7 +542,7 @@ OFCondition OverlapUtil::getSegmentsForLabelMapFrame(const Uint32 frameNumber, s
for (size_t n = 0; n < f_data->length; ++n)
{
Uint8 segmentNumber = f_data->pixData[n];
if ((segmentNumber == 0) || (segmentNumber > m_seg->getNumberOfSegments()))
if (segmentNumber > m_seg->getNumberOfSegments())
{
DCMSEG_ERROR("getSegmentsForFrame(): Segment number " << segmentNumber << " is out of range");
return EC_IllegalParameter;
Expand Down Expand Up @@ -605,32 +587,50 @@ OFCondition OverlapUtil::getSegmentForFrame(const Uint32 frameNumber, Uint32& se

OFCondition OverlapUtil::buildOverlapMatrix()
{
// Make 2 dimensional array matrix of Sint8 type for (segment numbers) X (segment numbers).
// Initialize with -1 (not checked yet)
m_segmentOverlapMatrix.clear();
m_segmentOverlapMatrix.resize(m_seg->getNumberOfSegments());
for (size_t i = 0; i < m_segmentOverlapMatrix.size(); ++i)
// TODO This corresponds to the list of segments in a binary or fractional
// segmentation IOD, we still need to get the actual list if this is a labelmap
const Uint32 numSegments = m_seg->getNumberOfSegments();
std::set<Uint32> segmentNumbers;
for (Uint32 i = 1; i <= numSegments; ++i)
{
m_segmentOverlapMatrix[i].resize(m_seg->getNumberOfSegments(), -1);
segmentNumbers.insert(i);
}

// Make 2 dimensional matrix of Sint8 type for (segment numbers) X (segment numbers).
// Initialize with -1 (not checked yet)
// Diagonal is always 0 (segment does not interfere/overlap with itself)
for (size_t i = 0; i < m_segmentOverlapMatrix.size(); ++i)
for (std::set<Uint32>::iterator it1 = segmentNumbers.begin();
it1 != segmentNumbers.end();
it1++)
{
m_segmentOverlapMatrix[i][i] = 0;
for (std::set<Uint32>::iterator it2 = segmentNumbers.begin();
it2 != segmentNumbers.end();
it2++)
{
if (*it1 == *it2)
{
m_segmentOverlapMatrix[*it1][*it2] = 0;
}
else
{
m_segmentOverlapMatrix[*it1][*it2] = -1;
}
}
}

// Go through all logical frame positions, and compare all segments at each position
size_t index1, index2;
index1 = index2 = 0;
for (size_t i = 0; i < m_segmentsByPosition.size(); ++i)
{
DCMSEG_DEBUG("getOverlappingSegments(): Comparing segments at logical frame position " << i);
// Compare all segments at this position
for (std::set<SegNumAndFrameNum>::iterator it = m_segmentsByPosition[i].begin();
for (SegmentsByPosition::value_type::iterator it = m_segmentsByPosition[i].begin();
it != m_segmentsByPosition[i].end();
++it)
{
index1++;
for (std::set<SegNumAndFrameNum>::iterator it2 = m_segmentsByPosition[i].begin();
for (SegmentsByPosition::value_type::iterator it2 = m_segmentsByPosition[i].begin();
it2 != m_segmentsByPosition[i].end();
++it2)
{
Expand All @@ -643,7 +643,7 @@ OFCondition OverlapUtil::buildOverlapMatrix()
{
// Check if we already have found an overlap on another logical frame, and if so, skip
Sint8 existing_result
= m_segmentOverlapMatrix[(*it).m_segmentNumber - 1][(*it2).m_segmentNumber - 1];
= m_segmentOverlapMatrix[(*it).m_segmentNumber][(*it2).m_segmentNumber];
if (existing_result == 1)
{
DCMSEG_DEBUG("getOverlappingSegments(): Skipping frame comparison on pos #"
Expand Down Expand Up @@ -671,21 +671,25 @@ OFCondition OverlapUtil::buildOverlapMatrix()
}

// Enter result into overlap matrix
m_segmentOverlapMatrix[(*it).m_segmentNumber - 1][(*it2).m_segmentNumber - 1] = overlap ? 1 : 0;
m_segmentOverlapMatrix[(*it2).m_segmentNumber - 1][(*it).m_segmentNumber - 1] = overlap ? 1 : 0;
m_segmentOverlapMatrix[(*it).m_segmentNumber][(*it2).m_segmentNumber] = overlap ? 1 : 0;
m_segmentOverlapMatrix[(*it2).m_segmentNumber][(*it).m_segmentNumber] = overlap ? 1 : 0;
}
}
}
}
// Since we don't compare all segments (since not all are showing up together on a single logical frame),
// we set all remaining entries that are still not initialized (-1) to 0 (no overlap)
for (size_t i = 0; i < m_segmentOverlapMatrix.size(); ++i)
for (std::set<Uint32>::iterator it1 = segmentNumbers.begin();
it1 != segmentNumbers.end();
it1++)
{
for (size_t j = 0; j < m_segmentOverlapMatrix[i].size(); ++j)
for (std::set<Uint32>::iterator it2 = segmentNumbers.begin();
it2 != segmentNumbers.end();
it2++)
{
if (m_segmentOverlapMatrix[i][j] == -1)
if (m_segmentOverlapMatrix[*it1][*it2] == -1)
{
m_segmentOverlapMatrix[i][j] = 0;
m_segmentOverlapMatrix[*it1][*it2] = 0;
}
}
}
Expand Down

0 comments on commit 4143036

Please sign in to comment.