/*=========================================================================

  Program:   Visualization Toolkit
  Module:    vtkGaussianSplatter.h

  Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
  All rights reserved.
  See Copyright.txt or http://www.kitware.com/Copyright.htm for details.

     This software is distributed WITHOUT ANY WARRANTY; without even
     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
     PURPOSE.  See the above copyright notice for more information.

=========================================================================*/
// .NAME vtkGaussianSplatter - splat points into a volume with an elliptical, Gaussian distribution
// .SECTION Description
// vtkGaussianSplatter is a filter that injects input points into a
// structured points (volume) dataset. As each point is injected, it "splats"
// or distributes values to nearby voxels. Data is distributed using an
// elliptical, Gaussian distribution function. The distribution function is
// modified using scalar values (expands distribution) or normals
// (creates ellipsoidal distribution rather than spherical).
//
// In general, the Gaussian distribution function f(x) around a given
// splat point p is given by
//
//     f(x) = ScaleFactor * exp( ExponentFactor*((r/Radius)**2) )
//
// where x is the current voxel sample point; r is the distance |x-p|
// ExponentFactor <= 0.0, and ScaleFactor can be multiplied by the scalar 
// value of the point p that is currently being splatted.
//
// If points normals are present (and NormalWarping is on), then the splat 
// function becomes elliptical (as compared to the spherical one described
// by the previous equation). The Gaussian distribution function then
// becomes:
// 
//     f(x) = ScaleFactor * 
//               exp( ExponentFactor*( ((rxy/E)**2 + z**2)/R**2) )
//
// where E is a user-defined eccentricity factor that controls the elliptical
// shape of the splat; z is the distance of the current voxel sample point
// along normal N; and rxy is the distance of x in the direction
// prependicular to N.
//
// This class is typically used to convert point-valued distributions into
// a volume representation. The volume is then usually iso-surfaced or
// volume rendered to generate a visualization. It can be used to create
// surfaces from point distributions, or to create structure (i.e.,
// topology) when none exists.

// .SECTION Caveats
// The input to this filter is any dataset type. This filter can be used 
// to resample any form of data, i.e., the input data need not be 
// unstructured. 
//
// Some voxels may never receive a contribution during the splatting process.
// The final value of these points can be specified with the "NullValue" 
// instance variable.

// .SECTION See Also
// vtkShepardMethod

#ifndef __vtkGaussianSplatter_h
#define __vtkGaussianSplatter_h

#include "vtkImageAlgorithm.h"

#define VTK_ACCUMULATION_MODE_MIN 0
#define VTK_ACCUMULATION_MODE_MAX 1
#define VTK_ACCUMULATION_MODE_SUM 2

class vtkDoubleArray;

class VTK_IMAGING_EXPORT vtkGaussianSplatter : public vtkImageAlgorithm 
{
public:
  vtkTypeMacro(vtkGaussianSplatter,vtkImageAlgorithm);
  void PrintSelf(ostream& os, vtkIndent indent);

  // Description:
  // Construct object with dimensions=(50,50,50); automatic computation of 
  // bounds; a splat radius of 0.1; an exponent factor of -5; and normal and 
  // scalar warping turned on.
  static vtkGaussianSplatter *New();

  // Description:
  // Set / get the dimensions of the sampling structured point set. Higher
  // values produce better results but are much slower.
  void SetSampleDimensions(int i, int j, int k);
  void SetSampleDimensions(int dim[3]);
  vtkGetVectorMacro(SampleDimensions,int,3);

  // Description:
  // Set / get the (xmin,xmax, ymin,ymax, zmin,zmax) bounding box in which
  // the sampling is performed. If any of the (min,max) bounds values are
  // min >= max, then the bounds will be computed automatically from the input
  // data. Otherwise, the user-specified bounds will be used.
  vtkSetVector6Macro(ModelBounds,double);
  vtkGetVectorMacro(ModelBounds,double,6);

  // Description:
  // Set / get the radius of propagation of the splat. This value is expressed
  // as a percentage of the length of the longest side of the sampling
  // volume. Smaller numbers greatly reduce execution time.
  vtkSetClampMacro(Radius,double,0.0,1.0);
  vtkGetMacro(Radius,double);

  // Description:
  // Multiply Gaussian splat distribution by this value. If ScalarWarping
  // is on, then the Scalar value will be multiplied by the ScaleFactor
  // times the Gaussian function.
  vtkSetClampMacro(ScaleFactor,double,0.0,VTK_DOUBLE_MAX);
  vtkGetMacro(ScaleFactor,double);

  // Description:
  // Set / get the sharpness of decay of the splats. This is the
  // exponent constant in the Gaussian equation. Normally this is
  // a negative value.
  vtkSetMacro(ExponentFactor,double);
  vtkGetMacro(ExponentFactor,double);

  // Description:
  // Turn on/off the generation of elliptical splats. If normal warping is
  // on, then the input normals affect the distribution of the splat. This
  // boolean is used in combination with the Eccentricity ivar.
  vtkSetMacro(NormalWarping,int);
  vtkGetMacro(NormalWarping,int);
  vtkBooleanMacro(NormalWarping,int);

  // Description:
  // Control the shape of elliptical splatting. Eccentricity is the ratio
  // of the major axis (aligned along normal) to the minor (axes) aligned
  // along other two axes. So Eccentricity > 1 creates needles with the
  // long axis in the direction of the normal; Eccentricity<1 creates
  // pancakes perpendicular to the normal vector.
  vtkSetClampMacro(Eccentricity,double,0.001,VTK_DOUBLE_MAX);
  vtkGetMacro(Eccentricity,double);

  // Description:
  // Turn on/off the scaling of splats by scalar value.
  vtkSetMacro(ScalarWarping,int);
  vtkGetMacro(ScalarWarping,int);
  vtkBooleanMacro(ScalarWarping,int);

  // Description:
  // Turn on/off the capping of the outer boundary of the volume
  // to a specified cap value. This can be used to close surfaces
  // (after iso-surfacing) and create other effects.
  vtkSetMacro(Capping,int);
  vtkGetMacro(Capping,int);
  vtkBooleanMacro(Capping,int);
  
  // Description:
  // Specify the cap value to use. (This instance variable only has effect 
  // if the ivar Capping is on.)
  vtkSetMacro(CapValue,double);
  vtkGetMacro(CapValue,double);

  // Description:
  // Specify the scalar accumulation mode. This mode expresses how scalar
  // values are combined when splats are overlapped. The Max mode acts
  // like a set union operation and is the most commonly used; the Min
  // mode acts like a set intersection, and the sum is just weird.
  vtkSetClampMacro(AccumulationMode,int,
                   VTK_ACCUMULATION_MODE_MIN,VTK_ACCUMULATION_MODE_SUM);
  vtkGetMacro(AccumulationMode,int);
  void SetAccumulationModeToMin()
    {this->SetAccumulationMode(VTK_ACCUMULATION_MODE_MIN);}
  void SetAccumulationModeToMax()
    {this->SetAccumulationMode(VTK_ACCUMULATION_MODE_MAX);}
  void SetAccumulationModeToSum()
    {this->SetAccumulationMode(VTK_ACCUMULATION_MODE_SUM);}
  const char *GetAccumulationModeAsString();

  // Description:
  // Set the Null value for output points not receiving a contribution from the
  // input points. (This is the initial value of the voxel samples.)
  vtkSetMacro(NullValue,double);
  vtkGetMacro(NullValue,double);

  // Description:
  // Compute the size of the sample bounding box automatically from the
  // input data. This is an internal helper function.
  void ComputeModelBounds(vtkDataSet *input, vtkImageData *output,
                          vtkInformation *outInfo);

protected:
  vtkGaussianSplatter();
  ~vtkGaussianSplatter() {};

  virtual int FillInputPortInformation(int port, vtkInformation* info);
  virtual int RequestInformation (vtkInformation *, 
                                  vtkInformationVector **, 
                                  vtkInformationVector *);
  virtual int RequestData(vtkInformation *, 
                          vtkInformationVector **, 
                          vtkInformationVector *);
  void Cap(vtkDoubleArray *s);

  int SampleDimensions[3]; // dimensions of volume to splat into
  double Radius; // maximum distance splat propagates (as fraction 0->1)
  double ExponentFactor; // scale exponent of gaussian function
  double ModelBounds[6]; // bounding box of splatting dimensions
  int NormalWarping; // on/off warping of splat via normal
  double Eccentricity;// elliptic distortion due to normals
  int ScalarWarping; // on/off warping of splat via scalar
  double ScaleFactor; // splat size influenced by scale factor
  int Capping; // Cap side of volume to close surfaces
  double CapValue; // value to use for capping
  int AccumulationMode; // how to combine scalar values

  double Gaussian(double x[3]);  
  double EccentricGaussian(double x[3]);  
  double ScalarSampling(double s) 
    {return this->ScaleFactor * s;}
  double PositionSampling(double) 
    {return this->ScaleFactor;}
  void SetScalar(int idx, double dist2, vtkDoubleArray *newScalars);

//BTX
private:
  double Radius2;
  double (vtkGaussianSplatter::*Sample)(double x[3]);
  double (vtkGaussianSplatter::*SampleFactor)(double s);
  char *Visited;
  double Eccentricity2;
  double *P;
  double *N;
  double S;
  double Origin[3];
  double Spacing[3];
  double SplatDistance[3];
  double NullValue;
//ETX

private:
  vtkGaussianSplatter(const vtkGaussianSplatter&);  // Not implemented.
  void operator=(const vtkGaussianSplatter&);  // Not implemented.
};

#endif