scalars to color mapping
Randy Heiland
heiland at ncsa.uiuc.edu
Fri Mar 24 13:34:32 EST 2000
Dingrong,
I'm surprised your code was working for small volumes even. The following line
in your code seems to prevent the multi-colored contour surfaces from
appearing:
aContour->ComputeScalarsOff();
I've attached a modified version of the pgm you sent. It generates a volume of
the Lorenz data and works for up to 200x200x200, generating 20 contours. Of
course, it probably depends on how much memory your machine has.
My question is: why does vtkKitwareContourFilter core dump?
--Randy
#define VTK_USE_PATENTED
#define LORENZ
#include "vtkMath.h"
#include "vtkShortArray.h"
#include "vtkRenderer.h"
#include "vtkRenderWindow.h"
#include "vtkRenderWindowInteractor.h"
#include "vtkSphereSource.h"
#include "vtkPolyDataMapper.h"
#include "vtkActor.h"
#include "vtkPlanes.h"
#include "vtkPoints.h"
#include "vtkNormals.h"
#include "vtkImplicitTextureCoords.h"
#include "vtkDataSetMapper.h"
#include "vtkStructuredPointsReader.h"
#include "vtkTexture.h"
#include "vtkDecimate.h"
#include "vtkOutlineFilter.h"
#include "vtkContourFilter.h"
#include "vtkKitwareContourFilter.h"
#include<stdlib.h>
#include<stdio.h>
#include<stream.h>
#include<fstream.h>
#include<math.h>
#define FALSE 0
#define TRUE 1
// Nasty Global Variables
// input parameters
float parLowThresh = 0.5;
float parHighThresh = 1.5;
float parColor[3] = {1.0, 0.0, 0.0};
float parRotation[3] = {0.0, 0.0, 0.0};
float parCamPos[3] = {0.0, 0.0, 0.0};
float parCamFoc[3] = {0.0, 0.0, 0.0};
int parCamera = FALSE;
float parZoom = 1.0;
int parSmooth = 100;
int parIterations = 9;
char outName[50];
char inName[50];
int outPPM = TRUE;
// ---- Lorenz stuff ------
float Pr = 10.0; // The Lorenz parameters
float b = 2.667;
float r = 28.0;
float x, y, z; // starting (and current) x, y, z
float h = 0.01; // integration step size
//int resolution=200; // slice resolution
//int resolution=20; // slice resolution
//int iter = 10000000; // number of iterations
int iter = 10000; // number of iterations
float xmin = -30.0; // x, y, z range for voxels
float xmax = 30.0;
float ymin = -30.0;
float ymax = 30.0;
float zmin = -10.0;
float zmax = 60.0;
int randomMode = 1;
float xIncr, yIncr, zIncr;
short *slice;
// ----------
int main (int argc, char * argv[])
{
if (argc < 3) {
cout <<"Usage: "<<argv[0]<<" <numContours> <domainRes>" <<endl;
cout <<" e.g., " <<argv[0] << " 25 30"<<endl;
cout <<" e.g., " <<argv[0] << " 20 200"<<endl;
exit(-1);
}
int numContours = atoi(argv[1]);
cout <<"numContours = "<< numContours <<endl;
int resolution = atoi(argv[2]);
cout <<"resolution = "<< resolution <<endl;
vtkRenderWindow *renWin;
vtkRenderer *ren1;
vtkRenderWindowInteractor *iren;
// create a window, renderer and interactor
renWin = vtkRenderWindow::New();
ren1 = vtkRenderer::New();
renWin->AddRenderer(ren1);
renWin->SetSize(500,500);
iren = vtkRenderWindowInteractor::New();
iren->SetRenderWindow(renWin);
#ifdef LORENZ
int i, j;
float xx, yy, zz;
short xxx, yyy, zzz;
int sliceSize;
short *s;
void options(int, char**);
int numPts, index;
// take a stab at an integration step size
xIncr = resolution / (xmax - xmin);
yIncr = resolution / (ymax - ymin);
zIncr = resolution / (zmax - zmin);
printf ("The Lorenz Attractor\n");
printf (" Pr = %f\n", Pr);
printf (" b = %f\n", b);
printf (" r = %f\n", r);
printf (" integration step size = %f\n", h);
printf (" slice resolution = %d\n", resolution);
printf (" # of iterations = %d\n", iter);
printf (" specified range:\n");
printf (" x: %f, %f\n", xmin, xmax);
printf (" y: %f, %f\n", ymin, ymax);
printf (" z: %f, %f\n", zmin, zmax);
x = vtkMath::Random(xmin,xmax);
y = vtkMath::Random(ymin,ymax);
z = vtkMath::Random(zmin,zmax);
printf (" starting at %f, %f, %f\n", x, y, z);
// allocate memory for the slices
sliceSize = resolution * resolution;
numPts = sliceSize * resolution;
vtkScalars *scalars = vtkScalars::New(VTK_SHORT);
s = ((vtkShortArray *)scalars->GetData())->WritePointer(0,numPts);
for (i=0; i < numPts; i++) s[i] = 0;
printf (" integrating...\n");
for (j = 0; j < iter; j++)
{
// integrate to next time step
xx = x + h * Pr * (y - x);
yy = y + h * (x * (r - z) - y);
zz = z + h * (x * y - (b * z));
x = xx; y = yy; z = zz;
// calculate voxel index
if (x < xmax && x > xmin && y < ymax && y > ymin && z < zmax && z > zmin)
{
xxx = (short) ((float)(xx - xmin) * xIncr);
yyy = (short) ((float)(yy - ymin) * yIncr);
zzz = (short) ((float)(zz - zmin) * zIncr);
index = xxx + yyy*resolution + zzz*sliceSize;
s[index] += 1;
}
}
vtkStructuredPoints *volume = vtkStructuredPoints::New();
volume->GetPointData()->SetScalars(scalars);
volume->SetDimensions(resolution,resolution,resolution);
volume->SetOrigin(xmin,ymin,zmin);
volume->SetSpacing((xmax-xmin)/resolution, (ymax-ymin)/resolution,
(zmax-zmin)/resolution);
#else
// try to read the data
vtkStructuredPointsReader *aReader;
aReader = vtkStructuredPointsReader::New();
aReader->SetFileName(argv[1]);
aReader->Update();
vtkStructuredPoints *aPoints;
aPoints = vtkStructuredPoints::New();
aPoints = aReader->GetOutput();
#endif
/*-----------
float max = 0;
// for (i = 0; i < aPoints->GetNumberOfPoints(); i++)
for (i = 0; i < volume->GetNumberOfPoints(); i++)
{
// if ( aPoints->GetPointData()->GetScalars()->GetScalar(i) > max)
if ( volume->GetPointData()->GetScalars()->GetScalar(i) > max)
{
// max = aPoints->GetPointData()->GetScalars()->GetScalar(i);
max = volume->GetPointData()->GetScalars()->GetScalar(i);
}
}
cout<<"Maximum value of the input data = "<<max<<endl;
-------------*/
float *srange = volume->GetPointData()->GetScalars()->GetRange();
cout<<"scalar range = "<<srange[0]<<", "<<srange[1]<<endl;
// Extract isosurface
vtkContourFilter *aContour = vtkContourFilter::New();
// hmmm, get a core dump when try to use this
// vtkKitwareContourFilter *aContour = vtkKitwareContourFilter::New();
// aContour->SetInput(aPoints);
aContour->SetInput(volume);
// aContour->ComputeScalarsOff();
/*----
aContour->ComputeGradientsOff();
aContour->ComputeNormalsOff();
---*/
aContour->UseScalarTreeOn();
float contourDelta = (srange[1] - srange[0]) / (numContours-1);
cout <<"contourDelta = "<< contourDelta <<endl;
int count = 0;
// for (i = 0; i <= max; i = i + 3)
for (float sval = srange[0]; sval <= srange[1]; sval = sval +
contourDelta)
{
aContour->SetValue(count, sval );
cout <<"--- setting aContour "<< count << ", value = "<< sval <<endl;
count = count + 1;
}
cout<<"Number of Contours = "<<aContour->GetNumberOfContours()<<endl;
aContour->SetNumberOfContours(count);
//derive my own lookup table for scalars color mapping.
vtkLookupTable *Lut = vtkLookupTable::New();
Lut->SetHueRange (0.0, 0.667);
Lut->SetSaturationRange (1.0, 1.0);
Lut->SetValueRange (1.0, 1.0);
Lut->SetAlphaRange( 1.0, 1.0);
Lut->SetNumberOfColors(256);
Lut->Build();
// map to graphics library
vtkPolyDataMapper *mapper;
mapper = vtkPolyDataMapper::New();
mapper->SetInput(aContour->GetOutput());
mapper->ScalarVisibilityOn();
//printf("setting scalar range = 0 %f\n",max);
// mapper->SetScalarRange(0, max);
mapper->SetScalarRange(srange);
// mapper->SetLookupTable(Lut);
cout<<"Finish mapper to graphics library."<<endl;
// actor coordinates geometry, properties, transformation
vtkActor *aObject;
aObject = vtkActor::New();
aObject->SetMapper(mapper);
aObject->GetProperty()->SetOpacity(0.5);
aObject->Update();
cout<<"Finish updating the actor."<<endl;
// Create box around object
// based on original dimensions of the image
vtkOutlineFilter *outline;
outline = vtkOutlineFilter::New();
// outline->SetInput(aReader->GetOutput());
outline->SetInput(volume);
vtkPolyDataMapper *outlineMapper;
outlineMapper = vtkPolyDataMapper::New();
outlineMapper->SetInput(outline->GetOutput());
vtkActor *outlineActor;
outlineActor = vtkActor::New();
outlineActor->SetMapper(outlineMapper);
outlineActor->GetProperty()->SetColor(0.0, 0.0, 1.0);
float * bounds;
bounds = (float *) malloc(sizeof(float)*6);
bounds = outlineActor->GetBounds();
float xc = (bounds[1]-bounds[0])/2.0;
float yc = (bounds[3]-bounds[2])/2.0;
float zc = (bounds[5]-bounds[4])/2.0;
cout<<"Origin xc, yc, zc = "<<xc<<" "<<yc<<" "<<zc<<endl;
// Set Origin of the objects at the same position
aObject->SetOrigin(xc,yc,zc);
outlineActor->SetOrigin(xc,yc,zc);
// Initial rotation
aObject->RotateWXYZ(parRotation[0],1,0,0);
aObject->RotateWXYZ(parRotation[1],0,1,0);
aObject->RotateWXYZ(parRotation[2],0,0,1);
outlineActor->RotateWXYZ(parRotation[0],1,0,0);
outlineActor->RotateWXYZ(parRotation[1],0,1,0);
outlineActor->RotateWXYZ(parRotation[2],0,0,1);
outlineActor->GetProperty()->SetOpacity(0.3);
//add actors to the renderer
ren1->AddActor(aObject);
ren1->AddActor(outlineActor);
ren1->SetBackground(1.0,1.0,1.0);
cout<<"Back ground color white"<<endl;
cout<<"Outline color = blue"<<endl;
cout<<"Ojbect color = Look up table"<<endl;
//if camera not specified use automatically created setup
if (parCamera==TRUE) {
ren1->GetActiveCamera()->SetPosition (parCamPos[0], parCamPos[1],
parCamPos[2]);
ren1->GetActiveCamera()->SetFocalPoint(parCamFoc[0], parCamFoc[1],
parCamFoc[2]);
ren1->GetActiveCamera()->SetViewUp(0,1,0);
}
// update the zoom
ren1->GetActiveCamera()->Zoom(parZoom);
// Render
renWin->Render();
cout<<"Render " <<endl;
/*-----
// Output PPM if asked
if (outPPM == TRUE) {
// char Name1[50];
float angle = 360.0/(float) parIterations;
for (int i=0; i<=parIterations; i++) {
cout<<"save "<<i<<"th scene. "<<endl;
cout<<"Type in an file name for PPM file: "<<endl;
cin>>outName;
renWin->SetFileName(outName);
renWin->SaveImageAsPPM();
aObject->RotateWXYZ((float) angle,0,1,0);
outlineActor->RotateWXYZ((float)angle,0,1,0);
//ren1->UpdateActors();
renWin->Render();
}
}
-----*/
cout<<"Begin mouse interaction. "<<endl;
// Begin mouse interaction
iren->Start();
return 0;
}
--------------------------------------------------------------------
This is the private VTK discussion list. Please keep messages on-topic.
Check the FAQ at: <http://public.kitware.com/cgi-bin/vtkfaq>
To UNSUBSCRIBE, send message body containing "unsubscribe vtkusers" to
<majordomo at public.kitware.com>. For help, send message body containing
"info vtkusers" to the same address.
--------------------------------------------------------------------
More information about the vtkusers
mailing list