Mercurial > hg > orthanc-stone

comparison Framework/Loaders/OrthancMultiframeVolumeLoader.cpp @ 1337:b1396be5aa27 broker

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Moved the fixed loaders back from the dead
author Benjamin Golinvaux <bgo@osimis.io>
date Fri, 03 Apr 2020 16:13:06 +0200
parents
children b1e6bef86955
comparison
equal deleted inserted replaced
1334:04055b6b9e2c 1337:b1396be5aa27
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-2020 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 "OrthancMultiframeVolumeLoader.h"
23
24 #include <Core/Endianness.h>
25 #include <Core/Toolbox.h>
26
27 namespace OrthancStone
28 {
29 class OrthancMultiframeVolumeLoader::LoadRTDoseGeometry : public LoaderStateMachine::State
30 {
31 private:
32 std::unique_ptr<Orthanc::DicomMap> dicom_;
33
34 public:
35 LoadRTDoseGeometry(OrthancMultiframeVolumeLoader& that,
36 Orthanc::DicomMap* dicom) :
37 State(that),
38 dicom_(dicom)
39 {
40 if (dicom == NULL)
41 {
42 throw Orthanc::OrthancException(Orthanc::ErrorCode_NullPointer);
43 }
44
45 }
46
47 virtual void Handle(const OrthancStone::OrthancRestApiCommand::SuccessMessage& message)
48 {
49 // Complete the DICOM tags with just-received "Grid Frame Offset Vector"
50 std::string s = Orthanc::Toolbox::StripSpaces(message.GetAnswer());
51 dicom_->SetValue(Orthanc::DICOM_TAG_GRID_FRAME_OFFSET_VECTOR, s, false);
52
53 GetLoader<OrthancMultiframeVolumeLoader>().SetGeometry(*dicom_);
54 }
55 };
56
57
58 static std::string GetSopClassUid(const Orthanc::DicomMap& dicom)
59 {
60 std::string s;
61 if (!dicom.LookupStringValue(s, Orthanc::DICOM_TAG_SOP_CLASS_UID, false))
62 {
63 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat,
64 "DICOM file without SOP class UID");
65 }
66 else
67 {
68 return s;
69 }
70 }
71
72
73 class OrthancMultiframeVolumeLoader::LoadGeometry : public State
74 {
75 public:
76 LoadGeometry(OrthancMultiframeVolumeLoader& that) :
77 State(that)
78 {
79 }
80
81 virtual void Handle(const OrthancStone::OrthancRestApiCommand::SuccessMessage& message)
82 {
83 OrthancMultiframeVolumeLoader& loader = GetLoader<OrthancMultiframeVolumeLoader>();
84
85 Json::Value body;
86 message.ParseJsonBody(body);
87
88 if (body.type() != Json::objectValue)
89 {
90 throw Orthanc::OrthancException(Orthanc::ErrorCode_NetworkProtocol);
91 }
92
93 std::unique_ptr<Orthanc::DicomMap> dicom(new Orthanc::DicomMap);
94 dicom->FromDicomAsJson(body);
95
96 if (OrthancStone::StringToSopClassUid(GetSopClassUid(*dicom)) == OrthancStone::SopClassUid_RTDose)
97 {
98 // Download the "Grid Frame Offset Vector" DICOM tag, that is
99 // mandatory for RT-DOSE, but is too long to be returned by default
100
101 std::unique_ptr<OrthancStone::OrthancRestApiCommand> command(new OrthancStone::OrthancRestApiCommand);
102 command->SetUri("/instances/" + loader.GetInstanceId() + "/content/" +
103 Orthanc::DICOM_TAG_GRID_FRAME_OFFSET_VECTOR.Format());
104 command->AcquirePayload(new LoadRTDoseGeometry(loader, dicom.release()));
105
106 Schedule(command.release());
107 }
108 else
109 {
110 loader.SetGeometry(*dicom);
111 }
112 }
113 };
114
115 class OrthancMultiframeVolumeLoader::LoadTransferSyntax : public State
116 {
117 public:
118 LoadTransferSyntax(OrthancMultiframeVolumeLoader& that) :
119 State(that)
120 {
121 }
122
123 virtual void Handle(const OrthancStone::OrthancRestApiCommand::SuccessMessage& message)
124 {
125 GetLoader<OrthancMultiframeVolumeLoader>().SetTransferSyntax(message.GetAnswer());
126 }
127 };
128
129 class OrthancMultiframeVolumeLoader::LoadUncompressedPixelData : public State
130 {
131 public:
132 LoadUncompressedPixelData(OrthancMultiframeVolumeLoader& that) :
133 State(that)
134 {
135 }
136
137 virtual void Handle(const OrthancStone::OrthancRestApiCommand::SuccessMessage& message)
138 {
139 GetLoader<OrthancMultiframeVolumeLoader>().SetUncompressedPixelData(message.GetAnswer());
140 }
141 };
142
143 const std::string& OrthancMultiframeVolumeLoader::GetInstanceId() const
144 {
145 if (IsActive())
146 {
147 return instanceId_;
148 }
149 else
150 {
151 LOG(ERROR) << "OrthancMultiframeVolumeLoader::GetInstanceId(): (!IsActive())";
152 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadSequenceOfCalls);
153 }
154 }
155
156 void OrthancMultiframeVolumeLoader::ScheduleFrameDownloads()
157 {
158 if (transferSyntaxUid_.empty() ||
159 !volume_->HasGeometry())
160 {
161 return;
162 }
163 /*
164 1.2.840.10008.1.2 Implicit VR Endian: Default Transfer Syntax for DICOM
165 1.2.840.10008.1.2.1 Explicit VR Little Endian
166 1.2.840.10008.1.2.2 Explicit VR Big Endian
167
168 See https://www.dicomlibrary.com/dicom/transfer-syntax/
169 */
170 if (transferSyntaxUid_ == "1.2.840.10008.1.2" ||
171 transferSyntaxUid_ == "1.2.840.10008.1.2.1" ||
172 transferSyntaxUid_ == "1.2.840.10008.1.2.2")
173 {
174 std::unique_ptr<OrthancStone::OrthancRestApiCommand> command(new OrthancStone::OrthancRestApiCommand);
175 command->SetHttpHeader("Accept-Encoding", "gzip");
176 command->SetUri("/instances/" + instanceId_ + "/content/" +
177 Orthanc::DICOM_TAG_PIXEL_DATA.Format() + "/0");
178 command->AcquirePayload(new LoadUncompressedPixelData(*this));
179 Schedule(command.release());
180 }
181 else
182 {
183 throw Orthanc::OrthancException(
184 Orthanc::ErrorCode_NotImplemented,
185 "No support for multiframe instances with transfer syntax: " + transferSyntaxUid_);
186 }
187 }
188
189 void OrthancMultiframeVolumeLoader::SetTransferSyntax(const std::string& transferSyntax)
190 {
191 transferSyntaxUid_ = Orthanc::Toolbox::StripSpaces(transferSyntax);
192 ScheduleFrameDownloads();
193 }
194
195 void OrthancMultiframeVolumeLoader::SetGeometry(const Orthanc::DicomMap& dicom)
196 {
197 OrthancStone::DicomInstanceParameters parameters(dicom);
198 volume_->SetDicomParameters(parameters);
199
200 Orthanc::PixelFormat format;
201 if (!parameters.GetImageInformation().ExtractPixelFormat(format, true))
202 {
203 throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented);
204 }
205
206 double spacingZ;
207 switch (parameters.GetSopClassUid())
208 {
209 case OrthancStone::SopClassUid_RTDose:
210 spacingZ = parameters.GetThickness();
211 break;
212
213 default:
214 throw Orthanc::OrthancException(
215 Orthanc::ErrorCode_NotImplemented,
216 "No support for multiframe instances with SOP class UID: " + GetSopClassUid(dicom));
217 }
218
219 const unsigned int width = parameters.GetImageInformation().GetWidth();
220 const unsigned int height = parameters.GetImageInformation().GetHeight();
221 const unsigned int depth = parameters.GetImageInformation().GetNumberOfFrames();
222
223 {
224 OrthancStone::VolumeImageGeometry geometry;
225 geometry.SetSizeInVoxels(width, height, depth);
226 geometry.SetAxialGeometry(parameters.GetGeometry());
227 geometry.SetVoxelDimensions(parameters.GetPixelSpacingX(),
228 parameters.GetPixelSpacingY(), spacingZ);
229 volume_->Initialize(geometry, format, true /* Do compute range */);
230 }
231
232 volume_->GetPixelData().Clear();
233
234 ScheduleFrameDownloads();
235
236
237
238 BroadcastMessage(OrthancStone::DicomVolumeImage::GeometryReadyMessage(*volume_));
239 }
240
241
242 ORTHANC_FORCE_INLINE
243 static void CopyPixel(uint32_t& target, const void* source)
244 {
245 // TODO - check alignement?
246 target = le32toh(*reinterpret_cast<const uint32_t*>(source));
247 }
248
249 ORTHANC_FORCE_INLINE
250 static void CopyPixel(uint16_t& target, const void* source)
251 {
252 // TODO - check alignement?
253 target = le16toh(*reinterpret_cast<const uint16_t*>(source));
254 }
255
256 ORTHANC_FORCE_INLINE
257 static void CopyPixel(int16_t& target, const void* source)
258 {
259 // byte swapping is the same for unsigned and signed integers
260 // (the sign bit is always stored with the MSByte)
261 uint16_t* targetUp = reinterpret_cast<uint16_t*>(&target);
262 CopyPixel(*targetUp, source);
263 }
264
265 template <typename T>
266 void OrthancMultiframeVolumeLoader::CopyPixelDataAndComputeDistribution(
267 const std::string& pixelData, std::map<T,uint64_t>& distribution)
268 {
269 OrthancStone::ImageBuffer3D& target = volume_->GetPixelData();
270
271 const unsigned int bpp = target.GetBytesPerPixel();
272 const unsigned int width = target.GetWidth();
273 const unsigned int height = target.GetHeight();
274 const unsigned int depth = target.GetDepth();
275
276 if (pixelData.size() != bpp * width * height * depth)
277 {
278 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat,
279 "The pixel data has not the proper size");
280 }
281
282 if (pixelData.empty())
283 {
284 return;
285 }
286
287 // first pass to initialize map
288 {
289 const uint8_t* source = reinterpret_cast<const uint8_t*>(pixelData.c_str());
290
291 for (unsigned int z = 0; z < depth; z++)
292 {
293 for (unsigned int y = 0; y < height; y++)
294 {
295 for (unsigned int x = 0; x < width; x++)
296 {
297 T value;
298 CopyPixel(value, source);
299 distribution[value] = 0;
300 source += bpp;
301 }
302 }
303 }
304 }
305
306 {
307 const uint8_t* source = reinterpret_cast<const uint8_t*>(pixelData.c_str());
308
309 for (unsigned int z = 0; z < depth; z++)
310 {
311 OrthancStone::ImageBuffer3D::SliceWriter writer(target, OrthancStone::VolumeProjection_Axial, z);
312
313 assert(writer.GetAccessor().GetWidth() == width &&
314 writer.GetAccessor().GetHeight() == height);
315
316 for (unsigned int y = 0; y < height; y++)
317 {
318 assert(sizeof(T) == Orthanc::GetBytesPerPixel(target.GetFormat()));
319
320 T* target = reinterpret_cast<T*>(writer.GetAccessor().GetRow(y));
321
322 for (unsigned int x = 0; x < width; x++)
323 {
324 CopyPixel(*target, source);
325
326 distribution[*target] += 1;
327
328 target++;
329 source += bpp;
330 }
331 }
332 }
333 }
334 }
335
336 template <typename T>
337 void OrthancMultiframeVolumeLoader::ComputeMinMaxWithOutlierRejection(
338 const std::map<T, uint64_t>& distribution)
339 {
340 if (distribution.size() == 0)
341 {
342 LOG(ERROR) << "ComputeMinMaxWithOutlierRejection -- Volume image empty.";
343 }
344 else
345 {
346 OrthancStone::ImageBuffer3D& target = volume_->GetPixelData();
347
348 const uint64_t width = target.GetWidth();
349 const uint64_t height = target.GetHeight();
350 const uint64_t depth = target.GetDepth();
351 const uint64_t voxelCount = width * height * depth;
352
353 // now that we have distribution[pixelValue] == numberOfPixelsWithValue
354 // compute number of values and check (assertion) that it is equal to
355 // width * height * depth
356 {
357 typename std::map<T, uint64_t>::const_iterator it = distribution.begin();
358 uint64_t totalCount = 0;
359 distributionRawMin_ = static_cast<float>(it->first);
360
361 while (it != distribution.end())
362 {
363 T pixelValue = it->first;
364 uint64_t count = it->second;
365 totalCount += count;
366 it++;
367 if (it == distribution.end())
368 distributionRawMax_ = static_cast<float>(pixelValue);
369 }
370 LOG(INFO) << "Volume image. First distribution value = "
371 << static_cast<float>(distributionRawMin_)
372 << " | Last distribution value = "
373 << static_cast<float>(distributionRawMax_);
374
375 if (totalCount != voxelCount)
376 {
377 LOG(ERROR) << "Internal error in dose distribution computation. TC ("
378 << totalCount << ") != VoxC (" << voxelCount;
379 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError);
380 }
381 }
382
383 // compute the number of voxels to reject at each end of the distribution
384 uint64_t endRejectionCount = static_cast<uint64_t>(
385 outliersHalfRejectionRate_ * voxelCount);
386
387 if (endRejectionCount > voxelCount)
388 {
389 LOG(ERROR) << "Internal error in dose distribution computation."
390 << " endRejectionCount = " << endRejectionCount
391 << " | voxelCount = " << voxelCount;
392 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError);
393 }
394
395 // this will contain the actual distribution minimum after outlier
396 // rejection
397 T resultMin = 0;
398
399 // then start from start and remove pixel values up to
400 // endRejectionCount voxels rejected
401 {
402 typename std::map<T, uint64_t>::const_iterator it = distribution.begin();
403
404 uint64_t currentCount = 0;
405
406 while (it != distribution.end())
407 {
408 T pixelValue = it->first;
409 uint64_t count = it->second;
410
411 // if this pixelValue crosses the rejection threshold, let's set it
412 // and exit the loop
413 if ((currentCount <= endRejectionCount) &&
414 (currentCount + count > endRejectionCount))
415 {
416 resultMin = pixelValue;
417 break;
418 }
419 else
420 {
421 currentCount += count;
422 }
423 // and continue walking along the distribution
424 it++;
425 }
426 }
427
428 // this will contain the actual distribution maximum after outlier
429 // rejection
430 T resultMax = 0;
431 // now start from END and remove pixel values up to
432 // endRejectionCount voxels rejected
433 {
434 typename std::map<T, uint64_t>::const_reverse_iterator it = distribution.rbegin();
435
436 uint64_t currentCount = 0;
437
438 while (it != distribution.rend())
439 {
440 T pixelValue = it->first;
441 uint64_t count = it->second;
442
443 if ((currentCount <= endRejectionCount) &&
444 (currentCount + count > endRejectionCount))
445 {
446 resultMax = pixelValue;
447 break;
448 }
449 else
450 {
451 currentCount += count;
452 }
453 // and continue walking along the distribution
454 it++;
455 }
456 }
457 if (resultMin > resultMax)
458 {
459 LOG(ERROR) << "Internal error in dose distribution computation! " <<
460 "resultMin (" << resultMin << ") > resultMax (" << resultMax << ")";
461 throw Orthanc::OrthancException(Orthanc::ErrorCode_InternalError);
462 }
463 computedDistributionMin_ = static_cast<float>(resultMin);
464 computedDistributionMax_ = static_cast<float>(resultMax);
465 }
466 }
467
468 template <typename T>
469 void OrthancMultiframeVolumeLoader::CopyPixelDataAndComputeMinMax(
470 const std::string& pixelData)
471 {
472 std::map<T, uint64_t> distribution;
473 CopyPixelDataAndComputeDistribution(pixelData, distribution);
474 ComputeMinMaxWithOutlierRejection(distribution);
475 }
476
477 void OrthancMultiframeVolumeLoader::SetUncompressedPixelData(const std::string& pixelData)
478 {
479 switch (volume_->GetPixelData().GetFormat())
480 {
481 case Orthanc::PixelFormat_Grayscale32:
482 CopyPixelDataAndComputeMinMax<uint32_t>(pixelData);
483 break;
484 case Orthanc::PixelFormat_Grayscale16:
485 CopyPixelDataAndComputeMinMax<uint16_t>(pixelData);
486 break;
487 case Orthanc::PixelFormat_SignedGrayscale16:
488 CopyPixelDataAndComputeMinMax<int16_t>(pixelData);
489 break;
490 default:
491 throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented);
492 }
493
494 volume_->IncrementRevision();
495
496 pixelDataLoaded_ = true;
497 BroadcastMessage(OrthancStone::DicomVolumeImage::ContentUpdatedMessage(*volume_));
498 }
499
500 bool OrthancMultiframeVolumeLoader::HasGeometry() const
501 {
502 return volume_->HasGeometry();
503 }
504
505 const OrthancStone::VolumeImageGeometry& OrthancMultiframeVolumeLoader::GetImageGeometry() const
506 {
507 return volume_->GetGeometry();
508 }
509
510 OrthancMultiframeVolumeLoader::OrthancMultiframeVolumeLoader(
511 OrthancStone::ILoadersContext& loadersContext,
512 boost::shared_ptr<OrthancStone::DicomVolumeImage> volume,
513 float outliersHalfRejectionRate)
514 : LoaderStateMachine(loadersContext)
515 , volume_(volume)
516 , pixelDataLoaded_(false)
517 , outliersHalfRejectionRate_(outliersHalfRejectionRate)
518 , distributionRawMin_(0)
519 , distributionRawMax_(0)
520 , computedDistributionMin_(0)
521 , computedDistributionMax_(0)
522 {
523 if (volume.get() == NULL)
524 {
525 throw Orthanc::OrthancException(Orthanc::ErrorCode_NullPointer);
526 }
527 }
528
529
530 boost::shared_ptr<OrthancMultiframeVolumeLoader>
531 OrthancMultiframeVolumeLoader::Create(
532 OrthancStone::ILoadersContext& loadersContext,
533 boost::shared_ptr<OrthancStone::DicomVolumeImage> volume,
534 float outliersHalfRejectionRate /*= 0.0005*/)
535 {
536 boost::shared_ptr<OrthancMultiframeVolumeLoader> obj(
537 new OrthancMultiframeVolumeLoader(
538 loadersContext,
539 volume,
540 outliersHalfRejectionRate));
541 obj->LoaderStateMachine::PostConstructor();
542 return obj;
543 }
544
545 OrthancMultiframeVolumeLoader::~OrthancMultiframeVolumeLoader()
546 {
547 LOG(TRACE) << "OrthancMultiframeVolumeLoader::~OrthancMultiframeVolumeLoader()";
548 }
549
550 void OrthancMultiframeVolumeLoader::GetDistributionMinMax
551 (float& minValue, float& maxValue) const
552 {
553 if (distributionRawMin_ == 0 && distributionRawMax_ == 0)
554 {
555 LOG(WARNING) << "GetDistributionMinMaxWithOutliersRejection called before computation!";
556 }
557 minValue = distributionRawMin_;
558 maxValue = distributionRawMax_;
559 }
560
561 void OrthancMultiframeVolumeLoader::GetDistributionMinMaxWithOutliersRejection
562 (float& minValue, float& maxValue) const
563 {
564 if (computedDistributionMin_ == 0 && computedDistributionMax_ == 0)
565 {
566 LOG(WARNING) << "GetDistributionMinMaxWithOutliersRejection called before computation!";
567 }
568 minValue = computedDistributionMin_;
569 maxValue = computedDistributionMax_;
570 }
571
572 void OrthancMultiframeVolumeLoader::LoadInstance(const std::string& instanceId)
573 {
574 Start();
575
576 instanceId_ = instanceId;
577
578 {
579 std::unique_ptr<OrthancStone::OrthancRestApiCommand> command(new OrthancStone::OrthancRestApiCommand);
580 command->SetHttpHeader("Accept-Encoding", "gzip");
581 command->SetUri("/instances/" + instanceId + "/tags");
582 command->AcquirePayload(new LoadGeometry(*this));
583 Schedule(command.release());
584 }
585
586 {
587 std::unique_ptr<OrthancStone::OrthancRestApiCommand> command(new OrthancStone::OrthancRestApiCommand);
588 command->SetUri("/instances/" + instanceId + "/metadata/TransferSyntax");
589 command->AcquirePayload(new LoadTransferSyntax(*this));
590 Schedule(command.release());
591 }
592 }
593 }