Mercurial > hg > orthanc-stone

comparison Framework/Toolbox/DicomInstanceParameters.cpp @ 746:d716bfb3e07c

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
reorganization
author Sebastien Jodogne <s.jodogne@gmail.com>
date Wed, 22 May 2019 12:48:57 +0200
parents
children f7c236894c1a
comparison
equal deleted inserted replaced
745:c44c1d2d3598 746:d716bfb3e07c
1 /**
2 * Stone of Orthanc
3 * Copyright (C) 2012-2016 Sebastien Jodogne, Medical Physics
4 * Department, University Hospital of Liege, Belgium
5 * Copyright (C) 2017-2019 Osimis S.A., Belgium
6 *
7 * This program is free software: you can redistribute it and/or
8 * modify it under the terms of the GNU Affero General Public License
9 * as published by the Free Software Foundation, either version 3 of
10 * the License, or (at your option) any later version.
11 *
12 * This program is distributed in the hope that it will be useful, but
13 * WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
15 * Affero General Public License for more details.
16 *
17 * You should have received a copy of the GNU Affero General Public License
18 * along with this program. If not, see <http://www.gnu.org/licenses/>.
19 **/
20
21
22 #include "DicomInstanceParameters.h"
23
24 #include "../Scene2D/ColorTextureSceneLayer.h"
25 #include "../Scene2D/FloatTextureSceneLayer.h"
26 #include "../Toolbox/GeometryToolbox.h"
27
28 #include <Core/Images/Image.h>
29 #include <Core/Images/ImageProcessing.h>
30 #include <Core/Logging.h>
31 #include <Core/OrthancException.h>
32 #include <Core/Toolbox.h>
33
34
35 namespace OrthancStone
36 {
37 void DicomInstanceParameters::Data::ComputeDoseOffsets(const Orthanc::DicomMap& dicom)
38 {
39 // http://dicom.nema.org/medical/Dicom/2016a/output/chtml/part03/sect_C.8.8.3.2.html
40
41 {
42 std::string increment;
43
44 if (dicom.CopyToString(increment, Orthanc::DICOM_TAG_FRAME_INCREMENT_POINTER, false))
45 {
46 Orthanc::Toolbox::ToUpperCase(increment);
47 if (increment != "3004,000C") // This is the "Grid Frame Offset Vector" tag
48 {
49 LOG(ERROR) << "RT-DOSE: Bad value for the \"FrameIncrementPointer\" tag";
50 return;
51 }
52 }
53 }
54
55 if (!LinearAlgebra::ParseVector(frameOffsets_, dicom, Orthanc::DICOM_TAG_GRID_FRAME_OFFSET_VECTOR) ||
56 frameOffsets_.size() < imageInformation_.GetNumberOfFrames())
57 {
58 LOG(ERROR) << "RT-DOSE: No information about the 3D location of some slice(s)";
59 frameOffsets_.clear();
60 }
61 else
62 {
63 if (frameOffsets_.size() >= 2)
64 {
65 thickness_ = frameOffsets_[1] - frameOffsets_[0];
66
67 if (thickness_ < 0)
68 {
69 thickness_ = -thickness_;
70 }
71 }
72 }
73 }
74
75
76 DicomInstanceParameters::Data::Data(const Orthanc::DicomMap& dicom) :
77 imageInformation_(dicom)
78 {
79 if (imageInformation_.GetNumberOfFrames() <= 0)
80 {
81 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
82 }
83
84 if (!dicom.CopyToString(studyInstanceUid_, Orthanc::DICOM_TAG_STUDY_INSTANCE_UID, false) ||
85 !dicom.CopyToString(seriesInstanceUid_, Orthanc::DICOM_TAG_SERIES_INSTANCE_UID, false) ||
86 !dicom.CopyToString(sopInstanceUid_, Orthanc::DICOM_TAG_SOP_INSTANCE_UID, false))
87 {
88 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
89 }
90
91 std::string s;
92 if (!dicom.CopyToString(s, Orthanc::DICOM_TAG_SOP_CLASS_UID, false))
93 {
94 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
95 }
96 else
97 {
98 sopClassUid_ = StringToSopClassUid(s);
99 }
100
101 if (!dicom.ParseDouble(thickness_, Orthanc::DICOM_TAG_SLICE_THICKNESS))
102 {
103 thickness_ = 100.0 * std::numeric_limits<double>::epsilon();
104 }
105
106 GeometryToolbox::GetPixelSpacing(pixelSpacingX_, pixelSpacingY_, dicom);
107
108 std::string position, orientation;
109 if (dicom.CopyToString(position, Orthanc::DICOM_TAG_IMAGE_POSITION_PATIENT, false) &&
110 dicom.CopyToString(orientation, Orthanc::DICOM_TAG_IMAGE_ORIENTATION_PATIENT, false))
111 {
112 geometry_ = CoordinateSystem3D(position, orientation);
113 }
114
115 if (sopClassUid_ == SopClassUid_RTDose)
116 {
117 ComputeDoseOffsets(dicom);
118 }
119
120 isColor_ = (imageInformation_.GetPhotometricInterpretation() != Orthanc::PhotometricInterpretation_Monochrome1 &&
121 imageInformation_.GetPhotometricInterpretation() != Orthanc::PhotometricInterpretation_Monochrome2);
122
123 double doseGridScaling;
124
125 if (dicom.ParseDouble(rescaleIntercept_, Orthanc::DICOM_TAG_RESCALE_INTERCEPT) &&
126 dicom.ParseDouble(rescaleSlope_, Orthanc::DICOM_TAG_RESCALE_SLOPE))
127 {
128 hasRescale_ = true;
129 }
130 else if (dicom.ParseDouble(doseGridScaling, Orthanc::DICOM_TAG_DOSE_GRID_SCALING))
131 {
132 hasRescale_ = true;
133 rescaleIntercept_ = 0;
134 rescaleSlope_ = doseGridScaling;
135 }
136 else
137 {
138 hasRescale_ = false;
139 }
140
141 Vector c, w;
142 if (LinearAlgebra::ParseVector(c, dicom, Orthanc::DICOM_TAG_WINDOW_CENTER) &&
143 LinearAlgebra::ParseVector(w, dicom, Orthanc::DICOM_TAG_WINDOW_WIDTH) &&
144 c.size() > 0 &&
145 w.size() > 0)
146 {
147 hasDefaultWindowing_ = true;
148 defaultWindowingCenter_ = static_cast<float>(c[0]);
149 defaultWindowingWidth_ = static_cast<float>(w[0]);
150 }
151 else
152 {
153 hasDefaultWindowing_ = false;
154 }
155
156 if (sopClassUid_ == SopClassUid_RTDose)
157 {
158 switch (imageInformation_.GetBitsStored())
159 {
160 case 16:
161 expectedPixelFormat_ = Orthanc::PixelFormat_Grayscale16;
162 break;
163
164 case 32:
165 expectedPixelFormat_ = Orthanc::PixelFormat_Grayscale32;
166 break;
167
168 default:
169 throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented);
170 }
171 }
172 else if (isColor_)
173 {
174 expectedPixelFormat_ = Orthanc::PixelFormat_RGB24;
175 }
176 else if (imageInformation_.IsSigned())
177 {
178 expectedPixelFormat_ = Orthanc::PixelFormat_SignedGrayscale16;
179 }
180 else
181 {
182 expectedPixelFormat_ = Orthanc::PixelFormat_Grayscale16;
183 }
184 }
185
186
187 CoordinateSystem3D DicomInstanceParameters::Data::GetFrameGeometry(unsigned int frame) const
188 {
189 if (frame == 0)
190 {
191 return geometry_;
192 }
193 else if (frame >= imageInformation_.GetNumberOfFrames())
194 {
195 throw Orthanc::OrthancException(Orthanc::ErrorCode_ParameterOutOfRange);
196 }
197 else if (sopClassUid_ == SopClassUid_RTDose)
198 {
199 if (frame >= frameOffsets_.size())
200 {
201 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError);
202 }
203
204 return CoordinateSystem3D(
205 geometry_.GetOrigin() + frameOffsets_[frame] * geometry_.GetNormal(),
206 geometry_.GetAxisX(),
207 geometry_.GetAxisY());
208 }
209 else
210 {
211 throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented);
212 }
213 }
214
215
216 bool DicomInstanceParameters::Data::IsPlaneWithinSlice(unsigned int frame,
217 const CoordinateSystem3D& plane) const
218 {
219 if (frame >= imageInformation_.GetNumberOfFrames())
220 {
221 throw Orthanc::OrthancException(Orthanc::ErrorCode_ParameterOutOfRange);
222 }
223
224 CoordinateSystem3D tmp = geometry_;
225
226 if (frame != 0)
227 {
228 tmp = GetFrameGeometry(frame);
229 }
230
231 double distance;
232
233 return (CoordinateSystem3D::GetDistance(distance, tmp, plane) &&
234 distance <= thickness_ / 2.0);
235 }
236
237
238 void DicomInstanceParameters::Data::ApplyRescale(Orthanc::ImageAccessor& image,
239 bool useDouble) const
240 {
241 if (image.GetFormat() != Orthanc::PixelFormat_Float32)
242 {
243 throw Orthanc::OrthancException(Orthanc::ErrorCode_IncompatibleImageFormat);
244 }
245
246 if (hasRescale_)
247 {
248 const unsigned int width = image.GetWidth();
249 const unsigned int height = image.GetHeight();
250
251 for (unsigned int y = 0; y < height; y++)
252 {
253 float* p = reinterpret_cast<float*>(image.GetRow(y));
254
255 if (useDouble)
256 {
257 // Slower, accurate implementation using double
258 for (unsigned int x = 0; x < width; x++, p++)
259 {
260 double value = static_cast<double>(*p);
261 *p = static_cast<float>(value * rescaleSlope_ + rescaleIntercept_);
262 }
263 }
264 else
265 {
266 // Fast, approximate implementation using float
267 for (unsigned int x = 0; x < width; x++, p++)
268 {
269 *p = (*p) * static_cast<float>(rescaleSlope_) + static_cast<float>(rescaleIntercept_);
270 }
271 }
272 }
273 }
274 }
275
276 double DicomInstanceParameters::GetRescaleIntercept() const
277 {
278 if (data_.hasRescale_)
279 {
280 return data_.rescaleIntercept_;
281 }
282 else
283 {
284 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadSequenceOfCalls);
285 }
286 }
287
288
289 double DicomInstanceParameters::GetRescaleSlope() const
290 {
291 if (data_.hasRescale_)
292 {
293 return data_.rescaleSlope_;
294 }
295 else
296 {
297 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadSequenceOfCalls);
298 }
299 }
300
301
302 float DicomInstanceParameters::GetDefaultWindowingCenter() const
303 {
304 if (data_.hasDefaultWindowing_)
305 {
306 return data_.defaultWindowingCenter_;
307 }
308 else
309 {
310 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadSequenceOfCalls);
311 }
312 }
313
314
315 float DicomInstanceParameters::GetDefaultWindowingWidth() const
316 {
317 if (data_.hasDefaultWindowing_)
318 {
319 return data_.defaultWindowingWidth_;
320 }
321 else
322 {
323 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadSequenceOfCalls);
324 }
325 }
326
327
328 TextureBaseSceneLayer* DicomInstanceParameters::CreateTexture(const Orthanc::ImageAccessor& pixelData) const
329 {
330 assert(sizeof(float) == 4);
331
332 Orthanc::PixelFormat sourceFormat = pixelData.GetFormat();
333
334 if (sourceFormat != GetExpectedPixelFormat())
335 {
336 throw Orthanc::OrthancException(Orthanc::ErrorCode_IncompatibleImageFormat);
337 }
338
339 if (sourceFormat == Orthanc::PixelFormat_RGB24)
340 {
341 // This is the case of a color image. No conversion has to be done.
342 return new ColorTextureSceneLayer(pixelData);
343 }
344 else
345 {
346 if (sourceFormat != Orthanc::PixelFormat_Grayscale16 &&
347 sourceFormat != Orthanc::PixelFormat_Grayscale32 &&
348 sourceFormat != Orthanc::PixelFormat_SignedGrayscale16)
349 {
350 throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented);
351 }
352
353 std::auto_ptr<FloatTextureSceneLayer> texture;
354
355 {
356 // This is the case of a grayscale frame. Convert it to Float32.
357 std::auto_ptr<Orthanc::Image> converted(new Orthanc::Image(Orthanc::PixelFormat_Float32,
358 pixelData.GetWidth(),
359 pixelData.GetHeight(),
360 false));
361 Orthanc::ImageProcessing::Convert(*converted, pixelData);
362
363 // Correct rescale slope/intercept if need be
364 data_.ApplyRescale(*converted, (sourceFormat == Orthanc::PixelFormat_Grayscale32));
365
366 texture.reset(new FloatTextureSceneLayer(*converted));
367 }
368
369 if (data_.hasDefaultWindowing_)
370 {
371 texture->SetCustomWindowing(data_.defaultWindowingCenter_,
372 data_.defaultWindowingWidth_);
373 }
374
375 return texture.release();
376 }
377 }
378 }