[vtkusers] x-ray effect in volume rendering
Divye
zombi2_84 at yahoo.com
Mon Oct 26 16:47:15 EDT 2009
Hey,
Were you able to try it successfully in java?
I am trying to do volume rendering of dicom images in java, and am runnning into lot of problems.
--- On Mon, 10/26/09, baliki2 at freemail.gr <baliki2 at freemail.gr> wrote:
From: baliki2 at freemail.gr <baliki2 at freemail.gr>
Subject: Re: [vtkusers] x-ray effect in volume rendering
To: vtkusers at vtk.org
Date: Monday, October 26, 2009, 6:59 AM
It works just great! I'll try it in Java, too.
Thanks a lot Mr. Gasteiger!
> Hi "Baliki",
>
> Below I have inserted a VTK-example program, which reads a raw volume
> dataset and performs normal volume rendering, X-ray and MIP projection. I
> use it for my students to illustrate these different approaches. The samples
> on each ray are interpolated with nearest neighbor (not the best, tri-linear
> is much better but expansive in computation time. It inFollowing steps you
> have to do:
> 1. Define "fname" (in main(...)) with the file name of you dataset
> 2 Specify your raw dataset for the vtkImageReader2 reader (In the example I
> use a volume version(!) of the foot dataset of the VTK distribution)
> 3. Turn on one of the render modes by uncomment one of the
> "currentRenderMode" line at the beginning of the example
>
> I have included some comments which hopefully helps you to understand each
> step.
>
> Best regards, Rocco
>
> //--------------------------------------------------------------------------
> ------------
> //
> // This example reads a volume dataset and then displays it.
> //
> #include "vtkRenderer.h"
> #include "vtkRenderWindow.h"
> #include "vtkRenderWindowInteractor.h"
> #include "vtkCamera.h"
> #include "vtkInteractorStyleSwitch.h"
> #include "vtkImageReader2.h"
> #include "vtkActor.h"
>
> #include "vtkColorTransferFunction.h"
> #include "vtkPiecewiseFunction.h"
> #include "vtkVolumeProperty.h"
> #include "vtkVolumeRayCastCompositeFunction.h"
> #include "vtkVolumeRayCastMapper.h"
> #include "vtkVolume.h"
>
>
> #define mode_MIP 0
> #define mode_XRay 1
> #define mode_Full 2
>
> static int currentRenderMode = mode_XRay;
>
> // Necessary for correct MIP and XRay scaling (0-1)
> static float globalMaximumScalarValue = 0;
> ////////////////////////////////////////////////////////////////
> // Definition and implementation of our own ray cast composite function
> ----------------------------------
> class MyVolumeRayCastCompositeFunction : public
> vtkVolumeRayCastCompositeFunction
> {
>
> public:
> static MyVolumeRayCastCompositeFunction *New()
> { return new MyVolumeRayCastCompositeFunction; }
>
> protected:
>
> ////////////////////////////////////////////////////////////////
> // Do not change this code
> //
> // This is called from RenderAnImage (in vtkDepthPARCMapper.cxx)
> // It uses the integer data type flag that is passed in
> // to determine what type of ray needs to be cast (which is handled
> // by a template function. It also uses the shading and
> // interpolation types to determine which template function
> // to call.
> void MyVolumeRayCastCompositeFunction::CastRay(
> vtkVolumeRayCastDynamicInfo *dynamicInfo,
> vtkVolumeRayCastStaticInfo *staticInfo )
> {
> void *data_ptr;
> data_ptr = staticInfo->ScalarDataPointer;
>
> switch ( staticInfo->ScalarDataType ) {
> case VTK_UNSIGNED_CHAR:
> vtkCastRay_OwnWork( (unsigned char *)data_ptr, dynamicInfo,
> staticInfo );
> break;
> case VTK_UNSIGNED_SHORT:
> vtkCastRay_OwnWork( (unsigned short *)data_ptr, dynamicInfo,
> staticInfo );
> break;
> default:
> vtkWarningMacro ( << "Unsigned char and unsigned short are the only
> supported datatypes for rendering" );
> break;
> }
> }
>
> // This is the template function that actually casts a ray and
> computes
> // the composite value. This version uses nearest neighbor
> interpolation, and performs shading.
> template <class T>
> void vtkCastRay_OwnWork(T *data_ptr, vtkVolumeRayCastDynamicInfo
> *dynamicInfo, vtkVolumeRayCastStaticInfo *staticInfo )
> {
> unsigned char *gmptr = NULL;
> float accum_red_intensity, accum_green_intensity,
> accum_blue_intensity, accum_opacity;
> float opacity;
> int xinc, yinc, zinc;
> int voxel[3];
> float ray_position[3];
> T *dptr;
> float *ScalarToOpacityTF;
> float *ScalarToColorTF;
> float *red_d_shade, *green_d_shade, *blue_d_shade;
> float *red_s_shade, *green_s_shade, *blue_s_shade;
> unsigned short *encoded_normals;
> int currentNormalIndex;
> float red_shaded_value, green_shaded_value,
> blue_shaded_value;
> int offset;
> int steps_this_ray;
> int scalar_value;
> float r, g, b;
> int num_steps;
> float *ray_start, *ray_increment;
> int shading = staticInfo->Shading;
>
>
> // For each ray some information are fetched which are
> needed for color and shading
> // computation.
>
> // Get ray parameters: Number of sample points, ray start
> position, and step size on the ray
> num_steps = dynamicInfo->NumberOfStepsToTake;
> ray_start = dynamicInfo->TransformedStart;
> ray_increment = dynamicInfo->TransformedIncrement;
>
> // Get diffuse shading table pointers. This table is
> preprocessed by VTK and stores
> // for different normal directions its diffuse portion of
> the phong illumination term.
> red_d_shade = staticInfo->RedDiffuseShadingTable;
> green_d_shade = staticInfo->GreenDiffuseShadingTable;
> blue_d_shade = staticInfo->BlueDiffuseShadingTable;
>
> // Get specular shading table pointers. This table is
> preprocessed by VTK and stores
> // for specular normal directions its diffuse portion of the
> phong illumination term.
> red_s_shade = staticInfo->RedSpecularShadingTable;
> green_s_shade = staticInfo->GreenSpecularShadingTable;
> blue_s_shade = staticInfo->BlueSpecularShadingTable;
>
> // Get a pointer to the encoded normals of this volume. This
> table is preprocessed by VTK and stores
> // for each voxel its normal.
> encoded_normals = staticInfo->EncodedNormals;
>
> // Get the scalar opacity transfer function which maps
> scalar input values to opacities
> ScalarToOpacityTF =
> staticInfo->Volume->GetCorrectedScalarOpacityArray();
>
> // Get the color transfer function which maps scalar input
> values to RGB values
> ScalarToColorTF = staticInfo->Volume->GetRGBArray();
>
> // store increments (distance between sample points in x, y,
> and z direction) in local variables
> xinc = staticInfo->DataIncrement[0];
> yinc = staticInfo->DataIncrement[1];
> zinc = staticInfo->DataIncrement[2];
> //cerr<<xinc<<", "<<yinc<<", "<<zinc<<", "<<endl; //constant
> for a data set
>
> // Initialize the ray position (in voxel space)
> ray_position[0] = ray_start[0];
> ray_position[1] = ray_start[1];
> ray_position[2] = ray_start[2];
>
>
> // Nearest neighbor interpolation - get nearest voxel (in
> voxel space) - you have to change this for tri-linear interpolation
> voxel[0] = vtkRoundFuncMacro( ray_position[0] );
> voxel[1] = vtkRoundFuncMacro( ray_position[1] );
> voxel[2] = vtkRoundFuncMacro( ray_position[2] );
>
> // So far we haven't accumulated anything
> accum_red_intensity = 0.0;
> accum_green_intensity = 0.0;
> accum_blue_intensity = 0.0;
> // accumulated opacity
> accum_opacity = 0.0;
>
> // added for MIP mode rendering
> float maxOpacity = 0.0;
> float maxValue = 0.0;
>
> currentNormalIndex = 0;
>
> // For each step along the ray do something (accumulate,
> find maximum opacity, ...)
> // "accum_opacity < 1.0" means early-ray-termination
> for ( steps_this_ray = 0; steps_this_ray < num_steps &&
> accum_opacity < 1.0; steps_this_ray++ ) {
>
> // get offset position of current voxel in linear
> list
> offset = voxel[2] * zinc + voxel[1] * yinc +
> voxel[0];
> // get pointer on current voxel in linear list
> dptr = data_ptr + offset;
>
> // get scalar value at current location
> scalar_value = (int) *(dptr);
>
> if (scalar_value > globalMaximumScalarValue)
> {
> globalMaximumScalarValue = scalar_value;
> }
>
> // get opacity value for this scalar_value out of
> the ScalarOpacityTransferFunction
> opacity = ScalarToOpacityTF[scalar_value];
>
> // If we have a valid opacity value, then compute
> the shading
> if ( opacity ) {
> switch (currentRenderMode)
> {
> case mode_MIP:
> // MIP - only maximal value or first
> local maximum above a threshold is visualized
> // check if current voxel opacity is
> above max opacity
> //if (opacity > maxOpacity)
> if (scalar_value > maxValue)
> {
> accum_red_intensity =
> opacity * ScalarToColorTF[(scalar_value) * 3 ];
> accum_green_intensity =
> opacity * ScalarToColorTF[(scalar_value) * 3 + 1];
> accum_blue_intensity =
> opacity * ScalarToColorTF[(scalar_value) * 3 + 2];
>
> //maxOpacity = opacity;
> maxValue = scalar_value;
> accum_opacity = 1.0-opacity;
>
> }
>
> break;
> case mode_XRay:
> // X-ray projection mode:
> Accumulation of the scalar values (or opacity values)
> accum_red_intensity =
> accum_green_intensity = accum_blue_intensity = accum_red_intensity +
> scalar_value;
>
> break;
> case mode_Full:
> // get r, g, and b value for current
> scalar value out of the ColorTransferFunction
> r = ScalarToColorTF[(scalar_value) *
> 3 ];
> g = ScalarToColorTF[(scalar_value) *
> 3 + 1];
> b = ScalarToColorTF[(scalar_value) *
> 3 + 2];
>
> red_shaded_value = opacity * r;
> green_shaded_value = opacity * g;
> blue_shaded_value = opacity * b;
>
>
> // Accumulate color
> // (the color of the current value
> is weighted by the remaining opacity of the voxels
> // in front of the current one)
> accum_red_intensity =
> red_shaded_value * (1.0-accum_opacity) + accum_red_intensity;
> accum_green_intensity =
> green_shaded_value * (1.0-accum_opacity) + accum_green_intensity;
> accum_blue_intensity =
> blue_shaded_value * (1.0-accum_opacity) + accum_blue_intensity;
>
> // Accumulate opacity
> accum_opacity =
> opacity*(1.0-accum_opacity) + accum_opacity;
>
> break;
> default:
> accum_red_intensity = 0.0;
> accum_green_intensity = 0.0;
> accum_blue_intensity = 0.0;
> accum_opacity = 1.0;
>
> break;
> }
>
> }
>
> // Increment our position and compute new voxel
> location
> ray_position[0] += ray_increment[0];
> ray_position[1] += ray_increment[1];
> ray_position[2] += ray_increment[2];
> voxel[0] = vtkRoundFuncMacro( ray_position[0] );
> voxel[1] = vtkRoundFuncMacro( ray_position[1] );
> voxel[2] = vtkRoundFuncMacro( ray_position[2] );
> }// endFor each step along the ray
> ----------------------------------------------------------------------------
> -----------
>
>
> // Now we stored in accum_{red,green,blue}_intensity the
> accumulated RGB values along a ray
> // added for XRay mode rendering
> if( currentRenderMode == mode_XRay)
> {
> // divide accumulated scalar values by the number of
> steps times maximum scalar value to get the average
> accum_red_intensity = accum_green_intensity =
> accum_blue_intensity = accum_red_intensity /
> float(num_steps*globalMaximumScalarValue);
> accum_opacity = 1.0-accum_red_intensity;
>
> }
>
> // Cap the accumulated intensities at 1.0
> if ( accum_red_intensity > 1.0 ) {
> accum_red_intensity = 1.0; }
> if ( accum_green_intensity > 1.0 ){ accum_green_intensity =
> 1.0;}
> if ( accum_blue_intensity > 1.0 ) {
> accum_blue_intensity = 1.0;}
>
> // Set the return pixel value.
> dynamicInfo->Color[0] = accum_red_intensity;
> dynamicInfo->Color[1] = accum_green_intensity;
> dynamicInfo->Color[2] = accum_blue_intensity;
> dynamicInfo->Color[3] = accum_opacity;
> dynamicInfo->NumberOfStepsTaken = steps_this_ray;
>
>
> }
>
> };
> //
> ----------------------------------------------------------------------------
> -
>
>
> ////////////////////////////////////////////////////////////////
> //
> int main (int argc, char **argv)
> {
> //check for arguments (e.g. a filename)
> const char* fname = "foot.raw";
>
> // Create the renderer, the render window, and the interactor. The
> renderer
> // draws into the render window, the interactor enables mouse- and
> // keyboard-based interaction with the data within the render
> window.
> vtkRenderer *renderer = vtkRenderer::New();
> //renderer->SetBackground( 0.2, 0.2, 0.4 );
> renderer->SetBackground( 1.0, 1.0, 1.0 );
> vtkCamera *camera = renderer->GetActiveCamera();
> camera->ParallelProjectionOff();
> camera->SetViewUp (0, 0, -1);
> camera->SetPosition (-1, 2, -0.5);
> vtkRenderWindow *renderWindow = vtkRenderWindow::New();
> renderWindow->SetSize( 640, 480 ); // Set the size of the render
> window (in pixel).
> renderWindow->AddRenderer(renderer);
> vtkRenderWindowInteractor *interactionRenderer =
> vtkRenderWindowInteractor::New();
> interactionRenderer->SetRenderWindow(renderWindow);
>
> reinterpret_cast<vtkInteractorStyleSwitch*>(interactionRenderer->GetInteract
> orStyle())->SetCurrentStyleToTrackballCamera();
>
> vtkImageReader2 *volumeDataset = vtkImageReader2::New();
> volumeDataset->SetFileName(fname);
> volumeDataset->SetFileDimensionality(3);
> volumeDataset->SetDataExtent(0, 127, 0, 127, 0, 63);
> volumeDataset->SetDataSpacing(1.60938, 1.60938, 1.95313);
> volumeDataset->SetDataOrigin(0.0, 0.0, 0.0);
> volumeDataset->SetDataScalarTypeToUnsignedShort();
> volumeDataset->SetDataByteOrderToLittleEndian();
> volumeDataset->UpdateWholeExtent();
>
>
>
> ////////////////////////////////////////////////////////////////////////////
> ///////////
> // Feel free to change this code
> // adjustments to the following code may be required and should be
> explored
>
> // Create transfer function - mapping scalar values [0-...] to COLOR
> [0-1, 0-1, 0-1]
> vtkColorTransferFunction *colorTransferFunction =
> vtkColorTransferFunction::New();
> colorTransferFunction->AddRGBPoint( 0.0, 0.0,0.0,0.0);
> colorTransferFunction->AddRGBPoint( 500.0, 0.9,0.5,0.3);
> colorTransferFunction->AddRGBPoint(1100.0, 0.8,0.8,0.6);
> colorTransferFunction->AddRGBPoint(1200.0, 0.6,0.6,0.6);
>
> vtkPiecewiseFunction *opacityTransferFunction =
> vtkPiecewiseFunction::New();
> opacityTransferFunction->AddPoint( 0, 0.0);
> //opacityTransferFunction->AddPoint( 980, 0.1);
> //opacityTransferFunction->AddPoint( 1055, 0.2);
> opacityTransferFunction->AddPoint( 1200, 0.05);
> opacityTransferFunction->AddPoint( 1300, 1.0);
>
> // The property describes how the data will look
> vtkVolumeProperty *volumeProperty = vtkVolumeProperty::New();
> volumeProperty->SetColor(colorTransferFunction);
> volumeProperty->SetScalarOpacity(opacityTransferFunction);
> //volumeProperty->ShadeOn(); // request 3d shading (german:
> beleuchtungsberechnung anfordern) //has to be implemented in
> vtkCastRay_OwnWork
>
> //vtkVolumeRayCastCompositeFunction *rayCompositeFunction =
> vtkVolumeRayCastCompositeFunction::New();
> MyVolumeRayCastCompositeFunction *rayCompositeFunction =
> MyVolumeRayCastCompositeFunction::New();
> // The mapper knows how to render the data
> vtkVolumeRayCastMapper *rayCastMapper =
> vtkVolumeRayCastMapper::New();
> rayCastMapper->SetInput( volumeDataset->GetOutput());
> rayCastMapper->SetVolumeRayCastFunction(rayCompositeFunction);
> //Feel free to change parameters for
> //rayCastMapper->AutoAdjustSampleDistancesOff(); //
> interactive framerates by lower quality (on/off)
> //rayCastMapper->SetSampleDistance(1.0);
> // adjust samples per ray rate if AutoAdjustSampleDistances is OFF
> //rayCastMapper->SetImageSampleDistance(1.0); //
> adjust rays per pixel rate (e.g. 0.5 => 4 rays for each pixel, 2.0 => 1 ray
> for 4 pixel) if AutoAdjustSampleDistances is OFF
> //rayCastMapper->SetMaximumImageSampleDistance(6.0); // max
> sample distance if AutoAdjustSampleDistances is ON
> //rayCastMapper->SetMinimumImageSampleDistance(1.0); // min
> sample distance if AutoAdjustSampleDistances is ON
>
> // The volume holds the mapper and the property and
> // can be used to position/orient the volume
> vtkVolume *volData= vtkVolume::New();
> volData->SetMapper(rayCastMapper);
> volData->SetProperty(volumeProperty);
>
> // Actors are added to the renderer.
> renderer->AddVolume(volData);
>
>
> // adjust camera for good initial view
> renderer->ResetCamera();
> renderer->GetActiveCamera()->Dolly(1);
> renderer->ResetCameraClippingRange();
>
> // Start the Interactor
> interactionRenderer->Initialize();
> interactionRenderer->Start();
>
> // It is important to delete all objects created previously to
> prevent
> // memory leaks. In this case, since the program is on its way to
> // exiting, it is not so important. But in applications it is
> // essential.
> renderer->Delete();
> renderWindow->Delete();
> interactionRenderer->Delete();
>
> return 0;
> }
>
>
> -----Original Message-----
> From: baliki2 at freemail.gr [mailto:baliki2 at freemail.gr]
> Sent: Friday, October 23, 2009 2:11 PM
> To: post at rocco-gasteiger.de
> Subject: Re: Re: [vtkusers] x-ray effect in volume rendering
>
>
> I don't really understand what you mean by "divide it by
> > numSteps*globalMaximumScalarValue at the end".
>
> What i would like to get is something like this:
> http://www.onlinetelemedicine.com/HTML/product/sam_images/X-Ray.jpg
> If you look carefully, the bones edges are more bright than the "inside"
> bone area.
> I've achieved a similar effect, but not this exact one:
> http://img12.imageshack.us/i/sendj.png/
> Please download it and zoom in to see it better.
>
> Thanks!
>
> > Hello,
> >
> > Yes, you can! ;-)
> >
> > It is quite simply. In your composite function accumulate all sampled
> scalar
> > values (or opacity values) on your ray and divide it by
> > numSteps*globalMaximumScalarValue at the end (to get the average between 0
> > and 1). I've implement it and it should work.
> >
> > Best regards, Rocco
> >
> > --------------------------------------------------
> > Dipl.-Ing. Rocco Gasteiger
> > Otto-von-Guericke University
> > Faculty of Computer Science
> > Department of Simulation and Graphics
> > Universit=E4tsplatz 2, 39106 Magdeburg, Germany
> >
> > Office: G29-223
> > Phone: +49 391 67 127 59
> > Fax:=A0=A0 +49 391 67 111 64
> > Website: http://wwwisg.cs.uni-magdeburg.de/cvcms/ =
> >
> > --------------------------------------------------
> >
> >
> >
> > -----Original Message-----
> > From: vtkusers-bounces at vtk.org [mailto:vtkusers-bounces at vtk.org] On Behalf
> > Of baliki2 at freemail.gr
> > Sent: Friday, October 23, 2009 11:30 AM
> > To: vtkusers at vtk.org
> > Subject: [vtkusers] x-ray effect in volume rendering
> >
> >
> > Is there a way that i can achieve an x-ray effect in volume rendering of
> > DICOM images?
> > I want to see the constructed volume as a 2D xray-like "image".
> >
> > I've tied the composite function combined with some tranfer and color
> > functions, but i didn't take the exact x-ray effect. The images i use are
> of
> > short type.
> >
> >
> > _______________________________________________
> > Powered by www.kitware.com
> >
> > Visit other Kitware open-source projects at
> > http://www.kitware.com/opensource/opensource.html
> >
> > Please keep messages on-topic and check the VTK FAQ at:
> > http://www.vtk.org/Wiki/VTK_FAQ
> >
> > Follow this link to subscribe/unsubscribe:
> > http://www.vtk.org/mailman/listinfo/vtkusers
> >
> > _______________________________________________
> > Powered by www.kitware.com
> >
> > Visit other Kitware open-source projects at
> http://www.kitware.com/opensour=
> > ce/opensource.html
> >
> > Please keep messages on-topic and check the VTK FAQ at:
> http://www.vtk.org/=
> > Wiki/VTK_FAQ
> >
> > Follow this link to subscribe/unsubscribe:
> > http://www.vtk.org/mailman/listinfo/vtkusers
>
> _______________________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at http://www.kitware.com/opensource/opensource.html
>
> Please keep messages on-topic and check the VTK FAQ at: http://www.vtk.org/Wiki/VTK_FAQ
>
> Follow this link to subscribe/unsubscribe:
> http://www.vtk.org/mailman/listinfo/vtkusers
_______________________________________________
Powered by www.kitware.com
Visit other Kitware open-source projects at http://www.kitware.com/opensource/opensource.html
Please keep messages on-topic and check the VTK FAQ at: http://www.vtk.org/Wiki/VTK_FAQ
Follow this link to subscribe/unsubscribe:
http://www.vtk.org/mailman/listinfo/vtkusers
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.vtk.org/pipermail/vtkusers/attachments/20091026/0910f8fb/attachment.htm>
More information about the vtkusers
mailing list