fw4spl
fwVtkWindowLevelLookupTable.cpp
1 /* ***** BEGIN LICENSE BLOCK *****
2  * FW4SPL - Copyright (C) IRCAD, 2004-2017.
3  * Distributed under the terms of the GNU Lesser General Public License (LGPL) as
4  * published by the Free Software Foundation.
5  * ****** END LICENSE BLOCK ****** */
6 
7 #include "fwRenderVTK/vtk/fwVtkWindowLevelLookupTable.hpp"
8 
9 #include <vtkBitArray.h>
10 #include <vtkMath.h>
11 #include <vtkObjectFactory.h>
12 
13 #include <cmath>
14 #include <iostream>
15 
16 vtkStandardNewMacro(fwVtkWindowLevelLookupTable);
17 
18 //----------------------------------------------------------------------------
19 fwVtkWindowLevelLookupTable::fwVtkWindowLevelLookupTable(int sze, int ext) :
20  vtkLookupTable(sze, ext)
21 {
22  this->Level = (this->TableRange[0] + this->TableRange[1])/2;
23  this->Window = (this->TableRange[1] - this->TableRange[0]);
24 
25  this->InverseVideo = 0;
26 
27  this->Clamping = 0;
28 
29  this->LeftClampValue[0] = 0.0;
30  this->LeftClampValue[1] = 0.0;
31  this->LeftClampValue[2] = 0.0;
32  this->LeftClampValue[3] = 0.0;
33 
34  this->RightClampValue[0] = 0.0;
35  this->RightClampValue[1] = 0.0;
36  this->RightClampValue[2] = 0.0;
37  this->RightClampValue[3] = 0.0;
38 
39  this->InvertTable = vtkUnsignedCharArray::New();
40  this->InvertTable->Register(this);
41  this->InvertTable->Delete();
42  this->InvertTable->SetNumberOfComponents(4);
43  this->InvertTable->Allocate(4*sze, 4*ext);
44 }
45 
46 fwVtkWindowLevelLookupTable::~fwVtkWindowLevelLookupTable()
47 {
48  this->InvertTable->UnRegister(this);
49  this->InvertTable = NULL;
50 }
51 
52 //------------------------------------------------------------------------------
53 
54 void fwVtkWindowLevelLookupTable::BuildInvert()
55 {
56  if (this->GetMTime() < this->InvertTime)
57  {
58  return;
59  }
60  if (this->Table->GetNumberOfTuples() < 1)
61  {
62  return;
63  }
64 
65  this->InvertTable->SetNumberOfTuples(this->NumberOfColors);
66 
67  unsigned char* tableRgba, * invertTableRgba2;
68 
69  int n = static_cast<int>(this->NumberOfColors-1);
70  for (int i = 0; i < this->NumberOfColors; i++)
71  {
72  tableRgba = this->Table->GetPointer(4*i);
73  invertTableRgba2 = this->InvertTable->WritePointer(4*(n-i), 4);
74 
75  invertTableRgba2[0] = tableRgba[0];
76  invertTableRgba2[1] = tableRgba[1];
77  invertTableRgba2[2] = tableRgba[2];
78  invertTableRgba2[3] = tableRgba[3];
79  }
80 
81  this->InvertTime.Modified();
82 
83 }
84 
85 //------------------------------------------------------------------------------
86 
87 unsigned char* fwVtkWindowLevelLookupTable::GetCurrentPointer(const vtkIdType id)
88 {
89  if(this->InverseVideo)
90  {
91  this->BuildInvert();
92  return this->InvertTable->GetPointer(4*id);
93  }
94  else
95  {
96  return this->Table->GetPointer(4*id);
97  }
98 }
99 
100 //----------------------------------------------------------------------------
101 void fwVtkWindowLevelLookupTable::PrintSelf(ostream& os, vtkIndent indent)
102 {
103  this->Superclass::PrintSelf(os, indent);
104 
105  os << indent << "Window: " << this->Window << "\n";
106  os << indent << "Level: " << this->Level << "\n";
107  os << indent << "InverseVideo: "
108  << (this->InverseVideo ? "On\n" : "Off\n");
109  os << indent << "LeftClampValue : ("
110  << this->LeftClampValue[0] << ", "
111  << this->LeftClampValue[1] << ", "
112  << this->LeftClampValue[2] << ", "
113  << this->LeftClampValue[3] << ")\n";
114  os << indent << "RightClampValue : ("
115  << this->RightClampValue[0] << ", "
116  << this->RightClampValue[1] << ", "
117  << this->RightClampValue[2] << ", "
118  << this->RightClampValue[3] << ")\n";
119 }
120 
121 //----------------------------------------------------------------------------
122 // Apply log to value, with appropriate constraints.
123 inline double vtkApplyLogScale(double v, const double range[2],
124  const double logRange[2])
125 {
126  // is the range set for negative numbers?
127  if (range[0] < 0)
128  {
129  if (v < 0)
130  {
131  v = log10(-static_cast<double>(v));
132  }
133  else if (range[0] > range[1])
134  {
135  v = logRange[0];
136  }
137  else
138  {
139  v = logRange[1];
140  }
141  }
142  else
143  {
144  if (v > 0)
145  {
146  v = log10(static_cast<double>(v));
147  }
148  else if (range[0] < range[1])
149  {
150  v = logRange[0];
151  }
152  else
153  {
154  v = logRange[1];
155  }
156  }
157  return v;
158 }
159 
160 //----------------------------------------------------------------------------
161 // Apply shift/scale to the scalar value v and do table lookup.
162 inline unsigned char* vtkLinearLookup(double v,
163  unsigned char* table,
164  double maxIndex,
165  double shift, double scale,
166  unsigned char* nanColor,
167  unsigned char* leftColor,
168  unsigned char* rightColor
169  )
170 {
171 
172  if (vtkMath::IsNan(v))
173  {
174  return nanColor;
175  }
176 
177  double findx = (v + shift)*scale;
178  if (findx < 0.0)
179  {
180  return leftColor;
181  }
182  else if (findx > maxIndex)
183  {
184  return rightColor;
185  }
186  return &table[4*static_cast<unsigned int>(findx)];
187  /* round
188  return &table[4*(unsigned int)(findx + 0.5f)];
189  */
190 }
191 
192 //----------------------------------------------------------------------------
193 // There is a little more to this than simply taking the log10 of the
194 // two range values: we do conversion of negative ranges to positive
195 // ranges, and conversion of zero to a 'very small number'
196 void fwVtkWindowLevelLookupTableLogRange(const double range[2], double logRange[2])
197 {
198  double rmin = range[0];
199  double rmax = range[1];
200 
201  if (rmin == 0)
202  {
203  rmin = 1.0e-6*(rmax - rmin);
204  if (rmax < 0)
205  {
206  rmin = -rmin;
207  }
208  }
209  if (rmax == 0)
210  {
211  rmax = 1.0e-6*(rmin - rmax);
212  if (rmin < 0)
213  {
214  rmax = -rmax;
215  }
216  }
217  if (rmin < 0 && rmax < 0)
218  {
219  logRange[0] = log10(-static_cast<double>(rmin));
220  logRange[1] = log10(-static_cast<double>(rmax));
221  }
222  else if (rmin > 0 && rmax > 0)
223  {
224  logRange[0] = log10(static_cast<double>(rmin));
225  logRange[1] = log10(static_cast<double>(rmax));
226  }
227 }
228 
229 //----------------------------------------------------------------------------
230 // accelerate the mapping by copying the data in 32-bit chunks instead
231 // of 8-bit chunks
232 template<class T>
233 void fwVtkWindowLevelLookupTableMapData(fwVtkWindowLevelLookupTable* self, T* input,
234  unsigned char* output, int length,
235  int inIncr, int outFormat)
236 {
237  int i = length;
238  double* range = self->GetTableRange();
239  double maxIndex = self->GetNumberOfColors() - 1;
240  double shift, scale;
241  unsigned char* table = self->GetCurrentPointer(0);
242  unsigned char* cptr;
243  double alpha;
244 
245  unsigned char nanColor[4];
246  unsigned char selfLeftColor[4];
247  unsigned char selfRightColor[4];
248 
249  for (int c = 0; c < 4; c++)
250  {
251  nanColor[c] = static_cast<unsigned char>(self->GetNanColor()[c]*255.0);
252  selfLeftColor[c] = static_cast<unsigned char>(self->GetLeftClampValue()[c]*255.0);
253  selfRightColor[c] = static_cast<unsigned char>(self->GetRightClampValue()[c]*255.0);
254  }
255 
256  unsigned char* leftColor = selfLeftColor;
257  unsigned char* rightColor = selfRightColor;
258 
259  if(self->GetClamping())
260  {
261  leftColor = &table[0];
262  rightColor = &table[(self->GetNumberOfColors()-1) * 4];
263  }
264 
265  if ( (alpha = self->GetAlpha()) >= 1.0 ) //no blending required
266  {
267  if (self->GetScale() == VTK_SCALE_LOG10)
268  {
269  double val;
270  double logRange[2];
271  fwVtkWindowLevelLookupTableLogRange(range, logRange);
272  shift = -logRange[0];
273  if (logRange[1] <= logRange[0])
274  {
275  scale = VTK_DOUBLE_MAX;
276  }
277  else
278  {
279  /* while this looks like the wrong scale, it is the correct scale
280  * taking into account the truncation to int that happens below. */
281  // scale = (maxIndex + 1)/(logRange[1] - logRange[0]);
282  // Fixed scale : finaly, no truncation happends below
283  scale = (maxIndex)/(logRange[1] - logRange[0]);
284  }
285  if (outFormat == VTK_RGBA)
286  {
287  while (--i >= 0)
288  {
289  val = vtkApplyLogScale(*input, range, logRange);
290  cptr = vtkLinearLookup(val, table, maxIndex, shift, scale, nanColor, leftColor, rightColor);
291  *output++ = *cptr++;
292  *output++ = *cptr++;
293  *output++ = *cptr++;
294  *output++ = *cptr++;
295  input += inIncr;
296  }
297  }
298  else if (outFormat == VTK_RGB)
299  {
300  while (--i >= 0)
301  {
302  val = vtkApplyLogScale(*input, range, logRange);
303  cptr = vtkLinearLookup(val, table, maxIndex, shift, scale, nanColor, leftColor, rightColor);
304  *output++ = *cptr++;
305  *output++ = *cptr++;
306  *output++ = *cptr++;
307  input += inIncr;
308  }
309  }
310  else if (outFormat == VTK_LUMINANCE_ALPHA)
311  {
312  while (--i >= 0)
313  {
314  val = vtkApplyLogScale(*input, range, logRange);
315  cptr = vtkLinearLookup(val, table, maxIndex, shift, scale, nanColor, leftColor, rightColor);
316  *output++ = static_cast<unsigned char>(cptr[0]*0.30 + cptr[1]*0.59 +
317  cptr[2]*0.11 + 0.5);
318  *output++ = cptr[3];
319  input += inIncr;
320  }
321  }
322  else // outFormat == VTK_LUMINANCE
323  {
324  while (--i >= 0)
325  {
326  val = vtkApplyLogScale(*input, range, logRange);
327  cptr = vtkLinearLookup(val, table, maxIndex, shift, scale, nanColor, leftColor, rightColor);
328  *output++ = static_cast<unsigned char>(cptr[0]*0.30 + cptr[1]*0.59 +
329  cptr[2]*0.11 + 0.5);
330  input += inIncr;
331  }
332  }
333  }//if log scale
334 
335  else //not log scale
336  {
337  shift = -range[0];
338  if (range[1] <= range[0])
339  {
340  scale = VTK_DOUBLE_MAX;
341  }
342  else
343  {
344  /* while this looks like the wrong scale, it is the correct scale
345  * taking into account the truncation to int that happens below. */
346  // scale = (maxIndex + 1)/(range[1] - range[0]);
347  // Fixed scale : finaly, no truncation happends below
348  scale = (maxIndex)/(range[1] - range[0]);
349  }
350 
351  if (outFormat == VTK_RGBA)
352  {
353  while (--i >= 0)
354  {
355  cptr = vtkLinearLookup(*input, table, maxIndex, shift, scale,
356  nanColor, leftColor, rightColor);
357  *output++ = *cptr++;
358  *output++ = *cptr++;
359  *output++ = *cptr++;
360  *output++ = *cptr++;
361  input += inIncr;
362  }
363  }
364  else if (outFormat == VTK_RGB)
365  {
366  while (--i >= 0)
367  {
368  cptr = vtkLinearLookup(*input, table, maxIndex, shift, scale,
369  nanColor, leftColor, rightColor);
370  *output++ = *cptr++;
371  *output++ = *cptr++;
372  *output++ = *cptr++;
373  input += inIncr;
374  }
375  }
376  else if (outFormat == VTK_LUMINANCE_ALPHA)
377  {
378  while (--i >= 0)
379  {
380  cptr = vtkLinearLookup(*input, table, maxIndex, shift, scale,
381  nanColor, leftColor, rightColor);
382  *output++ = static_cast<unsigned char>(cptr[0]*0.30 + cptr[1]*0.59 +
383  cptr[2]*0.11 + 0.5);
384  *output++ = cptr[3];
385  input += inIncr;
386  }
387  }
388  else // outFormat == VTK_LUMINANCE
389  {
390  while (--i >= 0)
391  {
392  cptr = vtkLinearLookup(*input, table, maxIndex, shift, scale,
393  nanColor, leftColor, rightColor);
394  *output++ = static_cast<unsigned char>(cptr[0]*0.30 + cptr[1]*0.59 +
395  cptr[2]*0.11 + 0.5);
396  input += inIncr;
397  }
398  }
399  }//if not log lookup
400  }//if blending not needed
401 
402  else //blend with the specified alpha
403  {
404  if (self->GetScale() == VTK_SCALE_LOG10)
405  {
406  double val;
407  double logRange[2];
408  fwVtkWindowLevelLookupTableLogRange(range, logRange);
409  shift = -logRange[0];
410  if (logRange[1] <= logRange[0])
411  {
412  scale = VTK_DOUBLE_MAX;
413  }
414  else
415  {
416  /* while this looks like the wrong scale, it is the correct scale
417  * taking into account the truncation to int that happens below. */
418  // scale = (maxIndex + 1)/(logRange[1] - logRange[0]);
419  // Fixed scale : finaly, no truncation happends below
420  scale = (maxIndex)/(logRange[1] - logRange[0]);
421  }
422  if (outFormat == VTK_RGBA)
423  {
424  while (--i >= 0)
425  {
426  val = vtkApplyLogScale(*input, range, logRange);
427  cptr = vtkLinearLookup(val, table, maxIndex, shift, scale, nanColor, leftColor, rightColor);
428  *output++ = *cptr++;
429  *output++ = *cptr++;
430  *output++ = *cptr++;
431  *output++ = static_cast<unsigned char>((*cptr)*alpha); cptr++;
432  input += inIncr;
433  }
434  }
435  else if (outFormat == VTK_RGB)
436  {
437  while (--i >= 0)
438  {
439  val = vtkApplyLogScale(*input, range, logRange);
440  cptr = vtkLinearLookup(val, table, maxIndex, shift, scale, nanColor, leftColor, rightColor);
441  *output++ = *cptr++;
442  *output++ = *cptr++;
443  *output++ = *cptr++;
444  input += inIncr;
445  }
446  }
447  else if (outFormat == VTK_LUMINANCE_ALPHA)
448  {
449  while (--i >= 0)
450  {
451  val = vtkApplyLogScale(*input, range, logRange);
452  cptr = vtkLinearLookup(val, table, maxIndex, shift, scale, nanColor, leftColor, rightColor);
453  *output++ = static_cast<unsigned char>(cptr[0]*0.30 + cptr[1]*0.59 +
454  cptr[2]*0.11 + 0.5);
455  *output++ = static_cast<unsigned char>(alpha*cptr[3]);
456  input += inIncr;
457  }
458  }
459  else // outFormat == VTK_LUMINANCE
460  {
461  while (--i >= 0)
462  {
463  val = vtkApplyLogScale(*input, range, logRange);
464  cptr = vtkLinearLookup(val, table, maxIndex, shift, scale, nanColor, leftColor, rightColor);
465  *output++ = static_cast<unsigned char>(cptr[0]*0.30 + cptr[1]*0.59 +
466  cptr[2]*0.11 + 0.5);
467  input += inIncr;
468  }
469  }
470  }//log scale with blending
471 
472  else //no log scale with blending
473  {
474  shift = -range[0];
475  if (range[1] <= range[0])
476  {
477  scale = VTK_DOUBLE_MAX;
478  }
479  else
480  {
481  /* while this looks like the wrong scale, it is the correct scale
482  * taking into account the truncation to int that happens below. */
483  // scale = (maxIndex + 1)/(range[1] - range[0]);
484  // Fixed scale : finaly, no truncation happends below
485  scale = (maxIndex)/(range[1] - range[0]);
486  }
487 
488  if (outFormat == VTK_RGBA)
489  {
490  while (--i >= 0)
491  {
492  cptr = vtkLinearLookup(*input, table, maxIndex, shift, scale,
493  nanColor, leftColor, rightColor);
494  *output++ = *cptr++;
495  *output++ = *cptr++;
496  *output++ = *cptr++;
497  *output++ = static_cast<unsigned char>((*cptr)*alpha); cptr++;
498  input += inIncr;
499  }
500  }
501  else if (outFormat == VTK_RGB)
502  {
503  while (--i >= 0)
504  {
505  cptr = vtkLinearLookup(*input, table, maxIndex, shift, scale,
506  nanColor, leftColor, rightColor);
507  *output++ = *cptr++;
508  *output++ = *cptr++;
509  *output++ = *cptr++;
510  input += inIncr;
511  }
512  }
513  else if (outFormat == VTK_LUMINANCE_ALPHA)
514  {
515  while (--i >= 0)
516  {
517  cptr = vtkLinearLookup(*input, table, maxIndex, shift, scale,
518  nanColor, leftColor, rightColor);
519  *output++ = static_cast<unsigned char>(cptr[0]*0.30 + cptr[1]*0.59 +
520  cptr[2]*0.11 + 0.5);
521  *output++ = static_cast<unsigned char>(cptr[3]*alpha);
522  input += inIncr;
523  }
524  }
525  else // outFormat == VTK_LUMINANCE
526  {
527  while (--i >= 0)
528  {
529  cptr = vtkLinearLookup(*input, table, maxIndex, shift, scale,
530  nanColor, leftColor, rightColor);
531  *output++ = static_cast<unsigned char>(cptr[0]*0.30 + cptr[1]*0.59 +
532  cptr[2]*0.11 + 0.5);
533  input += inIncr;
534  }
535  }
536  }//no log scale
537  }//alpha blending
538 
539 }
540 
541 //----------------------------------------------------------------------------
542 // Although this is a relatively expensive calculation,
543 // it is only done on the first render. Colors are cached
544 // for subsequent renders.
545 template<class T>
546 void fwVtkWindowLevelLookupTableMapMag(fwVtkWindowLevelLookupTable* self, T* input,
547  unsigned char* output, int length,
548  int inIncr, int outFormat)
549 {
550  double tmp, sum;
551  double* mag;
552  int i, j;
553 
554  mag = new double[length];
555  for (i = 0; i < length; ++i)
556  {
557  sum = 0;
558  for (j = 0; j < inIncr; ++j)
559  {
560  tmp = static_cast<double>(*input);
561  sum += (tmp * tmp);
562  ++input;
563  }
564  mag[i] = sqrt(sum);
565  }
566 
567  fwVtkWindowLevelLookupTableMapData(self, mag, output, length, 1, outFormat);
568 
569  delete [] mag;
570 }
571 
572 //----------------------------------------------------------------------------
573 void fwVtkWindowLevelLookupTable::MapScalarsThroughTable2(void* input,
574  unsigned char* output,
575  int inputDataType,
576  int numberOfValues,
577  int inputIncrement,
578  int outputFormat)
579 {
580  if (this->UseMagnitude && inputIncrement > 1)
581  {
582  switch (inputDataType)
583  {
584  vtkTemplateMacro(
585  fwVtkWindowLevelLookupTableMapMag(this, static_cast<VTK_TT*>(input), output,
586  numberOfValues, inputIncrement, outputFormat);
587  return
588  );
589  case VTK_BIT:
590  vtkErrorMacro("Cannot comput magnitude of bit array.");
591  break;
592  default:
593  vtkErrorMacro(<< "MapImageThroughTable: Unknown input ScalarType");
594  }
595  }
596 
597  switch (inputDataType)
598  {
599  case VTK_BIT:
600  {
601  vtkIdType i, id;
602  vtkBitArray* bitArray = vtkBitArray::New();
603  bitArray->SetVoidArray(input, numberOfValues, 1);
604  vtkUnsignedCharArray* newInput = vtkUnsignedCharArray::New();
605  newInput->SetNumberOfValues(numberOfValues);
606  for (id = i = 0; i < numberOfValues; i++, id += inputIncrement)
607  {
608  newInput->SetValue(i, static_cast<vtkUnsignedCharArray::ValueType>(bitArray->GetValue(id)));
609  }
610  fwVtkWindowLevelLookupTableMapData(this,
611  static_cast<unsigned char*>(newInput->GetPointer(0)),
612  output, numberOfValues,
613  inputIncrement, outputFormat);
614  newInput->Delete();
615  bitArray->Delete();
616  }
617  break;
618 
619  vtkTemplateMacro(
620  fwVtkWindowLevelLookupTableMapData(this, static_cast<VTK_TT*>(input), output,
621  numberOfValues, inputIncrement, outputFormat)
622  );
623  default:
624  vtkErrorMacro(<< "MapImageThroughTable: Unknown input ScalarType");
625  return;
626  }
627 }
628 
fwVtkWindowLevelLookupTable
Reinplementation of vtkWindowLevelLookupTable : add specific out-of-bounds colors.
Definition: fwVtkWindowLevelLookupTable.hpp:15