comparison OrthancStone/Sources/Toolbox/DicomInstanceParameters.cpp @ 2076:990f396484b1

fix rendering of RT-DOSE with negative GridFrameOffsetVector
author Sebastien Jodogne <s.jodogne@gmail.com>
date Tue, 11 Jul 2023 15:58:16 +0200
parents d84bdcbd8bf1
children 07964689cb0b
comparison
equal deleted inserted replaced
2075:d84bdcbd8bf1 2076:990f396484b1
112 { 112 {
113 height_ = 0; 113 height_ = 0;
114 } 114 }
115 115
116 116
117 bool sliceThicknessPresent = true; 117 if (dicom.ParseDouble(sliceThickness_, Orthanc::DICOM_TAG_SLICE_THICKNESS))
118 if (!dicom.ParseDouble(sliceThickness_, Orthanc::DICOM_TAG_SLICE_THICKNESS)) 118 {
119 { 119 hasSliceThickness_ = true;
120 }
121 else
122 {
123 hasSliceThickness_ = false;
124
120 if (numberOfFrames_ > 1) 125 if (numberOfFrames_ > 1)
121 { 126 {
122 LOG(INFO) << "The (non-madatory) slice thickness information is missing in a multiframe image"; 127 LOG(INFO) << "The (non-madatory) slice thickness information is missing in a multiframe image";
123 } 128 }
124
125 sliceThickness_ = 100.0 * std::numeric_limits<double>::epsilon();
126 sliceThicknessPresent = false;
127 } 129 }
128 130
129 hasPixelSpacing_ = GeometryToolbox::GetPixelSpacing(pixelSpacingX_, pixelSpacingY_, dicom); 131 hasPixelSpacing_ = GeometryToolbox::GetPixelSpacing(pixelSpacingX_, pixelSpacingY_, dicom);
130 132
131 std::string position, orientation; 133 std::string position, orientation;
137 139
138 // Must be AFTER setting "numberOfFrames_" 140 // Must be AFTER setting "numberOfFrames_"
139 if (numberOfFrames_ > 1) 141 if (numberOfFrames_ > 1)
140 { 142 {
141 ExtractFrameOffsets(frameOffsets_, dicom, numberOfFrames_); 143 ExtractFrameOffsets(frameOffsets_, dicom, numberOfFrames_);
142
143 // if the slice thickness is unknown, we try to infer it from the sequence of grid frame offsets
144 // this only works if:
145 // - the first offset is 0.0 (case (a) of http://dicom.nema.org/medical/Dicom/2017c/output/chtml/part03/sect_C.8.8.3.2.html)
146 // - the offsets are all equal, to some small tolerance
147 // - the offsets is positive (increasing throughout the frames)
148 if (!sliceThicknessPresent)
149 {
150 if (frameOffsets_.size() >= 2)
151 {
152 double sliceThickness = std::abs(frameOffsets_[1] - frameOffsets_[0]);
153
154 if (sliceThickness > 0)
155 {
156 bool sameSized = true;
157
158 for (size_t i = 2; i < frameOffsets_.size(); ++i)
159 {
160 double currentThickness = std::abs(frameOffsets_[i] - frameOffsets_[i-1]);
161 if (!LinearAlgebra::IsNear(sliceThickness, currentThickness))
162 {
163 LOG(ERROR) << "Unable to extract slice thickness from GridFrameOffsetVector (3004,000C) (reason: varying spacing)";
164 sameSized = false;
165 break;
166 }
167 }
168
169 if (sameSized)
170 {
171 sliceThickness_ = sliceThickness;
172 LOG(INFO) << "SliceThickness was not specified in the Dicom but was inferred from GridFrameOffsetVector (3004,000C).";
173 }
174 }
175 }
176 else
177 {
178 LOG(ERROR) << "Unable to extract slice thickness from GridFrameOffsetVector (3004,000C) (reason: GridFrameOffsetVector not present or too small)";
179 }
180 }
181 }
182 else
183 {
184 frameOffsets_.resize(0);
185 } 144 }
186 145
187 if (sopClassUid_ == SopClassUid_RTDose) 146 if (sopClassUid_ == SopClassUid_RTDose)
188 { 147 {
189 static const Orthanc::DicomTag DICOM_TAG_DOSE_UNITS(0x3004, 0x0002); 148 static const Orthanc::DicomTag DICOM_TAG_DOSE_UNITS(0x3004, 0x0002);
268 227
269 if (!dicom.HasTag(Orthanc::DICOM_TAG_INSTANCE_NUMBER) || 228 if (!dicom.HasTag(Orthanc::DICOM_TAG_INSTANCE_NUMBER) ||
270 !dicom.GetValue(Orthanc::DICOM_TAG_INSTANCE_NUMBER).ParseInteger32(instanceNumber_)) 229 !dicom.GetValue(Orthanc::DICOM_TAG_INSTANCE_NUMBER).ParseInteger32(instanceNumber_))
271 { 230 {
272 instanceNumber_ = 0; 231 instanceNumber_ = 0;
232 }
233 }
234
235
236 double DicomInstanceParameters::GetSliceThickness() const
237 {
238 if (data_.hasSliceThickness_)
239 {
240 return data_.sliceThickness_;
241 }
242 else
243 {
244 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadSequenceOfCalls);
273 } 245 }
274 } 246 }
275 247
276 248
277 const Orthanc::DicomImageInformation& DicomInstanceParameters::GetImageInformation() const 249 const Orthanc::DicomImageInformation& DicomInstanceParameters::GetImageInformation() const
713 685
714 return (value * scaling + offset); 686 return (value * scaling + offset);
715 } 687 }
716 688
717 689
718 bool DicomInstanceParameters::ComputeRegularSpacing(double& spacing) const 690 bool DicomInstanceParameters::ComputeFrameOffsetsSpacing(double& spacing) const
719 { 691 {
720 if (data_.frameOffsets_.size() == 0) // Not a RT-DOSE 692 if (data_.frameOffsets_.size() == 0) // Not a RT-DOSE
721 { 693 {
722 return false; 694 return false;
723 } 695 }
726 spacing = 1; // Edge case: RT-DOSE with one single frame 698 spacing = 1; // Edge case: RT-DOSE with one single frame
727 return true; 699 return true;
728 } 700 }
729 else 701 else
730 { 702 {
731 assert(data_.frameOffsets_.size() == GetNumberOfFrames()); 703 static double THRESHOLD = 0.001;
732 704
733 spacing = std::abs(data_.frameOffsets_[1] - data_.frameOffsets_[0]); 705 if (data_.frameOffsets_.size() != GetNumberOfFrames())
706 {
707 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError);
708 }
709
710 double a = data_.frameOffsets_[1] - data_.frameOffsets_[0];
734 711
735 for (size_t i = 1; i + 1 < data_.frameOffsets_.size(); i++) 712 for (size_t i = 1; i + 1 < data_.frameOffsets_.size(); i++)
736 { 713 {
737 double s = std::abs(data_.frameOffsets_[i + 1] - data_.frameOffsets_[i]); 714 double b = data_.frameOffsets_[i + 1] - data_.frameOffsets_[i];
738 if (!LinearAlgebra::IsNear(spacing, s, 0.001)) 715 if (!LinearAlgebra::IsNear(a, b, THRESHOLD))
739 { 716 {
717 LOG(ERROR) << "Unable to extract slice thickness from GridFrameOffsetVector (3004,000C) (reason: varying spacing)";
740 return false; 718 return false;
741 } 719 }
720 }
721
722 spacing = std::abs(a);
723
724 if (HasSliceThickness() &&
725 !LinearAlgebra::IsNear(spacing, GetSliceThickness(), THRESHOLD))
726 {
727 LOG(WARNING) << "SliceThickness and GridFrameOffsetVector (3004,000C) do not match";
742 } 728 }
743 729
744 return true; 730 return true;
745 } 731 }
746 } 732 }
814 10.0 * physicalDeltaY->asDouble()); 800 10.0 * physicalDeltaY->asDouble());
815 } 801 }
816 } 802 }
817 } 803 }
818 } 804 }
805
806
807 CoordinateSystem3D DicomInstanceParameters::GetMultiFrameGeometry() const
808 {
809 if (data_.frameOffsets_.empty())
810 {
811 return data_.geometry_;
812 }
813 else
814 {
815 assert(data_.frameOffsets_.size() == data_.numberOfFrames_);
816
817 double lowest = data_.frameOffsets_[0];
818 for (size_t i = 1; i < data_.frameOffsets_.size(); i++)
819 {
820 lowest = std::min(lowest, data_.frameOffsets_[i]);
821 }
822
823 return CoordinateSystem3D(
824 data_.geometry_.GetOrigin() + lowest * data_.geometry_.GetNormal(),
825 data_.geometry_.GetAxisX(),
826 data_.geometry_.GetAxisY());
827 }
828 }
829
830
831 bool DicomInstanceParameters::IsReversedFrameOffsets() const
832 {
833 if (data_.frameOffsets_.size() < 2)
834 {
835 return false;
836 }
837 else
838 {
839 return (data_.frameOffsets_[0] > data_.frameOffsets_[1]);
840 }
841 }
819 } 842 }