From 3d6401cf7a30c2c7db019784b8baee874590dea7 Mon Sep 17 00:00:00 2001 From: Shreeraj Jadhav Date: Tue, 16 Jul 2024 00:41:25 -0400 Subject: [PATCH] ENH: Add template for input ImageType for itk2Dicom --- include/dcmqi/ConverterBase.h | 3 +- include/dcmqi/Itk2DicomConverter.h | 3 +- libsrc/Itk2DicomConverter.cpp | 63 +++++++++++++++++++----------- 3 files changed, 45 insertions(+), 24 deletions(-) diff --git a/include/dcmqi/ConverterBase.h b/include/dcmqi/ConverterBase.h index bb4acedf..1485ca2c 100644 --- a/include/dcmqi/ConverterBase.h +++ b/include/dcmqi/ConverterBase.h @@ -298,8 +298,9 @@ namespace dcmqi { } // AF: I could not quickly figure out how to template this function over image type - suggestions are welcomed! + template static vector > getSliceMapForSegmentation2DerivationImage(const vector dcmDatasets, - const ShortImageType::ConstPointer &labelImage) { + const ImageSourceType& labelImage) { // Find mapping from the segmentation slice number to the derivation image // Assume that orientation of the segmentation is the same as the source series unsigned numLabelSlices = labelImage->GetLargestPossibleRegion().GetSize()[2]; diff --git a/include/dcmqi/Itk2DicomConverter.h b/include/dcmqi/Itk2DicomConverter.h index 0d671744..246a32da 100644 --- a/include/dcmqi/Itk2DicomConverter.h +++ b/include/dcmqi/Itk2DicomConverter.h @@ -62,8 +62,9 @@ namespace dcmqi { * added in the SharedFunctionalGroupsSequence without any geometry checks. * @return A pointer to the resulting DICOM Segmentation object. */ + template static DcmDataset* itkimage2dcmSegmentation(vector dcmDatasets, - vector segmentations, + vector> segmentations, const string &metaData, bool skipEmptySlices=true, bool useLabelIDAsSegmentNumber=false, diff --git a/libsrc/Itk2DicomConverter.cpp b/libsrc/Itk2DicomConverter.cpp index 15f3e400..b395668e 100644 --- a/libsrc/Itk2DicomConverter.cpp +++ b/libsrc/Itk2DicomConverter.cpp @@ -14,6 +14,7 @@ #endif #include +#include namespace dcmqi { @@ -25,14 +26,15 @@ namespace dcmqi { // ------------------------------------------------------------------------------------- + template DcmDataset* Itk2DicomConverter::itkimage2dcmSegmentation(vector dcmDatasets, - vector segmentations, + vector> segmentations, const string &metaData, bool skipEmptySlices, bool useLabelIDAsSegmentNumber, bool referencesGeometryCheck) { - ShortImageType::SizeType inputSize = segmentations[0]->GetBufferedRegion().GetSize(); + auto inputSize = segmentations[0]->GetBufferedRegion().GetSize(); JSONSegmentationMetaInformationHandler metaInfo(metaData.c_str()); metaInfo.read(); @@ -78,7 +80,7 @@ namespace dcmqi { // Shared FGs: PlaneOrientationPatientSequence { - ShortImageType::DirectionType labelDirMatrix = segmentations[0]->GetDirection(); + auto labelDirMatrix = segmentations[0]->GetDirection(); FGPlaneOrientationPatient *planor = FGPlaneOrientationPatient::createMinimal( @@ -96,7 +98,7 @@ namespace dcmqi { { FGPixelMeasures *pixmsr = new FGPixelMeasures(); - ShortImageType::SpacingType labelSpacing = segmentations[0]->GetSpacing(); + auto labelSpacing = segmentations[0]->GetSpacing(); ostringstream spacingSStream; spacingSStream << scientific << labelSpacing[1] << "\\" << labelSpacing[0]; CHECK_COND(pixmsr->setPixelSpacing(spacingSStream.str().c_str())); @@ -191,14 +193,15 @@ namespace dcmqi { // note that labels are the same in the input and output image produced // by this filter, see https://itk.org/Doxygen/html/classitk_1_1LabelImageToLabelMapFilter.html - LabelToLabelMapFilterType::Pointer l2lm = LabelToLabelMapFilterType::New(); + using LabelToLabelMapFilterType2 = itk::LabelImageToLabelMapFilter; + typename LabelToLabelMapFilterType2::Pointer l2lm = LabelToLabelMapFilterType2::New(); l2lm->SetInput(segmentations[segFileNumber]); l2lm->Update(); - typedef LabelToLabelMapFilterType::OutputImageType::LabelObjectType LabelType; - typedef itk::LabelStatisticsImageFilter LabelStatisticsType; + typedef typename LabelToLabelMapFilterType2::OutputImageType::LabelObjectType LabelType; + typedef itk::LabelStatisticsImageFilter LabelStatisticsType; - LabelStatisticsType::Pointer labelStats = LabelStatisticsType::New(); + typename LabelStatisticsType::Pointer labelStats = LabelStatisticsType::New(); cout << "Found " << l2lm->GetOutput()->GetNumberOfLabelObjects() << " label(s)" << endl; labelStats->SetInput(segmentations[segFileNumber]); @@ -208,21 +211,21 @@ namespace dcmqi { bool cropSegmentsBBox = false; if(cropSegmentsBBox){ cout << "WARNING: Crop operation enabled - WIP" << endl; - typedef itk::BinaryThresholdImageFilter ThresholdType; - ThresholdType::Pointer thresh = ThresholdType::New(); + typedef itk::BinaryThresholdImageFilter ThresholdType; + typename ThresholdType::Pointer thresh = ThresholdType::New(); thresh->SetInput(segmentations[segFileNumber]); thresh->SetLowerThreshold(1); thresh->SetLowerThreshold(100); thresh->SetInsideValue(1); thresh->Update(); - LabelStatisticsType::Pointer threshLabelStats = LabelStatisticsType::New(); + typename LabelStatisticsType::Pointer threshLabelStats = LabelStatisticsType::New(); threshLabelStats->SetInput(thresh->GetOutput()); threshLabelStats->SetLabelInput(thresh->GetOutput()); threshLabelStats->Update(); - LabelStatisticsType::BoundingBoxType threshBbox = threshLabelStats->GetBoundingBox(1); + auto threshBbox = threshLabelStats->GetBoundingBox(1); /* cout << "OVerall bounding box: " << threshBbox[0] << ", " << threshBbox[1] << threshBbox[2] << ", " << threshBbox[3] @@ -242,7 +245,7 @@ namespace dcmqi { cout << "Processing label " << label << endl; - LabelStatisticsType::BoundingBoxType bbox = labelStats->GetBoundingBox(label); + auto bbox = labelStats->GetBoundingBox(label); unsigned firstSlice, lastSlice; //bool skipEmptySlices = true; // TODO: what to do with that line? //bool skipEmptySlices = false; // TODO: what to do with that line? @@ -354,15 +357,15 @@ namespace dcmqi { // PerFrame FG: PlanePositionSequence { - ShortImageType::PointType sliceOriginPoint; - ShortImageType::IndexType sliceOriginIndex; + typename ImageSourceType::PointType sliceOriginPoint; + typename ImageSourceType::IndexType sliceOriginIndex; sliceOriginIndex.Fill(0); sliceOriginIndex[2] = sliceNumber; segmentations[segFileNumber]->TransformIndexToPhysicalPoint(sliceOriginIndex, sliceOriginPoint); ostringstream pppSStream; if(sliceNumber>0){ - ShortImageType::PointType prevOrigin; - ShortImageType::IndexType prevIndex; + typename ImageSourceType::PointType prevOrigin; + typename ImageSourceType::IndexType prevIndex; prevIndex.Fill(0); prevIndex[2] = sliceNumber-1; segmentations[segFileNumber]->TransformIndexToPhysicalPoint(prevIndex, prevOrigin); @@ -375,9 +378,9 @@ namespace dcmqi { /* Add frame that references this segment */ { - ShortImageType::RegionType sliceRegion; - ShortImageType::IndexType sliceIndex; - ShortImageType::SizeType sliceSize; + typename ImageSourceType::RegionType sliceRegion; + typename ImageSourceType::IndexType sliceIndex; + typename ImageSourceType::SizeType sliceSize; sliceIndex[0] = 0; sliceIndex[1] = 0; @@ -391,11 +394,11 @@ namespace dcmqi { sliceRegion.SetSize(sliceSize); unsigned framePixelCnt = 0; - itk::ImageRegionConstIteratorWithIndex sliceIterator(segmentations[segFileNumber], sliceRegion); + itk::ImageRegionConstIteratorWithIndex sliceIterator(segmentations[segFileNumber], sliceRegion); for(sliceIterator.GoToBegin();!sliceIterator.IsAtEnd();++sliceIterator,++framePixelCnt){ if(sliceIterator.Get() == label){ frameData[framePixelCnt] = 1; - ShortImageType::IndexType idx = sliceIterator.GetIndex(); + auto idx = sliceIterator.GetIndex(); } else frameData[framePixelCnt] = 0; } @@ -685,4 +688,20 @@ namespace dcmqi { return true; } + template DcmDataset* Itk2DicomConverter::itkimage2dcmSegmentation( + vector dcmDatasets, + vector segmentations, + const string& metaData, + bool skipEmptySlices, + bool useLabelIDAsSegmentNumber, + bool referencesGeometryCheck); + + using VectorImageAdapter = itk::VectorImageToImageAdaptor; + template DcmDataset* Itk2DicomConverter::itkimage2dcmSegmentation( + vector dcmDatasets, + vector segmentations, + const string& metaData, + bool skipEmptySlices, + bool useLabelIDAsSegmentNumber, + bool referencesGeometryCheck); }