Medical Imaging Interaction Toolkit  2018.4.99-a3d2e8fb
Medical Imaging Interaction Toolkit
vtkMitkOpenGLVolumeTextureMapper3D.cpp
Go to the documentation of this file.
1 /*============================================================================
2 
3 The Medical Imaging Interaction Toolkit (MITK)
4 
5 Copyright (c) German Cancer Research Center (DKFZ)
6 All rights reserved.
7 
8 Use of this source code is governed by a 3-clause BSD license that can be
9 found in the LICENSE file.
10 
11 ============================================================================*/
12 
13 #ifdef _OPENMP
14 #include <omp.h>
15 #endif
16 
17 #include "vtkWindows.h"
18 #include "mitkCommon.h"
19 #include "vtkMitkOpenGLVolumeTextureMapper3D.h"
20 
21 #define GPU_INFO MITK_INFO("mapper.vr")
22 #define GPU_WARN MITK_WARN("mapper.vr")
23 
24 #include "vtkCamera.h"
25 #include "vtkDataArray.h"
26 #include "vtkImageData.h"
27 #include "vtkLight.h"
28 #include "vtkLightCollection.h"
29 #include "vtkMath.h"
30 #include "vtkMatrix4x4.h"
31 #include "vtkObjectFactory.h"
32 #include "vtkOpenGLExtensionManager.h"
33 #include "vtkPlane.h"
34 #include "vtkPlaneCollection.h"
35 #include "vtkPointData.h"
36 #include "vtkRenderWindow.h"
37 #include "vtkRenderer.h"
38 #include "vtkTimerLog.h"
39 #include "vtkTransform.h"
40 #include "vtkVolumeProperty.h"
41 #include "vtkgl.h"
42 
43 #include "vtkOpenGLRenderWindow.h"
44 
45 #define myGL_COMPRESSED_RGB_S3TC_DXT1_EXT 0x83F0
46 #define myGL_COMPRESSED_LUMINANCE_ALPHA_LATC2_EXT 0x8C72
47 #define myGL_COMPRESSED_RGBA_S3TC_DXT5_EXT 0x83F3
48 
50 
51  //# We need some temporary variables
52  "TEMP index, normal, finalColor;\n"
53  "TEMP temp,temp1, temp2, temp3,temp4; \n"
54  "TEMP sampleColor;\n"
55  "TEMP ndotl, ndoth, ndotv; \n"
56  "TEMP lightInfo, lightResult;\n"
57 
58  //# We are going to use the first
59  //# texture coordinate
60  "ATTRIB tex0 = fragment.texcoord[0];\n"
61 
62  //# This is the lighting information
63  "PARAM lightDirection = program.local[0];\n"
64  "PARAM halfwayVector = program.local[1];\n"
65  "PARAM coefficient = program.local[2];\n"
66  "PARAM lightDiffColor = program.local[3]; \n"
67  "PARAM lightSpecColor = program.local[4]; \n"
68  "PARAM viewVector = program.local[5];\n"
69  "PARAM constants = program.local[6];\n"
70 
71  //# This is our output color
72  "OUTPUT out = result.color;\n"
73 
74  //# Look up the gradient direction
75  //# in the third volume
76  "TEX temp2, tex0, texture[0], 3D;\n"
77 
78  //# This normal is stored 0 to 1, change to -1 to 1
79  //# by multiplying by 2.0 then adding -1.0.
80  "MAD normal, temp2, constants.x, constants.y;\n"
81 
82  "DP3 temp4, normal, normal;\n"
83  "RSQ temp, temp4.x;\n"
84  "MUL normal, normal, temp;\n"
85 
86  //"RCP temp4,temp.x;\n"
87 
88  //"MUL temp2.w,temp2.w,temp4.x;\n"
89 
90  //"MUL_SAT temp2.w,temp2.w,6.0;\n"
91 
92  "TEX sampleColor, tex0, texture[1], 3D;\n"
93 
94  //# Take the dot product of the light
95  //# direction and the normal
96  "DP3 ndotl, normal, lightDirection;\n"
97 
98  //# Take the dot product of the halfway
99  //# vector and the normal
100  "DP3 ndoth, normal, halfwayVector;\n"
101 
102  "DP3 ndotv, normal, viewVector;\n"
103 
104  //# flip if necessary for two sided lighting
105  "MUL temp3, ndotl, constants.y; \n"
106  "CMP ndotl, ndotv, ndotl, temp3;\n"
107  "MUL temp3, ndoth, constants.y; \n"
108  "CMP ndoth, ndotv, ndoth, temp3;\n"
109 
110  //# put the pieces together for a LIT operation
111  "MOV lightInfo.x, ndotl.x; \n"
112  "MOV lightInfo.y, ndoth.x; \n"
113  "MOV lightInfo.w, coefficient.w; \n"
114 
115  //# compute the lighting
116  "LIT lightResult, lightInfo;\n"
117 
118  //# COLOR FIX
119  "MUL lightResult, lightResult, 4.0;\n"
120 
121  //# This is the ambient contribution
122  "MUL finalColor, coefficient.x, sampleColor;\n"
123 
124  //# This is the diffuse contribution
125  "MUL temp3, lightDiffColor, sampleColor;\n"
126  "MUL temp3, temp3, lightResult.y;\n"
127  "ADD finalColor, finalColor, temp3;\n"
128 
129  //# This is th specular contribution
130  "MUL temp3, lightSpecColor, lightResult.z; \n"
131 
132  //# Add specular into result so far, and replace
133  //# with the original alpha.
134  "ADD out, finalColor, temp3;\n"
135  "MOV out.w, temp2.w;\n"
136 
137  "END\n";
138 
140 
141  //# This is the fragment program for one
142  //# component data with shading
143 
144  //# We need some temporary variables
145  "TEMP index, normal, finalColor;\n"
146  "TEMP temp,temp1, temp2, temp3,temp4; \n"
147  "TEMP sampleColor;\n"
148  "TEMP ndotl, ndoth, ndotv; \n"
149  "TEMP lightInfo, lightResult;\n"
150 
151  //# We are going to use the first
152  //# texture coordinate
153  "ATTRIB tex0 = fragment.texcoord[0];\n"
154 
155  //# This is the lighting information
156  "PARAM lightDirection = program.local[0];\n"
157  "PARAM halfwayVector = program.local[1];\n"
158  "PARAM coefficient = program.local[2];\n"
159  "PARAM lightDiffColor = program.local[3]; \n"
160  "PARAM lightSpecColor = program.local[4]; \n"
161  "PARAM viewVector = program.local[5];\n"
162  "PARAM constants = program.local[6];\n"
163 
164  //# This is our output color
165  "OUTPUT out = result.color;\n"
166 
167  //# Look up the gradient direction
168  //# in the third volume
169  "TEX temp2, tex0, texture[0], 3D;\n"
170 
171  // Gradient Compution
172 
173  //# Look up the scalar value / gradient
174  //# magnitude in the first volume
175  //"TEX temp1, tex0, texture[0], 3D;\n"
176 
177  /*
178 
179  "ADD temp3,tex0,{-0.005,0,0};\n"
180  "TEX temp2,temp3, texture[0], 3D;\n"
181 
182  //"ADD temp3,tex0,{ 0.005,0,0};\n"
183  //"TEX temp1,temp3, texture[0], 3D;\n"
184 
185  "SUB normal.x,temp2.y,temp1.y;\n"
186 
187  "ADD temp3,tex0,{0,-0.005,0};\n"
188  "TEX temp2,temp3, texture[0], 3D;\n"
189 
190  //"ADD temp3,tex0,{0, 0.005,0};\n"
191  //"TEX temp1,temp3, texture[0], 3D;\n"
192 
193  "SUB normal.y,temp2.y,temp1.y;\n"
194 
195  "ADD temp3,tex0,{0,0,-0.005};\n"
196  "TEX temp2,temp3, texture[0], 3D;\n"
197 
198  //"ADD temp3,tex0,{0,0, 0.005};\n"
199  //"TEX temp1,temp3, texture[0], 3D;\n"
200 
201  "SUB normal.z,temp2.y,temp1.y;\n"
202 
203  */
204 
205  //"MOV normal,{1,1,1};\n"
206 
207  "MOV index.x,temp2.a;\n"
208 
209  //# This normal is stored 0 to 1, change to -1 to 1
210  //# by multiplying by 2.0 then adding -1.0.
211  "MAD normal, temp2, constants.x, constants.y;\n"
212 
213  //# Swizzle this to use (a,r) as texture
214  //# coordinates
215  //"SWZ index, temp1, a, r, 1, 1;\n"
216 
217  //# Use this coordinate to look up a
218  //# final color in the third texture
219  //# (this is a 2D texture)
220 
221  "DP3 temp4, normal, normal;\n"
222 
223  "RSQ temp, temp4.x;\n"
224 
225  "RCP temp4,temp.x;\n"
226 
227  "MUL normal, normal, temp;\n"
228 
229  "MOV index.y, temp4.x;\n"
230 
231  "TEX sampleColor, index, texture[1], 2D;\n"
232 
233  //"MUL sampleColor.w,sampleColor.w,temp4.x;\n"
234 
235  //# Take the dot product of the light
236  //# direction and the normal
237  "DP3 ndotl, normal, lightDirection;\n"
238 
239  //# Take the dot product of the halfway
240  //# vector and the normal
241  "DP3 ndoth, normal, halfwayVector;\n"
242 
243  "DP3 ndotv, normal, viewVector;\n"
244 
245  //# flip if necessary for two sided lighting
246  "MUL temp3, ndotl, constants.y; \n"
247  "CMP ndotl, ndotv, ndotl, temp3;\n"
248  "MUL temp3, ndoth, constants.y; \n"
249  "CMP ndoth, ndotv, ndoth, temp3;\n"
250 
251  //# put the pieces together for a LIT operation
252  "MOV lightInfo.x, ndotl.x; \n"
253  "MOV lightInfo.y, ndoth.x; \n"
254  "MOV lightInfo.w, coefficient.w; \n"
255 
256  //# compute the lighting
257  "LIT lightResult, lightInfo;\n"
258 
259  //# COLOR FIX
260  "MUL lightResult, lightResult, 4.0;\n"
261 
262  //# This is the ambient contribution
263  "MUL finalColor, coefficient.x, sampleColor;\n"
264 
265  //# This is the diffuse contribution
266  "MUL temp3, lightDiffColor, sampleColor;\n"
267  "MUL temp3, temp3, lightResult.y;\n"
268  "ADD finalColor, finalColor, temp3;\n"
269 
270  //# This is th specular contribution
271  "MUL temp3, lightSpecColor, lightResult.z; \n"
272 
273  //# Add specular into result so far, and replace
274  //# with the original alpha.
275  "ADD out, finalColor, temp3;\n"
276  "MOV out.w, sampleColor.w;\n"
277 
278  "END\n";
279 
280 //#ifndef VTK_IMPLEMENT_MESA_CXX
281 vtkStandardNewMacro(vtkMitkOpenGLVolumeTextureMapper3D);
282 //#endif
283 
284 vtkMitkOpenGLVolumeTextureMapper3D::vtkMitkOpenGLVolumeTextureMapper3D()
285 {
286  // GPU_INFO << "vtkMitkOpenGLVolumeTextureMapper3D";
287 
288  this->Initialized = 0;
289  this->Volume1Index = 0;
290  this->Volume2Index = 0;
291  this->Volume3Index = 0;
292  this->ColorLookupIndex = 0;
293  this->AlphaLookupIndex = 0;
294  this->RenderWindow = nullptr;
295  this->SupportsCompressedTexture = false;
296 
297  prgOneComponentShade = 0;
298  prgRGBAShade = 0;
299 }
300 
301 vtkMitkOpenGLVolumeTextureMapper3D::~vtkMitkOpenGLVolumeTextureMapper3D()
302 {
303  // GPU_INFO << "~vtkMitkOpenGLVolumeTextureMapper3D";
304  if (prgOneComponentShade)
305  vtkgl::DeleteProgramsARB(1, &prgOneComponentShade);
306 
307  if (prgRGBAShade)
308  vtkgl::DeleteProgramsARB(1, &prgRGBAShade);
309 }
310 
311 // Release the graphics resources used by this texture.
312 void vtkMitkOpenGLVolumeTextureMapper3D::ReleaseGraphicsResources(vtkWindow *renWin)
313 {
314  // GPU_INFO << "ReleaseGraphicsResources";
315 
316  if ((this->Volume1Index || this->Volume2Index || this->Volume3Index || this->ColorLookupIndex) && renWin)
317  {
318  static_cast<vtkRenderWindow *>(renWin)->MakeCurrent();
319 #ifdef GL_VERSION_1_1
320  // free any textures
321  this->DeleteTextureIndex(&this->Volume1Index);
322  this->DeleteTextureIndex(&this->Volume2Index);
323  this->DeleteTextureIndex(&this->Volume3Index);
324  this->DeleteTextureIndex(&this->ColorLookupIndex);
325  this->DeleteTextureIndex(&this->AlphaLookupIndex);
326 #endif
327  }
328  this->Volume1Index = 0;
329  this->Volume2Index = 0;
330  this->Volume3Index = 0;
331  this->ColorLookupIndex = 0;
332  this->RenderWindow = nullptr;
333  this->SupportsCompressedTexture = false;
334  this->SupportsNonPowerOfTwoTextures = false;
335 
336  this->Modified();
337 }
338 
339 // Release the graphics resources used by this texture.
340 void vtkMitkOpenGLVolumeTextureMapper3D::ReleaseGraphicsResources(mitk::BaseRenderer *renderer)
341 {
342  // GPU_INFO << "ReleaseGraphicsResources";
343 
344  vtkWindow *renWin = renderer->GetVtkRenderer()->GetRenderWindow();
345 
346  if ((this->Volume1Index || this->Volume2Index || this->Volume3Index || this->ColorLookupIndex) && renWin)
347  {
348  static_cast<vtkRenderWindow *>(renWin)->MakeCurrent();
349 #ifdef GL_VERSION_1_1
350  // free any textures
351  this->DeleteTextureIndex(&this->Volume1Index);
352  this->DeleteTextureIndex(&this->Volume2Index);
353  this->DeleteTextureIndex(&this->Volume3Index);
354  this->DeleteTextureIndex(&this->ColorLookupIndex);
355  this->DeleteTextureIndex(&this->AlphaLookupIndex);
356 #endif
357  }
358  this->Volume1Index = 0;
359  this->Volume2Index = 0;
360  this->Volume3Index = 0;
361  this->ColorLookupIndex = 0;
362  this->RenderWindow = nullptr;
363  this->SupportsCompressedTexture = false;
364  this->SupportsNonPowerOfTwoTextures = false;
365 
366  this->Modified();
367 }
368 
369 void vtkMitkOpenGLVolumeTextureMapper3D::Render(vtkRenderer *ren, vtkVolume *vol)
370 {
371  // GPU_INFO << "Render";
372 
373  ren->GetRenderWindow()->MakeCurrent();
374 
375  if (!this->Initialized)
376  {
377  // this->Initialize();
378  this->Initialize(ren);
379  }
380 
381  if (!this->RenderPossible)
382  {
383  vtkErrorMacro("required extensions not supported");
384  return;
385  }
386 
387  vtkMatrix4x4 *matrix = vtkMatrix4x4::New();
388  vtkPlaneCollection *clipPlanes;
389  vtkPlane *plane;
390  int numClipPlanes = 0;
391  double planeEquation[4];
392 
393  // build transformation
394  vol->GetMatrix(matrix);
395  matrix->Transpose();
396 
397  glPushAttrib(GL_ENABLE_BIT | GL_COLOR_BUFFER_BIT | GL_STENCIL_BUFFER_BIT | GL_DEPTH_BUFFER_BIT | GL_POLYGON_BIT |
398  GL_TEXTURE_BIT);
399 
400  int i;
401 
402  // Use the OpenGL clip planes
403  clipPlanes = this->ClippingPlanes;
404  if (clipPlanes)
405  {
406  numClipPlanes = clipPlanes->GetNumberOfItems();
407  if (numClipPlanes > 6)
408  {
409  vtkErrorMacro(<< "OpenGL guarantees only 6 additional clipping planes");
410  }
411 
412  for (i = 0; i < numClipPlanes; i++)
413  {
414  glEnable(static_cast<GLenum>(GL_CLIP_PLANE0 + i));
415 
416  plane = static_cast<vtkPlane *>(clipPlanes->GetItemAsObject(i));
417 
418  planeEquation[0] = plane->GetNormal()[0];
419  planeEquation[1] = plane->GetNormal()[1];
420  planeEquation[2] = plane->GetNormal()[2];
421  planeEquation[3] = -(planeEquation[0] * plane->GetOrigin()[0] + planeEquation[1] * plane->GetOrigin()[1] +
422  planeEquation[2] * plane->GetOrigin()[2]);
423  glClipPlane(static_cast<GLenum>(GL_CLIP_PLANE0 + i), planeEquation);
424  }
425  }
426 
427  // insert model transformation
428  glMatrixMode(GL_MODELVIEW);
429  glPushMatrix();
430  glMultMatrixd(matrix->Element[0]);
431 
432  glColor4f(1.0, 1.0, 1.0, 1.0);
433 
434  // Turn lighting off - the polygon textures already have illumination
435  glDisable(GL_LIGHTING);
436 
437  // vtkGraphicErrorMacro(ren->GetRenderWindow(),"Before actual render method");
438 
439  this->RenderFP(ren, vol);
440 
441  // pop transformation matrix
442  glMatrixMode(GL_MODELVIEW);
443  glPopMatrix();
444 
445  matrix->Delete();
446  glPopAttrib();
447 }
448 
449 void vtkMitkOpenGLVolumeTextureMapper3D::RenderFP(vtkRenderer *ren, vtkVolume *vol)
450 {
451  // GPU_INFO << "RenderFP";
452 
453  /*
454  glAlphaFunc (GL_GREATER, static_cast<GLclampf>(1.0/255.0));
455  glEnable (GL_ALPHA_TEST);
456  */
457 
458  glEnable(GL_BLEND);
459  glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
460 
461  int components = this->GetInput()->GetNumberOfScalarComponents();
462  switch (components)
463  {
464  case 1:
465  this->RenderOneIndependentShadeFP(ren, vol);
466  break;
467 
468  case 4:
469  this->RenderRGBAShadeFP(ren, vol);
470  break;
471  }
472 
473  vtkgl::ActiveTexture(vtkgl::TEXTURE2);
474  glDisable(GL_TEXTURE_2D);
475  glDisable(vtkgl::TEXTURE_3D);
476 
477  vtkgl::ActiveTexture(vtkgl::TEXTURE1);
478  glDisable(GL_TEXTURE_2D);
479  glDisable(vtkgl::TEXTURE_3D);
480 
481  vtkgl::ActiveTexture(vtkgl::TEXTURE0);
482  glDisable(GL_TEXTURE_2D);
483  glDisable(vtkgl::TEXTURE_3D);
484 
485  glDisable(GL_BLEND);
486 }
487 
488 void vtkMitkOpenGLVolumeTextureMapper3D::DeleteTextureIndex(GLuint *index)
489 {
490  // GPU_INFO << "DeleteTextureIndex";
491 
492  if (glIsTexture(*index))
493  {
494  GLuint tempIndex;
495  tempIndex = *index;
496  glDeleteTextures(1, &tempIndex);
497  *index = 0;
498  }
499 }
500 
501 void vtkMitkOpenGLVolumeTextureMapper3D::CreateTextureIndex(GLuint *index)
502 {
503  // GPU_INFO << "CreateTextureIndex";
504 
505  GLuint tempIndex = 0;
506  glGenTextures(1, &tempIndex);
507  *index = static_cast<long>(tempIndex);
508 }
509 
510 void vtkMitkOpenGLVolumeTextureMapper3D::RenderPolygons(vtkRenderer *ren, vtkVolume *vol, int stages[4])
511 {
512  // GPU_INFO << "RenderPolygons";
513 
514  vtkRenderWindow *renWin = ren->GetRenderWindow();
515 
516  if (renWin->CheckAbortStatus())
517  {
518  return;
519  }
520 
521  double bounds[27][6];
522  float distance2[27];
523 
524  int numIterations;
525  int i, j, k;
526 
527  // No cropping case - render the whole thing
528  if (!this->Cropping)
529  {
530  // Use the input data bounds - we'll take care of the volume's
531  // matrix during rendering
532  this->GetInput()->GetBounds(bounds[0]);
533  numIterations = 1;
534  }
535  // Simple cropping case - render the subvolume
536  else if (this->CroppingRegionFlags == 0x2000)
537  {
538  this->GetCroppingRegionPlanes(bounds[0]);
539  numIterations = 1;
540  }
541  // Complex cropping case - render each region in back-to-front order
542  else
543  {
544  // Get the camera position
545  double camPos[4];
546  ren->GetActiveCamera()->GetPosition(camPos);
547 
548  double volBounds[6];
549  this->GetInput()->GetBounds(volBounds);
550 
551  // Pass camera through inverse volume matrix
552  // so that we are in the same coordinate system
553  vtkMatrix4x4 *volMatrix = vtkMatrix4x4::New();
554  vol->GetMatrix(volMatrix);
555  camPos[3] = 1.0;
556  volMatrix->Invert();
557  volMatrix->MultiplyPoint(camPos, camPos);
558  volMatrix->Delete();
559  if (camPos[3])
560  {
561  camPos[0] /= camPos[3];
562  camPos[1] /= camPos[3];
563  camPos[2] /= camPos[3];
564  }
565 
566  // These are the region limits for x (first four), y (next four) and
567  // z (last four). The first region limit is the lower bound for
568  // that axis, the next two are the region planes along that axis, and
569  // the final one in the upper bound for that axis.
570  float limit[12];
571  for (i = 0; i < 3; i++)
572  {
573  limit[i * 4] = volBounds[i * 2];
574  limit[i * 4 + 1] = this->CroppingRegionPlanes[i * 2];
575  limit[i * 4 + 2] = this->CroppingRegionPlanes[i * 2 + 1];
576  limit[i * 4 + 3] = volBounds[i * 2 + 1];
577  }
578 
579  // For each of the 27 possible regions, find out if it is enabled,
580  // and if so, compute the bounds and the distance from the camera
581  // to the center of the region.
582  int numRegions = 0;
583  int region;
584  for (region = 0; region < 27; region++)
585  {
586  int regionFlag = 1 << region;
587 
588  if (this->CroppingRegionFlags & regionFlag)
589  {
590  // what is the coordinate in the 3x3x3 grid
591  int loc[3];
592  loc[0] = region % 3;
593  loc[1] = (region / 3) % 3;
594  loc[2] = (region / 9) % 3;
595 
596  // compute the bounds and center
597  float center[3];
598  for (i = 0; i < 3; i++)
599  {
600  bounds[numRegions][i * 2] = limit[4 * i + loc[i]];
601  bounds[numRegions][i * 2 + 1] = limit[4 * i + loc[i] + 1];
602  center[i] = (bounds[numRegions][i * 2] + bounds[numRegions][i * 2 + 1]) / 2.0;
603  }
604 
605  // compute the distance squared to the center
606  distance2[numRegions] = (camPos[0] - center[0]) * (camPos[0] - center[0]) +
607  (camPos[1] - center[1]) * (camPos[1] - center[1]) +
608  (camPos[2] - center[2]) * (camPos[2] - center[2]);
609 
610  // we've added one region
611  numRegions++;
612  }
613  }
614 
615  // Do a quick bubble sort on distance
616  for (i = 1; i < numRegions; i++)
617  {
618  for (j = i; j > 0 && distance2[j] > distance2[j - 1]; j--)
619  {
620  float tmpBounds[6];
621  float tmpDistance2;
622 
623  for (k = 0; k < 6; k++)
624  {
625  tmpBounds[k] = bounds[j][k];
626  }
627  tmpDistance2 = distance2[j];
628 
629  for (k = 0; k < 6; k++)
630  {
631  bounds[j][k] = bounds[j - 1][k];
632  }
633  distance2[j] = distance2[j - 1];
634 
635  for (k = 0; k < 6; k++)
636  {
637  bounds[j - 1][k] = tmpBounds[k];
638  }
639  distance2[j - 1] = tmpDistance2;
640  }
641  }
642 
643  numIterations = numRegions;
644  }
645 
646  // loop over all regions we need to render
647  for (int loop = 0; loop < numIterations; loop++)
648  {
649  // Compute the set of polygons for this region
650  // according to the bounds
651  this->ComputePolygons(ren, vol, bounds[loop]);
652 
653  // Loop over the polygons
654  for (i = 0; i < this->NumberOfPolygons; i++)
655  {
656  if (renWin->CheckAbortStatus())
657  {
658  return;
659  }
660 
661  float *ptr = this->PolygonBuffer + 36 * i;
662 
663  glBegin(GL_TRIANGLE_FAN);
664 
665  for (j = 0; j < 6; j++)
666  {
667  if (ptr[0] < 0.0)
668  {
669  break;
670  }
671 
672  for (k = 0; k < 4; k++)
673  {
674  if (stages[k])
675  {
676  vtkgl::MultiTexCoord3fv(vtkgl::TEXTURE0 + k, ptr);
677  }
678  }
679  glVertex3fv(ptr + 3);
680 
681  ptr += 6;
682  }
683  glEnd();
684  }
685  }
686 }
687 
688 // This method moves the scalars from the input volume into volume1 (and
689 // possibly volume2) which are the 3D texture maps used for rendering.
690 //
691 // In the case where our volume is a power of two, the copy is done
692 // directly. If we need to resample, then trilinear interpolation is used.
693 //
694 // A shift/scale is applied to the input scalar value to produce an 8 bit
695 // value for the texture volume.
696 //
697 // When the input data is one component, the scalar value is placed in the
698 // second component of the two component volume1. The first component is
699 // filled in later with the gradient magnitude.
700 //
701 // When the input data is two component non-independent, the first component
702 // of the input data is placed in the first component of volume1, and the
703 // second component of the input data is placed in the third component of
704 // volume1. Volume1 has three components - the second is filled in later with
705 // the gradient magnitude.
706 //
707 // When the input data is four component non-independent, the first three
708 // components of the input data are placed in volume1 (which has three
709 // components), and the fourth component is placed in the second component
710 // of volume2. The first component of volume2 is later filled in with the
711 // gradient magnitude.
712 
713 template <class T>
714 class ScalarGradientCompute
715 {
716  T *dataPtr;
717  unsigned char *tmpPtr;
718  unsigned char *tmpPtr2;
719  int sizeX;
720  int sizeY;
721  int sizeZ;
722  int sizeXY;
723  int sizeXm1;
724  int sizeYm1;
725  int sizeZm1;
726  int fullX;
727  int fullY;
728  int fullZ;
729  int fullXY;
730  int currentChunkStart;
731  int currentChunkEnd;
732 
733  int offZ;
734 
735  float offset;
736  float scale;
737 
738 public:
739  ScalarGradientCompute(T *_dataPtr,
740  unsigned char *_tmpPtr,
741  unsigned char *_tmpPtr2,
742  int _sizeX,
743  int _sizeY,
744  int _sizeZ,
745  int _fullX,
746  int _fullY,
747  int _fullZ,
748  float _offset,
749  float _scale)
750  {
751  dataPtr = _dataPtr;
752  tmpPtr = _tmpPtr;
753  tmpPtr2 = _tmpPtr2;
754  sizeX = _sizeX;
755  sizeY = _sizeY;
756  sizeZ = _sizeZ;
757  fullX = _fullX;
758  fullY = _fullY;
759  fullZ = _fullZ;
760  offset = _offset;
761  scale = _scale;
762 
763  sizeXY = sizeX * sizeY;
764  sizeXm1 = sizeX - 1;
765  sizeYm1 = sizeY - 1;
766  sizeZm1 = sizeZ - 1;
767 
768  fullXY = fullX * fullY;
769  currentChunkStart = 0;
770  currentChunkEnd = 0;
771  offZ = 0;
772  }
773 
774  inline float sample(int x, int y, int z) { return float(dataPtr[x + y * sizeX + z * sizeXY]); }
775  inline void fill(int x, int y, int z)
776  {
777  int doff = x + y * fullX + (z - offZ) * fullXY;
778 
779  tmpPtr[doff * 4 + 0] = 0;
780  tmpPtr[doff * 4 + 1] = 0;
781  tmpPtr[doff * 4 + 2] = 0;
782  tmpPtr[doff * 4 + 3] = 0;
783  /*
784 tmpPtr2[doff*3+0]= 0;
785 tmpPtr2[doff*3+1]= 0;
786 tmpPtr2[doff*3+2]= 0;
787 */
788  }
789 
790  inline int clamp(int x)
791  {
792  if (x < 0)
793  x = 0;
794  else if (x > 255)
795  x = 255;
796  return x;
797  }
798 
799  inline void write(int x, int y, int z, float grayValue, float gx, float gy, float gz)
800  {
801  /*
802  gx /= aspect[0];
803  gy /= aspect[1];
804  gz /= aspect[2];
805  */
806  // Compute the gradient magnitude
807 
808  int iGrayValue = static_cast<int>((grayValue + offset) * scale + 0.5f);
809 
810  gx *= scale;
811  gy *= scale;
812  gz *= scale;
813 
814  float t = sqrtf(gx * gx + gy * gy + gz * gz);
815 
816  if (t > 0.01f)
817  {
818  if (t < 2.0f)
819  {
820  float fac = 2.0f / t;
821  gx *= fac;
822  gy *= fac;
823  gz *= fac;
824  }
825  else if (t > 255.0f)
826  {
827  float fac = 255.0f / t;
828  gx *= fac;
829  gy *= fac;
830  gz *= fac;
831  }
832  }
833  else
834  {
835  gx = gy = gz = 0.0f;
836  }
837 
838  int nx = static_cast<int>(0.5f * gx + 127.5f);
839  int ny = static_cast<int>(0.5f * gy + 127.5f);
840  int nz = static_cast<int>(0.5f * gz + 127.5f);
841 
842  int doff = x + y * fullX + (z - offZ) * fullXY;
843 
844  // tmpPtr[doff*2+0]= 0;
845 
846  tmpPtr[doff * 4 + 0] = clamp(nx);
847  tmpPtr[doff * 4 + 1] = clamp(ny);
848  tmpPtr[doff * 4 + 2] = clamp(nz);
849  tmpPtr[doff * 4 + 3] = clamp(iGrayValue);
850 
851  /*
852  if( z == fullZ/2 )
853  if( y == fullY/2 )
854  MITK_INFO << x << " " << y << " " << z << " : " << iGrayValue << " : " << iGradient;
855  */
856  }
857 
858  inline void compute(int x, int y, int z)
859  {
860  float grayValue = sample(x, y, z);
861  float gx, gy, gz;
862 
863  gx = sample(x + 1, y, z) - sample(x - 1, y, z);
864  gy = sample(x, y + 1, z) - sample(x, y - 1, z);
865  gz = sample(x, y, z + 1) - sample(x, y, z - 1);
866 
867  write(x, y, z, grayValue, gx, gy, gz);
868  }
869 
870  inline void computeClamp(int x, int y, int z)
871  {
872  float grayValue = sample(x, y, z);
873  float gx, gy, gz;
874 
875  if (x == 0)
876  gx = 2.0f * (sample(x + 1, y, z) - grayValue);
877  else if (x == sizeXm1)
878  gx = 2.0f * (grayValue - sample(x - 1, y, z));
879  else
880  gx = sample(x + 1, y, z) - sample(x - 1, y, z);
881 
882  if (y == 0)
883  gy = 2.0f * (sample(x, y + 1, z) - grayValue);
884  else if (y == sizeYm1)
885  gy = 2.0f * (grayValue - sample(x, y - 1, z));
886  else
887  gy = sample(x, y + 1, z) - sample(x, y - 1, z);
888 
889  if (z == 0)
890  gz = 2.0f * (sample(x, y, z + 1) - grayValue);
891  else if (z == sizeZm1)
892  gz = 2.0f * (grayValue - sample(x, y, z - 1));
893  else
894  gz = sample(x, y, z + 1) - sample(x, y, z - 1);
895 
896  write(x, y, z, grayValue, gx, gy, gz);
897  }
898 
899  inline void compute1D(int y, int z)
900  {
901  int x;
902 
903  x = 0;
904  computeClamp(x, y, z);
905  x++;
906 
907  while (x < sizeX - 1)
908  {
909  compute(x, y, z);
910  x++;
911  }
912 
913  if (x < sizeX)
914  {
915  computeClamp(x, y, z);
916  x++;
917  }
918 
919  while (x < fullX)
920  {
921  fill(x, y, z);
922  x++;
923  }
924  }
925 
926  inline void fill1D(int y, int z)
927  {
928  int x;
929 
930  x = 0;
931  while (x < fullX)
932  {
933  fill(x, y, z);
934  x++;
935  }
936  }
937 
938  inline void computeClamp1D(int y, int z)
939  {
940  int x;
941 
942  x = 0;
943 
944  while (x < sizeX)
945  {
946  computeClamp(x, y, z);
947  x++;
948  }
949 
950  while (x < fullX)
951  {
952  fill(x, y, z);
953  x++;
954  }
955  }
956 
957  inline void computeClamp2D(int z)
958  {
959  int y;
960 
961  y = 0;
962 
963  while (y < sizeY)
964  {
965  computeClamp1D(y, z);
966  y++;
967  }
968 
969  while (y < fullY)
970  {
971  fill1D(y, z);
972  y++;
973  }
974  }
975 
976  inline void compute2D(int z)
977  {
978  int y;
979 
980  y = 0;
981  computeClamp1D(y, z);
982  y++;
983 
984  while (y < sizeY - 1)
985  {
986  compute1D(y, z);
987  y++;
988  }
989 
990  if (y < sizeY)
991  {
992  computeClamp1D(y, z);
993  y++;
994  }
995 
996  while (y < fullY)
997  {
998  fill1D(y, z);
999  y++;
1000  }
1001  }
1002 
1003  inline void fill2D(int z)
1004  {
1005  int y;
1006 
1007  y = 0;
1008  while (y < fullY)
1009  {
1010  fill1D(y, z);
1011  y++;
1012  }
1013  }
1014 
1015  inline void fillSlices(int currentChunkStart, int currentChunkEnd)
1016  {
1017  offZ = currentChunkStart;
1018 
1019 /*
1020  int num = omp_get_num_procs();
1021  MITK_INFO << "omp uses " << num << " processors";
1022 */
1023 
1024 #pragma omp parallel for
1025  for (int z = currentChunkStart; z <= currentChunkEnd; z++)
1026  {
1027  if (z == 0 || z == sizeZ - 1)
1028  computeClamp2D(z);
1029  else if (z >= sizeZ)
1030  fill2D(z);
1031  else
1032  compute2D(z);
1033  }
1034  }
1035 };
1036 
1037 template <class T>
1039  T *dataPtr, vtkMitkVolumeTextureMapper3D *me, float offset, float scale, GLuint volume1, GLuint /*volume2*/)
1040 {
1041  // T *inPtr;
1042  // unsigned char *outPtr, *outPtr2;
1043  // int i, j, k;
1044  // int idx;
1045 
1046  int inputDimensions[3];
1047  double inputSpacing[3];
1048  vtkImageData *input = me->GetInput();
1049 
1050  input->GetDimensions(inputDimensions);
1051  input->GetSpacing(inputSpacing);
1052 
1053  int outputDimensions[3];
1054  float outputSpacing[3];
1055  me->GetVolumeDimensions(outputDimensions);
1056  me->GetVolumeSpacing(outputSpacing);
1057 
1058  // int components = input->GetNumberOfScalarComponents();
1059 
1060  // double wx, wy, wz;
1061  // double fx, fy, fz;
1062  // int x, y, z;
1063 
1064  // double sampleRate[3];
1065  // sampleRate[0] = outputSpacing[0] / static_cast<double>(inputSpacing[0]);
1066  // sampleRate[1] = outputSpacing[1] / static_cast<double>(inputSpacing[1]);
1067  // sampleRate[2] = outputSpacing[2] / static_cast<double>(inputSpacing[2]);
1068 
1069  int fullX = outputDimensions[0];
1070  int fullY = outputDimensions[1];
1071  int fullZ = outputDimensions[2];
1072 
1073  int sizeX = inputDimensions[0];
1074  int sizeY = inputDimensions[1];
1075  int sizeZ = inputDimensions[2];
1076 
1077  int chunkSize = 64;
1078 
1079  if (fullZ < chunkSize)
1080  chunkSize = fullZ;
1081 
1082  int numChunks = (fullZ + (chunkSize - 1)) / chunkSize;
1083 
1084  // inPtr = dataPtr;
1085 
1086  unsigned char *tmpPtr = new unsigned char[fullX * fullY * chunkSize * 4];
1087  unsigned char *tmpPtr2 = 0; // new unsigned char[fullX*fullY*chunkSize*3];
1088 
1089  // For each Chunk
1090  {
1091  ScalarGradientCompute<T> sgc(dataPtr, tmpPtr, tmpPtr2, sizeX, sizeY, sizeZ, fullX, fullY, fullZ, offset, scale);
1092 
1093  int currentChunk = 0;
1094 
1095  while (currentChunk < numChunks)
1096  {
1097  int currentChunkStart = currentChunk * chunkSize;
1098  int currentChunkEnd = currentChunkStart + chunkSize - 1;
1099 
1100  if (currentChunkEnd > (fullZ - 1))
1101  currentChunkEnd = (fullZ - 1);
1102 
1103  int currentChunkSize = currentChunkEnd - currentChunkStart + 1;
1104 
1105  sgc.fillSlices(currentChunkStart, currentChunkEnd);
1106 
1107  glBindTexture(vtkgl::TEXTURE_3D, volume1);
1108  vtkgl::TexSubImage3D(vtkgl::TEXTURE_3D,
1109  0,
1110  0,
1111  0,
1112  currentChunkStart,
1113  fullX,
1114  fullY,
1115  currentChunkSize,
1116  GL_RGBA,
1117  GL_UNSIGNED_BYTE,
1118  tmpPtr);
1119  currentChunk++;
1120  }
1121  }
1122 
1123  delete[] tmpPtr;
1124 }
1125 
1126 class RGBACompute
1127 {
1128  unsigned char *dataPtr;
1129  unsigned char *tmpPtr;
1130  unsigned char *tmpPtr2;
1131  int sizeX;
1132  int sizeY;
1133  int sizeZ;
1134  int sizeXY;
1135  int sizeXm1;
1136  int sizeYm1;
1137  int sizeZm1;
1138  int fullX;
1139  int fullY;
1140  int fullZ;
1141  int fullXY;
1142  // int currentChunkStart;
1143  // int currentChunkEnd;
1144 
1145  int offZ;
1146 
1147 public:
1148  RGBACompute(unsigned char *_dataPtr,
1149  unsigned char *_tmpPtr,
1150  unsigned char *_tmpPtr2,
1151  int _sizeX,
1152  int _sizeY,
1153  int _sizeZ,
1154  int _fullX,
1155  int _fullY,
1156  int _fullZ)
1157  {
1158  dataPtr = _dataPtr;
1159  tmpPtr = _tmpPtr;
1160  tmpPtr2 = _tmpPtr2;
1161  sizeX = _sizeX;
1162  sizeY = _sizeY;
1163  sizeZ = _sizeZ;
1164  fullX = _fullX;
1165  fullY = _fullY;
1166  fullZ = _fullZ;
1167 
1168  sizeXY = sizeX * sizeY;
1169  sizeXm1 = sizeX - 1;
1170  sizeYm1 = sizeY - 1;
1171  sizeZm1 = sizeZ - 1;
1172 
1173  fullXY = fullX * fullY;
1174  offZ = 0;
1175  }
1176 
1177  inline int sample(int x, int y, int z) { return dataPtr[(x + y * sizeX + z * sizeXY) * 4 + 3]; }
1178  inline void fill(int x, int y, int z)
1179  {
1180  int doff = x + y * fullX + (z - offZ) * fullXY;
1181 
1182  tmpPtr[doff * 4 + 0] = 0;
1183  tmpPtr[doff * 4 + 1] = 0;
1184  tmpPtr[doff * 4 + 2] = 0;
1185  tmpPtr[doff * 4 + 3] = 0;
1186 
1187  tmpPtr2[doff * 3 + 0] = 0;
1188  tmpPtr2[doff * 3 + 1] = 0;
1189  tmpPtr2[doff * 3 + 2] = 0;
1190  }
1191 
1192  inline int clamp(int x)
1193  {
1194  if (x < 0)
1195  x = 0;
1196  else if (x > 255)
1197  x = 255;
1198  return x;
1199  }
1200 
1201  inline void write(int x, int y, int z, int iGrayValue, int gx, int gy, int gz)
1202  {
1203  /*
1204  gx /= aspect[0];
1205  gy /= aspect[1];
1206  gz /= aspect[2];
1207  */
1208  int nx = static_cast<int>(0.5f * gx + 127.5f);
1209  int ny = static_cast<int>(0.5f * gy + 127.5f);
1210  int nz = static_cast<int>(0.5f * gz + 127.5f);
1211 
1212  int doff = x + y * fullX + (z - offZ) * fullXY;
1213 
1214  // tmpPtr[doff*2+0]= 0;
1215 
1216  tmpPtr[doff * 4 + 0] = clamp(nx);
1217  tmpPtr[doff * 4 + 1] = clamp(ny);
1218  tmpPtr[doff * 4 + 2] = clamp(nz);
1219  tmpPtr[doff * 4 + 3] = clamp(iGrayValue);
1220 
1221  int soff = x + y * sizeX + z * sizeXY;
1222 
1223  tmpPtr2[doff * 3 + 0] = dataPtr[soff * 4 + 0];
1224  tmpPtr2[doff * 3 + 1] = dataPtr[soff * 4 + 1];
1225  tmpPtr2[doff * 3 + 2] = dataPtr[soff * 4 + 2];
1226 
1227  /*
1228  if( z == fullZ/2 )
1229  if( y == fullY/2 )
1230  MITK_INFO << x << " " << y << " " << z << " : " << iGrayValue << " : " << iGradient;
1231  */
1232  }
1233 
1234  inline void compute(int x, int y, int z)
1235  {
1236  int grayValue = sample(x, y, z);
1237  int gx, gy, gz;
1238 
1239  gx = sample(x + 1, y, z) - sample(x - 1, y, z);
1240  gy = sample(x, y + 1, z) - sample(x, y - 1, z);
1241  gz = sample(x, y, z + 1) - sample(x, y, z - 1);
1242 
1243  write(x, y, z, grayValue, gx, gy, gz);
1244  }
1245 
1246  inline void computeClamp(int x, int y, int z)
1247  {
1248  int grayValue = sample(x, y, z);
1249  int gx, gy, gz;
1250 
1251  if (x == 0)
1252  gx = 2 * (sample(x + 1, y, z) - grayValue);
1253  else if (x == sizeXm1)
1254  gx = 2 * (grayValue - sample(x - 1, y, z));
1255  else
1256  gx = sample(x + 1, y, z) - sample(x - 1, y, z);
1257 
1258  if (y == 0)
1259  gy = 2 * (sample(x, y + 1, z) - grayValue);
1260  else if (y == sizeYm1)
1261  gy = 2 * (grayValue - sample(x, y - 1, z));
1262  else
1263  gy = sample(x, y + 1, z) - sample(x, y - 1, z);
1264 
1265  if (z == 0)
1266  gz = 2 * (sample(x, y, z + 1) - grayValue);
1267  else if (z == sizeZm1)
1268  gz = 2 * (grayValue - sample(x, y, z - 1));
1269  else
1270  gz = sample(x, y, z + 1) - sample(x, y, z - 1);
1271 
1272  write(x, y, z, grayValue, gx, gy, gz);
1273  }
1274 
1275  inline void compute1D(int y, int z)
1276  {
1277  int x = 0;
1278 
1279  computeClamp(x, y, z);
1280  x++;
1281 
1282  while (x < sizeX - 1)
1283  {
1284  compute(x, y, z);
1285  x++;
1286  }
1287 
1288  if (x < sizeX)
1289  {
1290  computeClamp(x, y, z);
1291  x++;
1292  }
1293 
1294  while (x < fullX)
1295  {
1296  fill(x, y, z);
1297  x++;
1298  }
1299  }
1300 
1301  inline void fill1D(int y, int z)
1302  {
1303  int x = 0;
1304 
1305  while (x < fullX)
1306  {
1307  fill(x, y, z);
1308  x++;
1309  }
1310  }
1311 
1312  inline void computeClamp1D(int y, int z)
1313  {
1314  int x = 0;
1315 
1316  while (x < sizeX)
1317  {
1318  computeClamp(x, y, z);
1319  x++;
1320  }
1321 
1322  while (x < fullX)
1323  {
1324  fill(x, y, z);
1325  x++;
1326  }
1327  }
1328 
1329  inline void computeClamp2D(int z)
1330  {
1331  int y = 0;
1332 
1333  while (y < sizeY)
1334  {
1335  computeClamp1D(y, z);
1336  y++;
1337  }
1338 
1339  while (y < fullY)
1340  {
1341  fill1D(y, z);
1342  y++;
1343  }
1344  }
1345 
1346  inline void compute2D(int z)
1347  {
1348  int y = 0;
1349 
1350  computeClamp1D(y, z);
1351  y++;
1352 
1353  while (y < sizeY - 1)
1354  {
1355  compute1D(y, z);
1356  y++;
1357  }
1358 
1359  if (y < sizeY)
1360  {
1361  computeClamp1D(y, z);
1362  y++;
1363  }
1364 
1365  while (y < fullY)
1366  {
1367  fill1D(y, z);
1368  y++;
1369  }
1370  }
1371 
1372  inline void fill2D(int z)
1373  {
1374  int y = 0;
1375 
1376  while (y < fullY)
1377  {
1378  fill1D(y, z);
1379  y++;
1380  }
1381  }
1382 
1383  inline void fillSlices(int currentChunkStart, int currentChunkEnd)
1384  {
1385  offZ = currentChunkStart;
1386 
1387 #pragma omp parallel for
1388  for (int z = currentChunkStart; z <= currentChunkEnd; z++)
1389  {
1390  if (z == 0 || z == sizeZ - 1)
1391  computeClamp2D(z);
1392  else if (z >= sizeZ)
1393  fill2D(z);
1394  else
1395  compute2D(z);
1396  }
1397  }
1398 };
1399 
1400 void vtkVolumeTextureMapper3DComputeRGBA(unsigned char *dataPtr,
1401  vtkMitkVolumeTextureMapper3D *me,
1402  GLuint volume1,
1403  GLuint volume2)
1404 {
1405  // unsigned char *inPtr;
1406  // unsigned char *outPtr, *outPtr2;
1407  // int i, j, k;
1408  // int idx;
1409 
1410  int inputDimensions[3];
1411  double inputSpacing[3];
1412  vtkImageData *input = me->GetInput();
1413 
1414  input->GetDimensions(inputDimensions);
1415  input->GetSpacing(inputSpacing);
1416 
1417  int outputDimensions[3];
1418  float outputSpacing[3];
1419  me->GetVolumeDimensions(outputDimensions);
1420  me->GetVolumeSpacing(outputSpacing);
1421 
1422  int components = input->GetNumberOfScalarComponents();
1423 
1424  MITK_INFO << "components are " << components;
1425 
1426  // double wx, wy, wz;
1427  // double fx, fy, fz;
1428  // int x, y, z;
1429 
1430  // double sampleRate[3];
1431  // sampleRate[0] = outputSpacing[0] / static_cast<double>(inputSpacing[0]);
1432  // sampleRate[1] = outputSpacing[1] / static_cast<double>(inputSpacing[1]);
1433  // sampleRate[2] = outputSpacing[2] / static_cast<double>(inputSpacing[2]);
1434 
1435  int fullX = outputDimensions[0];
1436  int fullY = outputDimensions[1];
1437  int fullZ = outputDimensions[2];
1438 
1439  int sizeX = inputDimensions[0];
1440  int sizeY = inputDimensions[1];
1441  int sizeZ = inputDimensions[2];
1442 
1443  int chunkSize = 64;
1444 
1445  if (fullZ < chunkSize)
1446  chunkSize = fullZ;
1447 
1448  int numChunks = (fullZ + (chunkSize - 1)) / chunkSize;
1449 
1450  // inPtr = dataPtr;
1451 
1452  unsigned char *tmpPtr = new unsigned char[fullX * fullY * chunkSize * 4];
1453  unsigned char *tmpPtr2 = new unsigned char[fullX * fullY * chunkSize * 3];
1454 
1455  // For each Chunk
1456  {
1457  RGBACompute sgc(dataPtr, tmpPtr, tmpPtr2, sizeX, sizeY, sizeZ, fullX, fullY, fullZ);
1458 
1459  int currentChunk = 0;
1460 
1461  while (currentChunk < numChunks)
1462  {
1463  // MITK_INFO << "processing chunk " << currentChunk;
1464 
1465  int currentChunkStart = currentChunk * chunkSize;
1466  int currentChunkEnd = currentChunkStart + chunkSize - 1;
1467 
1468  if (currentChunkEnd > (fullZ - 1))
1469  currentChunkEnd = (fullZ - 1);
1470 
1471  int currentChunkSize = currentChunkEnd - currentChunkStart + 1;
1472 
1473  sgc.fillSlices(currentChunkStart, currentChunkEnd);
1474 
1475  glBindTexture(vtkgl::TEXTURE_3D, volume1);
1476  vtkgl::TexSubImage3D(vtkgl::TEXTURE_3D,
1477  0,
1478  0,
1479  0,
1480  currentChunkStart,
1481  fullX,
1482  fullY,
1483  currentChunkSize,
1484  GL_RGBA,
1485  GL_UNSIGNED_BYTE,
1486  tmpPtr);
1487 
1488  glBindTexture(vtkgl::TEXTURE_3D, volume2);
1489  vtkgl::TexSubImage3D(vtkgl::TEXTURE_3D,
1490  0,
1491  0,
1492  0,
1493  currentChunkStart,
1494  fullX,
1495  fullY,
1496  currentChunkSize,
1497  GL_RGB,
1498  GL_UNSIGNED_BYTE,
1499  tmpPtr2);
1500 
1501  currentChunk++;
1502  }
1503  }
1504 
1505  delete[] tmpPtr;
1506  delete[] tmpPtr2;
1507 }
1508 
1509 //-----------------------------------------------------------------------------
1510 void vtkMitkOpenGLVolumeTextureMapper3D::ComputeVolumeDimensions()
1511 {
1512  // Get the image data
1513  vtkImageData *input = this->GetInput();
1514 
1515  // How big does the Volume need to be?
1516  int dim[3];
1517  input->GetDimensions(dim);
1518 
1519  int powerOfTwoDim[3];
1520 
1521  if (this->SupportsNonPowerOfTwoTextures)
1522  {
1523  for (int i = 0; i < 3; i++)
1524  powerOfTwoDim[i] = (dim[i] + 1) & ~1;
1525 
1526  // MITK_INFO << "using non-power-two even textures (" <<
1527  // (1.0-double(dim[0]*dim[1]*dim[2])/double(powerOfTwoDim[0]*powerOfTwoDim[1]*powerOfTwoDim[2])) * 100.0 << "%
1528  // memory wasted)";
1529  }
1530  else
1531  {
1532  for (int i = 0; i < 3; i++)
1533  {
1534  powerOfTwoDim[i] = 4;
1535  while (powerOfTwoDim[i] < dim[i])
1536  powerOfTwoDim[i] *= 2;
1537  }
1538 
1539  MITK_WARN << "using power-two textures ("
1540  << (1.0 -
1541  double(dim[0] * dim[1] * dim[2]) / double(powerOfTwoDim[0] * powerOfTwoDim[1] * powerOfTwoDim[2])) *
1542  100.0
1543  << "% memory wasted)";
1544  }
1545 
1546  // Save the volume size
1547  this->VolumeDimensions[0] = powerOfTwoDim[0];
1548  this->VolumeDimensions[1] = powerOfTwoDim[1];
1549  this->VolumeDimensions[2] = powerOfTwoDim[2];
1550 
1551  // What is the spacing?
1552  double spacing[3];
1553  input->GetSpacing(spacing);
1554 
1555  // Compute the new spacing
1556  this->VolumeSpacing[0] = (dim[0] - 1.01) * spacing[0] / static_cast<double>(this->VolumeDimensions[0] - 1);
1557  this->VolumeSpacing[1] = (dim[1] - 1.01) * spacing[1] / static_cast<double>(this->VolumeDimensions[1] - 1);
1558  this->VolumeSpacing[2] = ((dim[2]) - 1.01) * spacing[2] / static_cast<double>(this->VolumeDimensions[2] - 1);
1559 }
1560 
1561 //-----------------------------------------------------------------------------
1562 bool vtkMitkOpenGLVolumeTextureMapper3D::UpdateVolumes(vtkVolume *vtkNotUsed(vol))
1563 {
1564  // Get the image data
1565  vtkImageData *input = this->GetInput();
1566  // input->Update(); //VTK6_TODO
1567 
1568  bool needUpdate = false;
1569 
1570  // Has the volume changed in some way?
1571  if (this->SavedTextureInput != input || this->SavedTextureMTime.GetMTime() < input->GetMTime())
1572  needUpdate = true;
1573 
1574  // Do we have any volume on the gpu already?
1575  if (!this->Volume1Index)
1576  needUpdate = true;
1577 
1578  if (!needUpdate)
1579  return true;
1580 
1581  ComputeVolumeDimensions();
1582 
1583  int components = input->GetNumberOfScalarComponents();
1584 
1585  // Find the scalar range
1586  double scalarRange[2];
1587  input->GetPointData()->GetScalars()->GetRange(scalarRange, components - 1);
1588 
1589  // Is the difference between max and min less than 4096? If so, and if
1590  // the data is not of float or double type, use a simple offset mapping.
1591  // If the difference between max and min is 4096 or greater, or the data
1592  // is of type float or double, we must use an offset / scaling mapping.
1593  // In this case, the array size will be 4096 - we need to figure out the
1594  // offset and scale factor.
1595  float offset;
1596  float scale;
1597 
1598  int arraySizeNeeded;
1599 
1600  int scalarType = input->GetScalarType();
1601 
1602  if (scalarType == VTK_FLOAT || scalarType == VTK_DOUBLE || scalarRange[1] - scalarRange[0] > 255)
1603  {
1604  arraySizeNeeded = 256;
1605  offset = -scalarRange[0];
1606  scale = 255.0 / (scalarRange[1] - scalarRange[0]);
1607  }
1608  else
1609  {
1610  arraySizeNeeded = static_cast<int>(scalarRange[1] - scalarRange[0] + 1);
1611  offset = -scalarRange[0];
1612  scale = 1.0;
1613  }
1614 
1615  this->ColorTableSize = arraySizeNeeded;
1616  this->ColorTableOffset = offset;
1617  this->ColorTableScale = scale;
1618 
1619  // Allocating volume on gpu
1620  {
1621  // Deleting old textures
1622  this->DeleteTextureIndex(&this->Volume1Index);
1623  this->DeleteTextureIndex(&this->Volume2Index);
1624  this->DeleteTextureIndex(&this->Volume3Index);
1625 
1626  this->CreateTextureIndex(&this->Volume1Index);
1627  // this->CreateTextureIndex(&this->Volume2Index);
1628 
1629  int dim[3];
1630  this->GetVolumeDimensions(dim);
1631 
1632  vtkgl::ActiveTexture(vtkgl::TEXTURE0);
1633 
1634  MITK_INFO << "allocating volume on gpu";
1635 
1636  GLint gradientScalarTextureFormat = GL_RGBA8;
1637 
1638  if (this->UseCompressedTexture && SupportsCompressedTexture)
1639  gradientScalarTextureFormat = myGL_COMPRESSED_RGBA_S3TC_DXT5_EXT;
1640 
1641  glBindTexture(vtkgl::TEXTURE_3D, this->Volume1Index);
1642  vtkgl::TexImage3D(
1643  vtkgl::TEXTURE_3D, 0, gradientScalarTextureFormat, dim[0], dim[1], dim[2], 0, GL_RGBA, GL_UNSIGNED_BYTE, 0);
1644  this->Setup3DTextureParameters(true);
1645  }
1646 
1647  // Transfer the input volume to the RGBA volume
1648  void *dataPtr = input->GetScalarPointer();
1649 
1650  switch (scalarType)
1651  {
1652  vtkTemplateMacro(vtkVolumeTextureMapper3DComputeScalars(
1653  static_cast<VTK_TT *>(dataPtr), this, offset, scale, this->Volume1Index, this->Volume2Index));
1654  }
1655 
1656  this->SavedTextureInput = input;
1657  this->SavedTextureMTime.Modified();
1658 
1659  return true;
1660 }
1661 
1662 //-----------------------------------------------------------------------------
1663 bool vtkMitkOpenGLVolumeTextureMapper3D::UpdateVolumesRGBA(vtkVolume *vtkNotUsed(vol))
1664 {
1665  // Get the image data
1666  vtkImageData *input = this->GetInput();
1667  // input->Update(); //VTK6_TODO
1668 
1669  bool needUpdate = false;
1670 
1671  // Has the volume changed in some way?
1672  if (this->SavedTextureInput != input || this->SavedTextureMTime.GetMTime() < input->GetMTime())
1673  needUpdate = true;
1674 
1675  // Do we have any volume on the gpu already?
1676  if (!this->Volume1Index)
1677  needUpdate = true;
1678 
1679  if (!needUpdate)
1680  return true;
1681 
1682  MITK_INFO << "updating rgba volume";
1683 
1684  ComputeVolumeDimensions();
1685 
1686  // Allocating volume on gpu
1687  {
1688  // Deleting old textures
1689  this->DeleteTextureIndex(&this->Volume1Index);
1690  this->DeleteTextureIndex(&this->Volume2Index);
1691  this->DeleteTextureIndex(&this->Volume3Index);
1692 
1693  this->CreateTextureIndex(&this->Volume1Index);
1694  this->CreateTextureIndex(&this->Volume2Index);
1695 
1696  int dim[3];
1697  this->GetVolumeDimensions(dim);
1698 
1699  MITK_INFO << "allocating volume on gpu";
1700 
1701  GLint gradientScalarTextureFormat = GL_RGBA8;
1702  GLint colorTextureFormat = GL_RGB8;
1703 
1704  if (this->UseCompressedTexture && SupportsCompressedTexture)
1705  {
1706  gradientScalarTextureFormat = myGL_COMPRESSED_RGBA_S3TC_DXT5_EXT;
1707  colorTextureFormat = myGL_COMPRESSED_RGB_S3TC_DXT1_EXT;
1708  }
1709 
1710  vtkgl::ActiveTexture(vtkgl::TEXTURE0);
1711  glBindTexture(vtkgl::TEXTURE_3D, this->Volume1Index);
1712  vtkgl::TexImage3D(
1713  vtkgl::TEXTURE_3D, 0, gradientScalarTextureFormat, dim[0], dim[1], dim[2], 0, GL_RGBA, GL_UNSIGNED_BYTE, 0);
1714  this->Setup3DTextureParameters(true);
1715 
1716  glBindTexture(vtkgl::TEXTURE_3D, this->Volume2Index);
1717  vtkgl::TexImage3D(vtkgl::TEXTURE_3D, 0, colorTextureFormat, dim[0], dim[1], dim[2], 0, GL_RGB, GL_UNSIGNED_BYTE, 0);
1718  this->Setup3DTextureParameters(true);
1719  }
1720 
1721  // Transfer the input volume to the RGBA volume
1722  unsigned char *dataPtr = (unsigned char *)input->GetScalarPointer();
1723  vtkVolumeTextureMapper3DComputeRGBA(dataPtr, this, this->Volume1Index, this->Volume2Index);
1724 
1725  this->SavedTextureInput = input;
1726  this->SavedTextureMTime.Modified();
1727 
1728  return true;
1729 }
1730 
1731 void vtkMitkOpenGLVolumeTextureMapper3D::Setup3DTextureParameters(bool linear)
1732 {
1733  // GPU_INFO << "Setup3DTextureParameters";
1734 
1735  if (linear)
1736  {
1737  glTexParameterf(vtkgl::TEXTURE_3D, GL_TEXTURE_MIN_FILTER, GL_LINEAR);
1738  glTexParameterf(vtkgl::TEXTURE_3D, GL_TEXTURE_MAG_FILTER, GL_LINEAR);
1739  }
1740  else
1741  {
1742  glTexParameterf(vtkgl::TEXTURE_3D, GL_TEXTURE_MIN_FILTER, GL_NEAREST);
1743  glTexParameterf(vtkgl::TEXTURE_3D, GL_TEXTURE_MAG_FILTER, GL_NEAREST);
1744  }
1745 
1746  glTexParameterf(vtkgl::TEXTURE_3D, GL_TEXTURE_WRAP_S, GL_CLAMP);
1747  glTexParameterf(vtkgl::TEXTURE_3D, GL_TEXTURE_WRAP_T, GL_CLAMP);
1748 }
1749 
1750 void vtkMitkOpenGLVolumeTextureMapper3D::SetupOneIndependentTextures(vtkRenderer *vtkNotUsed(ren), vtkVolume *vol)
1751 {
1752  // Update the volume containing the 2 byte scalar / gradient magnitude
1753  this->UpdateVolumes(vol);
1754 
1755  // Update the dependent 2D color table mapping scalar value and
1756  // gradient magnitude to RGBA
1757  if (this->UpdateColorLookup(vol) || !this->ColorLookupIndex)
1758  {
1759  this->DeleteTextureIndex(&this->ColorLookupIndex);
1760  this->DeleteTextureIndex(&this->AlphaLookupIndex);
1761 
1762  this->CreateTextureIndex(&this->ColorLookupIndex);
1763 
1764  vtkgl::ActiveTexture(vtkgl::TEXTURE1);
1765  glBindTexture(GL_TEXTURE_2D, this->ColorLookupIndex);
1766 
1767  glTexParameterf(GL_TEXTURE_2D, GL_TEXTURE_MIN_FILTER, GL_NEAREST);
1768  glTexParameterf(GL_TEXTURE_2D, GL_TEXTURE_MAG_FILTER, GL_NEAREST);
1769  glTexParameterf(GL_TEXTURE_2D, GL_TEXTURE_WRAP_S, GL_CLAMP);
1770  glTexParameterf(GL_TEXTURE_2D, GL_TEXTURE_WRAP_T, GL_CLAMP);
1771 
1772  // MITK_INFO << "uploading transferfunction";
1773 
1774  GLint colorLookupTextureFormat = GL_RGBA8;
1775 
1776  if (this->UseCompressedTexture && SupportsCompressedTexture)
1777  colorLookupTextureFormat = myGL_COMPRESSED_RGBA_S3TC_DXT5_EXT;
1778 
1779  glTexImage2D(GL_TEXTURE_2D, 0, colorLookupTextureFormat, 256, 256, 0, GL_RGBA, GL_UNSIGNED_BYTE, this->ColorLookup);
1780  }
1781 }
1782 
1783 void vtkMitkOpenGLVolumeTextureMapper3D::SetupRGBATextures(vtkRenderer *vtkNotUsed(ren), vtkVolume *vol)
1784 {
1785  MITK_INFO << "SetupFourDependentTextures";
1786 
1787  this->UpdateVolumesRGBA(vol);
1788 
1789  /*
1790 vtkgl::ActiveTexture( vtkgl::TEXTURE0 );
1791 glDisable( GL_TEXTURE_2D );
1792 glEnable( vtkgl::TEXTURE_3D );
1793 
1794 vtkgl::ActiveTexture( vtkgl::TEXTURE1 );
1795 glDisable( GL_TEXTURE_2D );
1796 glEnable( vtkgl::TEXTURE_3D );
1797 
1798 vtkgl::ActiveTexture( vtkgl::TEXTURE2 );
1799 glDisable( GL_TEXTURE_2D );
1800 glEnable( vtkgl::TEXTURE_3D );
1801 
1802 // Update the volume containing the 3 byte scalars / gradient magnitude
1803 if ( this->UpdateVolumes( vol ) || !this->Volume1Index ||
1804  !this->Volume2Index || !this->Volume3Index )
1805  {
1806  int dim[3];
1807  this->GetVolumeDimensions(dim);
1808 
1809  vtkgl::ActiveTexture( vtkgl::TEXTURE0 );
1810  glBindTexture(vtkgl::TEXTURE_3D,0);
1811  this->DeleteTextureIndex(&this->Volume1Index);
1812  this->CreateTextureIndex(&this->Volume1Index);
1813  glBindTexture(vtkgl::TEXTURE_3D, this->Volume1Index);
1814  vtkgl::TexImage3D(vtkgl::TEXTURE_3D,0,this->InternalRGB,dim[0],dim[1],
1815  dim[2],0,GL_RGB,GL_UNSIGNED_BYTE,this->Volume1);
1816 
1817  vtkgl::ActiveTexture( vtkgl::TEXTURE1 );
1818  glBindTexture(vtkgl::TEXTURE_3D,0);
1819  this->DeleteTextureIndex(&this->Volume2Index);
1820  this->CreateTextureIndex(&this->Volume2Index);
1821  glBindTexture(vtkgl::TEXTURE_3D, this->Volume2Index);
1822  vtkgl::TexImage3D(vtkgl::TEXTURE_3D,0,this->InternalLA,dim[0],dim[1],
1823  dim[2],0,GL_LUMINANCE_ALPHA,GL_UNSIGNED_BYTE,
1824  this->Volume2);
1825 
1826  vtkgl::ActiveTexture( vtkgl::TEXTURE2 );
1827  glBindTexture(vtkgl::TEXTURE_3D,0);
1828  this->DeleteTextureIndex(&this->Volume3Index);
1829  this->CreateTextureIndex(&this->Volume3Index);
1830  glBindTexture(vtkgl::TEXTURE_3D, this->Volume3Index);
1831  vtkgl::TexImage3D(vtkgl::TEXTURE_3D,0,this->InternalRGB,dim[0],dim[1],
1832  dim[2],0,GL_RGB,GL_UNSIGNED_BYTE,this->Volume3);
1833  }
1834 
1835 vtkgl::ActiveTexture( vtkgl::TEXTURE0 );
1836 glBindTexture(vtkgl::TEXTURE_3D, this->Volume1Index);
1837 this->Setup3DTextureParameters( true );
1838 
1839 vtkgl::ActiveTexture( vtkgl::TEXTURE1 );
1840 glBindTexture(vtkgl::TEXTURE_3D, this->Volume2Index);
1841 this->Setup3DTextureParameters( true );
1842 
1843 vtkgl::ActiveTexture( vtkgl::TEXTURE2 );
1844 glBindTexture(vtkgl::TEXTURE_3D_EXT, this->Volume3Index);
1845 this->Setup3DTextureParameters( true );
1846 
1847 vtkgl::ActiveTexture( vtkgl::TEXTURE3 );
1848 glEnable( GL_TEXTURE_2D );
1849 glDisable( vtkgl::TEXTURE_3D );
1850 
1851 // Update the dependent 2D table mapping scalar value and
1852 // gradient magnitude to opacity
1853 if ( this->UpdateColorLookup( vol ) || !this->AlphaLookupIndex )
1854  {
1855  this->DeleteTextureIndex(&this->ColorLookupIndex);
1856 
1857  vtkgl::ActiveTexture( vtkgl::TEXTURE3 );
1858  glBindTexture(GL_TEXTURE_2D,0);
1859  this->DeleteTextureIndex(&this->AlphaLookupIndex);
1860  this->CreateTextureIndex(&this->AlphaLookupIndex);
1861  glBindTexture(GL_TEXTURE_2D, this->AlphaLookupIndex);
1862 
1863  glTexParameterf( GL_TEXTURE_2D, GL_TEXTURE_MIN_FILTER, GL_NEAREST);
1864  glTexParameterf( GL_TEXTURE_2D, GL_TEXTURE_MAG_FILTER, GL_NEAREST );
1865  glTexParameterf( GL_TEXTURE_2D, GL_TEXTURE_WRAP_S, GL_CLAMP );
1866  glTexParameterf( GL_TEXTURE_2D, GL_TEXTURE_WRAP_T, GL_CLAMP );
1867 
1868  //MITK_INFO << "uploading transferfunction";
1869 
1870  glTexImage2D(GL_TEXTURE_2D,0,this->InternalAlpha, 256, 256, 0,
1871  GL_ALPHA, GL_UNSIGNED_BYTE, this->AlphaLookup );
1872  }
1873 
1874 vtkgl::ActiveTexture( vtkgl::TEXTURE3 );
1875 glBindTexture(GL_TEXTURE_2D, this->AlphaLookupIndex);
1876 
1877 */
1878 }
1879 
1880 void vtkMitkOpenGLVolumeTextureMapper3D::RenderOneIndependentShadeFP(vtkRenderer *ren, vtkVolume *vol)
1881 {
1882  // GPU_INFO << "RenderOneIndependentShadeFP";
1883 
1884  this->SetupOneIndependentTextures(ren, vol);
1885 
1886  glEnable(vtkgl::FRAGMENT_PROGRAM_ARB);
1887 
1888  vtkgl::BindProgramARB(vtkgl::FRAGMENT_PROGRAM_ARB, prgOneComponentShade);
1889 
1890  this->SetupProgramLocalsForShadingFP(ren, vol);
1891 
1892  // Bind Textures
1893  {
1894  vtkgl::ActiveTexture(vtkgl::TEXTURE0);
1895  glDisable(GL_TEXTURE_2D);
1896  glEnable(vtkgl::TEXTURE_3D);
1897  glBindTexture(vtkgl::TEXTURE_3D, this->Volume1Index);
1898 
1899  vtkgl::ActiveTexture(vtkgl::TEXTURE1);
1900  glEnable(GL_TEXTURE_2D);
1901  glDisable(vtkgl::TEXTURE_3D);
1902  glBindTexture(GL_TEXTURE_2D, this->ColorLookupIndex);
1903 
1904  vtkgl::ActiveTexture(vtkgl::TEXTURE2);
1905  glDisable(GL_TEXTURE_2D);
1906  glEnable(vtkgl::TEXTURE_3D);
1907  glBindTexture(vtkgl::TEXTURE_3D, this->Volume2Index);
1908  }
1909 
1910  int stages[4] = {1, 1, 1, 0};
1911  this->RenderPolygons(ren, vol, stages);
1912 
1913  glDisable(vtkgl::FRAGMENT_PROGRAM_ARB);
1914 }
1915 
1916 void vtkMitkOpenGLVolumeTextureMapper3D::RenderRGBAShadeFP(vtkRenderer *ren, vtkVolume *vol)
1917 {
1918  this->SetupRGBATextures(ren, vol);
1919 
1920  glEnable(vtkgl::FRAGMENT_PROGRAM_ARB);
1921 
1922  vtkgl::BindProgramARB(vtkgl::FRAGMENT_PROGRAM_ARB, prgRGBAShade);
1923 
1924  this->SetupProgramLocalsForShadingFP(ren, vol);
1925 
1926  // Bind Textures
1927  {
1928  vtkgl::ActiveTexture(vtkgl::TEXTURE0);
1929  glDisable(GL_TEXTURE_2D);
1930  glEnable(vtkgl::TEXTURE_3D);
1931  glBindTexture(vtkgl::TEXTURE_3D, this->Volume1Index);
1932 
1933  vtkgl::ActiveTexture(vtkgl::TEXTURE1);
1934  glDisable(GL_TEXTURE_2D);
1935  glEnable(vtkgl::TEXTURE_3D);
1936  glBindTexture(vtkgl::TEXTURE_3D, this->Volume2Index);
1937  }
1938 
1939  int stages[4] = {1, 1, 1, 0};
1940  this->RenderPolygons(ren, vol, stages);
1941 
1942  glDisable(vtkgl::FRAGMENT_PROGRAM_ARB);
1943 }
1944 
1945 void vtkMitkOpenGLVolumeTextureMapper3D::GetLightInformation(vtkRenderer *ren,
1946  vtkVolume *vol,
1947  GLfloat lightDirection[2][4],
1948  GLfloat lightDiffuseColor[2][4],
1949  GLfloat lightSpecularColor[2][4],
1950  GLfloat halfwayVector[2][4],
1951  GLfloat ambientColor[4])
1952 {
1953  // GPU_INFO << "GetLightInformation";
1954 
1955  float ambient = vol->GetProperty()->GetAmbient();
1956  float diffuse = vol->GetProperty()->GetDiffuse();
1957  float specular = vol->GetProperty()->GetSpecular();
1958 
1959  vtkTransform *volumeTransform = vtkTransform::New();
1960 
1961  volumeTransform->SetMatrix(vol->GetMatrix());
1962  volumeTransform->Inverse();
1963 
1964  vtkLightCollection *lights = ren->GetLights();
1965  lights->InitTraversal();
1966 
1967  vtkLight *light[2];
1968  light[0] = lights->GetNextItem();
1969  light[1] = lights->GetNextItem();
1970 
1971  int lightIndex = 0;
1972 
1973  double cameraPosition[3];
1974  double cameraFocalPoint[3];
1975 
1976  ren->GetActiveCamera()->GetPosition(cameraPosition);
1977  ren->GetActiveCamera()->GetFocalPoint(cameraFocalPoint);
1978 
1979  double viewDirection[3];
1980 
1981  volumeTransform->TransformPoint(cameraPosition, cameraPosition);
1982  volumeTransform->TransformPoint(cameraFocalPoint, cameraFocalPoint);
1983 
1984  viewDirection[0] = cameraFocalPoint[0] - cameraPosition[0];
1985  viewDirection[1] = cameraFocalPoint[1] - cameraPosition[1];
1986  viewDirection[2] = cameraFocalPoint[2] - cameraPosition[2];
1987 
1988  vtkMath::Normalize(viewDirection);
1989 
1990  ambientColor[0] = 0.0;
1991  ambientColor[1] = 0.0;
1992  ambientColor[2] = 0.0;
1993  ambientColor[3] = 0.0;
1994 
1995  for (lightIndex = 0; lightIndex < 2; lightIndex++)
1996  {
1997  float dir[3] = {0, 0, 0};
1998  float half[3] = {0, 0, 0};
1999 
2000  if (light[lightIndex] == nullptr || light[lightIndex]->GetSwitch() == 0)
2001  {
2002  lightDiffuseColor[lightIndex][0] = 0.0;
2003  lightDiffuseColor[lightIndex][1] = 0.0;
2004  lightDiffuseColor[lightIndex][2] = 0.0;
2005  lightDiffuseColor[lightIndex][3] = 0.0;
2006 
2007  lightSpecularColor[lightIndex][0] = 0.0;
2008  lightSpecularColor[lightIndex][1] = 0.0;
2009  lightSpecularColor[lightIndex][2] = 0.0;
2010  lightSpecularColor[lightIndex][3] = 0.0;
2011  }
2012  else
2013  {
2014  float lightIntensity = light[lightIndex]->GetIntensity();
2015  double lightColor[3];
2016 
2017  light[lightIndex]->GetDiffuseColor(lightColor);
2018 
2019  double lightPosition[3];
2020  double lightFocalPoint[3];
2021  light[lightIndex]->GetTransformedPosition(lightPosition);
2022  light[lightIndex]->GetTransformedFocalPoint(lightFocalPoint);
2023 
2024  volumeTransform->TransformPoint(lightPosition, lightPosition);
2025  volumeTransform->TransformPoint(lightFocalPoint, lightFocalPoint);
2026 
2027  dir[0] = lightPosition[0] - lightFocalPoint[0];
2028  dir[1] = lightPosition[1] - lightFocalPoint[1];
2029  dir[2] = lightPosition[2] - lightFocalPoint[2];
2030 
2031  vtkMath::Normalize(dir);
2032 
2033  lightDiffuseColor[lightIndex][0] = lightColor[0] * diffuse * lightIntensity;
2034  lightDiffuseColor[lightIndex][1] = lightColor[1] * diffuse * lightIntensity;
2035  lightDiffuseColor[lightIndex][2] = lightColor[2] * diffuse * lightIntensity;
2036  lightDiffuseColor[lightIndex][3] = 1.0;
2037 
2038  lightSpecularColor[lightIndex][0] = lightColor[0] * specular * lightIntensity;
2039  lightSpecularColor[lightIndex][1] = lightColor[1] * specular * lightIntensity;
2040  lightSpecularColor[lightIndex][2] = lightColor[2] * specular * lightIntensity;
2041  lightSpecularColor[lightIndex][3] = 0.0;
2042 
2043  half[0] = dir[0] - viewDirection[0];
2044  half[1] = dir[1] - viewDirection[1];
2045  half[2] = dir[2] - viewDirection[2];
2046 
2047  vtkMath::Normalize(half);
2048 
2049  ambientColor[0] += ambient * lightColor[0];
2050  ambientColor[1] += ambient * lightColor[1];
2051  ambientColor[2] += ambient * lightColor[2];
2052  }
2053 
2054  lightDirection[lightIndex][0] = (dir[0] + 1.0) / 2.0;
2055  lightDirection[lightIndex][1] = (dir[1] + 1.0) / 2.0;
2056  lightDirection[lightIndex][2] = (dir[2] + 1.0) / 2.0;
2057  lightDirection[lightIndex][3] = 0.0;
2058 
2059  halfwayVector[lightIndex][0] = (half[0] + 1.0) / 2.0;
2060  halfwayVector[lightIndex][1] = (half[1] + 1.0) / 2.0;
2061  halfwayVector[lightIndex][2] = (half[2] + 1.0) / 2.0;
2062  halfwayVector[lightIndex][3] = 0.0;
2063  }
2064 
2065  volumeTransform->Delete();
2066 }
2067 
2068 void vtkMitkOpenGLVolumeTextureMapper3D::SetupProgramLocalsForShadingFP(vtkRenderer *ren, vtkVolume *vol)
2069 {
2070  // GPU_INFO << "SetupProgramLocalsForShadingFP";
2071 
2072  GLfloat lightDirection[2][4];
2073  GLfloat lightDiffuseColor[2][4];
2074  GLfloat lightSpecularColor[2][4];
2075  GLfloat halfwayVector[2][4];
2076  GLfloat ambientColor[4];
2077 
2078  float ambient = vol->GetProperty()->GetAmbient();
2079  float diffuse = vol->GetProperty()->GetDiffuse();
2080  float specular = vol->GetProperty()->GetSpecular();
2081  float specularPower = vol->GetProperty()->GetSpecularPower();
2082 
2083  vtkTransform *volumeTransform = vtkTransform::New();
2084 
2085  volumeTransform->SetMatrix(vol->GetMatrix());
2086  volumeTransform->Inverse();
2087 
2088  vtkLightCollection *lights = ren->GetLights();
2089  lights->InitTraversal();
2090 
2091  vtkLight *light[2];
2092  light[0] = lights->GetNextItem();
2093  light[1] = lights->GetNextItem();
2094 
2095  int lightIndex = 0;
2096 
2097  double cameraPosition[3];
2098  double cameraFocalPoint[3];
2099 
2100  ren->GetActiveCamera()->GetPosition(cameraPosition);
2101  ren->GetActiveCamera()->GetFocalPoint(cameraFocalPoint);
2102 
2103  volumeTransform->TransformPoint(cameraPosition, cameraPosition);
2104  volumeTransform->TransformPoint(cameraFocalPoint, cameraFocalPoint);
2105 
2106  double viewDirection[4];
2107 
2108  viewDirection[0] = cameraFocalPoint[0] - cameraPosition[0];
2109  viewDirection[1] = cameraFocalPoint[1] - cameraPosition[1];
2110  viewDirection[2] = cameraFocalPoint[2] - cameraPosition[2];
2111  viewDirection[3] = 0.0;
2112 
2113  vtkMath::Normalize(viewDirection);
2114 
2115  ambientColor[0] = 0.0;
2116  ambientColor[1] = 0.0;
2117  ambientColor[2] = 0.0;
2118  ambientColor[3] = 0.0;
2119 
2120  for (lightIndex = 0; lightIndex < 2; lightIndex++)
2121  {
2122  float dir[3] = {0, 0, 0};
2123  float half[3] = {0, 0, 0};
2124 
2125  if (light[lightIndex] == nullptr || light[lightIndex]->GetSwitch() == 0)
2126  {
2127  lightDiffuseColor[lightIndex][0] = 0.0;
2128  lightDiffuseColor[lightIndex][1] = 0.0;
2129  lightDiffuseColor[lightIndex][2] = 0.0;
2130  lightDiffuseColor[lightIndex][3] = 0.0;
2131 
2132  lightSpecularColor[lightIndex][0] = 0.0;
2133  lightSpecularColor[lightIndex][1] = 0.0;
2134  lightSpecularColor[lightIndex][2] = 0.0;
2135  lightSpecularColor[lightIndex][3] = 0.0;
2136  }
2137  else
2138  {
2139  float lightIntensity = light[lightIndex]->GetIntensity();
2140  double lightColor[3];
2141 
2142  light[lightIndex]->GetDiffuseColor(lightColor);
2143 
2144  double lightPosition[3];
2145  double lightFocalPoint[3];
2146  light[lightIndex]->GetTransformedPosition(lightPosition);
2147  light[lightIndex]->GetTransformedFocalPoint(lightFocalPoint);
2148 
2149  volumeTransform->TransformPoint(lightPosition, lightPosition);
2150  volumeTransform->TransformPoint(lightFocalPoint, lightFocalPoint);
2151 
2152  dir[0] = lightPosition[0] - lightFocalPoint[0];
2153  dir[1] = lightPosition[1] - lightFocalPoint[1];
2154  dir[2] = lightPosition[2] - lightFocalPoint[2];
2155 
2156  vtkMath::Normalize(dir);
2157 
2158  lightDiffuseColor[lightIndex][0] = lightColor[0] * diffuse * lightIntensity;
2159  lightDiffuseColor[lightIndex][1] = lightColor[1] * diffuse * lightIntensity;
2160  lightDiffuseColor[lightIndex][2] = lightColor[2] * diffuse * lightIntensity;
2161  lightDiffuseColor[lightIndex][3] = 0.0;
2162 
2163  lightSpecularColor[lightIndex][0] = lightColor[0] * specular * lightIntensity;
2164  lightSpecularColor[lightIndex][1] = lightColor[1] * specular * lightIntensity;
2165  lightSpecularColor[lightIndex][2] = lightColor[2] * specular * lightIntensity;
2166  lightSpecularColor[lightIndex][3] = 0.0;
2167 
2168  half[0] = dir[0] - viewDirection[0];
2169  half[1] = dir[1] - viewDirection[1];
2170  half[2] = dir[2] - viewDirection[2];
2171 
2172  vtkMath::Normalize(half);
2173 
2174  ambientColor[0] += ambient * lightColor[0];
2175  ambientColor[1] += ambient * lightColor[1];
2176  ambientColor[2] += ambient * lightColor[2];
2177  }
2178 
2179  lightDirection[lightIndex][0] = dir[0];
2180  lightDirection[lightIndex][1] = dir[1];
2181  lightDirection[lightIndex][2] = dir[2];
2182  lightDirection[lightIndex][3] = 0.0;
2183 
2184  halfwayVector[lightIndex][0] = half[0];
2185  halfwayVector[lightIndex][1] = half[1];
2186  halfwayVector[lightIndex][2] = half[2];
2187  halfwayVector[lightIndex][3] = 0.0;
2188  }
2189 
2190  volumeTransform->Delete();
2191 
2192  vtkgl::ProgramLocalParameter4fARB(vtkgl::FRAGMENT_PROGRAM_ARB,
2193  0,
2194  lightDirection[0][0],
2195  lightDirection[0][1],
2196  lightDirection[0][2],
2197  lightDirection[0][3]);
2198 
2199  vtkgl::ProgramLocalParameter4fARB(
2200  vtkgl::FRAGMENT_PROGRAM_ARB, 1, halfwayVector[0][0], halfwayVector[0][1], halfwayVector[0][2], halfwayVector[0][3]);
2201 
2202  vtkgl::ProgramLocalParameter4fARB(vtkgl::FRAGMENT_PROGRAM_ARB, 2, ambient, diffuse, specular, specularPower);
2203 
2204  vtkgl::ProgramLocalParameter4fARB(vtkgl::FRAGMENT_PROGRAM_ARB,
2205  3,
2206  lightDiffuseColor[0][0],
2207  lightDiffuseColor[0][1],
2208  lightDiffuseColor[0][2],
2209  lightDiffuseColor[0][3]);
2210 
2211  vtkgl::ProgramLocalParameter4fARB(vtkgl::FRAGMENT_PROGRAM_ARB,
2212  4,
2213  lightSpecularColor[0][0],
2214  lightSpecularColor[0][1],
2215  lightSpecularColor[0][2],
2216  lightSpecularColor[0][3]);
2217 
2218  vtkgl::ProgramLocalParameter4fARB(
2219  vtkgl::FRAGMENT_PROGRAM_ARB, 5, viewDirection[0], viewDirection[1], viewDirection[2], viewDirection[3]);
2220 
2221  vtkgl::ProgramLocalParameter4fARB(vtkgl::FRAGMENT_PROGRAM_ARB, 6, 2.0, -1.0, 0.0, 0.0);
2222 }
2223 
2224 int vtkMitkOpenGLVolumeTextureMapper3D::IsRenderSupported(vtkRenderer *renderer, vtkVolumeProperty * /*property*/)
2225 {
2226  // GPU_INFO << "IsRenderSupported";
2227 
2228  if (!this->Initialized)
2229  {
2230  // this->Initialize();
2231  this->Initialize(renderer);
2232  }
2233 
2234  if (!this->RenderPossible)
2235  {
2236  MITK_WARN << "vtkMitkOpenGLVolumeTextureMapper3D::IsRenderSupported Rendering not possible";
2237  return 0;
2238  }
2239 
2240  if (!this->GetInput())
2241  {
2242  MITK_WARN << "vtkMitkOpenGLVolumeTextureMapper3D::IsRenderSupported No input available";
2243  return 0;
2244  }
2245 
2246  return 1;
2247 }
2248 
2249 void vtkMitkOpenGLVolumeTextureMapper3D::Initialize(vtkRenderer *renderer)
2250 {
2251  // GPU_INFO << "Initialize";
2252 
2253  this->Initialized = 1;
2254  // vtkOpenGLExtensionManager * extensions = vtkOpenGLExtensionManager::New();
2255  // extensions->SetRenderWindow(nullptr); // set render window to the current one.
2256  vtkOpenGLExtensionManager *extensions =
2257  static_cast<vtkOpenGLRenderWindow *>(renderer->GetRenderWindow())->GetExtensionManager();
2258 
2259  int supports_texture3D = extensions->ExtensionSupported("GL_VERSION_1_2");
2260  if (supports_texture3D)
2261  {
2262  extensions->LoadExtension("GL_VERSION_1_2");
2263  }
2264  else
2265  {
2266  supports_texture3D = extensions->ExtensionSupported("GL_EXT_texture3D");
2267  if (supports_texture3D)
2268  {
2269  extensions->LoadCorePromotedExtension("GL_EXT_texture3D");
2270  }
2271  }
2272 
2273  int supports_multitexture = extensions->ExtensionSupported("GL_VERSION_1_3");
2274  if (supports_multitexture)
2275  {
2276  extensions->LoadExtension("GL_VERSION_1_3");
2277  }
2278  else
2279  {
2280  supports_multitexture = extensions->ExtensionSupported("GL_ARB_multitexture");
2281  if (supports_multitexture)
2282  {
2283  extensions->LoadCorePromotedExtension("GL_ARB_multitexture");
2284  }
2285  }
2286 
2287  this->SupportsCompressedTexture = extensions->ExtensionSupported("GL_VERSION_1_3") == 1;
2288  if (!this->SupportsCompressedTexture)
2289  {
2290  this->SupportsCompressedTexture = extensions->ExtensionSupported("GL_ARB_texture_compression") == 1;
2291  if (this->SupportsCompressedTexture)
2292  {
2293  extensions->LoadCorePromotedExtension("GL_ARB_texture_compression");
2294  }
2295  }
2296  // GPU_INFO(this->SupportsCompressedTexture) << "supporting compressed textures";
2297 
2298  this->SupportsNonPowerOfTwoTextures = extensions->ExtensionSupported("GL_VERSION_2_0") ||
2299  extensions->ExtensionSupported("GL_ARB_texture_non_power_of_two");
2300 
2301  // GPU_INFO << "np2: " << (this->SupportsNonPowerOfTwoTextures?1:0);
2302 
2303  int supports_GL_ARB_fragment_program = extensions->ExtensionSupported("GL_ARB_fragment_program");
2304  if (supports_GL_ARB_fragment_program)
2305  {
2306  extensions->LoadExtension("GL_ARB_fragment_program");
2307  }
2308 
2309  int supports_GL_ARB_vertex_program = extensions->ExtensionSupported("GL_ARB_vertex_program");
2310  if (supports_GL_ARB_vertex_program)
2311  {
2312  extensions->LoadExtension("GL_ARB_vertex_program");
2313  }
2314 
2315  RenderPossible = 0;
2316 
2317  if (supports_texture3D && supports_multitexture && supports_GL_ARB_fragment_program &&
2318  supports_GL_ARB_vertex_program && vtkgl::TexImage3D && vtkgl::ActiveTexture && vtkgl::MultiTexCoord3fv &&
2319  vtkgl::GenProgramsARB && vtkgl::DeleteProgramsARB && vtkgl::BindProgramARB && vtkgl::ProgramStringARB &&
2320  vtkgl::ProgramLocalParameter4fARB)
2321  {
2322  RenderPossible = 1;
2323  }
2324  else
2325  {
2326  std::string errString =
2327  "no gpu-acceleration possible cause following extensions/methods are missing or unsupported:";
2328  if (!supports_texture3D)
2329  errString += " EXT_TEXTURE3D";
2330  if (!supports_multitexture)
2331  errString += " EXT_MULTITEXTURE";
2332  if (!supports_GL_ARB_fragment_program)
2333  errString += " ARB_FRAGMENT_PROGRAM";
2334  if (!supports_GL_ARB_vertex_program)
2335  errString += " ARB_VERTEX_PROGRAM";
2336  if (!vtkgl::TexImage3D)
2337  errString += " glTexImage3D";
2338  if (!vtkgl::ActiveTexture)
2339  errString += " glActiveTexture";
2340  if (!vtkgl::MultiTexCoord3fv)
2341  errString += " glMultiTexCoord3fv";
2342  if (!vtkgl::GenProgramsARB)
2343  errString += " glGenProgramsARB";
2344  if (!vtkgl::DeleteProgramsARB)
2345  errString += " glDeleteProgramsARB";
2346  if (!vtkgl::BindProgramARB)
2347  errString += " glBindProgramARB";
2348  if (!vtkgl::ProgramStringARB)
2349  errString += " glProgramStringARB";
2350  if (!vtkgl::ProgramLocalParameter4fARB)
2351  errString += " glProgramLocalParameter4fARB";
2352  GPU_WARN << errString;
2353  };
2354 
2355  if (RenderPossible)
2356  {
2357  vtkgl::GenProgramsARB(1, &prgOneComponentShade);
2358  vtkgl::BindProgramARB(vtkgl::FRAGMENT_PROGRAM_ARB, prgOneComponentShade);
2359  vtkgl::ProgramStringARB(vtkgl::FRAGMENT_PROGRAM_ARB,
2360  vtkgl::PROGRAM_FORMAT_ASCII_ARB,
2361  static_cast<GLsizei>(strlen(vtkMitkVolumeTextureMapper3D_OneComponentShadeFP)),
2363 
2364  vtkgl::GenProgramsARB(1, &prgRGBAShade);
2365  vtkgl::BindProgramARB(vtkgl::FRAGMENT_PROGRAM_ARB, prgRGBAShade);
2366  vtkgl::ProgramStringARB(vtkgl::FRAGMENT_PROGRAM_ARB,
2367  vtkgl::PROGRAM_FORMAT_ASCII_ARB,
2368  static_cast<GLsizei>(strlen(vtkMitkVolumeTextureMapper3D_FourDependentShadeFP)),
2370  }
2371 }
2372 
2373 // ----------------------------------------------------------------------------
2374 // Print the vtkMitkOpenGLVolumeTextureMapper3D
2375 void vtkMitkOpenGLVolumeTextureMapper3D::PrintSelf(ostream &os, vtkIndent indent)
2376 {
2377  // vtkOpenGLExtensionManager * extensions = vtkOpenGLExtensionManager::New();
2378  // extensions->SetRenderWindow(nullptr); // set render window to current render window
2379 
2380  os << indent << "Initialized " << this->Initialized << endl;
2381  /* if ( this->Initialized )
2382  {
2383  os << indent << "Supports GL_VERSION_1_2:"
2384  << extensions->ExtensionSupported( "GL_VERSION_1_2" ) << endl;
2385  os << indent << "Supports GL_EXT_texture3D:"
2386  << extensions->ExtensionSupported( "GL_EXT_texture3D" ) << endl;
2387  os << indent << "Supports GL_VERSION_1_3:"
2388  << extensions->ExtensionSupported( "GL_VERSION_1_3" ) << endl;
2389  os << indent << "Supports GL_ARB_multitexture: "
2390  << extensions->ExtensionSupported( "GL_ARB_multitexture" ) << endl;
2391  os << indent << "Supports GL_NV_texture_shader2: "
2392  << extensions->ExtensionSupported( "GL_NV_texture_shader2" ) << endl;
2393  os << indent << "Supports GL_NV_register_combiners2: "
2394  << extensions->ExtensionSupported( "GL_NV_register_combiners2" )
2395  << endl;
2396  os << indent << "Supports GL_ATI_fragment_shader: "
2397  << extensions->ExtensionSupported( "GL_ATI_fragment_shader" ) << endl;
2398  os << indent << "Supports GL_ARB_fragment_program: "
2399  << extensions->ExtensionSupported( "GL_ARB_fragment_program" ) << endl;
2400  os << indent << "Supports GL_ARB_texture_compression: "
2401  << extensions->ExtensionSupported( "GL_ARB_texture_compression" )
2402  << endl;
2403  os << indent << "Supports GL_VERSION_2_0:"
2404  << extensions->ExtensionSupported( "GL_VERSION_2_0" )
2405  << endl;
2406  os << indent << "Supports GL_ARB_texture_non_power_of_two:"
2407  << extensions->ExtensionSupported( "GL_ARB_texture_non_power_of_two" )
2408  << endl;
2409  }
2410  extensions->Delete();
2411  */
2412 
2413  if (this->RenderWindow != 0)
2414  {
2415  vtkOpenGLExtensionManager *extensions =
2416  static_cast<vtkOpenGLRenderWindow *>(this->RenderWindow)->GetExtensionManager();
2417 
2418  if (this->Initialized)
2419  {
2420  os << indent << "Supports GL_VERSION_1_2:" << extensions->ExtensionSupported("GL_VERSION_1_2") << endl;
2421  os << indent << "Supports GL_EXT_texture3D:" << extensions->ExtensionSupported("GL_EXT_texture3D") << endl;
2422  os << indent << "Supports GL_VERSION_1_3:" << extensions->ExtensionSupported("GL_VERSION_1_3") << endl;
2423  os << indent << "Supports GL_ARB_multitexture: " << extensions->ExtensionSupported("GL_ARB_multitexture") << endl;
2424  os << indent << "Supports GL_NV_texture_shader2: " << extensions->ExtensionSupported("GL_NV_texture_shader2")
2425  << endl;
2426  os << indent
2427  << "Supports GL_NV_register_combiners2: " << extensions->ExtensionSupported("GL_NV_register_combiners2")
2428  << endl;
2429  os << indent << "Supports GL_ATI_fragment_shader: " << extensions->ExtensionSupported("GL_ATI_fragment_shader")
2430  << endl;
2431  os << indent << "Supports GL_ARB_fragment_program: " << extensions->ExtensionSupported("GL_ARB_fragment_program")
2432  << endl;
2433  os << indent
2434  << "Supports GL_ARB_texture_compression: " << extensions->ExtensionSupported("GL_ARB_texture_compression")
2435  << endl;
2436  os << indent << "Supports GL_VERSION_2_0:" << extensions->ExtensionSupported("GL_VERSION_2_0") << endl;
2437  os << indent << "Supports GL_ARB_texture_non_power_of_two:"
2438  << extensions->ExtensionSupported("GL_ARB_texture_non_power_of_two") << endl;
2439  }
2440  }
2441 
2442  this->Superclass::PrintSelf(os, indent);
2443 }
const char * vtkMitkVolumeTextureMapper3D_FourDependentShadeFP
float k(1.0)
#define MITK_INFO
Definition: mitkLogMacros.h:18
vtkRenderer * GetVtkRenderer() const
const char * vtkMitkVolumeTextureMapper3D_OneComponentShadeFP
Organizes the rendering process.
void vtkVolumeTextureMapper3DComputeScalars(T *dataPtr, vtkMitkVolumeTextureMapper3D *me, float offset, float scale, GLuint volume1, GLuint)
#define myGL_COMPRESSED_RGB_S3TC_DXT1_EXT
#define myGL_COMPRESSED_RGBA_S3TC_DXT5_EXT
static Vector3D offset
#define MITK_WARN
Definition: mitkLogMacros.h:19
vtkStandardNewMacro(vtkMitkOpenGLVolumeTextureMapper3D)
void Normalize(itk::Image< TPixel, VImageDimension > *itkImage, mitk::Image::Pointer im2, mitk::Image::Pointer mask1, std::string output)
Definition: CLBrainMask.cpp:36
void vtkVolumeTextureMapper3DComputeRGBA(unsigned char *dataPtr, vtkMitkVolumeTextureMapper3D *me, GLuint volume1, GLuint volume2)