|
32
|
1 /**
|
|
|
2 * SPDX-FileCopyrightText: 2023-2024 Sebastien Jodogne, UCLouvain, Belgium
|
|
|
3 * SPDX-License-Identifier: GPL-3.0-or-later
|
|
|
4 */
|
|
|
5
|
|
|
6 /**
|
|
|
7 * STL plugin for Orthanc
|
|
|
8 * Copyright (C) 2023-2024 Sebastien Jodogne, UCLouvain, Belgium
|
|
|
9 *
|
|
|
10 * This program is free software: you can redistribute it and/or
|
|
|
11 * modify it under the terms of the GNU General Public License as
|
|
|
12 * published by the Free Software Foundation, either version 3 of the
|
|
|
13 * License, or (at your option) any later version.
|
|
|
14 *
|
|
|
15 * This program is distributed in the hope that it will be useful, but
|
|
|
16 * WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
|
17 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
|
|
|
18 * General Public License for more details.
|
|
|
19 *
|
|
|
20 * You should have received a copy of the GNU General Public License
|
|
|
21 * along with this program. If not, see <http://www.gnu.org/licenses/>.
|
|
|
22 **/
|
|
|
23
|
|
|
24
|
|
|
25 #include "VTKToolbox.h"
|
|
|
26
|
|
|
27 #include <ChunkedBuffer.h>
|
|
|
28 #include <Compression/GzipCompressor.h>
|
|
|
29 #include <Endianness.h>
|
|
|
30 #include <OrthancException.h>
|
|
|
31 #include <Toolbox.h>
|
|
|
32
|
|
33
|
33 #include <vtkImageConstantPad.h>
|
|
32
|
34 #include <vtkImageResize.h>
|
|
|
35 #include <vtkMarchingCubes.h>
|
|
33
|
36 #include <vtkNew.h>
|
|
|
37 #include <vtkPolyDataNormals.h>
|
|
32
|
38 #include <vtkSmoothPolyDataFilter.h>
|
|
33
|
39 #include <vtkTriangle.h>
|
|
32
|
40
|
|
|
41 #include <nifti1_io.h>
|
|
|
42
|
|
|
43
|
|
|
44 namespace VTKToolbox
|
|
|
45 {
|
|
|
46 static void WriteFloat(Orthanc::ChunkedBuffer& buffer,
|
|
|
47 float value,
|
|
|
48 Orthanc::Endianness endianness)
|
|
|
49 {
|
|
|
50 switch (endianness)
|
|
|
51 {
|
|
|
52 case Orthanc::Endianness_Little:
|
|
|
53 buffer.AddChunk(&value, sizeof(float));
|
|
|
54 break;
|
|
|
55
|
|
|
56 case Orthanc::Endianness_Big:
|
|
|
57 {
|
|
|
58 uint8_t tmp[4];
|
|
|
59 tmp[0] = reinterpret_cast<uint8_t*>(&value) [3];
|
|
|
60 tmp[1] = reinterpret_cast<uint8_t*>(&value) [2];
|
|
|
61 tmp[2] = reinterpret_cast<uint8_t*>(&value) [1];
|
|
|
62 tmp[3] = reinterpret_cast<uint8_t*>(&value) [0];
|
|
|
63 buffer.AddChunk(&tmp, sizeof(float));
|
|
|
64 break;
|
|
|
65 }
|
|
|
66
|
|
|
67 default:
|
|
|
68 throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented);
|
|
|
69 }
|
|
|
70 }
|
|
|
71
|
|
|
72
|
|
|
73 static void WriteInteger(Orthanc::ChunkedBuffer& buffer,
|
|
|
74 uint32_t value,
|
|
|
75 Orthanc::Endianness endianness)
|
|
|
76 {
|
|
|
77 switch (endianness)
|
|
|
78 {
|
|
|
79 case Orthanc::Endianness_Little:
|
|
|
80 buffer.AddChunk(&value, sizeof(uint32_t));
|
|
|
81 break;
|
|
|
82
|
|
|
83 case Orthanc::Endianness_Big:
|
|
|
84 {
|
|
|
85 uint8_t tmp[4];
|
|
|
86 tmp[0] = reinterpret_cast<uint8_t*>(&value) [3];
|
|
|
87 tmp[1] = reinterpret_cast<uint8_t*>(&value) [2];
|
|
|
88 tmp[2] = reinterpret_cast<uint8_t*>(&value) [1];
|
|
|
89 tmp[3] = reinterpret_cast<uint8_t*>(&value) [0];
|
|
|
90 buffer.AddChunk(&tmp, sizeof(uint32_t));
|
|
|
91 break;
|
|
|
92 }
|
|
|
93
|
|
|
94 default:
|
|
|
95 throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented);
|
|
|
96 }
|
|
|
97 }
|
|
|
98
|
|
|
99
|
|
|
100 void EncodeSTL(std::string& target /* out */,
|
|
|
101 vtkPolyData& mesh /* in */)
|
|
|
102 {
|
|
|
103 const Orthanc::Endianness endianness = Orthanc::Toolbox::DetectEndianness();
|
|
|
104
|
|
|
105 Orthanc::ChunkedBuffer buffer;
|
|
|
106
|
|
|
107 uint8_t header[80];
|
|
|
108 memset(header, 0, sizeof(header));
|
|
|
109 buffer.AddChunk(header, sizeof(header)); // This doesn't depend on endianness
|
|
|
110
|
|
|
111 WriteInteger(buffer, mesh.GetNumberOfCells(), endianness);
|
|
|
112
|
|
|
113 for (vtkIdType i = 0; i < mesh.GetNumberOfCells(); i++)
|
|
|
114 {
|
|
|
115 vtkCell* cell = mesh.GetCell(i);
|
|
|
116 vtkTriangle* triangle = dynamic_cast<vtkTriangle*>(cell);
|
|
|
117
|
|
|
118 double p0[3];
|
|
|
119 double p1[3];
|
|
|
120 double p2[3];
|
|
|
121 triangle->GetPoints()->GetPoint(0, p0);
|
|
|
122 triangle->GetPoints()->GetPoint(1, p1);
|
|
|
123 triangle->GetPoints()->GetPoint(2, p2);
|
|
|
124
|
|
|
125 double normal[3];
|
|
|
126 vtkTriangle::ComputeNormal(p0, p1, p2, normal);
|
|
|
127
|
|
|
128 WriteFloat(buffer, normal[0], endianness);
|
|
|
129 WriteFloat(buffer, normal[1], endianness);
|
|
|
130 WriteFloat(buffer, normal[2], endianness);
|
|
|
131 WriteFloat(buffer, p0[0], endianness);
|
|
|
132 WriteFloat(buffer, p0[1], endianness);
|
|
|
133 WriteFloat(buffer, p0[2], endianness);
|
|
|
134 WriteFloat(buffer, p1[0], endianness);
|
|
|
135 WriteFloat(buffer, p1[1], endianness);
|
|
|
136 WriteFloat(buffer, p1[2], endianness);
|
|
|
137 WriteFloat(buffer, p2[0], endianness);
|
|
|
138 WriteFloat(buffer, p2[1], endianness);
|
|
|
139 WriteFloat(buffer, p2[2], endianness);
|
|
|
140
|
|
|
141 uint16_t a = 0;
|
|
|
142 buffer.AddChunk(&a, sizeof(a)); // This doesn't depend on endianness
|
|
|
143 }
|
|
|
144
|
|
|
145 buffer.Flatten(target);
|
|
|
146 }
|
|
|
147
|
|
|
148
|
|
|
149 bool EncodeVolume(std::string& stl,
|
|
|
150 vtkImageData* volume,
|
|
|
151 unsigned int resolution,
|
|
|
152 bool smooth)
|
|
|
153 {
|
|
|
154 if (volume == NULL)
|
|
|
155 {
|
|
|
156 throw Orthanc::OrthancException(Orthanc::ErrorCode_NullPointer);
|
|
|
157 }
|
|
|
158
|
|
|
159 vtkNew<vtkImageResize> resize;
|
|
|
160 resize->SetOutputDimensions(resolution, resolution, resolution);
|
|
|
161 resize->SetInputData(volume);
|
|
|
162 resize->Update();
|
|
|
163
|
|
|
164 vtkNew<vtkImageConstantPad> padding;
|
|
|
165 padding->SetConstant(0);
|
|
|
166 padding->SetOutputNumberOfScalarComponents(1);
|
|
|
167 padding->SetOutputWholeExtent(-1, resolution, -1, resolution, -1, resolution);
|
|
|
168 padding->SetInputData(resize->GetOutput());
|
|
|
169 padding->Update();
|
|
|
170
|
|
|
171 double range[2];
|
|
|
172 padding->GetOutput()->GetScalarRange(range);
|
|
|
173
|
|
|
174 const double isoValue = (range[0] + range[1]) / 2.0;
|
|
|
175
|
|
|
176 vtkNew<vtkMarchingCubes> surface;
|
|
|
177 surface->SetInputData(padding->GetOutput());
|
|
|
178 surface->ComputeNormalsOn();
|
|
|
179 surface->SetValue(0, isoValue);
|
|
|
180 surface->Update();
|
|
|
181
|
|
|
182 if (smooth)
|
|
|
183 {
|
|
|
184 vtkNew<vtkSmoothPolyDataFilter> smoothFilter;
|
|
|
185 // Apply volume smoothing
|
|
|
186 // https://examples.vtk.org/site/Cxx/PolyData/SmoothPolyDataFilter/
|
|
|
187 smoothFilter->SetInputConnection(surface->GetOutputPort());
|
|
|
188 smoothFilter->SetNumberOfIterations(15);
|
|
|
189 smoothFilter->SetRelaxationFactor(0.1);
|
|
|
190 smoothFilter->FeatureEdgeSmoothingOff();
|
|
|
191 smoothFilter->BoundarySmoothingOn();
|
|
|
192 smoothFilter->Update();
|
|
|
193
|
|
|
194 vtkNew<vtkPolyDataNormals> normalGenerator;
|
|
|
195 normalGenerator->SetInputConnection(smoothFilter->GetOutputPort());
|
|
|
196 normalGenerator->ComputePointNormalsOn();
|
|
|
197 normalGenerator->ComputeCellNormalsOn();
|
|
|
198 normalGenerator->Update();
|
|
|
199
|
|
|
200 EncodeSTL(stl, *normalGenerator->GetOutput());
|
|
|
201 }
|
|
|
202 else
|
|
|
203 {
|
|
|
204 EncodeSTL(stl, *surface->GetOutput());
|
|
|
205 }
|
|
|
206
|
|
|
207 return true;
|
|
|
208 }
|
|
|
209
|
|
|
210
|
|
|
211 namespace
|
|
|
212 {
|
|
|
213 class NiftiHeader : public boost::noncopyable
|
|
|
214 {
|
|
|
215 private:
|
|
|
216 nifti_image* image_;
|
|
|
217
|
|
|
218 public:
|
|
40
|
219 explicit NiftiHeader(const std::string& nifti)
|
|
32
|
220 {
|
|
|
221 nifti_1_header header;
|
|
|
222 if (nifti.size() < sizeof(header))
|
|
|
223 {
|
|
|
224 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
|
|
|
225 }
|
|
|
226
|
|
|
227 memcpy(&header, nifti.c_str(), sizeof(header));
|
|
|
228 if (!nifti_hdr_looks_good(&header))
|
|
|
229 {
|
|
|
230 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
|
|
|
231 }
|
|
|
232
|
|
|
233 image_ = nifti_convert_nhdr2nim(header, "dummy_filename");
|
|
|
234 if (image_ == NULL)
|
|
|
235 {
|
|
|
236 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat);
|
|
|
237 }
|
|
|
238 }
|
|
|
239
|
|
|
240 ~NiftiHeader()
|
|
|
241 {
|
|
|
242 nifti_image_free(image_);
|
|
|
243 }
|
|
|
244
|
|
|
245 const nifti_image& GetInfo() const
|
|
|
246 {
|
|
|
247 assert(image_ != NULL);
|
|
|
248 return *image_;
|
|
|
249 }
|
|
|
250 };
|
|
|
251 }
|
|
|
252
|
|
|
253
|
|
|
254 void LoadNifti(vtkImageData* volume /* out */,
|
|
|
255 std::string& nifti /* in */)
|
|
|
256 {
|
|
|
257 if (volume == NULL)
|
|
|
258 {
|
|
|
259 throw Orthanc::OrthancException(Orthanc::ErrorCode_NullPointer);
|
|
|
260 }
|
|
|
261
|
|
|
262 const uint8_t* p = reinterpret_cast<const uint8_t*>(nifti.c_str());
|
|
|
263
|
|
|
264 if (nifti.size() >= 2 &&
|
|
|
265 p[0] == 0x1f &&
|
|
|
266 p[1] == 0x8b)
|
|
|
267 {
|
|
|
268 Orthanc::GzipCompressor compressor;
|
|
|
269 std::string uncompressed;
|
|
|
270 Orthanc::IBufferCompressor::Uncompress(uncompressed, compressor, nifti);
|
|
|
271 nifti.swap(uncompressed);
|
|
|
272 }
|
|
|
273
|
|
|
274 NiftiHeader header(nifti);
|
|
|
275
|
|
|
276 if (header.GetInfo().ndim != 3)
|
|
|
277 {
|
|
|
278 throw Orthanc::OrthancException(Orthanc::ErrorCode_BadFileFormat,
|
|
|
279 "Only 3D NIfTI volumes are allowed");
|
|
|
280 }
|
|
|
281
|
|
|
282 size_t itemSize;
|
|
|
283 int vtkType;
|
|
|
284
|
|
|
285 switch (header.GetInfo().datatype)
|
|
|
286 {
|
|
|
287 case DT_UNSIGNED_CHAR:
|
|
|
288 itemSize = 1;
|
|
|
289 vtkType = VTK_UNSIGNED_CHAR;
|
|
|
290 break;
|
|
|
291
|
|
|
292 case DT_FLOAT:
|
|
|
293 itemSize = sizeof(float);
|
|
|
294 vtkType = VTK_FLOAT;
|
|
|
295 break;
|
|
|
296
|
|
|
297 case DT_DOUBLE:
|
|
|
298 itemSize = sizeof(double);
|
|
|
299 vtkType = VTK_DOUBLE;
|
|
|
300 break;
|
|
|
301
|
|
|
302 default:
|
|
|
303 throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented);
|
|
|
304 }
|
|
|
305
|
|
|
306 assert(static_cast<int>(header.GetInfo().nvox) == header.GetInfo().nx * header.GetInfo().ny * header.GetInfo().nz);
|
|
|
307
|
|
|
308 const size_t pixelDataOffset = sizeof(nifti_1_header) + 4 /* extension */;
|
|
|
309
|
|
|
310 if (nifti.size() != pixelDataOffset + header.GetInfo().nvox * itemSize)
|
|
|
311 {
|
|
|
312 throw Orthanc::OrthancException(Orthanc::ErrorCode_CorruptedFile);
|
|
|
313 }
|
|
|
314
|
|
|
315 volume->SetDimensions(header.GetInfo().nx, header.GetInfo().ny, header.GetInfo().nz);
|
|
|
316 volume->AllocateScalars(vtkType, 1);
|
|
|
317 volume->SetSpacing(header.GetInfo().dx, header.GetInfo().dy, header.GetInfo().dz);
|
|
|
318 memcpy(volume->GetScalarPointer(), &nifti[pixelDataOffset], header.GetInfo().nvox * itemSize);
|
|
|
319 }
|
|
|
320 }
|
