Commit 56961fdb authored by Will Schroeder's avatar Will Schroeder
Browse files

Initial revision

parent 984d0509
/*=========================================================================
Program: Visualization Library
Module: Splatter.hh
Language: C++
Date: $Date$
Version: $Revision$
Description:
---------------------------------------------------------------------------
This file is part of the Visualization Library. No part of this file
or its contents may be copied, reproduced or altered in any way
without the express written consent of the authors.
Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen 1993, 1994
=========================================================================*/
//
// Splat points using Gaussian distribution
//
#ifndef __vlGaussianSplatter_h
#define __vlGaussianSplatter_h
#include "DS2SPtsF.hh"
class vlGaussianSplatter : public vlDataSetToStructuredPointsFilter
{
public:
vlGaussianSplatter();
~vlGaussianSplatter() {};
char *GetClassName() {return "vlGaussianSplatter";};
void PrintSelf(ostream& os, vlIndent indent);
void ComputeModelBounds();
void SetSampleDimensions(int i, int j, int k);
void SetSampleDimensions(int dim[3]);
vlGetVectorMacro(SampleDimensions,int);
vlSetClampMacro(Radius,float,0.0,1.0);
vlGetMacro(Radius,float);
vlSetClampMacro(ScaleFactor,float,0.0,LARGE_FLOAT);
vlGetMacro(ScaleFactor,float);
vlSetMacro(ExponentFactor,float);
vlGetMacro(ExponentFactor,float);
vlSetClampMacro(Eccentricity,float,0.001,LARGE_FLOAT);
vlGetMacro(Eccentricity,float);
void SetModelBounds(float *bounds);
void SetModelBounds(float xmin, float xmax, float ymin, float ymax, float zmin, float zmax);
vlGetVectorMacro(ModelBounds,float);
vlSetMacro(NormalWarping,int);
vlGetMacro(NormalWarping,int);
vlBooleanMacro(NormalWarping,int);
vlSetMacro(ScalarWarping,int);
vlGetMacro(ScalarWarping,int);
vlBooleanMacro(ScalarWarping,int);
vlSetMacro(Capping,int);
vlGetMacro(Capping,int);
vlBooleanMacro(Capping,int);
vlSetMacro(CapValue,float);
vlGetMacro(CapValue,float);
protected:
void Execute();
void Cap(vlFloatScalars *s);
int SampleDimensions[3]; // dimensions of volume to splat into
float Radius; // maximum distance splat propagates (as fraction 0->1)
float ExponentFactor; // scale exponent of gaussian function
float ModelBounds[6]; // bounding box of splatting dimensions
int NormalWarping; // on/off warping of splat via normal
float Eccentricity;// elliptic distortion due to normals
int ScalarWarping; // on/off warping of splat via scalar
float ScaleFactor; // splat size influenced by scale factor
int Capping; // Cap side of volume to close surfaces
float CapValue; // value to use for capping
// recursive propagation of splat
void SplitIJK(int i, int idir, int j, int jdir, int k, int kdir);
void SplitIJ(int i, int idir, int j, int jdir, int k);
void SplitIK(int i, int idir, int j, int k, int kdir);
void SplitJK(int i, int j, int jdir, int k, int kdir);
void SplitI(int i, int idir, int j, int k);
void SplitJ(int i, int j, int jdir, int k);
void SplitK(int i, int j, int k, int kdir);
float Gaussian(float x[3]);
float EccentricGaussian(float x[3]);
float ScalarSampling(float s) {return this->ScaleFactor * s;};
float PositionSampling(float s) {return this->ScaleFactor;};
void SetScalar(int idx, float dist2);
};
#endif
/*=========================================================================
Program: Visualization Library
Module: Splatter.cc
Language: C++
Date: $Date$
Version: $Revision$
Description:
---------------------------------------------------------------------------
This file is part of the Visualization Library. No part of this file
or its contents may be copied, reproduced or altered in any way
without the express written consent of the authors.
Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen 1993, 1994
=========================================================================*/
#include <math.h>
#include "Splatter.hh"
#include "FScalars.hh"
vlGaussianSplatter::vlGaussianSplatter()
{
this->SampleDimensions[0] = 50;
this->SampleDimensions[1] = 50;
this->SampleDimensions[2] = 50;
this->Radius = 0.1;
this->ExponentFactor = -5.0;
this->ModelBounds[0] = 0.0;
this->ModelBounds[1] = 0.0;
this->ModelBounds[2] = 0.0;
this->ModelBounds[3] = 0.0;
this->ModelBounds[4] = 0.0;
this->ModelBounds[5] = 0.0;
this->NormalWarping = 1;
this->Eccentricity = 2.5;
this->ScalarWarping = 1;
this->ScaleFactor = 1.0;
this->Capping = 1;
this->CapValue = LARGE_FLOAT;
}
void vlGaussianSplatter::PrintSelf(ostream& os, vlIndent indent)
{
if (this->ShouldIPrint(vlGaussianSplatter::GetClassName()))
{
vlDataSetToStructuredPointsFilter::PrintSelf(os,indent);
os << indent << "Sample Dimensions: (" << this->SampleDimensions[0] << ", "
<< this->SampleDimensions[1] << ", "
<< this->SampleDimensions[2] << ")\n";
os << indent << "Radius: " << this->Radius << "\n";
os << indent << "Exponent Factor: " << this->ExponentFactor << "\n";
os << indent << "ModelBounds: \n";
os << indent << " Xmin,Xmax: (" << this->ModelBounds[0] << ", " << this->ModelBounds[1] << ")\n";
os << indent << " Ymin,Ymax: (" << this->ModelBounds[2] << ", " << this->ModelBounds[3] << ")\n";
os << indent << " Zmin,Zmax: (" << this->ModelBounds[4] << ", " << this->ModelBounds[5] << ")\n";
os << indent << "Normal Warping: " << (this->NormalWarping ? "On\n" : "Off\n");
os << indent << "Eccentricity: " << this->Eccentricity << "\n";
os << indent << "Scalar Warping: " << (this->ScalarWarping ? "On\n" : "Off\n");
os << indent << "Scale Factor: " << this->ScaleFactor << "\n";
os << indent << "Capping: " << (this->Capping ? "On\n" : "Off\n");
os << indent << "Cap Value: " << this->CapValue << "\n";
}
}
void vlGaussianSplatter::SetModelBounds(float *bounds)
{
vlGaussianSplatter::SetModelBounds(bounds[0], bounds[1], bounds[2], bounds[3], bounds[4], bounds[5]);
}
void vlGaussianSplatter::SetModelBounds(float xmin, float xmax, float ymin, float ymax, float zmin, float zmax)
{
if (this->ModelBounds[0] != xmin || this->ModelBounds[1] != xmax ||
this->ModelBounds[2] != ymin || this->ModelBounds[3] != ymax ||
this->ModelBounds[4] != zmin || this->ModelBounds[5] != zmax )
{
float length;
this->Modified();
this->ModelBounds[0] = xmin;
this->ModelBounds[1] = xmax;
this->ModelBounds[2] = ymin;
this->ModelBounds[3] = ymax;
this->ModelBounds[4] = zmin;
this->ModelBounds[5] = zmax;
this->Origin[0] = xmin;
this->Origin[1] = ymin;
this->Origin[2] = zmin;
if ( (length = xmin - xmax) == 0.0 ) length = 1.0;
this->AspectRatio[0] = 1.0;
this->AspectRatio[1] = (ymax - ymin) / length;
this->AspectRatio[2] = (zmax - zmin) / length;
}
}
//
// Static variables aid recursion
//
static vlFloatScalars *NewScalars;
static float Radius2;
static int Dim[3];
static float (vlGaussianSplatter::*Sample)(float x[3]);
static float (vlGaussianSplatter::*SampleFactor)(float s);
static char *Visited;
static float Eccentricity2;
static float *P, *N, S;
void vlGaussianSplatter::Execute()
{
int numSplatPts, numPts;
int ptId, i, j, k;
vlPointData *pd;
vlNormals *inNormals;
vlScalars *inScalars;
int loc[3], ip, jp, kp, idir, jdir, kdir;
vlDebugMacro(<< "Splatting data");
//
// Initialize self; create output objects
//
this->Initialize();
//
// Make sure points are available
//
if ( (numSplatPts=this->Input->GetNumberOfPoints()) < 1 )
{
vlErrorMacro(<<"No points to splat!");
return;
}
//
// Compute the radius of influence of the points. If an
// automatically generated bounding box has been generated, increase
// its size slightly to acoomodate the radius of influence.
//
Eccentricity2 = this->Eccentricity * this->Eccentricity;
numPts = this->SampleDimensions[0] * this->SampleDimensions[1]
* this->SampleDimensions[2];
NewScalars = new vlFloatScalars(numPts);
for (i=0; i<numPts; i++) NewScalars->SetScalar(i,0.0);
Visited = new char[numPts];
for (i=0; i < numPts; i++) Visited[i] = 0;
this->SetDimensions(this->GetSampleDimensions());
this->ComputeModelBounds();
//
// Set up proper function pointers
//
pd = this->Input->GetPointData();
if ( this->NormalWarping && (inNormals=pd->GetNormals()) != NULL )
Sample = this->Gaussian;
else
Sample = this->EccentricGaussian;
if ( this->ScalarWarping && (inScalars=pd->GetScalars()) != NULL )
SampleFactor = this->ScalarSampling;
else
SampleFactor = this->PositionSampling;
//
// Traverse all points - injecting into volume.
// For each input point, determine which cell it is in. Then start
// the recursive distribution of sampling function.
//
for (ptId=0; ptId < this->Input->GetNumberOfPoints(); ptId++)
{
P = this->Input->GetPoint(ptId);
if ( inNormals != NULL ) N = inNormals->GetNormal(ptId);
if ( inScalars != NULL ) S = inScalars->GetScalar(ptId);
if ( ! (ptId % 5000) && ptId > 0 )
vlDebugMacro(<< "Vertex #" << ptId);
for (i=0; i<3; i++)
loc[i] = (P[i] - this->Origin[i]) / this->AspectRatio[i];
//
// For each of the eight corners of the cell, need to evaluate sample
// function and then begin recursive distribution
//
for (i=0; i < 2; i++)
{
for (j=0; j < 2; j++)
{
for (k=0; k < 2; k++)
{
if ( (ip=loc[0]+i) >= 0 && ip < this->SampleDimensions[0] &&
(jp=loc[1]+j) >= 0 && jp < this->SampleDimensions[1] &&
(kp=loc[2]+k) >= 0 && kp < this->SampleDimensions[2] )
{
idir = (i==0 ? -1 : 1);
jdir = (j==0 ? -1 : 1);
kdir = (k==0 ? -1 : 1);
SplitIJK (ip, idir, jp, jdir, kp, kdir);
}
}
}
}
}
//
// If capping is turned on, set the distances of the outside of the volume
// to the CapValue.
//
if ( this->Capping )
{
this->Cap(NewScalars);
}
//
// Update self
//
delete [] Visited;
this->PointData.SetScalars(NewScalars);
}
void vlGaussianSplatter::ComputeModelBounds()
{
float *bounds, maxDist;
int i, adjustBounds=0;
// compute model bounds if not set previously
if ( this->ModelBounds[0] >= this->ModelBounds[1] ||
this->ModelBounds[2] >= this->ModelBounds[3] ||
this->ModelBounds[4] >= this->ModelBounds[5] )
{
adjustBounds = 1;
bounds = this->Input->GetBounds();
}
else
{
bounds = this->ModelBounds;
}
for (maxDist=0.0, i=0; i<3; i++)
if ( (bounds[2*i+1] - bounds[2*i]) > maxDist )
maxDist = bounds[2*i+1] - bounds[2*i];
maxDist *= this->Radius;
// adjust bounds so model fits strictly inside (only if not set previously)
if ( adjustBounds )
{
for (i=0; i<3; i++)
{
this->ModelBounds[2*i] = bounds[2*i] - maxDist;
this->ModelBounds[2*i+1] = bounds[2*i+1] + maxDist;
}
}
// Set volume origin and aspect ratio
for (i=0; i<3; i++)
{
this->Origin[i] = this->ModelBounds[2*i];
this->AspectRatio[i] = (this->ModelBounds[2*i+1] - this->ModelBounds[2*i])
/ (this->SampleDimensions[i] - 1);
}
}
void vlGaussianSplatter::SetSampleDimensions(int i, int j, int k)
{
int dim[3];
dim[0] = i;
dim[1] = j;
dim[2] = k;
this->SetSampleDimensions(dim);
}
void vlGaussianSplatter::SetSampleDimensions(int dim[3])
{
vlDebugMacro(<< " setting SampleDimensions to (" << dim[0] << "," << dim[1] << "," << dim[2] << ")");
if ( dim[0] != this->SampleDimensions[0] || dim[1] != SampleDimensions[1] ||
dim[2] != SampleDimensions[2] )
{
if ( dim[0]<1 || dim[1]<1 || dim[2]<1 )
{
vlErrorMacro (<< "Bad Sample Dimensions, retaining previous values");
return;
}
for (int dataDim=0, i=0; i<3 ; i++) if (dim[i] > 1) dataDim++;
if ( dataDim < 3 )
{
vlErrorMacro(<<"Sample dimensions must define a volume!");
return;
}
for ( i=0; i<3; i++) this->SampleDimensions[i] = dim[i];
this->Modified();
}
}
void vlGaussianSplatter::Cap(vlFloatScalars *s)
{
int i,j,k;
int idx;
int d01=this->SampleDimensions[0]*this->SampleDimensions[1];
// i-j planes
k = 0;
for (j=0; j<this->SampleDimensions[1]; j++)
for (i=0; i<this->SampleDimensions[0]; i++)
s->SetScalar(i+j*this->SampleDimensions[1], this->CapValue);
k = this->SampleDimensions[2] - 1;
idx = k*d01;
for (j=0; j<this->SampleDimensions[1]; j++)
for (i=0; i<this->SampleDimensions[0]; i++)
s->SetScalar(idx+i+j*this->SampleDimensions[1], this->CapValue);
// j-k planes
i = 0;
for (k=0; k<this->SampleDimensions[2]; k++)
for (j=0; j<this->SampleDimensions[1]; j++)
s->SetScalar(j*this->SampleDimensions[0]+k*d01, this->CapValue);
i = this->SampleDimensions[0] - 1;
for (k=0; k<this->SampleDimensions[2]; k++)
for (j=0; j<this->SampleDimensions[1]; j++)
s->SetScalar(i+j*this->SampleDimensions[0]+k*d01, this->CapValue);
// i-k planes
j = 0;
for (k=0; k<this->SampleDimensions[2]; k++)
for (i=0; i<this->SampleDimensions[0]; i++)
s->SetScalar(i+k*d01, this->CapValue);
j = this->SampleDimensions[1] - 1;
idx = j*this->SampleDimensions[0];
for (k=0; k<this->SampleDimensions[2]; k++)
for (i=0; i<this->SampleDimensions[0]; i++)
s->SetScalar(idx+i+k*d01, this->CapValue);
}
void vlGaussianSplatter::SplitIJK (int i, int idir, int j, int jdir,
int k, int kdir)
{
int idx, ip, jp, kp;
float cx[3], dist2;
cx[0] = this->Origin[0] + this->AspectRatio[0]*i;
cx[1] = this->Origin[1] + this->AspectRatio[1]*j;
cx[2] = this->Origin[2] + this->AspectRatio[2]*k;
if ( (dist2=Sample(cx)) <= Radius2 )
{
idx = i + j*this->SampleDimensions[0] +
k*this->SampleDimensions[0]*this->SampleDimensions[1];
this->SetScalar(idx,dist2);
//
// Continue sampling on opposite cell vertex, cell walls that emanate
// from this vertex, and cell edges that emanate from this vertex.
//
ip = i + idir;
jp = j + jdir;
kp = k + kdir;
if ( ip >= 0 && ip < this->SampleDimensions[0] &&
jp >= 0 && jp < this->SampleDimensions[1] &&
kp >= 0 && kp < this->SampleDimensions[2] )
SplitIJK (ip, idir, jp, jdir, kp, kdir);
//
// Don't forget cell walls emanating from this vertex
//
if ( ip >= 0 && ip < this->SampleDimensions[0] &&
jp >= 0 && jp < this->SampleDimensions[1] )
SplitIJ (ip, idir, jp, jdir, k);
if ( jp >= 0 && jp < this->SampleDimensions[1] &&
kp >= 0 && kp < this->SampleDimensions[2] )
SplitJK (i, jp, jdir, kp, kdir);
if ( ip >= 0 && ip < this->SampleDimensions[0] &&
kp >= 0 && kp < this->SampleDimensions[2] )
SplitIK (ip, idir, j, kp, kdir);
//
// Don't forget the cell edges emanating from this vertex
//
if ( ip >= 0 && ip < this->SampleDimensions[0] )
SplitI (ip, idir, j, k);
if ( jp >= 0 && jp < this->SampleDimensions[1] )
SplitJ (i, jp, jdir, k);
if ( kp >= 0 && kp < this->SampleDimensions[2] )
SplitK (i, j, kp, kdir);
}
}
void vlGaussianSplatter::SplitIJ (int i, int idir, int j, int jdir, int k)
{
int idx, ip, jp, kp;
float cx[3], dist2;
void SplitI(), SplitJ();
cx[0] = this->Origin[0] + this->AspectRatio[0]*i;
cx[1] = this->Origin[1] + this->AspectRatio[1]*j;
cx[2] = this->Origin[2] + this->AspectRatio[2]*k;
if ( (dist2=Sample(cx)) <= Radius2 )
{
idx = i + j*this->SampleDimensions[0] +
k*this->SampleDimensions[0]*this->SampleDimensions[1];
this->SetScalar(idx,dist2);
//
// Continue sampling on opposite cell vertex and cell edges that
// emanate from this vertex.
//
ip = i + idir;
jp = j + jdir;
if ( ip >= 0 && ip < this->SampleDimensions[0] &&
jp >= 0 && jp < this->SampleDimensions[1] )
SplitIJ (ip, idir, jp, jdir, k);
//
// Don't forget the cell edges emanating from this vertex
//
if ( ip >= 0 && ip < this->SampleDimensions[0] )
SplitI (ip, idir, j, k);
if ( jp >= 0 && jp < this->SampleDimensions[1] )
SplitJ (i, jp, jdir, k);
}
}
void vlGaussianSplatter::SplitJK (int i, int j, int jdir, int k, int kdir)
{
int idx, ip, jp, kp;
float cx[3], dist2;
void SplitJ(), SplitK();
cx[0] = this->Origin[0] + this->AspectRatio[0]*i;
cx[1] = this->Origin[1] + this->AspectRatio[1]*j;
cx[2] = this->Origin[2] + this->AspectRatio[2]*k;
if ( (dist2=Sample(cx)) <= Radius2 )
{
idx = i + j*this->SampleDimensions[0] +
k*this->SampleDimensions[0]*this->SampleDimensions[1];
this->SetScalar(idx,dist2);
//
// Continue sampling on opposite cell vertex and cell edges that
// emanate from this vertex.
//
jp = j + jdir;
kp = k + kdir;
if ( jp >= 0 && jp < this->SampleDimensions[1] &&
kp >= 0 && kp < this->SampleDimensions[2] )
SplitJK (i, jp, jdir, kp, kdir);
//
// Don't forget the cell edges emanating from this vertex
//
if ( jp >= 0 && jp < this->SampleDimensions[1] )
SplitJ (i, jp, jdir, k);
if ( kp >= 0 && kp < this->SampleDimensions[2] )
SplitK (i, j, kp, kdir);
}
}
void vlGaussianSplatter::SplitIK (int i, int idir, int j, int k, int kdir)
{
int idx, ip, jp, kp;
float cx[3], dist2;
void SplitI(), SplitK();
cx[0] = this->Origin[0] + this->AspectRatio[0]*i;
cx[1] = this->Origin[1] + this->AspectRatio[1]*j;
cx[2] = this->Origin[2] + this->AspectRatio[2]*k;
if ( (dist2=Sample(cx)) <= Radius2 )
{