1 /*=========================================================================
2
3 Program: GDCM (Grassroots DICOM). A DICOM library
4
5 Copyright (c) 2006-2011 Mathieu Malaterre
6 All rights reserved.
7 See Copyright.txt or http://gdcm.sourceforge.net/Copyright.html for details.
8
9 This software is distributed WITHOUT ANY WARRANTY; without even
10 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
11 PURPOSE. See the above copyright notice for more information.
12
13 =========================================================================*/
14 #include "gdcmImageHelper.h"
15 #include "gdcmMediaStorage.h"
16 #include "gdcmFile.h"
17 #include "gdcmDataSet.h"
18 #include "gdcmDataElement.h"
19 #include "gdcmItem.h"
20 #include "gdcmSequenceOfItems.h"
21 #include "gdcmGlobal.h"
22 #include "gdcmDictEntry.h"
23 #include "gdcmDicts.h"
24 #include "gdcmAttribute.h"
25 #include "gdcmImage.h"
26 #include "gdcmDirectionCosines.h"
27 #include "gdcmSegmentedPaletteColorLookupTable.h"
28 #include "gdcmByteValue.h"
29
30 #include <math.h> // fabs
31
32 /* TODO:
33 *
34 * (0028,9145) SQ (Sequence with undefined length #=1) # u/l, 1 PixelValueTransformationSequence
35 * (fffe,e000) na (Item with undefined length #=3) # u/l, 1 Item
36 * (0028,1052) DS [0.00000] # 8, 1 RescaleIntercept
37 * (0028,1053) DS [1.00000] # 8, 1 RescaleSlope
38 * (0028,1054) LO [US] # 2, 1 RescaleType
39 * (fffe,e00d) na (ItemDelimitationItem) # 0, 0 ItemDelimitationItem
40 * (fffe,e0dd) na (SequenceDelimitationItem) # 0, 0 SequenceDelimitationItem
41 *
42 * Same goes for window level...
43 */
44
45 namespace gdcm
46 {
47
48 bool ImageHelper::ForceRescaleInterceptSlope = false;
49 bool ImageHelper::PMSRescaleInterceptSlope = true;
50 bool ImageHelper::ForcePixelSpacing = false;
51
GetOriginValueFromSequence(const DataSet & ds,const Tag & tfgs,std::vector<double> & ori)52 static bool GetOriginValueFromSequence(const DataSet& ds, const Tag& tfgs, std::vector<double> &ori)
53 {
54 if( !ds.FindDataElement( tfgs ) ) return false;
55 //const SequenceOfItems * sqi = ds.GetDataElement( tfgs ).GetSequenceOfItems();
56 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( tfgs ).GetValueAsSQ();
57 if( !(sqi && sqi->GetNumberOfItems() > 0) ) return false;
58 // Get first item:
59 const Item &item = sqi->GetItem(1);
60 const DataSet & subds = item.GetNestedDataSet();
61 // Plane position Sequence
62 const Tag tpms(0x0020,0x9113);
63 if( !subds.FindDataElement(tpms) ) return false;
64 //const SequenceOfItems * sqi2 = subds.GetDataElement( tpms ).GetSequenceOfItems();
65 SmartPointer<SequenceOfItems> sqi2 = subds.GetDataElement( tpms ).GetValueAsSQ();
66 if( sqi2 && !sqi2->IsEmpty() )
67 {
68 const Item &item2 = sqi2->GetItem(1);
69 const DataSet & subds2 = item2.GetNestedDataSet();
70 //
71 const Tag tps(0x0020,0x0032);
72 if( !subds2.FindDataElement(tps) ) return false;
73 const DataElement &de = subds2.GetDataElement( tps );
74 //assert( bv );
75 Attribute<0x0020,0x0032> at;
76 at.SetFromDataElement( de );
77 //at.Print( std::cout );
78 ori.push_back( at.GetValue(0) );
79 ori.push_back( at.GetValue(1) );
80 ori.push_back( at.GetValue(2) );
81
82 return true;
83 }
84 return false;
85 }
86
GetDirectionCosinesValueFromSequence(const DataSet & ds,const Tag & tfgs,std::vector<double> & dircos)87 static bool GetDirectionCosinesValueFromSequence(const DataSet& ds, const Tag& tfgs, std::vector<double> &dircos)
88 {
89 if( !ds.FindDataElement( tfgs ) ) return false;
90 //const SequenceOfItems * sqi = ds.GetDataElement( tfgs ).GetSequenceOfItems();
91 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( tfgs ).GetValueAsSQ();
92 if( !(sqi && sqi->GetNumberOfItems() > 0) ) return false;
93 // Get first item:
94 const Item &item = sqi->GetItem(1);
95 const DataSet & subds = item.GetNestedDataSet();
96 // Plane position Sequence
97 const Tag tpms(0x0020,0x9116);
98 if( !subds.FindDataElement(tpms) ) return false;
99 //const SequenceOfItems * sqi2 = subds.GetDataElement( tpms ).GetSequenceOfItems();
100 SmartPointer<SequenceOfItems> sqi2 = subds.GetDataElement( tpms ).GetValueAsSQ();
101 if( !(sqi2 && sqi2->GetNumberOfItems()) ) return false;
102 assert( sqi2 && sqi2->GetNumberOfItems() );
103 // Take it from the first item
104 const Item &item2 = sqi2->GetItem(1);
105 const DataSet & subds2 = item2.GetNestedDataSet();
106 //
107 const Tag tps(0x0020,0x0037);
108 if( !subds2.FindDataElement(tps) ) return false;
109 const DataElement &de = subds2.GetDataElement( tps );
110 //assert( bv );
111 Attribute<0x0020,0x0037> at;
112 at.SetFromDataElement( de );
113 dircos.push_back( at.GetValue(0) );
114 dircos.push_back( at.GetValue(1) );
115 dircos.push_back( at.GetValue(2) );
116 dircos.push_back( at.GetValue(3) );
117 dircos.push_back( at.GetValue(4) );
118 dircos.push_back( at.GetValue(5) );
119
120 return true;
121 }
122
GetInterceptSlopeValueFromSequence(const DataSet & ds,const Tag & tfgs,std::vector<double> & intslope)123 static bool GetInterceptSlopeValueFromSequence(const DataSet& ds, const Tag& tfgs, std::vector<double> &intslope)
124 {
125 if( !ds.FindDataElement( tfgs ) ) return false;
126 //const SequenceOfItems * sqi = ds.GetDataElement( tfgs ).GetSequenceOfItems();
127 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( tfgs ).GetValueAsSQ();
128 if( !(sqi && sqi->GetNumberOfItems() > 0) ) return false;
129 // Get first item:
130 const Item &item = sqi->GetItem(1);
131 const DataSet & subds = item.GetNestedDataSet();
132 // (0028,9145) SQ (Sequence with undefined length) # u/l,1 Pixel Value Transformation Sequence
133 const Tag tpms(0x0028,0x9145);
134 if( !subds.FindDataElement(tpms) ) return false;
135 //const SequenceOfItems * sqi2 = subds.GetDataElement( tpms ).GetSequenceOfItems();
136 SmartPointer<SequenceOfItems> sqi2 = subds.GetDataElement( tpms ).GetValueAsSQ();
137 assert( sqi2 );
138 const Item &item2 = sqi2->GetItem(1);
139 const DataSet & subds2 = item2.GetNestedDataSet();
140 {
141 // (0028,1052) DS [0] # 2,1 Rescale Intercept
142 const Tag tps(0x0028,0x1052);
143 if( !subds2.FindDataElement(tps) ) return false;
144 const DataElement &de = subds2.GetDataElement( tps );
145 //assert( bv );
146 Attribute<0x0028,0x1052> at;
147 at.SetFromDataElement( de );
148 //at.Print( std::cout );
149 intslope.push_back( at.GetValue() );
150 }
151 {
152 // (0028,1053) DS [5.65470085470085] # 16,1 Rescale Slope
153 const Tag tps(0x0028,0x1053);
154 if( !subds2.FindDataElement(tps) ) return false;
155 const DataElement &de = subds2.GetDataElement( tps );
156 //assert( bv );
157 Attribute<0x0028,0x1053> at;
158 at.SetFromDataElement( de );
159 //at.Print( std::cout );
160 intslope.push_back( at.GetValue() );
161 }
162
163 assert( intslope.size() == 2 );
164 return true;
165 }
166
ComputeZSpacingFromIPP(const DataSet & ds,double & zspacing)167 static bool ComputeZSpacingFromIPP(const DataSet &ds, double &zspacing)
168 {
169 // first we need to get the direction cosines:
170 const Tag t1(0x5200,0x9229);
171 const Tag t2(0x5200,0x9230);
172 std::vector<double> cosines;
173 // For some reason TOSHIBA-EnhancedCT.dcm is storing the direction cosines in the per-frame section
174 // and not the shared one... oh well
175 bool b1 = GetDirectionCosinesValueFromSequence(ds, t1, cosines)
176 || GetDirectionCosinesValueFromSequence(ds,t2,cosines);
177 if(!b1)
178 {
179 cosines.resize( 6 );
180 bool b2 = ImageHelper::GetDirectionCosinesFromDataSet(ds, cosines);
181 if( b2 )
182 {
183 gdcmWarningMacro( "Image Orientation (Patient) cannot be stored here!. Continuing" );
184 b1 = b2;
185 }
186 else
187 {
188 gdcmErrorMacro( "Image Orientation (Patient) was not found" );
189 cosines[0] = 1;
190 cosines[1] = 0;
191 cosines[2] = 0;
192 cosines[3] = 0;
193 cosines[4] = 1;
194 cosines[5] = 0;
195 }
196 }
197 assert( b1 && cosines.size() == 6 ); // yeah we really need that
198
199 const Tag tfgs(0x5200,0x9230);
200 if( !ds.FindDataElement( tfgs ) ) return false;
201 //const SequenceOfItems * sqi = ds.GetDataElement( tfgs ).GetSequenceOfItems();
202 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( tfgs ).GetValueAsSQ();
203 if( !sqi ) return false;
204 assert( sqi );
205 double normal[3];
206 DirectionCosines dc( &cosines[0] );
207 dc.Cross( normal );
208
209 // For each item
210 SequenceOfItems::SizeType nitems = sqi->GetNumberOfItems();
211 if( nitems > 1 ) {
212 std::vector<double> dircos_subds2; dircos_subds2.resize(6);
213 std::vector<double> distances;
214 for(SequenceOfItems::SizeType i0 = 1; i0 <= nitems; ++i0)
215 {
216 const Item &item = sqi->GetItem(i0);
217 const DataSet & subds = item.GetNestedDataSet();
218 // (0020,9113) SQ (Sequence with undefined length #=1) # u/l, 1 PlanePositionSequence
219 const Tag tpms(0x0020,0x9113);
220 if( !subds.FindDataElement(tpms) ) return false;
221 //const SequenceOfItems * sqi2 = subds.GetDataElement( tpms ).GetSequenceOfItems();
222 SmartPointer<SequenceOfItems> sqi2 = subds.GetDataElement( tpms ).GetValueAsSQ();
223 assert( sqi2 );
224 const Item &item2 = sqi2->GetItem(1);
225 const DataSet & subds2 = item2.GetNestedDataSet();
226 // Check Image Orientation (Patient)
227 if( ImageHelper::GetDirectionCosinesFromDataSet(subds2, dircos_subds2) )
228 {
229 assert( dircos_subds2 == cosines );
230 }
231 // (0020,0032) DS [-82.5\-82.5\1153.75] # 20, 3 ImagePositionPatient
232 const Tag tps(0x0020,0x0032);
233 if( !subds2.FindDataElement(tps) ) return false;
234 const DataElement &de = subds2.GetDataElement( tps );
235 Attribute<0x0020,0x0032> ipp;
236 ipp.SetFromDataElement(de);
237 double dist = 0;
238 for (int i = 0; i < 3; ++i) dist += normal[i]*ipp[i];
239 distances.push_back( dist );
240 }
241 assert( distances.size() == nitems );
242 double meanspacing = 0;
243 double prev = distances[0];
244 for(unsigned int i = 1; i < nitems; ++i)
245 {
246 const double current = distances[i] - prev;
247 meanspacing += current;
248 prev = distances[i];
249 }
250 bool timeseries = false;
251 if( nitems > 1 )
252 {
253 meanspacing /= (double)(nitems - 1);
254 if( meanspacing == 0.0 )
255 {
256 // Could be a time series. Assume time spacing of 1. for now:
257 gdcmDebugMacro( "Assuming time series for Z-spacing" );
258 meanspacing = 1.0;
259 timeseries = true;
260 }
261 }
262
263 zspacing = meanspacing;
264 if( nitems > 1 )
265 assert( zspacing != 0.0 ); // technically this should not happen
266
267 if( !timeseries )
268 {
269 // Check spacing is consistent:
270 const double ZTolerance = 1e-3; // ??? FIXME
271 prev = distances[0];
272 for(unsigned int i = 1; i < nitems; ++i)
273 {
274 const double current = distances[i] - prev;
275 if( fabs(current - zspacing) > ZTolerance )
276 {
277 // For now simply gives up
278 gdcmErrorMacro( "This Enhanced Multiframe is not supported for now. Sorry" );
279 return false;
280 }
281 prev = distances[i];
282 }
283 }
284 } else {
285 // single slice, this is not an error to not find the zspacing in this case.
286 zspacing = 1.0;
287 const Tag tfgs0(0x5200,0x9229);
288 if( !ds.FindDataElement( tfgs0 ) ) return true;
289 SmartPointer<SequenceOfItems> sqi0 = ds.GetDataElement( tfgs0 ).GetValueAsSQ();
290 if( !(sqi0 && sqi0->GetNumberOfItems() > 0) ) return true;
291 // Get first item:
292 const Item &item = sqi0->GetItem(1);
293 const DataSet & subds = item.GetNestedDataSet();
294 // <entry group="0028" element="9110" vr="SQ" vm="1" name="Pixel Measures Sequence"/>
295 const Tag tpms(0x0028,0x9110);
296 if( !subds.FindDataElement(tpms) ) return true;
297 //const SequenceOfItems * sqi2 = subds.GetDataElement( tpms ).GetSequenceOfItems();
298 SmartPointer<SequenceOfItems> sqi2 = subds.GetDataElement( tpms ).GetValueAsSQ();
299 assert( sqi2 );
300 const Item &item2 = sqi2->GetItem(1);
301 const DataSet & subds2 = item2.GetNestedDataSet();
302 // <entry group="0028" element="0030" vr="DS" vm="2" name="Pixel Spacing"/>
303 const Tag tps(0x0018,0x0088);
304 if( !subds2.FindDataElement(tps) ) return true;
305 const DataElement &de = subds2.GetDataElement( tps );
306 Attribute<0x0018,0x0088> at;
307 at.SetFromDataElement( de );
308 zspacing = at.GetValue();
309 }
310 return true;
311 }
312
313 // EnhancedMRImageStorage & EnhancedCTImageStorage
GetSpacingValueFromSequence(const DataSet & ds,const Tag & tfgs,std::vector<double> & sp)314 static bool GetSpacingValueFromSequence(const DataSet& ds, const Tag& tfgs, std::vector<double> &sp)
315 {
316 // (0028,9110) SQ (Sequence with undefined length #=1) # u/l, 1 PixelMeasuresSequence
317 // (fffe,e000) na (Item with undefined length #=2) # u/l, 1 Item
318 // (0018,0050) DS [0.5] # 4, 1 SliceThickness
319 // (0028,0030) DS [0.322\0.322] # 12, 2 PixelSpacing
320 // <entry group="5200" element="9229" vr="SQ" vm="1" name="Shared Functional Groups Sequence"/>
321 //const Tag tfgs(0x5200,0x9229);
322 //const Tag tfgs(0x5200,0x9230);
323 //assert( ds.FindDataElement( tfgs ) );
324 if( !ds.FindDataElement( tfgs ) ) return false;
325 //const SequenceOfItems * sqi = ds.GetDataElement( tfgs ).GetSequenceOfItems();
326 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( tfgs ).GetValueAsSQ();
327 if( !(sqi && sqi->GetNumberOfItems() > 0) ) return false;
328 // Get first item:
329 const Item &item = sqi->GetItem(1);
330 const DataSet & subds = item.GetNestedDataSet();
331 // <entry group="0028" element="9110" vr="SQ" vm="1" name="Pixel Measures Sequence"/>
332 const Tag tpms(0x0028,0x9110);
333 if( !subds.FindDataElement(tpms) ) return false;
334 //const SequenceOfItems * sqi2 = subds.GetDataElement( tpms ).GetSequenceOfItems();
335 SmartPointer<SequenceOfItems> sqi2 = subds.GetDataElement( tpms ).GetValueAsSQ();
336 assert( sqi2 );
337 const Item &item2 = sqi2->GetItem(1);
338 const DataSet & subds2 = item2.GetNestedDataSet();
339 // <entry group="0028" element="0030" vr="DS" vm="2" name="Pixel Spacing"/>
340 const Tag tps(0x0028,0x0030);
341 if( !subds2.FindDataElement(tps) ) return false;
342 const DataElement &de = subds2.GetDataElement( tps );
343 //assert( bv );
344 Attribute<0x0028,0x0030> at;
345 at.SetFromDataElement( de );
346 //at.Print( std::cout );
347 sp.push_back( at.GetValue(1) );
348 sp.push_back( at.GetValue(0) );
349
350 // BUG ! Check for instance:
351 // gdcmData/BRTUM001.dcm
352 // Slice Thickness is 5.0 while the Zspacing should be 6.0 !
353 #if 0
354 // Do the 3rd dimension zspacing:
355 // <entry group="0018" element="0050" vr="DS" vm="1" name="Slice Thickness"/>
356 const Tag tst(0x0018,0x0050);
357 if( !subds2.FindDataElement(tst) ) return false;
358 const DataElement &de2 = subds2.GetDataElement( tst );
359 Attribute<0x0018,0x0050> at2;
360 at2.SetFromDataElement( de2 );
361 //at2.Print( std::cout );
362 sp.push_back( at2.GetValue(0) );
363 #endif
364 double zspacing;
365 bool b = ComputeZSpacingFromIPP(ds, zspacing);
366 if( !b ) return false;
367
368 sp.push_back( zspacing );
369
370 return true;
371 }
372
InsertOrReplaceSQ(DataSet & ds,const Tag & tag)373 static SmartPointer<SequenceOfItems> InsertOrReplaceSQ( DataSet & ds, const Tag &tag )
374 {
375 SmartPointer<SequenceOfItems> sqi;
376 if( !ds.FindDataElement( tag ) )
377 {
378 sqi = new SequenceOfItems;
379 DataElement de( tag );
380 de.SetVR( VR::SQ );
381 de.SetValue( *sqi );
382 assert( de.GetVL().IsUndefined() );
383 de.SetVLToUndefined();
384 ds.Insert( de );
385 }
386 sqi = ds.GetDataElement( tag ).GetValueAsSQ();
387 sqi->SetLengthToUndefined();
388 DataElement de_dup = ds.GetDataElement( tag );
389 de_dup.SetValue( *sqi );
390 ds.Replace( de_dup );
391 return sqi;
392 }
393
394
395 // UltrasoundMultiframeImageStorage
GetUltraSoundSpacingValueFromSequence(const DataSet & ds,std::vector<double> & sp)396 static bool GetUltraSoundSpacingValueFromSequence(const DataSet& ds, std::vector<double> &sp)
397 {
398 /*
399 (0018,6011) SQ (Sequence with explicit length #=1) # 196, 1 SequenceOfUltrasoundRegions
400 (fffe,e000) na (Item with explicit length #=15) # 188, 1 Item
401 (0018,6012) US 1 # 2, 1 RegionSpatialFormat
402 (0018,6014) US 1 # 2, 1 RegionDataType
403 (0018,6016) UL 0 # 4, 1 RegionFlags
404 (0018,6018) UL 0 # 4, 1 RegionLocationMinX0
405 (0018,601a) UL 0 # 4, 1 RegionLocationMinY0
406 (0018,601c) UL 479 # 4, 1 RegionLocationMaxX1
407 (0018,601e) UL 479 # 4, 1 RegionLocationMaxY1
408 (0018,6020) SL 0 # 4, 1 ReferencePixelX0
409 (0018,6022) SL 0 # 4, 1 ReferencePixelY0
410 (0018,6024) US 3 # 2, 1 PhysicalUnitsXDirection
411 (0018,6026) US 3 # 2, 1 PhysicalUnitsYDirection
412 (0018,6028) FD 0 # 8, 1 ReferencePixelPhysicalValueX
413 (0018,602a) FD 0 # 8, 1 ReferencePixelPhysicalValueY
414 (0018,602c) FD 0.002 # 8, 1 PhysicalDeltaX
415 (0018,602e) FD 0.002 # 8, 1 PhysicalDeltaY
416 (fffe,e00d) na (ItemDelimitationItem for re-encoding) # 0, 0 ItemDelimitationItem
417 (fffe,e0dd) na (SequenceDelimitationItem for re-encod.) # 0, 0 SequenceDelimitationItem
418 */
419 const Tag tsqusreg(0x0018,0x6011);
420 if( !ds.FindDataElement( tsqusreg ) ) return false;
421 //const SequenceOfItems * sqi = ds.GetDataElement( tsqusreg ).GetSequenceOfItems();
422 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( tsqusreg ).GetValueAsSQ();
423 if( !sqi ) return false;
424 assert( sqi );
425 // Get first item:
426 const Item &item = sqi->GetItem(1);
427 const DataSet & subds = item.GetNestedDataSet();
428 // (0018,602c) FD 0.002 # 8, 1 PhysicalDeltaX
429 // (0018,602e) FD 0.002 # 8, 1 PhysicalDeltaY
430 Attribute<0x0018,0x602c> at1;
431 Attribute<0x0018,0x602e> at2;
432 const DataElement &de1 = subds.GetDataElement( at1.GetTag() );
433 at1.SetFromDataElement( de1 );
434 assert( at1.GetNumberOfValues() == 1 );
435 const DataElement &de2 = subds.GetDataElement( at2.GetTag() );
436 at2.SetFromDataElement( de2 );
437 assert( at2.GetNumberOfValues() == 1 );
438 sp.push_back( at1.GetValue() );
439 sp.push_back( at2.GetValue() );
440
441 return true;
442 }
SetUltraSoundSpacingValueFromSequence(DataSet & ds,std::vector<double> const & spacing)443 static void SetUltraSoundSpacingValueFromSequence(DataSet& ds, std::vector<double> const &spacing)
444 {
445 /*
446 (0018,6011) SQ (Sequence with explicit length #=1) # 196, 1 SequenceOfUltrasoundRegions
447 (fffe,e000) na (Item with explicit length #=15) # 188, 1 Item
448 (0018,6012) US 1 # 2, 1 RegionSpatialFormat
449 (0018,6014) US 1 # 2, 1 RegionDataType
450 (0018,6016) UL 0 # 4, 1 RegionFlags
451 (0018,6018) UL 0 # 4, 1 RegionLocationMinX0
452 (0018,601a) UL 0 # 4, 1 RegionLocationMinY0
453 (0018,601c) UL 479 # 4, 1 RegionLocationMaxX1
454 (0018,601e) UL 479 # 4, 1 RegionLocationMaxY1
455 (0018,6020) SL 0 # 4, 1 ReferencePixelX0
456 (0018,6022) SL 0 # 4, 1 ReferencePixelY0
457 (0018,6024) US 3 # 2, 1 PhysicalUnitsXDirection
458 (0018,6026) US 3 # 2, 1 PhysicalUnitsYDirection
459 (0018,6028) FD 0 # 8, 1 ReferencePixelPhysicalValueX
460 (0018,602a) FD 0 # 8, 1 ReferencePixelPhysicalValueY
461 (0018,602c) FD 0.002 # 8, 1 PhysicalDeltaX
462 (0018,602e) FD 0.002 # 8, 1 PhysicalDeltaY
463 (fffe,e00d) na (ItemDelimitationItem for re-encoding) # 0, 0 ItemDelimitationItem
464 (fffe,e0dd) na (SequenceDelimitationItem for re-encod.) # 0, 0 SequenceDelimitationItem
465 */
466 const Tag tsqusreg(0x0018,0x6011);
467 SmartPointer<SequenceOfItems> sqi = InsertOrReplaceSQ( ds, tsqusreg);
468 if( !sqi->GetNumberOfItems() )
469 {
470 Item item; //( Tag(0xfffe,0xe000) );
471 sqi->AddItem( item );
472 }
473 Item &item1 = sqi->GetItem(1);
474 item1.SetVLToUndefined();
475 DataSet &subds = item1.GetNestedDataSet();
476
477 // (0018,602c) FD 0.002 # 8, 1 PhysicalDeltaX
478 // (0018,602e) FD 0.002 # 8, 1 PhysicalDeltaY
479 Attribute<0x0018,0x602c> at1;
480 at1.SetValue( spacing[0] );
481 Attribute<0x0018,0x602e> at2;
482 at2.SetValue( spacing[1] );
483 subds.Replace( at1.GetAsDataElement() );
484 subds.Replace( at2.GetAsDataElement() );
485
486 }
487
488
489 /* Enhanced stuff looks like:
490
491 (0020,9113) SQ (Sequence with undefined length #=1) # u/l, 1 PlanePositionSequence
492 (fffe,e000) na (Item with undefined length #=1) # u/l, 1 Item
493 (0020,0032) DS [73.5100815890831\-129.65028828174\189.777023529388] # 50, 3 ImagePositionPatient
494 (fffe,e00d) na (ItemDelimitationItem) # 0, 0 ItemDelimitationItem
495 (fffe,e0dd) na (SequenceDelimitationItem) # 0, 0 SequenceDelimitationItem
496 (0020,9116) SQ (Sequence with undefined length #=1) # u/l, 1 PlaneOrientationSequence
497 (fffe,e000) na (Item with undefined length #=1) # u/l, 1 Item
498 (0020,0037) DS [0.01604138687252\0.99942564964294\-0.0298495516180\0.0060454937629... # 102, 6 ImageOrientationPatient
499 (fffe,e00d) na (ItemDelimitationItem) # 0, 0 ItemDelimitationItem
500 (fffe,e0dd) na (SequenceDelimitationItem) # 0, 0 SequenceDelimitationItem
501 (0028,9110) SQ (Sequence with undefined length #=1) # u/l, 1 PixelMeasuresSequence
502 (fffe,e000) na (Item with undefined length #=2) # u/l, 1 Item
503 (0018,0050) DS [1] # 2, 1 SliceThickness
504 (0028,0030) DS [0.83333331346511\0.83333331346511] # 34, 2 PixelSpacing
505 (fffe,e00d) na (ItemDelimitationItem) # 0, 0 ItemDelimitationItem
506 (fffe,e0dd) na (SequenceDelimitationItem) # 0, 0 SequenceDelimitationItem
507 */
508
GetOriginValue(File const & f)509 std::vector<double> ImageHelper::GetOriginValue(File const & f)
510 {
511 std::vector<double> ori;
512 MediaStorage ms;
513 ms.SetFromFile(f);
514 const DataSet& ds = f.GetDataSet();
515
516 if( ms == MediaStorage::EnhancedCTImageStorage
517 || ms == MediaStorage::EnhancedMRImageStorage
518 || ms == MediaStorage::EnhancedMRColorImageStorage
519 || ms == MediaStorage::EnhancedPETImageStorage
520 || ms == MediaStorage::OphthalmicTomographyImageStorage
521 || ms == MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage
522 || ms == MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage
523 || ms == MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage
524 || ms == MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage
525 || ms == MediaStorage::XRay3DAngiographicImageStorage
526 || ms == MediaStorage::XRay3DCraniofacialImageStorage
527 || ms == MediaStorage::SegmentationStorage
528 || ms == MediaStorage::IVOCTForProcessing
529 || ms == MediaStorage::IVOCTForPresentation
530 || ms == MediaStorage::BreastTomosynthesisImageStorage
531 || ms == MediaStorage::BreastProjectionXRayImageStorageForPresentation
532 || ms == MediaStorage::BreastProjectionXRayImageStorageForProcessing
533 || ms == MediaStorage::LegacyConvertedEnhancedMRImageStorage
534 || ms == MediaStorage::LegacyConvertedEnhancedCTImageStorage
535 || ms == MediaStorage::LegacyConvertedEnhancedPETImageStorage )
536 {
537 const Tag t1(0x5200,0x9229);
538 const Tag t2(0x5200,0x9230);
539 if( GetOriginValueFromSequence(ds,t1, ori)
540 || GetOriginValueFromSequence(ds, t2, ori) )
541 {
542 assert( ori.size() == 3 );
543 return ori;
544 }
545 ori.resize( 3 );
546 gdcmWarningMacro( "Could not find Origin" );
547 return ori;
548 }
549 if( ms == MediaStorage::NuclearMedicineImageStorage )
550 {
551 const Tag t1(0x0054,0x0022);
552 if( ds.FindDataElement( t1 ) )
553 {
554 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( t1 ).GetValueAsSQ();
555 if( sqi && sqi->GetNumberOfItems() >= 1)
556 {
557 // Get first item:
558 const Item &item = sqi->GetItem(1);
559 const DataSet & subds = item.GetNestedDataSet();
560 const Tag timagepositionpatient(0x0020, 0x0032);
561 assert( subds.FindDataElement( timagepositionpatient ) );
562 Attribute<0x0020,0x0032> at = {{0,0,0}}; // default value if empty
563 at.SetFromDataSet( subds );
564 ori.resize( at.GetNumberOfValues() );
565 for( unsigned int i = 0; i < at.GetNumberOfValues(); ++i )
566 {
567 ori[i] = at.GetValue(i);
568 }
569 return ori;
570 }
571 }
572 ori.resize( 3 );
573 gdcmWarningMacro( "Could not find Origin" );
574 return ori;
575 }
576 ori.resize( 3 );
577
578 // else
579 const Tag timagepositionpatient(0x0020, 0x0032);
580 if( ms != MediaStorage::SecondaryCaptureImageStorage && ds.FindDataElement( timagepositionpatient ) )
581 {
582 const DataElement& de = ds.GetDataElement( timagepositionpatient );
583 Attribute<0x0020,0x0032> at = {{0,0,0}}; // default value if empty
584 at.SetFromDataElement( de );
585 for( unsigned int i = 0; i < at.GetNumberOfValues(); ++i )
586 {
587 ori[i] = at.GetValue(i);
588 }
589 }
590 else
591 {
592 ori[0] = 0;
593 ori[1] = 0;
594 ori[2] = 0;
595 }
596 assert( ori.size() == 3 );
597 return ori;
598 }
599
GetDirectionCosinesFromDataSet(DataSet const & ds,std::vector<double> & dircos)600 bool ImageHelper::GetDirectionCosinesFromDataSet(DataSet const & ds, std::vector<double> & dircos)
601 {
602 // \precondition: this dataset is not a secondary capture !
603 // else
604 const Tag timageorientationpatient(0x0020, 0x0037);
605 if( ds.FindDataElement( timageorientationpatient ) /*&& !ds.GetDataElement( timageorientationpatient ).IsEmpty()*/ )
606 {
607 const DataElement& de = ds.GetDataElement( timageorientationpatient );
608 Attribute<0x0020,0x0037> at = {{1,0,0,0,1,0}}; // default value if empty
609 at.SetFromDataElement( de );
610 for( unsigned int i = 0; i < at.GetNumberOfValues(); ++i )
611 {
612 dircos[i] = at.GetValue(i);
613 }
614 DirectionCosines dc( &dircos[0] );
615 if( !dc.IsValid() )
616 {
617 dc.Normalize();
618 if( dc.IsValid() )
619 {
620 gdcmWarningMacro( "DirectionCosines was not normalized. Fixed" );
621 const double * p = dc;
622 dircos = std::vector<double>(p, p + 6);
623 //return dircos;
624 }
625 else
626 {
627 // PAPYRUS_CR_InvalidIOP.dcm
628 gdcmWarningMacro( "Could not get DirectionCosines. Will be set to unit vector." );
629 //dircos[0] = 1;
630 //dircos[1] = 0;
631 //dircos[2] = 0;
632 //dircos[3] = 0;
633 //dircos[4] = 1;
634 //dircos[5] = 0;
635 return false;
636 }
637 }
638 return true;
639 }
640 return false;
641 }
642
GetDirectionCosinesValue(File const & f)643 std::vector<double> ImageHelper::GetDirectionCosinesValue(File const & f)
644 {
645 std::vector<double> dircos;
646 MediaStorage ms;
647 ms.SetFromFile(f);
648 const DataSet& ds = f.GetDataSet();
649
650 if( ms == MediaStorage::EnhancedCTImageStorage
651 || ms == MediaStorage::EnhancedMRImageStorage
652 || ms == MediaStorage::EnhancedMRColorImageStorage
653 || ms == MediaStorage::EnhancedPETImageStorage
654 || ms == MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage
655 || ms == MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage
656 || ms == MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage
657 || ms == MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage
658 || ms == MediaStorage::XRay3DAngiographicImageStorage
659 || ms == MediaStorage::XRay3DCraniofacialImageStorage
660 || ms == MediaStorage::SegmentationStorage
661 || ms == MediaStorage::IVOCTForPresentation
662 || ms == MediaStorage::IVOCTForProcessing
663 || ms == MediaStorage::BreastTomosynthesisImageStorage
664 || ms == MediaStorage::BreastProjectionXRayImageStorageForPresentation
665 || ms == MediaStorage::BreastProjectionXRayImageStorageForProcessing
666 || ms == MediaStorage::LegacyConvertedEnhancedMRImageStorage
667 || ms == MediaStorage::LegacyConvertedEnhancedCTImageStorage
668 || ms == MediaStorage::LegacyConvertedEnhancedPETImageStorage )
669 {
670 const Tag t1(0x5200,0x9229);
671 const Tag t2(0x5200,0x9230);
672 if( GetDirectionCosinesValueFromSequence(ds,t1, dircos)
673 || GetDirectionCosinesValueFromSequence(ds, t2, dircos) )
674 {
675 assert( dircos.size() == 6 );
676 return dircos;
677 }
678 else
679 {
680 dircos.resize( 6 );
681 bool b2 = ImageHelper::GetDirectionCosinesFromDataSet(ds, dircos);
682 if( b2 )
683 {
684 gdcmWarningMacro( "Image Orientation (Patient) cannot be stored here!. Continuing" );
685 }
686 else
687 {
688 gdcmErrorMacro( "Image Orientation (Patient) was not found" );
689 dircos[0] = 1;
690 dircos[1] = 0;
691 dircos[2] = 0;
692 dircos[3] = 0;
693 dircos[4] = 1;
694 dircos[5] = 0;
695 }
696 return dircos;
697 }
698 }
699 if( ms == MediaStorage::NuclearMedicineImageStorage )
700 {
701 const Tag t1(0x0054,0x0022);
702 if( ds.FindDataElement( t1 ) )
703 {
704 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( t1 ).GetValueAsSQ();
705 if( sqi && sqi->GetNumberOfItems() >= 1 )
706 {
707 // Get first item:
708 const Item &item = sqi->GetItem(1);
709 const DataSet & subds = item.GetNestedDataSet();
710
711 dircos.resize( 6 );
712 bool b2 = ImageHelper::GetDirectionCosinesFromDataSet(subds, dircos);
713 if( b2 )
714 {
715 }
716 else
717 {
718 gdcmErrorMacro( "Image Orientation (Patient) was not found" );
719 dircos[0] = 1;
720 dircos[1] = 0;
721 dircos[2] = 0;
722 dircos[3] = 0;
723 dircos[4] = 1;
724 dircos[5] = 0;
725 }
726 return dircos;
727 }
728 }
729 }
730
731 dircos.resize( 6 );
732 if( ms == MediaStorage::SecondaryCaptureImageStorage || !GetDirectionCosinesFromDataSet(ds, dircos) )
733 {
734 dircos[0] = 1;
735 dircos[1] = 0;
736 dircos[2] = 0;
737 dircos[3] = 0;
738 dircos[4] = 1;
739 dircos[5] = 0;
740 }
741
742 assert( dircos.size() == 6 );
743 return dircos;
744 }
745
SetForceRescaleInterceptSlope(bool b)746 void ImageHelper::SetForceRescaleInterceptSlope(bool b)
747 {
748 ForceRescaleInterceptSlope = b;
749 }
750
GetForceRescaleInterceptSlope()751 bool ImageHelper::GetForceRescaleInterceptSlope()
752 {
753 return ForceRescaleInterceptSlope;
754 }
755
SetPMSRescaleInterceptSlope(bool b)756 void ImageHelper::SetPMSRescaleInterceptSlope(bool b)
757 {
758 PMSRescaleInterceptSlope = b;
759 }
760
GetPMSRescaleInterceptSlope()761 bool ImageHelper::GetPMSRescaleInterceptSlope()
762 {
763 return PMSRescaleInterceptSlope;
764 }
765
SetForcePixelSpacing(bool b)766 void ImageHelper::SetForcePixelSpacing(bool b)
767 {
768 ForcePixelSpacing = b;
769 }
770
GetForcePixelSpacing()771 bool ImageHelper::GetForcePixelSpacing()
772 {
773 return ForcePixelSpacing;
774 }
775
GetRescaleInterceptSlopeValueFromDataSet(const DataSet & ds,std::vector<double> & interceptslope)776 bool GetRescaleInterceptSlopeValueFromDataSet(const DataSet& ds, std::vector<double> & interceptslope)
777 {
778 Attribute<0x0028,0x1052> at1;
779 bool intercept = ds.FindDataElement(at1.GetTag());
780 if( intercept )
781 {
782 if( !ds.GetDataElement(at1.GetTag()).IsEmpty() )
783 {
784 at1.SetFromDataElement( ds.GetDataElement(at1.GetTag()) );
785 interceptslope[0] = at1.GetValue();
786 }
787 }
788 Attribute<0x0028,0x1053> at2;
789 bool slope = ds.FindDataElement(at2.GetTag());
790 if ( slope )
791 {
792 if( !ds.GetDataElement(at2.GetTag()).IsEmpty() )
793 {
794 at2.SetFromDataElement( ds.GetDataElement(at2.GetTag()) );
795 interceptslope[1] = at2.GetValue();
796 if( interceptslope[1] == 0 )
797 {
798 // come' on ! WTF
799 gdcmDebugMacro( "Cannot have slope == 0. Defaulting to 1.0 instead" );
800 interceptslope[1] = 1;
801 }
802 }
803 }
804 return intercept || slope;
805 }
806
807
808 /// This function returns pixel information about an image from its dataset
809 /// That includes samples per pixel and bit depth (in that order)
810 /// Returns a PixelFormat
GetPixelFormatValue(const File & f)811 PixelFormat ImageHelper::GetPixelFormatValue(const File& f)
812 {
813 // D 0028|0011 [US] [Columns] [512]
814 //[10/20/10 9:05:07 AM] Mathieu Malaterre:
815 PixelFormat pf;
816 const DataSet& ds = f.GetDataSet();
817 // D 0028|0100 [US] [Bits Allocated] [16]
818 {
819 //const DataElement& de = ds.GetDataElement( Tag(0x0028, 0x0100) );
820 Attribute<0x0028,0x0100> at = { 0 };
821 at.SetFromDataSet( ds );
822 pf.SetBitsAllocated( at.GetValue() );
823 }
824 // D 0028|0101 [US] [Bits Stored] [12]
825 {
826 //const DataElement& de = ds.GetDataElement( Tag(0x0028, 0x0101) );
827 Attribute<0x0028,0x0101> at = { 0 };
828 at.SetFromDataSet( ds );
829 pf.SetBitsStored( at.GetValue() );
830 }
831 // D 0028|0102 [US] [High Bit] [11]
832 {
833 //const DataElement& de = ds.GetDataElement( Tag(0x0028, 0x0102) );
834 Attribute<0x0028,0x0102> at = { 0 };
835 at.SetFromDataSet( ds );
836 pf.SetHighBit( at.GetValue() );
837 }
838 // D 0028|0103 [US] [Pixel Representation] [0]
839 {
840 //const DataElement& de = ds.GetDataElement( Tag(0x0028, 0x0103) );
841 Attribute<0x0028,0x0103> at = { 0 };
842 at.SetFromDataSet( ds );
843 pf.SetPixelRepresentation( at.GetValue() );
844 }
845 // (0028,0002) US 1 # 2, 1 SamplesPerPixel
846 {
847 //if( ds.FindDataElement( Tag(0x0028, 0x0002) ) )
848 {
849 //const DataElement& de = ds.GetDataElement( Tag(0x0028, 0x0002) );
850 Attribute<0x0028,0x0002> at = { 1 };
851 at.SetFromDataSet( ds );
852 pf.SetSamplesPerPixel( at.GetValue() );
853 }
854 // else pf will default to 1...
855 }
856 return pf;
857
858 }
859 /// This function checks tags (0x0028, 0x0010) and (0x0028, 0x0011) for the
860 /// rows and columns of the image in pixels (as opposed to actual distances).
861 /// Also checks 0054, 0081 for the number of z slices for a 3D image
862 /// If that tag is not present, default the z dimension to 1
GetDimensionsValue(const File & f)863 std::vector<unsigned int> ImageHelper::GetDimensionsValue(const File& f)
864 {
865 DataSet const & ds = f.GetDataSet();
866
867 MediaStorage ms;
868 ms.SetFromFile(f);
869 std::vector<unsigned int> theReturn(3);
870 #if 0
871 if( ms == MediaStorage::VLWholeSlideMicroscopyImageStorage )
872 {
873 {
874 Attribute<0x0048,0x0006> at = { 0 };
875 at.SetFromDataSet( ds );
876 theReturn[0] = at.GetValue();
877 }
878 {
879 Attribute<0x0048,0x0007> at = { 0 };
880 at.SetFromDataSet( ds );
881 theReturn[1] = at.GetValue();
882 }
883 theReturn[2] = 1;
884 }
885 else
886 #endif
887 {
888 {
889 Attribute<0x0028,0x0011> at = { 0 }; // Columns
890 at.SetFromDataSet( ds );
891 theReturn[0] = at.GetValue();
892 }
893 {
894 Attribute<0x0028,0x0010> at = { 0 }; // Rows
895 at.SetFromDataSet( ds );
896 theReturn[1] = at.GetValue();
897 }
898 {
899 Attribute<0x0028,0x0008> at = { 0 };
900 at.SetFromDataSet( ds );
901 int numberofframes = at.GetValue();
902 theReturn[2] = 1;
903 if( numberofframes > 1 )
904 {
905 theReturn[2] = at.GetValue();
906 }
907 }
908 // ACR-NEMA legacy
909 {
910 Attribute<0x0028,0x0005> at = { 0 };
911 if( ds.FindDataElement( at.GetTag() ) )
912 {
913 const DataElement &de = ds.GetDataElement( at.GetTag() );
914 // SIEMENS_MAGNETOM-12-MONO2-Uncompressed.dcm picks VR::SS instead...
915 if( at.GetVR().Compatible( de.GetVR() ) )
916 {
917 at.SetFromDataSet( ds );
918 int imagedimensions = at.GetValue();
919 if( imagedimensions == 3 )
920 {
921 Attribute<0x0028,0x0012> at2 = { 0 };
922 at2.SetFromDataSet( ds );
923 theReturn[2] = at2.GetValue();
924 }
925 }
926 else
927 {
928 gdcmWarningMacro( "Sorry cannot read attribute (wrong VR): " << at.GetTag() );
929 }
930 }
931 }
932 }
933
934 return theReturn;
935 }
936
SetDimensionsValue(File & f,const Pixmap & img)937 void ImageHelper::SetDimensionsValue(File& f, const Pixmap & img)
938 {
939 const unsigned int *dims = img.GetDimensions();
940 MediaStorage ms;
941 ms.SetFromFile(f);
942 DataSet& ds = f.GetDataSet();
943 assert( MediaStorage::IsImage( ms ) );
944 {
945 Attribute<0x0028,0x0010> rows;
946 rows.SetValue( (uint16_t)dims[1] );
947 ds.Replace( rows.GetAsDataElement() );
948 Attribute<0x0028,0x0011> columns;
949 columns.SetValue( (uint16_t)dims[0] );
950 ds.Replace( columns.GetAsDataElement() );
951 Attribute<0x0028,0x0008> numframes = { 0 };
952 numframes.SetValue( dims[2] );
953 if( img.GetNumberOfDimensions() == 3 && dims[2] >= 1 )
954 {
955 if( ms.MediaStorage::GetModalityDimension() > 2 )
956 ds.Replace( numframes.GetAsDataElement() );
957 else if( ms.MediaStorage::GetModalityDimension() == 2 && dims[2] == 1 )
958 ds.Remove( numframes.GetTag() );
959 else
960 {
961 gdcmErrorMacro( "MediaStorage does not allow 3rd dimension. But value is: " << dims[2] );
962 gdcmAssertAlwaysMacro( "Could not set third dimension" );
963 }
964 }
965 else if( img.GetNumberOfDimensions() == 2 && dims[2] == 1 )
966 {
967 // This is a MF instances, need to set Number of Frame to 1 when Required
968 if (ms.MediaStorage::GetModalityDimension() > 2)
969 {
970 // Only include Multi-Frame when required (not Conditional):
971 if( ms == MediaStorage::XRayAngiographicImageStorage // A.14.3 XA Image IOD Module Table: Multi-frame C.7.6.6 C - Required if pixel data is Multi - frame Cine data
972 )
973 {
974 ds.Remove(numframes.GetTag());
975 }
976 else
977 {
978 ds.Replace(numframes.GetAsDataElement());
979 }
980 }
981 }
982 else // cleanup
983 ds.Remove( numframes.GetTag() );
984 }
985 // cleanup pass:
986 if( ms == MediaStorage::EnhancedCTImageStorage
987 || ms == MediaStorage::EnhancedMRImageStorage
988 || ms == MediaStorage::EnhancedMRColorImageStorage
989 || ms == MediaStorage::EnhancedPETImageStorage
990 || ms == MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage
991 || ms == MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage
992 || ms == MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage
993 || ms == MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage
994 || ms == MediaStorage::XRay3DAngiographicImageStorage
995 || ms == MediaStorage::XRay3DCraniofacialImageStorage
996 || ms == MediaStorage::SegmentationStorage
997 || ms == MediaStorage::IVOCTForProcessing
998 || ms == MediaStorage::IVOCTForPresentation
999 || ms == MediaStorage::BreastTomosynthesisImageStorage
1000 || ms == MediaStorage::BreastProjectionXRayImageStorageForPresentation
1001 || ms == MediaStorage::BreastProjectionXRayImageStorageForProcessing
1002 || ms == MediaStorage::LegacyConvertedEnhancedMRImageStorage
1003 || ms == MediaStorage::LegacyConvertedEnhancedCTImageStorage
1004 || ms == MediaStorage::LegacyConvertedEnhancedPETImageStorage )
1005 {
1006 const Tag tfgs(0x5200,0x9230);
1007 if( ds.FindDataElement( tfgs ) )
1008 {
1009 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( tfgs ).GetValueAsSQ();
1010 assert( sqi );
1011 sqi->SetLengthToUndefined();
1012 sqi->SetNumberOfItems( dims[2] );
1013 }
1014 }
1015
1016 }
1017
GetRescaleInterceptSlopeValue(File const & f)1018 std::vector<double> ImageHelper::GetRescaleInterceptSlopeValue(File const & f)
1019 {
1020 std::vector<double> interceptslope;
1021 MediaStorage ms;
1022 ms.SetFromFile(f);
1023 const DataSet& ds = f.GetDataSet();
1024
1025 if( ms == MediaStorage::EnhancedCTImageStorage
1026 || ms == MediaStorage::EnhancedMRImageStorage
1027 /*|| ms == MediaStorage::EnhancedMRColorImageStorage*/
1028 || ms == MediaStorage::EnhancedPETImageStorage
1029 || ms == MediaStorage::XRay3DAngiographicImageStorage
1030 || ms == MediaStorage::XRay3DCraniofacialImageStorage
1031 || ms == MediaStorage::SegmentationStorage
1032 || ms == MediaStorage::BreastTomosynthesisImageStorage
1033 || ms == MediaStorage::BreastProjectionXRayImageStorageForPresentation
1034 || ms == MediaStorage::BreastProjectionXRayImageStorageForProcessing
1035 || ms == MediaStorage::LegacyConvertedEnhancedMRImageStorage
1036 || ms == MediaStorage::LegacyConvertedEnhancedCTImageStorage
1037 || ms == MediaStorage::LegacyConvertedEnhancedPETImageStorage )
1038 {
1039 const Tag t1(0x5200,0x9229);
1040 const Tag t2(0x5200,0x9230);
1041 if( GetInterceptSlopeValueFromSequence(ds,t1, interceptslope)
1042 || GetInterceptSlopeValueFromSequence(ds,t2, interceptslope) )
1043 {
1044 assert( interceptslope.size() == 2 );
1045 return interceptslope;
1046 }
1047
1048 // Workaround to get intercept/slope from some Philips
1049 // XRay3DAngiographic files
1050 if (ms == MediaStorage::XRay3DAngiographicImageStorage && ForceRescaleInterceptSlope)
1051 {
1052 const Tag t3(0x0018,0x9530);
1053 if(ds.FindDataElement( t3 ))
1054 {
1055 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( t3 ).GetValueAsSQ();
1056 if(sqi && sqi->GetNumberOfItems() > 0)
1057 {
1058 const Item &item = sqi->GetItem(1);
1059 const DataSet & subds = item.GetNestedDataSet();
1060 const Tag tpi(0x0028,0x1052);
1061 const Tag tps(0x0028,0x1053);
1062 if( subds.FindDataElement(tps) && subds.FindDataElement(tpi))
1063 {
1064 const DataElement &dei = subds.GetDataElement( tpi );
1065 Attribute<0x0028,0x1052> ati;
1066 ati.SetFromDataElement( dei );
1067 interceptslope.push_back( ati.GetValue() );
1068 const DataElement &des = subds.GetDataElement( tps );
1069 Attribute<0x0028,0x1053> ats;
1070 ats.SetFromDataElement( des );
1071 interceptslope.push_back( ats.GetValue() );
1072 return interceptslope;
1073 }
1074 }
1075 }
1076 }
1077
1078 //else
1079 // {
1080 // interceptslope.resize( 2 );
1081 // interceptslope[0] = 0;
1082 // interceptslope[1] = 1;
1083 // bool b = GetRescaleInterceptSlopeValueFromDataSet(ds, interceptslope);
1084 // gdcmAssertMacro( b ); (void)b;
1085 // return interceptslope;
1086 // }
1087 }
1088
1089 // else
1090 interceptslope.resize( 2 );
1091 interceptslope[0] = 0;
1092 interceptslope[1] = 1;
1093 if( ms == MediaStorage::CTImageStorage
1094 || ms == MediaStorage::ComputedRadiographyImageStorage
1095 || ms == MediaStorage::PETImageStorage
1096 || ms == MediaStorage::SecondaryCaptureImageStorage
1097 || ms == MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage
1098 || ms == MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage
1099 || ms == MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage
1100 || ms == MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage
1101 || ForceRescaleInterceptSlope
1102 )
1103 {
1104 bool b = GetRescaleInterceptSlopeValueFromDataSet(ds, interceptslope);
1105 if( !b )
1106 {
1107 gdcmDebugMacro( "No Modality LUT found (Rescale Intercept/Slope)" );
1108 }
1109 }
1110 else if ( ms == MediaStorage::MRImageStorage )
1111 {
1112 #if 0
1113 const Tag trwvms(0x0040,0x9096); // Real World Value Mapping Sequence
1114 if( ds.FindDataElement( trwvms ) )
1115 {
1116 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( trwvms ).GetValueAsSQ();
1117 if( sqi )
1118 {
1119 const Tag trwvlutd(0x0040,0x9212); // Real World Value LUT Data
1120 if( ds.FindDataElement( trwvlutd ) )
1121 {
1122 gdcmAssertAlwaysMacro(0); // Not supported !
1123 }
1124 // don't know how to handle multiples:
1125 gdcmAssertAlwaysMacro( sqi->GetNumberOfItems() == 1 );
1126 const Item &item = sqi->GetItem(1);
1127 const DataSet & subds = item.GetNestedDataSet();
1128 //const Tag trwvi(0x0040,0x9224); // Real World Value Intercept
1129 //const Tag trwvs(0x0040,0x9225); // Real World Value Slope
1130 Attribute<0x0040,0x9224> at1 = {0};
1131 at1.SetFromDataSet( subds );
1132 Attribute<0x0040,0x9225> at2 = {1};
1133 at2.SetFromDataSet( subds );
1134 interceptslope[0] = at1.GetValue();
1135 interceptslope[1] = at2.GetValue();
1136 }
1137 }
1138 #else
1139 // See the long thread at:
1140 // https://groups.google.com/d/msg/comp.protocols.dicom/M4kdqcrs50Y/_TSx0EjtAQAJ
1141 // in particular this paper:
1142 // Errors in Quantitative Image Analysis due to Platform-Dependent Image Scaling
1143 // http://www.ncbi.nlm.nih.gov/pmc/articles/PMC3998685/
1144 const PrivateTag tpriv_rescaleintercept( 0x2005,0x09,"Philips MR Imaging DD 005" );
1145 const PrivateTag tpriv_rescaleslope( 0x2005,0x0a,"Philips MR Imaging DD 005" );
1146 if( ds.FindDataElement( tpriv_rescaleintercept ) && ds.FindDataElement( tpriv_rescaleslope ) )
1147 {
1148 // The following will work out of the box for Philips whether or not
1149 // "Combine MR Rescaling" was set:
1150 // PMS DICOM CS states that Modality LUT for MR Image Storage is to be
1151 // used for image processing. VOI LUT are always recomputed, so output
1152 // from GDCM may not look right for display (sorry!)
1153 const DataElement &priv_rescaleintercept = ds.GetDataElement( tpriv_rescaleintercept );
1154 const DataElement &priv_rescaleslope = ds.GetDataElement( tpriv_rescaleslope );
1155 Element<VR::DS,VM::VM1> el_ri = {{ 0 }};
1156 el_ri.SetFromDataElement( priv_rescaleintercept );
1157 Element<VR::DS,VM::VM1> el_rs = {{ 1 }};
1158 el_rs.SetFromDataElement( priv_rescaleslope );
1159 if( PMSRescaleInterceptSlope )
1160 {
1161 interceptslope[0] = el_ri.GetValue();
1162 interceptslope[1] = el_rs.GetValue();
1163 if( interceptslope[1] == 0 )
1164 interceptslope[1] = 1;
1165 gdcmWarningMacro( "PMS Modality LUT loaded for MR Image Storage: [" << interceptslope[0] << "," << interceptslope[1] << "]" );
1166 }
1167 }
1168 else
1169 {
1170 std::vector<double> dummy(2);
1171 if( GetRescaleInterceptSlopeValueFromDataSet(ds, dummy) )
1172 {
1173 // for everyone else, read your DCS, and set: ForceRescaleInterceptSlope = true if needed
1174 gdcmDebugMacro( "Modality LUT unused for MR Image Storage: [" << dummy[0] << "," << dummy[1] << "]" );
1175 }
1176 }
1177 #endif
1178 }
1179 else if (
1180 ms == MediaStorage::RTDoseStorage
1181 )
1182 {
1183 // TODO. Should I check FrameIncrementPointer ? (0028,0009) AT (3004,000c)
1184 Attribute<0x3004,0x000e> gridscaling = { 0 };
1185 gridscaling.SetFromDataSet( ds );
1186 interceptslope[0] = 0;
1187 interceptslope[1] = gridscaling.GetValue();
1188 if( interceptslope[1] == 0 )
1189 {
1190 // come' on ! WTF
1191 gdcmWarningMacro( "Cannot have slope == 0. Defaulting to 1.0 instead" );
1192 interceptslope[1] = 1;
1193 }
1194 }
1195
1196 // \post condition slope can never be 0:
1197 assert( interceptslope[1] != 0. );
1198 return interceptslope;
1199 }
1200
GetSpacingTagFromMediaStorage(MediaStorage const & ms)1201 Tag ImageHelper::GetSpacingTagFromMediaStorage(MediaStorage const &ms)
1202 {
1203 Tag t;
1204
1205 // gdcmData/MR00010001.dcm => GeneralElectricMagneticResonanceImageStorage
1206 // (0018,0088) DS [] # 4, 1 SpacingBetweenSlices
1207 // (0028,0030) DS [ 0.8593750000\0.8593750000] # 26, 2 PixelSpacing
1208 switch(ms)
1209 {
1210 case MediaStorage::EnhancedUSVolumeStorage:
1211 // Enhanced stuff are handled elsewhere... look carefully :)
1212 //case MediaStorage::EnhancedMRImageStorage:
1213 //case MediaStorage::EnhancedCTImageStorage:
1214 //case MediaStorage::XRay3DAngiographicImageStorage
1215 //case MediaStorage::XRay3DCraniofacialImageStorage
1216 // gdcmWarningMacro( "Enhanced image are not currently supported. Spacing will be wrong" );
1217 case MediaStorage::CTImageStorage:
1218 case MediaStorage::MRImageStorage:
1219 case MediaStorage::RTDoseStorage:
1220 case MediaStorage::NuclearMedicineImageStorage:
1221 case MediaStorage::PETImageStorage:
1222 case MediaStorage::GeneralElectricMagneticResonanceImageStorage:
1223 case MediaStorage::PhilipsPrivateMRSyntheticImageStorage:
1224 case MediaStorage::VLPhotographicImageStorage: // VL Image IOD
1225 case MediaStorage::VLMicroscopicImageStorage:
1226 case MediaStorage::IVOCTForProcessing:
1227 case MediaStorage::IVOCTForPresentation:
1228 // (0028,0030) DS [2.0\2.0] # 8, 2 PixelSpacing
1229 t = Tag(0x0028,0x0030);
1230 break;
1231 case MediaStorage::ComputedRadiographyImageStorage: // See pixelspacingtestimages/DISCIMG/IMAGES/CRIMAGE
1232 case MediaStorage::DigitalXRayImageStorageForPresentation:
1233 case MediaStorage::DigitalXRayImageStorageForProcessing:
1234 case MediaStorage::DigitalMammographyImageStorageForPresentation:
1235 case MediaStorage::DigitalMammographyImageStorageForProcessing:
1236 case MediaStorage::DigitalIntraoralXrayImageStorageForPresentation:
1237 case MediaStorage::DigitalIntraoralXRayImageStorageForProcessing:
1238 case MediaStorage::XRayAngiographicImageStorage:
1239 case MediaStorage::XRayRadiofluoroscopingImageStorage:
1240 case MediaStorage::XRayAngiographicBiPlaneImageStorageRetired:
1241 case MediaStorage::FujiPrivateMammoCRImageStorage:
1242 // (0018,1164) DS [0.5\0.5] # 8, 2 ImagerPixelSpacing
1243 t = Tag(0x0018,0x1164);
1244 break;
1245 case MediaStorage::RTImageStorage: // gdcmData/SYNGORTImage.dcm
1246 t = Tag(0x3002,0x0011); // ImagePlanePixelSpacing
1247 break;
1248 case MediaStorage::SecondaryCaptureImageStorage:
1249 case MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage:
1250 case MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage:
1251 case MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage:
1252 case MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage:
1253 // See PS 3.3-2008. Table C.8-25 SC IMAGE MODULE ATTRIBUTES
1254 // and Table C.8-25b SC MULTI-FRAME IMAGE MODULE ATTRIBUTES
1255 t = Tag(0x0018,0x2010);
1256 break;
1257 case MediaStorage::HardcopyGrayscaleImageStorage:
1258 case MediaStorage::HardcopyColorImageStorage:
1259 t = Tag(0x0018,0x2010); // Nominal Scanned Pixel Spacing
1260 break;
1261 case MediaStorage::GEPrivate3DModelStorage: // FIXME FIXME !!!
1262 case MediaStorage::Philips3D:
1263 case MediaStorage::VideoEndoscopicImageStorage:
1264 gdcmWarningMacro( "FIXME" );
1265 t = Tag(0xffff,0xffff);
1266 break;
1267 case MediaStorage::UltrasoundMultiFrameImageStorage:
1268 // gdcmData/US-MONO2-8-8x-execho.dcm
1269 // this should be handled somewhere else
1270 //assert(0);
1271 gdcmWarningMacro( "FIXME" );
1272 t = Tag(0xffff,0xffff);
1273 break;
1274 case MediaStorage::UltrasoundImageStorage: // ??
1275 case MediaStorage::UltrasoundImageStorageRetired:
1276 case MediaStorage::UltrasoundMultiFrameImageStorageRetired:
1277 // (0028,0034) IS [4\3] # 4, 2 PixelAspectRatio
1278 t = Tag(0xffff,0xffff); // FIXME
1279 t = Tag(0x0028,0x0034); // FIXME
1280 break;
1281 default:
1282 gdcmDebugMacro( "Do not handle: " << ms );
1283 t = Tag(0xffff,0xffff);
1284 break;
1285 }
1286
1287 // should only override unless Modality set it already
1288 // basically only Secondary Capture should reach that point
1289 if( ForcePixelSpacing && t == Tag(0xffff,0xffff) )
1290 {
1291 t = Tag(0x0028,0x0030);
1292 }
1293
1294 return t;
1295 }
1296
GetZSpacingTagFromMediaStorage(MediaStorage const & ms)1297 Tag ImageHelper::GetZSpacingTagFromMediaStorage(MediaStorage const &ms)
1298 {
1299 Tag t;
1300
1301 switch(ms)
1302 {
1303 case MediaStorage::EnhancedUSVolumeStorage:
1304 case MediaStorage::MRImageStorage:
1305 case MediaStorage::NuclearMedicineImageStorage: // gdcmData/Nm.dcm
1306 case MediaStorage::GeneralElectricMagneticResonanceImageStorage:
1307 // (0018,0088) DS [3] # 2, 1 SpacingBetweenSlices
1308 t = Tag(0x0018,0x0088);
1309 break;
1310 // No spacing AFAIK for those:
1311 /*
1312 $ dciodvfy gdcmData/D_CLUNIE_CT1_JPLL.dcm
1313 CTImage
1314 Warning - Attribute is not present in standard DICOM IOD - (0x0018,0x0088) DS Spacing Between Slices
1315 Warning - Dicom dataset contains attributes not present in standard DICOM IOD - this is a Standard Extended SOP Class
1316 */
1317 case MediaStorage::PETImageStorage: // ??
1318 case MediaStorage::CTImageStorage:
1319 case MediaStorage::RTImageStorage:
1320 // ImagerPixelSpacing section:
1321 case MediaStorage::ComputedRadiographyImageStorage: // ??
1322 case MediaStorage::DigitalXRayImageStorageForPresentation:
1323 case MediaStorage::DigitalXRayImageStorageForProcessing:
1324 case MediaStorage::DigitalMammographyImageStorageForPresentation:
1325 case MediaStorage::DigitalMammographyImageStorageForProcessing:
1326 case MediaStorage::DigitalIntraoralXrayImageStorageForPresentation:
1327 case MediaStorage::DigitalIntraoralXRayImageStorageForProcessing:
1328 case MediaStorage::XRayAngiographicImageStorage:
1329 case MediaStorage::XRayRadiofluoroscopingImageStorage:
1330 case MediaStorage::XRayAngiographicBiPlaneImageStorageRetired:
1331 // US:
1332 case MediaStorage::UltrasoundImageStorage:
1333 case MediaStorage::UltrasoundMultiFrameImageStorage:
1334 case MediaStorage::UltrasoundImageStorageRetired:
1335 case MediaStorage::UltrasoundMultiFrameImageStorageRetired:
1336 // SC:
1337 case MediaStorage::SecondaryCaptureImageStorage:
1338 case MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage:
1339 case MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage:
1340 case MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage:
1341 case MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage:
1342 case MediaStorage::HardcopyGrayscaleImageStorage:
1343 t = Tag(0xffff,0xffff);
1344 break;
1345 case MediaStorage::RTDoseStorage: // gdcmData/BogugsItemAndSequenceLengthCorrected.dcm
1346 t = Tag(0x3004,0x000c);
1347 break;
1348 case MediaStorage::GEPrivate3DModelStorage:
1349 case MediaStorage::Philips3D:
1350 case MediaStorage::VideoEndoscopicImageStorage:
1351 gdcmWarningMacro( "FIXME" );
1352 t = Tag(0xffff,0xffff);
1353 break;
1354 default:
1355 gdcmDebugMacro( "Do not handle Z spacing for: " << ms );
1356 t = Tag(0xffff,0xffff);
1357 break;
1358 }
1359
1360 if( ForcePixelSpacing && t == Tag(0xffff,0xffff) )
1361 {
1362 t = Tag(0x0018,0x0088);
1363 }
1364 return t;
1365 }
1366
GetSpacingValue(File const & f)1367 std::vector<double> ImageHelper::GetSpacingValue(File const & f)
1368 {
1369 std::vector<double> sp;
1370 sp.reserve(3);
1371 MediaStorage ms;
1372 ms.SetFromFile(f);
1373 const DataSet& ds = f.GetDataSet();
1374
1375 if( ms == MediaStorage::EnhancedCTImageStorage
1376 || ms == MediaStorage::EnhancedMRImageStorage
1377 || ms == MediaStorage::EnhancedMRColorImageStorage
1378 || ms == MediaStorage::EnhancedPETImageStorage
1379 || ms == MediaStorage::OphthalmicTomographyImageStorage
1380 || ms == MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage
1381 || ms == MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage
1382 || ms == MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage
1383 || ms == MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage
1384 || ms == MediaStorage::XRay3DAngiographicImageStorage
1385 || ms == MediaStorage::XRay3DCraniofacialImageStorage
1386 || ms == MediaStorage::SegmentationStorage
1387 || ms == MediaStorage::IVOCTForProcessing
1388 || ms == MediaStorage::IVOCTForPresentation
1389 || ms == MediaStorage::BreastTomosynthesisImageStorage
1390 || ms == MediaStorage::BreastProjectionXRayImageStorageForPresentation
1391 || ms == MediaStorage::BreastProjectionXRayImageStorageForProcessing
1392 || ms == MediaStorage::LegacyConvertedEnhancedMRImageStorage
1393 || ms == MediaStorage::LegacyConvertedEnhancedCTImageStorage
1394 || ms == MediaStorage::LegacyConvertedEnhancedPETImageStorage)
1395 {
1396 // <entry group="5200" element="9230" vr="SQ" vm="1" name="Per-frame Functional Groups Sequence"/>
1397 const Tag t1(0x5200,0x9229);
1398 const Tag t2(0x5200,0x9230);
1399 if( GetSpacingValueFromSequence(ds,t1, sp)
1400 || GetSpacingValueFromSequence(ds, t2, sp) )
1401 {
1402 assert( sp.size() == 3 );
1403 return sp;
1404 }
1405 // Else.
1406 // How do I send an error ?
1407 sp.resize( 3 ); // FIXME !!
1408 sp[0] = 1.;
1409 sp[1] = 1.;
1410 sp[2] = 1.;
1411 gdcmWarningMacro( "Could not find Spacing" );
1412 return sp;
1413 }
1414 else if( ms == MediaStorage::UltrasoundMultiFrameImageStorage || ms == MediaStorage::UltrasoundImageStorage )
1415 {
1416 if( GetUltraSoundSpacingValueFromSequence(ds, sp) )
1417 {
1418 // 3rd dimension is too difficult to handle for now...
1419 // (0018,1065) DS [0\ 957\ 990\ 990\1023\1023\ 990\ 990\1023\ 957\1023\1023\1023\ 990... # 562,113 FrameTimeVector
1420 sp.push_back( 1.0 );
1421 return sp;
1422 }
1423 else
1424 {
1425 // TODO this one is easy:
1426 // (0028,0009) AT (0018,1063) # 4, 1 FrameIncrementPointer
1427 // -> (0018,1063) DS [76.000000] # 10, 1 FrameTime
1428
1429 gdcmWarningMacro( "No spacing value found" );
1430 sp.push_back( 1.0 );
1431 sp.push_back( 1.0 );
1432 sp.push_back( 1.0 );
1433 return sp;
1434 }
1435 }
1436
1437 Tag spacingtag = GetSpacingTagFromMediaStorage(ms);
1438 if( spacingtag != Tag(0xffff,0xffff) && ds.FindDataElement( spacingtag ) && !ds.GetDataElement( spacingtag ).IsEmpty() )
1439 {
1440 const DataElement& de = ds.GetDataElement( spacingtag );
1441 const Global &g = GlobalInstance;
1442 const Dicts &dicts = g.GetDicts();
1443 const DictEntry &entry = dicts.GetDictEntry(de.GetTag());
1444 const VR & vr = entry.GetVR();
1445 assert( vr.Compatible( de.GetVR() ) );
1446 switch(vr)
1447 {
1448 case VR::DS:
1449 {
1450 Element<VR::DS,VM::VM1_n> el;
1451 std::stringstream ss;
1452 const ByteValue *bv = de.GetByteValue();
1453 assert( bv );
1454 std::string s = std::string( bv->GetPointer(), bv->GetLength() );
1455 ss.str( s );
1456 // Stupid file: CT-MONO2-8-abdo.dcm
1457 // The spacing is something like that: [0.2\0\0.200000]
1458 // I would need to throw an exception that VM is not compatible
1459 el.SetLength( entry.GetVM().GetLength() * entry.GetVR().GetSizeof() );
1460 std::string::size_type found = s.find('\\');
1461 if( found != std::string::npos )
1462 {
1463 el.Read( ss );
1464 assert( el.GetLength() == 2 );
1465 for(unsigned int i = 0; i < el.GetLength(); ++i)
1466 {
1467 if( el.GetValue(i) )
1468 {
1469 sp.push_back( el.GetValue(i) );
1470 }
1471 else
1472 {
1473 gdcmWarningMacro( "Cannot have a spacing of 0" );
1474 sp.push_back( 1.0 );
1475 }
1476 }
1477 std::swap( sp[0], sp[1]);
1478 }
1479 else
1480 {
1481 double singleval;
1482 ss >> singleval;
1483 if( singleval == 0.0 )
1484 {
1485 singleval = 1.0;
1486 }
1487 sp.push_back( singleval );
1488 sp.push_back( singleval );
1489 }
1490 assert( sp.size() == (unsigned int)entry.GetVM() );
1491 }
1492 break;
1493 case VR::IS:
1494 {
1495 Element<VR::IS,VM::VM1_n> el;
1496 std::stringstream ss;
1497 const ByteValue *bv = de.GetByteValue();
1498 assert( bv );
1499 std::string s = std::string( bv->GetPointer(), bv->GetLength() );
1500 ss.str( s );
1501 el.SetLength( entry.GetVM().GetLength() * entry.GetVR().GetSizeof() );
1502 el.Read( ss );
1503 for(unsigned int i = 0; i < el.GetLength(); ++i)
1504 {
1505 if( el.GetValue(i) )
1506 sp.push_back( el.GetValue(i) );
1507 else
1508 {
1509 gdcmWarningMacro( "Cannot have a spacing of 0" );
1510 sp.push_back( 1.0 );
1511 }
1512 }
1513 std::swap( sp[0], sp[1]);
1514 assert( sp.size() == (unsigned int)entry.GetVM() );
1515 }
1516 break;
1517 default:
1518 assert(0);
1519 break;
1520 }
1521 }
1522 else
1523 {
1524 sp.push_back( 1.0 );
1525 sp.push_back( 1.0 );
1526 }
1527 assert( sp.size() == 2 );
1528 // Make sure multiframe:
1529 std::vector<unsigned int> dims = ImageHelper::GetDimensionsValue( f );
1530
1531 // Do Z:
1532 Tag zspacingtag = ImageHelper::GetZSpacingTagFromMediaStorage(ms);
1533 if( zspacingtag != Tag(0xffff,0xffff) && ds.FindDataElement( zspacingtag ) )
1534 {
1535 const DataElement& de = ds.GetDataElement( zspacingtag );
1536 if( de.IsEmpty() )
1537 {
1538 sp.push_back( 1.0 );
1539 }
1540 else
1541 {
1542 const Global &g = GlobalInstance;
1543 const Dicts &dicts = g.GetDicts();
1544 const DictEntry &entry = dicts.GetDictEntry(de.GetTag());
1545 const VR & vr = entry.GetVR();
1546 assert( de.GetVR() == vr || de.GetVR() == VR::INVALID || de.GetVR() == VR::UN );
1547 if( entry.GetVM() == VM::VM1 )
1548 {
1549 switch(vr)
1550 {
1551 case VR::DS:
1552 {
1553 Element<VR::DS,VM::VM1_n> el;
1554 std::stringstream ss;
1555 const ByteValue *bv = de.GetByteValue();
1556 assert( bv );
1557 std::string s = std::string( bv->GetPointer(), bv->GetLength() );
1558 ss.str( s );
1559 el.SetLength( entry.GetVM().GetLength() * entry.GetVR().GetSizeof() );
1560 el.Read( ss );
1561 for(unsigned int i = 0; i < el.GetLength(); ++i)
1562 {
1563 const double value = el.GetValue(i);
1564 sp.push_back( value );
1565 }
1566 //assert( sp.size() == entry.GetVM() );
1567 }
1568 break;
1569 default:
1570 assert(0);
1571 break;
1572 }
1573 }
1574 else
1575 {
1576 assert( entry.GetVM() == VM::VM2_n );
1577 assert( vr == VR::DS );
1578 Attribute<0x28,0x8> numberoframes;
1579 const DataElement& de1 = ds.GetDataElement( numberoframes.GetTag() );
1580 numberoframes.SetFromDataElement( de1 );
1581 Attribute<0x3004,0x000c> gridframeoffsetvector;
1582 const DataElement& de2 = ds.GetDataElement( gridframeoffsetvector.GetTag() );
1583 gridframeoffsetvector.SetFromDataElement( de2 );
1584 double v1 = gridframeoffsetvector[0];
1585 double v2 = gridframeoffsetvector[1];
1586 // FIXME. I should check consistency
1587 sp.push_back( v2 - v1 );
1588 }
1589 }
1590 }
1591 else if( ds.FindDataElement( Tag(0x0028,0x0009) ) ) // Frame Increment Pointer
1592 {
1593 const DataElement& de = ds.GetDataElement( Tag(0x0028,0x0009) );
1594 Attribute<0x0028,0x0009,VR::AT,VM::VM1> at;
1595 at.SetFromDataElement( de );
1596 assert( ds.FindDataElement( at.GetTag() ) );
1597 if( ds.FindDataElement( at.GetValue() ) )
1598 {
1599 /*
1600 $ dcmdump D_CLUNIE_NM1_JPLL.dcm" | grep 0028,0009
1601 (0028,0009) AT (0054,0010)\(0054,0020) # 8, 2 FrameIncrementPointer
1602 */
1603 const DataElement& de2 = ds.GetDataElement( at.GetValue() );
1604 if( at.GetValue() == Tag(0x0018,0x1063) && at.GetNumberOfValues() == 1 )
1605 {
1606 Attribute<0x0018,0x1063> at2;
1607 at2.SetFromDataElement( de2 );
1608 if( dims[2] > 1 )
1609 {
1610 sp.push_back( at2.GetValue() );
1611 }
1612 else
1613 {
1614 if( at2.GetValue() != 0. )
1615 {
1616 gdcmErrorMacro( "Number of Frame should be equal to 0" );
1617 sp.push_back( 0.0 );
1618 }
1619 else
1620 {
1621 sp.push_back( 1.0 );
1622 }
1623 }
1624 }
1625 else
1626 {
1627 gdcmWarningMacro( "Don't know how to handle spacing for: " << de );
1628 sp.push_back( 1.0 );
1629 }
1630 }
1631 else
1632 {
1633 gdcmErrorMacro( "Tag: " << at.GetTag() << " was found to point to missing"
1634 "Tag: " << at.GetValue() << " default to 1.0." );
1635 sp.push_back( 1.0 );
1636 }
1637 }
1638 else
1639 {
1640 sp.push_back( 1.0 );
1641 }
1642
1643 assert( sp.size() == 3 );
1644 assert( sp[0] != 0. );
1645 assert( sp[1] != 0. );
1646 //if( ms != MediaStorage::MRImageStorage )
1647 // assert( sp[2] != 0. );
1648 return sp;
1649 }
1650
SetSpacingValue(DataSet & ds,const std::vector<double> & spacing)1651 void ImageHelper::SetSpacingValue(DataSet & ds, const std::vector<double> & spacing)
1652 {
1653 MediaStorage ms;
1654 ms.SetFromDataSet(ds);
1655 assert( MediaStorage::IsImage( ms ) );
1656 if( ms == MediaStorage::SecondaryCaptureImageStorage )
1657 {
1658 Tag pixelspacing(0x0028,0x0030);
1659 Tag imagerpixelspacing(0x0018,0x1164);
1660 Tag spacingbetweenslice(0x0018,0x0088);
1661 //ds.Remove( pixelspacing );
1662 //ds.Remove( imagerpixelspacing );
1663 //ds.Remove( spacingbetweenslice );
1664 //return;
1665 }
1666 assert( spacing.size() == 3 );
1667
1668 if( ms == MediaStorage::EnhancedCTImageStorage
1669 || ms == MediaStorage::EnhancedMRImageStorage
1670 || ms == MediaStorage::EnhancedMRColorImageStorage
1671 || ms == MediaStorage::EnhancedPETImageStorage
1672 || ms == MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage
1673 || ms == MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage
1674 || ms == MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage
1675 || ms == MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage
1676 || ms == MediaStorage::XRay3DAngiographicImageStorage
1677 || ms == MediaStorage::XRay3DCraniofacialImageStorage
1678 || ms == MediaStorage::SegmentationStorage
1679 || ms == MediaStorage::IVOCTForPresentation
1680 || ms == MediaStorage::IVOCTForProcessing
1681 || ms == MediaStorage::BreastTomosynthesisImageStorage
1682 || ms == MediaStorage::BreastProjectionXRayImageStorageForPresentation
1683 || ms == MediaStorage::BreastProjectionXRayImageStorageForProcessing
1684 || ms == MediaStorage::LegacyConvertedEnhancedMRImageStorage
1685 || ms == MediaStorage::LegacyConvertedEnhancedCTImageStorage
1686 || ms == MediaStorage::LegacyConvertedEnhancedPETImageStorage )
1687 {
1688 /*
1689 (0028,9110) SQ (Sequence with undefined length #=1) # u/l, 1 PixelMeasuresSequence
1690 (fffe,e000) na (Item with undefined length #=2) # u/l, 1 Item
1691 (0018,0050) DS [5.00000] # 8, 1 SliceThickness
1692 (0028,0030) DS [0.820312\0.820312] # 18, 2 PixelSpacing
1693 (fffe,e00d) na (ItemDelimitationItem) # 0, 0 ItemDelimitationItem
1694 (fffe,e0dd) na (SequenceDelimitationItem) # 0, 0 SequenceDelimitationItem
1695 */
1696 {
1697 const Tag tfgs(0x5200,0x9229);
1698 SmartPointer<SequenceOfItems> sqi = InsertOrReplaceSQ( ds, tfgs );
1699 if( !sqi->GetNumberOfItems() )
1700 {
1701 Item item; //( Tag(0xfffe,0xe000) );
1702 sqi->AddItem( item );
1703 }
1704 Item &item1 = sqi->GetItem(1);
1705 item1.SetVLToUndefined();
1706 DataSet &subds = item1.GetNestedDataSet();
1707 const Tag tpms(0x0028,0x9110);
1708 sqi = InsertOrReplaceSQ( subds, tpms );
1709 if( !sqi->GetNumberOfItems() )
1710 {
1711 Item item; //( Tag(0xfffe,0xe000) );
1712 sqi->AddItem( item );
1713 }
1714 Item &item2 = sqi->GetItem(1);
1715 item2.SetVLToUndefined();
1716 DataSet &subds2 = item2.GetNestedDataSet();
1717
1718 // <entry group="0028" element="9110" vr="SQ" vm="1" name="Pixel Measures Sequence"/>
1719 // do not set spacing between slices since GDCM always recompute it from the IOP/IPP
1720 Attribute<0x0018,0x0088> at2;
1721 at2.SetValue( fabs(spacing[2]) );
1722 Attribute<0x0028,0x0030> at1;
1723 at1.SetValue( spacing[1], 0 );
1724 at1.SetValue( spacing[0], 1 );
1725 subds2.Replace( at1.GetAsDataElement() );
1726 subds2.Replace( at2.GetAsDataElement() );
1727 }
1728 // cleanup per-frame
1729 {
1730 const Tag tfgs(0x5200,0x9230);
1731 if( ds.FindDataElement( tfgs ) )
1732 {
1733 SmartPointer<SequenceOfItems> sqi = InsertOrReplaceSQ( ds, tfgs );
1734 SequenceOfItems::SizeType nitems = sqi->GetNumberOfItems();
1735 for(SequenceOfItems::SizeType i0 = 1; i0 <= nitems; ++i0)
1736 {
1737 Item &item = sqi->GetItem(i0);
1738 item.SetVLToUndefined();
1739 DataSet & subds = item.GetNestedDataSet();
1740 const Tag tpms(0x0028,0x9110);
1741 subds.Remove(tpms);
1742 }
1743 }
1744 }
1745 // cleanup root (famous MR -> EMR case)
1746 {
1747 const Tag t1(0x0018,0x0088);
1748 ds.Remove(t1);
1749 const Tag t2(0x0028,0x0030);
1750 ds.Remove(t2);
1751 }
1752
1753 return;
1754 }
1755 else if( ms == MediaStorage::UltrasoundMultiFrameImageStorage || ms == MediaStorage::UltrasoundImageStorage )
1756 {
1757 SetUltraSoundSpacingValueFromSequence(ds, spacing);
1758 return;
1759 }
1760
1761
1762 Tag spacingtag = GetSpacingTagFromMediaStorage(ms);
1763 Tag zspacingtag = GetZSpacingTagFromMediaStorage(ms);
1764 //std::vector<Tag> spacingtags;
1765 //spacingtags.push_back( spacingtag );
1766 //spacingtags.push_back( zspacingtag );
1767 {
1768 const Tag ¤tspacing = spacingtag;
1769 if( currentspacing != Tag(0xffff,0xffff) )
1770 {
1771 DataElement de( currentspacing );
1772 const Global &g = GlobalInstance;
1773 const Dicts &dicts = g.GetDicts();
1774 const DictEntry &entry = dicts.GetDictEntry(de.GetTag());
1775 const VR & vr = entry.GetVR();
1776 const VM & vm = entry.GetVM(); (void)vm;
1777 assert( de.GetVR() == vr || de.GetVR() == VR::INVALID );
1778 switch(vr)
1779 {
1780 case VR::DS:
1781 {
1782 Element<VR::DS,VM::VM1_n> el;
1783 el.SetLength( entry.GetVM().GetLength() * vr.GetSizeof() );
1784 assert( entry.GetVM() == VM::VM2 );
1785 for( unsigned int i = 0; i < entry.GetVM().GetLength(); ++i)
1786 {
1787 el.SetValue( spacing[i], i );
1788 }
1789 el.SetValue( spacing[1], 0 );
1790 el.SetValue( spacing[0], 1 );
1791 //assert( el.GetValue(0) == spacing[0] && el.GetValue(1) == spacing[1] );
1792 std::stringstream os;
1793 el.Write( os );
1794 de.SetVR( VR::DS );
1795 if( os.str().size() % 2 ) os << " ";
1796 VL::Type osStrSize = (VL::Type)os.str().size();
1797 de.SetByteValue( os.str().c_str(),osStrSize );
1798 ds.Replace( de );
1799 }
1800 break;
1801 case VR::IS:
1802 {
1803 Element<VR::IS,VM::VM1_n> el;
1804 el.SetLength( entry.GetVM().GetLength() * vr.GetSizeof() );
1805 assert( entry.GetVM() == VM::VM2 );
1806 for( unsigned int i = 0; i < entry.GetVM().GetLength(); ++i)
1807 {
1808 el.SetValue( (int)spacing[i], i );
1809 }
1810 //assert( el.GetValue(0) == spacing[0] && el.GetValue(1) == spacing[1] );
1811 std::stringstream os;
1812 el.Write( os );
1813 de.SetVR( VR::IS );
1814 if( os.str().size() % 2 ) os << " ";
1815 VL::Type osStrSize = (VL::Type)os.str().size();
1816 de.SetByteValue( os.str().c_str(), osStrSize );
1817 ds.Replace( de );
1818 }
1819 break;
1820 default:
1821 assert(0);
1822 }
1823 }
1824 }
1825 {
1826 const Tag ¤tspacing = zspacingtag;
1827 if( currentspacing != Tag(0xffff,0xffff) )
1828 {
1829 DataElement de( currentspacing );
1830 const Global &g = GlobalInstance;
1831 const Dicts &dicts = g.GetDicts();
1832 const DictEntry &entry = dicts.GetDictEntry(de.GetTag());
1833 const VR & vr = entry.GetVR();
1834 const VM & vm = entry.GetVM(); (void)vm;
1835 assert( de.GetVR() == vr || de.GetVR() == VR::INVALID );
1836 if( entry.GetVM() == VM::VM2_n )
1837 {
1838 assert( vr == VR::DS );
1839 assert( de.GetTag() == Tag(0x3004,0x000c) );
1840 Attribute<0x28,0x8> numberoframes;
1841 // Make we are multiframes:
1842 if( ds.FindDataElement( numberoframes.GetTag() ) )
1843 {
1844 const DataElement& de1 = ds.GetDataElement( numberoframes.GetTag() );
1845 numberoframes.SetFromDataElement( de1 );
1846
1847 Element<VR::DS,VM::VM2_n> el;
1848 el.SetLength( numberoframes.GetValue() * vr.GetSizeof() );
1849 assert( entry.GetVM() == VM::VM2_n );
1850 double spacing_start = 0;
1851 assert( 0 < numberoframes.GetValue() );
1852 for( int i = 0; i < numberoframes.GetValue(); ++i)
1853 {
1854 el.SetValue( spacing_start, i );
1855 spacing_start += spacing[2];
1856 }
1857 //assert( el.GetValue(0) == spacing[0] && el.GetValue(1) == spacing[1] );
1858 std::stringstream os;
1859 el.Write( os );
1860 de.SetVR( VR::DS );
1861 if( os.str().size() % 2 ) os << " ";
1862 VL::Type osStrSize = (VL::Type)os.str().size();
1863 de.SetByteValue( os.str().c_str(), osStrSize );
1864 ds.Replace( de );
1865 }
1866 }
1867 else
1868 {
1869 switch(vr)
1870 {
1871 case VR::DS:
1872 {
1873 Element<VR::DS,VM::VM1_n> el;
1874 el.SetLength( entry.GetVM().GetLength() * vr.GetSizeof() );
1875 assert( entry.GetVM() == VM::VM1 );
1876 for( unsigned int i = 0; i < entry.GetVM().GetLength(); ++i)
1877 {
1878 el.SetValue( spacing[i+2], i );
1879 }
1880 //assert( el.GetValue(0) == spacing[0] && el.GetValue(1) == spacing[1] );
1881 std::stringstream os;
1882 el.Write( os );
1883 de.SetVR( VR::DS );
1884 if( os.str().size() % 2 ) os << " ";
1885 VL::Type osStrSize = (VL::Type)os.str().size();
1886 de.SetByteValue( os.str().c_str(), osStrSize );
1887 ds.Replace( de );
1888 }
1889 break;
1890 default:
1891 assert(0);
1892 }
1893 }
1894 }
1895 }
1896
1897 }
1898
SetDataElementInSQAsItemNumber(DataSet & ds,DataElement const & de,Tag const & sqtag,unsigned int itemidx)1899 static void SetDataElementInSQAsItemNumber(DataSet & ds, DataElement const & de, Tag const & sqtag, unsigned int itemidx)
1900 {
1901 const Tag tfgs = sqtag; //(0x5200,0x9230);
1902 SmartPointer<SequenceOfItems> sqi = InsertOrReplaceSQ( ds, tfgs );
1903 if( sqi->GetNumberOfItems() < itemidx )
1904 {
1905 Item item; //( Tag(0xfffe,0xe000) );
1906 sqi->AddItem( item );
1907 }
1908 Item &item1 = sqi->GetItem(itemidx);
1909 item1.SetVLToUndefined();
1910 DataSet &subds = item1.GetNestedDataSet();
1911 const Tag tpms(0x0020,0x9113);
1912 sqi = InsertOrReplaceSQ( subds, tpms );
1913 if( !sqi->GetNumberOfItems() )
1914 {
1915 Item item; //( Tag(0xfffe,0xe000) );
1916 sqi->AddItem( item );
1917 }
1918 Item &item2 = sqi->GetItem(1);
1919 item2.SetVLToUndefined();
1920 DataSet &subds2 = item2.GetNestedDataSet();
1921
1922 //Attribute<0x0020,0x0032> ipp = {{0,0,0}}; // default value
1923 //ipp.SetValue( origin[0], 0);
1924 //ipp.SetValue( origin[1], 1);
1925 //ipp.SetValue( origin[2], 2);
1926
1927 subds2.Replace( de );
1928 }
1929
SetOriginValue(DataSet & ds,const Image & image)1930 void ImageHelper::SetOriginValue(DataSet & ds, const Image & image)
1931 {
1932 const double *origin = image.GetOrigin();
1933 //assert( origin.size() == 3 );
1934 MediaStorage ms;
1935 ms.SetFromDataSet(ds);
1936 assert( MediaStorage::IsImage( ms ) );
1937
1938 if( ms == MediaStorage::SecondaryCaptureImageStorage )
1939 {
1940 // https://sourceforge.net/p/gdcm/bugs/322/
1941 // default behavior is simply to pass
1942 return;
1943 }
1944
1945 // FIXME Hardcoded
1946 if( ms != MediaStorage::CTImageStorage
1947 && ms != MediaStorage::MRImageStorage
1948 && ms != MediaStorage::RTDoseStorage
1949 && ms != MediaStorage::PETImageStorage
1950 //&& ms != MediaStorage::ComputedRadiographyImageStorage
1951 && ms != MediaStorage::SegmentationStorage
1952 && ms != MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage
1953 && ms != MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage
1954 && ms != MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage
1955 && ms != MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage
1956 && ms != MediaStorage::XRay3DAngiographicImageStorage
1957 && ms != MediaStorage::XRay3DCraniofacialImageStorage
1958 && ms != MediaStorage::EnhancedMRImageStorage
1959 && ms != MediaStorage::EnhancedMRColorImageStorage
1960 && ms != MediaStorage::EnhancedPETImageStorage
1961 && ms != MediaStorage::EnhancedCTImageStorage
1962 && ms != MediaStorage::IVOCTForPresentation
1963 && ms != MediaStorage::IVOCTForProcessing
1964 && ms != MediaStorage::BreastTomosynthesisImageStorage
1965 && ms != MediaStorage::BreastProjectionXRayImageStorageForPresentation
1966 && ms != MediaStorage::BreastProjectionXRayImageStorageForProcessing
1967 && ms != MediaStorage::LegacyConvertedEnhancedMRImageStorage
1968 && ms != MediaStorage::LegacyConvertedEnhancedCTImageStorage
1969 && ms != MediaStorage::LegacyConvertedEnhancedPETImageStorage )
1970 {
1971 // FIXME: should I remove the ipp tag ???
1972 return;
1973 }
1974
1975 if( ms == MediaStorage::EnhancedCTImageStorage
1976 || ms == MediaStorage::EnhancedMRImageStorage
1977 || ms == MediaStorage::EnhancedMRColorImageStorage
1978 || ms == MediaStorage::EnhancedPETImageStorage
1979 || ms == MediaStorage::XRay3DAngiographicImageStorage
1980 || ms == MediaStorage::XRay3DCraniofacialImageStorage
1981 || ms == MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage
1982 || ms == MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage
1983 || ms == MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage
1984 || ms == MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage
1985 || ms == MediaStorage::SegmentationStorage
1986 || ms == MediaStorage::IVOCTForPresentation
1987 || ms == MediaStorage::IVOCTForProcessing
1988 || ms == MediaStorage::BreastTomosynthesisImageStorage
1989 || ms == MediaStorage::BreastProjectionXRayImageStorageForPresentation
1990 || ms == MediaStorage::BreastProjectionXRayImageStorageForProcessing
1991 || ms == MediaStorage::LegacyConvertedEnhancedMRImageStorage
1992 || ms == MediaStorage::LegacyConvertedEnhancedCTImageStorage
1993 || ms == MediaStorage::LegacyConvertedEnhancedPETImageStorage)
1994 {
1995 /*
1996 (0020,9113) SQ (Sequence with undefined length #=1) # u/l, 1 PlanePositionSequence
1997 (fffe,e000) na (Item with undefined length #=1) # u/l, 1 Item
1998 (0020,0032) DS [40.0000\-105.000\105.000] # 24, 3 ImagePositionPatient
1999 (fffe,e00d) na (ItemDelimitationItem) # 0, 0 ItemDelimitationItem
2000 (fffe,e0dd) na (SequenceDelimitationItem) # 0, 0 SequenceDelimitationItem
2001 */
2002
2003 const Tag tfgs(0x5200,0x9230);
2004
2005 Attribute<0x0020,0x0032> ipp = {{0,0,0}}; // default value
2006 double zspacing = image.GetSpacing(2);
2007 unsigned int dimz = image.GetDimension(2);
2008 const double *cosines = image.GetDirectionCosines();
2009 DirectionCosines dc( cosines );
2010
2011 double normal[3];
2012 dc.Cross( normal );
2013
2014 for(unsigned int i = 0; i < dimz; ++i )
2015 {
2016 double new_origin[3];
2017 for (int j = 0; j < 3; j++)
2018 {
2019 // the n'th slice is n * z-spacing aloung the IOP-derived
2020 // z-axis
2021 new_origin[j] = origin[j] + normal[j] * i * zspacing;
2022 }
2023
2024 ipp.SetValue( new_origin[0], 0);
2025 ipp.SetValue( new_origin[1], 1);
2026 ipp.SetValue( new_origin[2], 2);
2027 SetDataElementInSQAsItemNumber(ds, ipp.GetAsDataElement(), tfgs, i+1);
2028 }
2029 // cleanup the sharedgroup:
2030 {
2031 const Tag tfgs0(0x5200,0x9229);
2032 if( ds.FindDataElement( tfgs0 ) )
2033 {
2034 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( tfgs0 ).GetValueAsSQ();
2035 assert( sqi );
2036 SequenceOfItems::SizeType nitems = sqi->GetNumberOfItems();
2037 for(SequenceOfItems::SizeType i0 = 1; i0 <= nitems; ++i0)
2038 {
2039 // Get first item:
2040 Item &item = sqi->GetItem(i0);
2041 item.SetVLToUndefined();
2042 DataSet & subds = item.GetNestedDataSet();
2043 const Tag tpms(0x0020,0x9113);
2044 subds.Remove(tpms);
2045 }
2046 }
2047 }
2048 // Cleanup root level:
2049 {
2050 const Tag tiop(0x0020,0x0032);
2051 ds.Remove(tiop);
2052 }
2053
2054
2055 // C.7.6.6.1.2 Frame Increment Pointer
2056 // (0028,0009) AT (0018,2005) # 4,1-n Frame Increment Pointer
2057 if( ms == MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage
2058 || ms == MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage
2059 || ms == MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage
2060 || ms == MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage )
2061 {
2062 if( dimz > 1 ) {
2063 Attribute<0x0028,0x0009> fip;
2064 fip.SetNumberOfValues( 1 );
2065 fip.SetValue( tfgs );
2066 ds.Replace( fip.GetAsDataElement() );
2067 }
2068 }
2069
2070 return;
2071 }
2072
2073 // Image Position (Patient)
2074 Attribute<0x0020,0x0032> ipp = {{0,0,0}}; // default value
2075 ipp.SetValue( origin[0], 0);
2076 ipp.SetValue( origin[1], 1);
2077 ipp.SetValue( origin[2], 2);
2078
2079 ds.Replace( ipp.GetAsDataElement() );
2080 }
2081
SetDirectionCosinesValue(DataSet & ds,const std::vector<double> & dircos)2082 void ImageHelper::SetDirectionCosinesValue(DataSet & ds, const std::vector<double> & dircos)
2083 {
2084 MediaStorage ms;
2085 ms.SetFromDataSet(ds);
2086 assert( MediaStorage::IsImage( ms ) );
2087
2088 if( ms == MediaStorage::SecondaryCaptureImageStorage )
2089 {
2090 // https://sourceforge.net/p/gdcm/bugs/322/
2091 // default behavior is simply to pass
2092 return;
2093 }
2094
2095 // FIXME Hardcoded
2096 if( ms != MediaStorage::CTImageStorage
2097 && ms != MediaStorage::MRImageStorage
2098 && ms != MediaStorage::RTDoseStorage
2099 && ms != MediaStorage::PETImageStorage
2100 //&& ms != MediaStorage::ComputedRadiographyImageStorage
2101 && ms != MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage
2102 && ms != MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage
2103 && ms != MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage
2104 && ms != MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage
2105 && ms != MediaStorage::SegmentationStorage
2106 && ms != MediaStorage::XRay3DAngiographicImageStorage
2107 && ms != MediaStorage::XRay3DCraniofacialImageStorage
2108 && ms != MediaStorage::EnhancedMRImageStorage
2109 && ms != MediaStorage::EnhancedMRColorImageStorage
2110 && ms != MediaStorage::EnhancedPETImageStorage
2111 && ms != MediaStorage::EnhancedCTImageStorage
2112 && ms != MediaStorage::IVOCTForPresentation
2113 && ms != MediaStorage::IVOCTForProcessing
2114 && ms != MediaStorage::BreastTomosynthesisImageStorage
2115 && ms != MediaStorage::BreastProjectionXRayImageStorageForPresentation
2116 && ms != MediaStorage::BreastProjectionXRayImageStorageForProcessing
2117 && ms != MediaStorage::LegacyConvertedEnhancedMRImageStorage
2118 && ms != MediaStorage::LegacyConvertedEnhancedCTImageStorage
2119 && ms != MediaStorage::LegacyConvertedEnhancedPETImageStorage )
2120 {
2121 // FIXME: should I remove the iop tag ???
2122 return;
2123 }
2124
2125 // Image Orientation (Patient)
2126 Attribute<0x0020,0x0037> iop = {{1,0,0,0,1,0}}; // default value
2127
2128 assert( dircos.size() == 6 );
2129 DirectionCosines dc( &dircos[0] );
2130 if( !dc.IsValid() )
2131 {
2132 gdcmWarningMacro( "Direction Cosines are not valid. Using default value (1\\0\\0\\0\\1\\0)" );
2133 }
2134 else
2135 {
2136 iop.SetValue( dircos[0], 0);
2137 iop.SetValue( dircos[1], 1);
2138 iop.SetValue( dircos[2], 2);
2139 iop.SetValue( dircos[3], 3);
2140 iop.SetValue( dircos[4], 4);
2141 iop.SetValue( dircos[5], 5);
2142 }
2143
2144 if( ms == MediaStorage::EnhancedCTImageStorage
2145 || ms == MediaStorage::EnhancedMRImageStorage
2146 || ms == MediaStorage::EnhancedMRColorImageStorage
2147 || ms == MediaStorage::EnhancedPETImageStorage
2148 || ms == MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage
2149 || ms == MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage
2150 || ms == MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage
2151 || ms == MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage
2152 || ms == MediaStorage::XRay3DAngiographicImageStorage
2153 || ms == MediaStorage::XRay3DCraniofacialImageStorage
2154 || ms == MediaStorage::SegmentationStorage
2155 || ms == MediaStorage::IVOCTForPresentation
2156 || ms == MediaStorage::IVOCTForProcessing
2157 || ms == MediaStorage::BreastTomosynthesisImageStorage
2158 || ms == MediaStorage::BreastProjectionXRayImageStorageForPresentation
2159 || ms == MediaStorage::BreastProjectionXRayImageStorageForProcessing
2160 || ms == MediaStorage::LegacyConvertedEnhancedMRImageStorage
2161 || ms == MediaStorage::LegacyConvertedEnhancedCTImageStorage
2162 || ms == MediaStorage::LegacyConvertedEnhancedPETImageStorage )
2163 {
2164 /*
2165 (0020,9116) SQ (Sequence with undefined length #=1) # u/l, 1 PlaneOrientationSequence
2166 (fffe,e000) na (Item with undefined length #=1) # u/l, 1 Item
2167 (0020,0037) DS [0.00000\1.00000\0.00000\0.00000\0.00000\-1.00000] # 48, 6 ImageOrientationPatient
2168 (fffe,e00d) na (ItemDelimitationItem) # 0, 0 ItemDelimitationItem
2169 (fffe,e0dd) na (SequenceDelimitationItem) # 0, 0 SequenceDelimitationItem
2170 */
2171 {
2172 const Tag tfgs(0x5200,0x9229);
2173 SmartPointer<SequenceOfItems> sqi = InsertOrReplaceSQ( ds, tfgs );
2174 if( !sqi->GetNumberOfItems() )
2175 {
2176 Item item; //( Tag(0xfffe,0xe000) );
2177 sqi->AddItem( item );
2178 }
2179 Item &item1 = sqi->GetItem(1);
2180 item1.SetVLToUndefined();
2181 DataSet &subds = item1.GetNestedDataSet();
2182 const Tag tpms(0x0020,0x9116);
2183 sqi = InsertOrReplaceSQ( subds, tpms );
2184 if( !sqi->GetNumberOfItems() )
2185 {
2186 Item item; //( Tag(0xfffe,0xe000) );
2187 sqi->AddItem( item );
2188 }
2189 Item &item2 = sqi->GetItem(1);
2190 item2.SetVLToUndefined();
2191 DataSet &subds2 = item2.GetNestedDataSet();
2192
2193 subds2.Replace( iop.GetAsDataElement() );
2194 }
2195 // cleanup per-frame
2196 {
2197 const Tag tfgs(0x5200,0x9230);
2198 if( ds.FindDataElement( tfgs ) )
2199 {
2200 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( tfgs ).GetValueAsSQ();
2201 assert( sqi );
2202 SequenceOfItems::SizeType nitems = sqi->GetNumberOfItems();
2203 for(SequenceOfItems::SizeType i0 = 1; i0 <= nitems; ++i0)
2204 {
2205 // Get first item:
2206 Item &item = sqi->GetItem(i0);
2207 item.SetVLToUndefined();
2208 DataSet & subds = item.GetNestedDataSet();
2209 const Tag tpms(0x0020,0x9116);
2210 subds.Remove(tpms);
2211 }
2212 }
2213 }
2214 // Cleanup root level:
2215 {
2216 const Tag tiop(0x0020,0x0037);
2217 ds.Remove(tiop);
2218 }
2219
2220 return;
2221 }
2222
2223 ds.Replace( iop.GetAsDataElement() );
2224 }
2225
SetRescaleInterceptSlopeValue(File & f,const Image & img)2226 void ImageHelper::SetRescaleInterceptSlopeValue(File & f, const Image & img)
2227 {
2228 MediaStorage ms;
2229 // SetFromFile is required here, SetFromDataSet is not enough for all cases
2230 ms.SetFromFile(f);
2231 assert( MediaStorage::IsImage( ms ) );
2232 DataSet &ds = f.GetDataSet();
2233
2234 // FIXME Hardcoded
2235 if( ms != MediaStorage::CTImageStorage
2236 && ms != MediaStorage::ComputedRadiographyImageStorage
2237 && ms != MediaStorage::MRImageStorage // FIXME !
2238 && ms != MediaStorage::PETImageStorage
2239 && ms != MediaStorage::RTDoseStorage
2240 && ms != MediaStorage::SecondaryCaptureImageStorage
2241 && ms != MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage
2242 && ms != MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage
2243 && ms != MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage
2244 && ms != MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage
2245 && ms != MediaStorage::EnhancedMRImageStorage
2246 /*&& ms != MediaStorage::EnhancedMRColorImageStorage*/
2247 && ms != MediaStorage::EnhancedCTImageStorage
2248 && ms != MediaStorage::EnhancedPETImageStorage
2249 && ms != MediaStorage::XRay3DAngiographicImageStorage
2250 && ms != MediaStorage::XRay3DCraniofacialImageStorage
2251 && ms != MediaStorage::SegmentationStorage
2252 && ms != MediaStorage::IVOCTForPresentation
2253 && ms != MediaStorage::IVOCTForProcessing
2254 && ms != MediaStorage::BreastTomosynthesisImageStorage
2255 && ms != MediaStorage::BreastProjectionXRayImageStorageForPresentation
2256 && ms != MediaStorage::BreastProjectionXRayImageStorageForProcessing
2257 && ms != MediaStorage::LegacyConvertedEnhancedMRImageStorage
2258 && ms != MediaStorage::LegacyConvertedEnhancedCTImageStorage
2259 && ms != MediaStorage::LegacyConvertedEnhancedPETImageStorage )
2260 {
2261 if( img.GetIntercept() != 0. || img.GetSlope() != 1. )
2262 {
2263 gdcmErrorMacro( "Cannot have non-identity intercept/slope values." );
2264 throw "Impossible"; // Please report
2265 }
2266 return;
2267 }
2268
2269 if( ms == MediaStorage::SegmentationStorage ) return; // seg storage cannot have rescale slope
2270 if( ms == MediaStorage::EnhancedCTImageStorage
2271 || ms == MediaStorage::EnhancedMRImageStorage
2272 /*|| ms == MediaStorage::EnhancedMRColorImageStorage*/
2273 || ms == MediaStorage::EnhancedPETImageStorage
2274 || ms == MediaStorage::XRay3DAngiographicImageStorage
2275 || ms == MediaStorage::XRay3DCraniofacialImageStorage
2276 || ms == MediaStorage::BreastTomosynthesisImageStorage
2277 || ms == MediaStorage::BreastProjectionXRayImageStorageForPresentation
2278 || ms == MediaStorage::BreastProjectionXRayImageStorageForProcessing
2279 || ms == MediaStorage::LegacyConvertedEnhancedMRImageStorage
2280 || ms == MediaStorage::LegacyConvertedEnhancedCTImageStorage
2281 || ms == MediaStorage::LegacyConvertedEnhancedPETImageStorage )
2282 {
2283 /*
2284 (0020,9116) SQ (Sequence with undefined length #=1) # u/l, 1 PlaneOrientationSequence
2285 (fffe,e000) na (Item with undefined length #=1) # u/l, 1 Item
2286 (0020,0037) DS [0.00000\1.00000\0.00000\0.00000\0.00000\-1.00000] # 48, 6 ImageOrientationPatient
2287 (fffe,e00d) na (ItemDelimitationItem) # 0, 0 ItemDelimitationItem
2288 (fffe,e0dd) na (SequenceDelimitationItem) # 0, 0 SequenceDelimitationItem
2289 */
2290 {
2291 const Tag tfgs(0x5200,0x9229);
2292 SmartPointer<SequenceOfItems> sqi = InsertOrReplaceSQ( ds, tfgs );
2293 if( !sqi->GetNumberOfItems() )
2294 {
2295 Item item; //( Tag(0xfffe,0xe000) );
2296 sqi->AddItem( item );
2297 }
2298 Item &item1 = sqi->GetItem(1);
2299 item1.SetVLToUndefined();
2300 DataSet &subds = item1.GetNestedDataSet();
2301 const Tag tpms(0x0028,0x9145);
2302 sqi = InsertOrReplaceSQ( subds, tpms );
2303 if( !sqi->GetNumberOfItems() )
2304 {
2305 Item item; //( Tag(0xfffe,0xe000) );
2306 sqi->AddItem( item );
2307 }
2308 Item &item2 = sqi->GetItem(1);
2309 item2.SetVLToUndefined();
2310 DataSet &subds2 = item2.GetNestedDataSet();
2311
2312 Attribute<0x0028,0x1052> at1;
2313 at1.SetValue( img.GetIntercept() );
2314 subds2.Replace( at1.GetAsDataElement() );
2315 Attribute<0x0028,0x1053> at2;
2316 at2.SetValue( img.GetSlope() );
2317 subds2.Replace( at2.GetAsDataElement() );
2318 }
2319
2320 // cleanup per-frame
2321 {
2322 const Tag tfgs(0x5200,0x9230);
2323 if( ds.FindDataElement( tfgs ) )
2324 {
2325 SmartPointer<SequenceOfItems> sqi = ds.GetDataElement( tfgs ).GetValueAsSQ();
2326 assert( sqi );
2327 SequenceOfItems::SizeType nitems = sqi->GetNumberOfItems();
2328 for(SequenceOfItems::SizeType i0 = 1; i0 <= nitems; ++i0)
2329 {
2330 // Get first item:
2331 Item &item = sqi->GetItem(i0);
2332 DataSet & subds = item.GetNestedDataSet();
2333 // (0028,9145) SQ (Sequence with undefined length) # u/l,1 Pixel Value Transformation Sequence
2334 const Tag tpms(0x0028,0x9145);
2335 subds.Remove(tpms);
2336 }
2337 }
2338 }
2339 // cleanup root (famous MR -> EMR case)
2340 {
2341 const Tag t1(0x0028,0x1052);
2342 ds.Remove(t1);
2343 const Tag t2(0x0028,0x1053);
2344 ds.Remove(t2);
2345 const Tag t3(0x0028,0x1053);
2346 ds.Remove(t3);
2347 }
2348 return;
2349 }
2350
2351 if( ms == MediaStorage::RTDoseStorage )
2352 {
2353 if( img.GetIntercept() != 0 )
2354 {
2355 gdcmErrorMacro( "Cannot have an intercept value for RTDOSE, only Scaling allowed" );
2356 return;
2357 }
2358 Attribute<0x3004,0x00e> at2;
2359 at2.SetValue( img.GetSlope() );
2360 ds.Replace( at2.GetAsDataElement() );
2361
2362 Attribute<0x0028,0x0009> framePointer;
2363 framePointer.SetNumberOfValues(1);
2364 framePointer.SetValue( Tag(0x3004,0x000C) );
2365 ds.Replace( framePointer.GetAsDataElement() );
2366
2367 return;
2368 }
2369
2370 if( ms == MediaStorage::MRImageStorage )
2371 {
2372 #if 0
2373 /*
2374 * http://dicom.nema.org/medical/dicom/current/output/chtml/part03/sect_C.7.6.16.2.html#table_C.7.6.16-12b
2375 (0040,9096) SQ (Sequence with undefined length) # u/l,1 Real World Value Mapping Sequence
2376 (fffe,e000) na (Item with defined length)
2377 (0028,3003) LO [Grey Scale LUT] # 14,1 LUT Explanation
2378 (0040,08ea) SQ (Sequence with undefined length) # u/l,1 Measurement Units Code Sequence
2379 (fffe,e000) na (Item with defined length)
2380 (0008,0100) SH [mm2/s ] # 6,1 Code Value
2381 (0008,0102) SH [UCUM] # 4,1 Coding Scheme Designator
2382 (0008,0103) SH [1.4 ] # 4,1 Coding Scheme Version
2383 (0008,0104) LO [mm2/s ] # 6,1 Code Meaning
2384 (fffe,e0dd)
2385 (0040,9210) SH [GE_GREY ] # 8,1 LUT Label
2386 (0040,9211) US 4904 # 2,1 Real World Value Last Value Mapped
2387 (0040,9216) US 359 # 2,1 Real World Value First Value Mapped
2388 (0040,9224) FD 0 # 8,1 Real World Value Intercept
2389 (0040,9225) FD 1e-06 # 8,1 Real World Value Slope
2390 */
2391 if( img.GetIntercept() != 0.0 || img.GetSlope() != 1.0 )
2392 {
2393 SmartPointer<SequenceOfItems> sq = new SequenceOfItems;
2394 Item it;
2395 DataSet & subds = it.GetNestedDataSet();
2396 Attribute<0x0040,0x9224> at1 = {0};
2397 at1.SetValue( img.GetIntercept() );
2398 Attribute<0x0040,0x9225> at2 = {1};
2399 at2.SetValue( img.GetSlope() );
2400 subds.Insert( at1.GetAsDataElement() );
2401 subds.Insert( at2.GetAsDataElement() );
2402 sq->AddItem( it );
2403 const Tag trwvms(0x0040,0x9096); // Real World Value Mapping Sequence
2404 DataElement de( trwvms );
2405 de.SetVR( VR::SQ );
2406 de.SetValue(*sq);
2407 ds.Replace( de );
2408 }
2409
2410 ds.Remove( Tag(0x28,0x1052) );
2411 ds.Remove( Tag(0x28,0x1053) );
2412 ds.Remove( Tag(0x28,0x1054) );
2413 #else
2414 //if( img.GetIntercept() != 0.0 || img.GetSlope() != 1.0 )
2415 {
2416 if( ForceRescaleInterceptSlope )
2417 {
2418 gdcmDebugMacro( "Forcing MR Image Storage / Modality LUT: [" << img.GetIntercept() << "," << img.GetSlope() );
2419 Attribute<0x0028,0x1052> at1;
2420 at1.SetValue( img.GetIntercept() );
2421 ds.Replace( at1.GetAsDataElement() );
2422 Attribute<0x0028,0x1053> at2;
2423 at2.SetValue( img.GetSlope() );
2424 ds.Replace( at2.GetAsDataElement() );
2425
2426 Attribute<0x0028,0x1054> at3; // Rescale Type
2427 at3.SetValue( "US" ); // Compatible with Enhanced MR Image Storage
2428 ds.Replace( at3.GetAsDataElement() );
2429 }
2430 }
2431 #endif
2432 }
2433 else
2434 {
2435 Attribute<0x0028,0x1052> at1;
2436 at1.SetValue( img.GetIntercept() );
2437 ds.Replace( at1.GetAsDataElement() );
2438 Attribute<0x0028,0x1053> at2;
2439 at2.SetValue( img.GetSlope() );
2440 ds.Replace( at2.GetAsDataElement() );
2441
2442 Attribute<0x0028,0x1054> at3; // Rescale Type
2443 at3.SetValue( "US" ); // FIXME
2444 if( ms == MediaStorage::SecondaryCaptureImageStorage )
2445 {
2446 // As per 3-2009, US is the only valid enumerated value:
2447 ds.Replace( at3.GetAsDataElement() );
2448 }
2449 else if( ms == MediaStorage::PETImageStorage )
2450 {
2451 // not there anyway:
2452 ds.Remove( at3.GetTag() );
2453 }
2454 else
2455 {
2456 // In case user decide to override the default:
2457 ds.ReplaceEmpty( at3.GetAsDataElement() );
2458 }
2459 }
2460 }
2461
GetRealWorldValueMappingContent(File const & f,RealWorldValueMappingContent & ret)2462 bool ImageHelper::GetRealWorldValueMappingContent(File const & f, RealWorldValueMappingContent & ret)
2463 {
2464 MediaStorage ms;
2465 ms.SetFromFile(f);
2466 const DataSet& ds = f.GetDataSet();
2467
2468 if( ms == MediaStorage::MRImageStorage )
2469 {
2470 const Tag trwvms(0x0040,0x9096); // Real World Value Mapping Sequence
2471 if( ds.FindDataElement( trwvms ) )
2472 {
2473 SmartPointer<SequenceOfItems> sqi0 = ds.GetDataElement( trwvms ).GetValueAsSQ();
2474 if( sqi0 )
2475 {
2476 const Tag trwvlutd(0x0040,0x9212); // Real World Value LUT Data
2477 if( ds.FindDataElement( trwvlutd ) )
2478 {
2479 gdcmAssertAlwaysMacro(0); // Not supported !
2480 }
2481 // don't know how to handle multiples:
2482 gdcmAssertAlwaysMacro( sqi0->GetNumberOfItems() == 1 );
2483 const Item &item0 = sqi0->GetItem(1);
2484 const DataSet & subds0 = item0.GetNestedDataSet();
2485 //const Tag trwvi(0x0040,0x9224); // Real World Value Intercept
2486 //const Tag trwvs(0x0040,0x9225); // Real World Value Slope
2487 {
2488 Attribute<0x0040,0x9224> at1 = {0};
2489 at1.SetFromDataSet( subds0 );
2490 Attribute<0x0040,0x9225> at2 = {1};
2491 at2.SetFromDataSet( subds0 );
2492 ret.RealWorldValueIntercept = at1.GetValue();
2493 ret.RealWorldValueSlope = at2.GetValue();
2494 }
2495 const Tag tmucs(0x0040,0x08ea); // Measurement Units Code Sequence
2496 if( subds0.FindDataElement( tmucs ) )
2497 {
2498 SmartPointer<SequenceOfItems> sqi = subds0.GetDataElement( tmucs ).GetValueAsSQ();
2499 if( sqi )
2500 {
2501 gdcmAssertAlwaysMacro( sqi->GetNumberOfItems() == 1 );
2502 const Item &item = sqi->GetItem(1);
2503 const DataSet & subds = item.GetNestedDataSet();
2504 Attribute<0x0008,0x0100> at1;
2505 at1.SetFromDataSet( subds );
2506 Attribute<0x0008,0x0104> at2;
2507 at2.SetFromDataSet( subds );
2508 ret.CodeValue = at1.GetValue().Trim();
2509 ret.CodeMeaning = at2.GetValue().Trim();
2510 }
2511 }
2512 }
2513 return true;
2514 }
2515 }
2516 return false;
2517 }
2518
ComputeSpacingFromImagePositionPatient(const std::vector<double> & imageposition,std::vector<double> & spacing)2519 bool ImageHelper::ComputeSpacingFromImagePositionPatient(const std::vector<double> & imageposition, std::vector<double> & spacing)
2520 {
2521 if( imageposition.size() % 3 != 0 )
2522 {
2523 return false;
2524 }
2525 std::vector<double>::const_iterator it = imageposition.begin();
2526 //const double x0 = *it++;
2527 //const double y0 = *it++;
2528 //const double z0 = *it++;
2529 spacing[0] = spacing[1] = spacing[2] = 0.;
2530 for( ; it != imageposition.end(); ++it)
2531 {
2532 const double x = *it++;
2533 const double y = *it++;
2534 const double z = *it;
2535 spacing[0] += x;
2536 spacing[1] += y;
2537 spacing[2] += z;
2538 }
2539 size_t n = imageposition.size() / 3;
2540 spacing[0] /= (double)n;
2541 spacing[1] /= (double)n;
2542 spacing[2] /= (double)n;
2543
2544 return true;
2545 }
2546
2547
2548 //functions to get more information from a file
2549 //useful for the stream image reader, which fills in necessary image information
2550 //distinctly from the reader-style data input
2551 //code is borrowed from gdcmPixmapReader::ReadImage(MediaStorage const &ms)
GetPhotometricInterpretationValue(File const & f)2552 PhotometricInterpretation ImageHelper::GetPhotometricInterpretationValue(File const& f)
2553 {
2554 // 5. Photometric Interpretation
2555 // D 0028|0004 [CS] [Photometric Interpretation] [MONOCHROME2 ]
2556 PixelFormat pf = GetPixelFormatValue(f);
2557 const Tag tphotometricinterpretation(0x0028, 0x0004);
2558 const ByteValue *photometricinterpretation =
2559 ImageHelper::GetPointerFromElement(tphotometricinterpretation, f);
2560 PhotometricInterpretation pi = PhotometricInterpretation::UNKNOWN;
2561 if( photometricinterpretation )
2562 {
2563 std::string photometricinterpretation_str(
2564 photometricinterpretation->GetPointer(),
2565 photometricinterpretation->GetLength() );
2566 pi = PhotometricInterpretation::GetPIType( photometricinterpretation_str.c_str() );
2567 }
2568 else
2569 {
2570 if( pf.GetSamplesPerPixel() == 1 )
2571 {
2572 gdcmWarningMacro( "No PhotometricInterpretation found, default to MONOCHROME2" );
2573 pi = PhotometricInterpretation::MONOCHROME2;
2574 }
2575 else if( pf.GetSamplesPerPixel() == 3 )
2576 {
2577 gdcmWarningMacro( "No PhotometricInterpretation found, default to RGB" );
2578 pi = PhotometricInterpretation::RGB;
2579 }
2580 else if( pf.GetSamplesPerPixel() == 4 )
2581 {
2582 gdcmWarningMacro( "No PhotometricInterpretation found, default to ARGB" );
2583 pi = PhotometricInterpretation::ARGB;
2584 }
2585 }
2586
2587 bool isacrnema = false;
2588 DataSet ds = f.GetDataSet();
2589 const Tag trecognitioncode(0x0008,0x0010);
2590 if( ds.FindDataElement( trecognitioncode ) && !ds.GetDataElement( trecognitioncode ).IsEmpty() )
2591 {
2592 // PHILIPS_Gyroscan-12-MONO2-Jpeg_Lossless.dcm
2593 // PHILIPS_Gyroscan-12-Jpeg_Extended_Process_2_4.dcm
2594 gdcmDebugMacro( "Mixture of ACR NEMA and DICOM file" );
2595 isacrnema = true;
2596 }
2597 if( !pf.GetSamplesPerPixel() || (pi.GetSamplesPerPixel() != pf.GetSamplesPerPixel()) )
2598 {
2599 if( pi != PhotometricInterpretation::UNKNOWN )
2600 {
2601 pf.SetSamplesPerPixel( pi.GetSamplesPerPixel() );
2602 }
2603 else if ( isacrnema )
2604 {
2605 assert ( pf.GetSamplesPerPixel() == 0 );
2606 assert ( pi == PhotometricInterpretation::UNKNOWN );
2607 pf.SetSamplesPerPixel( 1 );
2608 pi = PhotometricInterpretation::MONOCHROME2;
2609 }
2610 else
2611 {
2612 gdcmWarningMacro( "Cannot recognize image type. Does not looks like ACR-NEMA and is missing both Sample Per Pixel AND PhotometricInterpretation. Please report" );
2613 }
2614 }
2615 return pi;
2616 }
2617 //returns the configuration of colors in a plane, either RGB RGB RGB or RRR GGG BBB
2618 //code is borrowed from gdcmPixmapReader::ReadImage(MediaStorage const &ms)
GetPlanarConfigurationValue(const File & f)2619 unsigned int ImageHelper::GetPlanarConfigurationValue(const File& f)
2620 {
2621 // 4. Planar Configuration
2622 // D 0028|0006 [US] [Planar Configuration] [1]
2623 const Tag planarconfiguration = Tag(0x0028, 0x0006);
2624 PixelFormat pf = GetPixelFormatValue(f);
2625 unsigned int pc = 0;
2626 // FIXME: Whatif planaconfiguration is send in a grayscale image... it would be empty...
2627 // well hopefully :(
2628 DataSet const & ds = f.GetDataSet();
2629 if( ds.FindDataElement( planarconfiguration ) && !ds.GetDataElement( planarconfiguration ).IsEmpty() )
2630 {
2631 const DataElement& de = ds.GetDataElement( planarconfiguration );
2632 Attribute<0x0028,0x0006> at = { 0 };
2633 at.SetFromDataElement( de );
2634
2635 //unsigned int pc = ReadUSFromTag( planarconfiguration, ss, conversion );
2636 pc = at.GetValue();
2637 if( pc && pf.GetSamplesPerPixel() != 3 )
2638 {
2639 gdcmDebugMacro( "Cannot have PlanarConfiguration=1, when Sample Per Pixel != 3" );
2640 pc = 0;
2641 }
2642 }
2643 return pc;
2644 }
2645
2646 //returns the lookup table of an image file
GetLUT(File const & f)2647 SmartPointer<LookupTable> ImageHelper::GetLUT(File const& f)
2648 {
2649 DataSet const & ds = f.GetDataSet();
2650 PixelFormat pf = GetPixelFormatValue(f);
2651 PhotometricInterpretation pi = GetPhotometricInterpretationValue(f);
2652 // Do the Palette Color:
2653 // 1. Modality LUT Sequence
2654 bool modlut = ds.FindDataElement(Tag(0x0028,0x3000) );
2655 if( modlut )
2656 {
2657 gdcmWarningMacro( "Modality LUT (0028,3000) are not handled. Image will not be displayed properly" );
2658 }
2659 // 2. LUTData (0028,3006)
2660 // technically I do not need to warn about LUTData since either modality lut XOR VOI LUT need to
2661 // be sent to require a LUT Data...
2662 bool lutdata = ds.FindDataElement(Tag(0x0028,0x3006) );
2663 if( lutdata )
2664 {
2665 gdcmWarningMacro( "LUT Data (0028,3006) are not handled. Image will not be displayed properly" );
2666 }
2667 // 3. VOILUTSequence (0028,3010)
2668 bool voilut = ds.FindDataElement(Tag(0x0028,0x3010) );
2669 if( voilut )
2670 {
2671 gdcmWarningMacro( "VOI LUT (0028,3010) are not handled. Image will not be displayed properly" );
2672 }
2673 // (0028,0120) US 32767 # 2, 1 PixelPaddingValue
2674 bool pixelpaddingvalue = ds.FindDataElement(Tag(0x0028,0x0120));
2675
2676 // PS 3.3 - 2008 / C.7.5.1.1.2 Pixel Padding Value and Pixel Padding Range Limit
2677 if(pixelpaddingvalue)
2678 {
2679 // Technically if Pixel Padding Value is 0 on MONOCHROME2 image, then appearance should be fine...
2680 bool vizissue = true;
2681 if( pf.GetPixelRepresentation() == 0 )
2682 {
2683 Element<VR::US,VM::VM1> ppv;
2684 if( !ds.GetDataElement(Tag(0x0028,0x0120) ).IsEmpty() )
2685 {
2686 ppv.SetFromDataElement( ds.GetDataElement(Tag(0x0028,0x0120)) ); //.GetValue() );
2687 if( pi == PhotometricInterpretation::MONOCHROME2 && ppv.GetValue() == 0 )
2688 {
2689 vizissue = false;
2690 }
2691 }
2692 }
2693 else if( pf.GetPixelRepresentation() == 1 )
2694 {
2695 gdcmDebugMacro( "TODO" );
2696 }
2697 // test if there is any viz issue:
2698 if( vizissue )
2699 {
2700 gdcmDebugMacro( "Pixel Padding Value (0028,0120) is not handled. Image will not be displayed properly" );
2701 }
2702 }
2703 SmartPointer<LookupTable> lut = new LookupTable;
2704 const Tag testseglut(0x0028, (0x1221 + 0));
2705 if( ds.FindDataElement( testseglut ) )
2706 {
2707 lut = new SegmentedPaletteColorLookupTable;
2708 }
2709 //SmartPointer<SegmentedPaletteColorLookupTable> lut = new SegmentedPaletteColorLookupTable;
2710 lut->Allocate( pf.GetBitsAllocated() );
2711
2712 // for each red, green, blue:
2713 for(int i=0; i<3; ++i)
2714 {
2715 // (0028,1101) US 0\0\16
2716 // (0028,1102) US 0\0\16
2717 // (0028,1103) US 0\0\16
2718 const Tag tdescriptor(0x0028, (uint16_t)(0x1101 + i));
2719 //const Tag tdescriptor(0x0028, 0x3002);
2720 Element<VR::US,VM::VM3> el_us3 = {{ 0, 0, 0}};
2721 // Now pass the byte array to a DICOMizer:
2722 el_us3.SetFromDataElement( ds[tdescriptor] ); //.GetValue() );
2723 lut->InitializeLUT( LookupTable::LookupTableType(i),
2724 el_us3[0], el_us3[1], el_us3[2] );
2725
2726 // (0028,1201) OW
2727 // (0028,1202) OW
2728 // (0028,1203) OW
2729 const Tag tlut(0x0028, (uint16_t)(0x1201 + i));
2730 //const Tag tlut(0x0028, 0x3006);
2731
2732 // Segmented LUT
2733 // (0028,1221) OW
2734 // (0028,1222) OW
2735 // (0028,1223) OW
2736 const Tag seglut(0x0028, (uint16_t)(0x1221 + i));
2737 if( ds.FindDataElement( tlut ) )
2738 {
2739 const ByteValue *lut_raw = ds.GetDataElement( tlut ).GetByteValue();
2740 if( lut_raw )
2741 {
2742 // LookupTableType::RED == 0
2743 lut->SetLUT( LookupTable::LookupTableType(i),
2744 (const unsigned char*)lut_raw->GetPointer(), lut_raw->GetLength() );
2745 //assert( pf.GetBitsAllocated() == el_us3.GetValue(2) );
2746 }
2747 else
2748 {
2749 lut->Clear();
2750 }
2751
2752 unsigned long check =
2753 (el_us3.GetValue(0) ? el_us3.GetValue(0) : 65536)
2754 * el_us3.GetValue(2) / 8;
2755 assert( !lut->Initialized() || check == lut_raw->GetLength() ); (void)check;
2756 }
2757 else if( ds.FindDataElement( seglut ) )
2758 {
2759 const ByteValue *lut_raw = ds.GetDataElement( seglut ).GetByteValue();
2760 if( lut_raw )
2761 {
2762 lut->SetLUT( LookupTable::LookupTableType(i),
2763 (const unsigned char*)lut_raw->GetPointer(), lut_raw->GetLength() );
2764 //assert( pf.GetBitsAllocated() == el_us3.GetValue(2) );
2765 }
2766 else
2767 {
2768 lut->Clear();
2769 }
2770
2771 //unsigned long check =
2772 // (el_us3.GetValue(0) ? el_us3.GetValue(0) : 65536)
2773 // * el_us3.GetValue(2) / 8;
2774 //assert( check == lut_raw->GetLength() ); (void)check;
2775 }
2776 else
2777 {
2778 assert(0);
2779 }
2780 }
2781 if( ! lut->Initialized() ) {
2782 gdcmDebugMacro("LUT was uninitialized!");
2783 }
2784 return lut;
2785 }
2786
2787
GetPointerFromElement(Tag const & tag,const File & inF)2788 const ByteValue* ImageHelper::GetPointerFromElement(Tag const &tag, const File& inF)
2789 {
2790 const DataSet &ds = inF.GetDataSet();
2791 if( ds.FindDataElement( tag ) )
2792 {
2793 const DataElement &de = ds.GetDataElement( tag );
2794 return de.GetByteValue();
2795 }
2796 return nullptr;
2797 }
2798
ComputeMediaStorageFromModality(const char * modality,unsigned int dimension,PixelFormat const & pixeltype,PhotometricInterpretation const & pi,double intercept,double slope)2799 MediaStorage ImageHelper::ComputeMediaStorageFromModality(const char *modality,
2800 unsigned int dimension, PixelFormat const & pixeltype,
2801 PhotometricInterpretation const & pi,
2802 double intercept , double slope
2803 )
2804 {
2805 // FIXME: Planar Configuration (0028,0006) shall not be present
2806 MediaStorage ms = MediaStorage::SecondaryCaptureImageStorage;
2807 ms.GuessFromModality(modality, dimension );
2808
2809 // refine for SC family
2810 if( dimension != 2 &&
2811 (ms == MediaStorage::SecondaryCaptureImageStorage // dim 2
2812 || ms == MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage ) // dim 3
2813 )
2814 {
2815 // A.8.3.4 Multi-frame Grayscale Byte SC Image IOD Content Constraints
2816 /*
2817 - Samples per Pixel (0028,0002) shall be 1
2818 - Photometric Interpretation (0028,0004) shall be MONOCHROME2
2819 - Bits Allocated (0028,0100) shall be 8
2820 - Bits Stored (0028,0101) shall be 8
2821 - High Bit (0028,0102) shall be 7
2822 - Pixel Representation (0028,0103) shall be 0
2823 - Planar Configuration (0028,0006) shall not be present
2824 */
2825 if( dimension == 3 &&
2826 pixeltype.GetSamplesPerPixel() == 1 &&
2827 pi == PhotometricInterpretation::MONOCHROME2 &&
2828 pixeltype.GetBitsAllocated() == 8 &&
2829 pixeltype.GetBitsStored() == 8 &&
2830 pixeltype.GetHighBit() == 7 &&
2831 pixeltype.GetPixelRepresentation() == 0
2832 )
2833 {
2834 ms = MediaStorage::MultiframeGrayscaleByteSecondaryCaptureImageStorage;
2835 if( intercept != 0 || slope != 1 )
2836 {
2837 // A.8.3.4 Multi-frame Grayscale Byte SC Image IOD Content Constraints
2838 gdcmErrorMacro( "Cannot have shift/scale" );
2839 return MediaStorage::MS_END;
2840 }
2841 }
2842 else if( dimension == 3 &&
2843 pixeltype.GetSamplesPerPixel() == 1 &&
2844 pi == PhotometricInterpretation::MONOCHROME2 &&
2845 pixeltype.GetBitsAllocated() == 1 &&
2846 pixeltype.GetBitsStored() == 1 &&
2847 pixeltype.GetHighBit() == 0 &&
2848 pixeltype.GetPixelRepresentation() == 0
2849 )
2850 {
2851 ms = MediaStorage::MultiframeSingleBitSecondaryCaptureImageStorage;
2852 // FIXME: GDCM does not handle bit packing...
2853 if( intercept != 0 || slope != 1 )
2854 {
2855 gdcmDebugMacro( "Cannot have shift/scale" );
2856 return MediaStorage::MS_END;
2857 }
2858 }
2859 else if( dimension == 3 &&
2860 pixeltype.GetSamplesPerPixel() == 1 &&
2861 pi == PhotometricInterpretation::MONOCHROME2 &&
2862 pixeltype.GetBitsAllocated() == 16 &&
2863 pixeltype.GetBitsStored() <= 16 && pixeltype.GetBitsStored() >= 9 &&
2864 pixeltype.GetHighBit() == pixeltype.GetBitsStored() - 1 &&
2865 pixeltype.GetPixelRepresentation() == 0
2866 )
2867 {
2868 ms = MediaStorage::MultiframeGrayscaleWordSecondaryCaptureImageStorage;
2869 // A.8.4.4 Multi-frame Grayscale Word SC Image IOD Content Constraints
2870 // Rescale Slope and Rescale Intercept are not constrained in this IOD to
2871 // any particular values. E.g., they may be used to recover floating
2872 // point values scaled to the integer range of the stored pixel values,
2873 // Rescale Slope may be less than one, e.g., a Rescale Slope of 1.0/65535
2874 // would allow represent floating point values from 0 to 1.0.
2875 }
2876 else if( dimension == 3 && /* A.8.5.4 Multi-frame True Color SC Image IOD Content Constraints */
2877 pixeltype.GetSamplesPerPixel() == 3 &&
2878 ( pi == PhotometricInterpretation::RGB
2879 || pi == PhotometricInterpretation::YBR_RCT
2880 || pi == PhotometricInterpretation::YBR_ICT
2881 || pi == PhotometricInterpretation::YBR_PARTIAL_420
2882 || pi == PhotometricInterpretation::YBR_FULL_422 ) &&
2883 pixeltype.GetBitsAllocated() == 8 &&
2884 pixeltype.GetBitsStored() == 8 &&
2885 pixeltype.GetHighBit() == 7 &&
2886 pixeltype.GetPixelRepresentation() == 0
2887 )
2888 {
2889 ms = MediaStorage::MultiframeTrueColorSecondaryCaptureImageStorage;
2890 if( intercept != 0 || slope != 1 )
2891 {
2892 gdcmDebugMacro( "Cannot have shift/scale" );
2893 return MediaStorage::MS_END;
2894 }
2895 }
2896 else
2897 {
2898 gdcmDebugMacro( "Cannot handle Multi Frame image in SecondaryCaptureImageStorage" );
2899 return MediaStorage::MS_END;
2900 }
2901 }
2902 // check MR Image Storage
2903 if( ms == MediaStorage::MRImageStorage )
2904 {
2905 if( intercept != 0.0 || slope != 1.0 )
2906 {
2907 if( !ForceRescaleInterceptSlope )
2908 {
2909 // hopefully this is not a lame choice:
2910 ms = MediaStorage::EnhancedMRImageStorage;
2911 }
2912 }
2913 }
2914 // check 'DX' / dim == 2:
2915 if( ms == MediaStorage::DigitalXRayImageStorageForPresentation )
2916 {
2917 if( intercept != 0.0 || slope != 1.0 )
2918 {
2919 // hopefully this is not a lame choice:
2920 ms = MediaStorage::XRay3DCraniofacialImageStorage;
2921 }
2922 }
2923 return ms;
2924 }
2925
2926 } // end namespace gdcm
2927