VTK
Static Public Member Functions | List of all members
vtkMomentsHelper Class Reference

rotation invariant pattern detetction More...

#include <vtkMomentsHelper.h>

Static Public Member Functions

static std::vector< vtkMomentsTensororthonormalizeMoments (int dimension, std::vector< vtkMomentsTensor > moments, double radius)
 The monomial basis is not orthonormal. More...
 
static std::vector< vtkMomentsTensorallMoments (int dimension, int order, int fieldRank, double radius, double center[3], vtkImageData *stencil, std::string nameOfPointData)
 This function computes the moments at a given location and radius the moments are the projections of the function to the monomial basis they are evaluated using a numerical integration over uniformly sampled 2D or 3D space. More...
 
static double getVolume (vtkCell *cell, vtkDataSet *source)
 This function approximates the volume of a cell for unstructured grids. More...
 
static std::vector< vtkMomentsTensorallMomentsOrigRes (int dimension, int order, int fieldRank, double radius, double center[3], vtkDataSet *dataset, std::string nameOfPointData)
 This function computes the moments at a given location and radius the moments are the projections of the function to the monomial basis they are evaluated using a numerical integration over the original unstructured dataset. More...
 
static std::vector< vtkMomentsTensorallMomentsOrigResImageData (int dimension, int order, int fieldRank, double radius, int ptID, vtkImageData *dataset, std::string nameOfPointData)
 This function computes the moments at a given location and radius the moments are the projections of the function to the monomial basis they are evaluated using a numerical integration over the original dataset if it is structured data. More...
 
static double translationFactorAnalytic (double radius, int dimension, size_t p, size_t q, size_t r)
 This function computes the factor that needs to be removed for the translational normalization it corresponds to the moment of the function identical to one for the lowest orders, we know the analytic solution. More...
 
static double translationFactor (double radius, size_t p, size_t q, size_t r, vtkImageData *stencil)
 This function computes the factor that needs to be removed for the translational normalization it corresponds to the moment of the function identical to one if we do not know the analytic solution, we evaluate it numerically. More...
 
static void BuildStencil (vtkImageData *stencil, double radius, int numberOfIntegrationSteps, int dimension, vtkDataSet *source, std::string nameOfPointData)
 This function generates the stencil, which contains the locations at which the dataset is evaluated for the integration. More...
 
static bool CenterStencil (double *center, vtkDataSet *source, vtkImageData *stencil, int numberOfIntegrationSteps, std::string nameOfPointData)
 This function moves the stencil to the current location, where the integration is supposed o be performed. More...
 
static std::vector< size_t > getTensorIndicesFromFieldIndex (size_t index, int dimension, int order, int fieldRank)
 Inverse function to getFieldIndexFromTensorIndices The output contains a vector with the tensor indices that describe the basis function that belongs to the given output array. More...
 
static std::string getTensorIndicesFromFieldIndexAsString (size_t index, int dimension, int order, int fieldRank)
 The output contains the tensor indices that describe the basis function that belongs to the given output array as string. More...
 
static size_t getFieldIndexFromTensorIndices (size_t radiusIndex, std::vector< size_t > indices, int dimension, int fieldRank, int numberOfBasisFunctions)
 Inverse function to getTensorIndicesFromFieldIndex given a vector with tensor indices and a radius, this function returns the index in the output of this algorithm that corresponds tothis basis function. More...
 
static bool isCloseToEdge (int dimension, int ptId, double radius, vtkImageData *field)
 the function returns true if the point lies within radius of the boundary of the dataset More...
 
static bool isEdge (int dimension, int ptId, vtkImageData *field)
 the function returns true if the point lies on the boundary of the dataset More...
 
static vtkIdType getArrayIndex (std::vector< int > coord, std::vector< int > dimensions)
 Calculates the index of the coordinate in the 1D array that is treated as dimensions[0] x dimensions[1] x dimensions[2] matrix. More...
 
static std::vector< intgetCoord (vtkIdType index, std::vector< int > dimensions)
 Calculates the coordinate as if we are in a dimensions[0] x dimensions[1] x dimensions[2] matrix based on the index of the 1D array. More...
 
static vtkImageDatatranslateToOrigin (vtkImageData *data)
 Translates the data to the origin (0, 0, 0) More...
 
static vtkImageDatapadField (vtkImageData *field, vtkImageData *kernel, int dimension, std::string nameOfPointData)
 Pad the field to a square where the size is max(field->Dimensions) + max(kernel->Dimensions) More...
 
static vtkImageDatapadKernel (vtkImageData *kernel, vtkImageData *paddedField)
 Pad the kernel to the same size as paddedField The center of the kernel is the origin of the final output and the rest is wrapped accordingly. More...
 

Detailed Description

rotation invariant pattern detetction

vtkMomentsHelper is a helper class that contains functions that will be used by more than one algorithm in the moments module the theory and the algorithm is described in Roxana Bujack and Hans Hagen: "Moment Invariants for Multi-Dimensional Data" http://www.informatik.uni-leipzig.de/~bujack/2017TensorDagstuhl.pdf

Thanks:
Developed by Roxana Bujack at Los Alamos National Laboratory.

Definition at line 87 of file vtkMomentsHelper.h.

Member Function Documentation

◆ orthonormalizeMoments()

static std::vector<vtkMomentsTensor> vtkMomentsHelper::orthonormalizeMoments ( int  dimension,
std::vector< vtkMomentsTensor moments,
double  radius 
)
static

The monomial basis is not orthonormal.

We need this function for the reconstruction of the function from the moments. This function uses Gram Schmidt

Parameters
dimension2D or 3D
momentsthe moments at a point
radiusthe corresponding integration radius
Returns
the orthonormal moments

◆ allMoments()

static std::vector<vtkMomentsTensor> vtkMomentsHelper::allMoments ( int  dimension,
int  order,
int  fieldRank,
double  radius,
double  center[3],
vtkImageData stencil,
std::string  nameOfPointData 
)
static

This function computes the moments at a given location and radius the moments are the projections of the function to the monomial basis they are evaluated using a numerical integration over uniformly sampled 2D or 3D space.

Parameters
dimension2D or 3D
orderthe maximal order up to which the moments are computed
fieldRank0 for scalar, 1 for vector and 2 for matrix
radiusthe integration radius at which the moments are computed
centerlocation where the moments are computed
stencilcontains the locations at which the dataset is evaluated for the integration
nameOfPointDatathe name of the array in the point data of which the momens are computed.
Returns
the moments

◆ getVolume()

static double vtkMomentsHelper::getVolume ( vtkCell cell,
vtkDataSet source 
)
static

This function approximates the volume of a cell for unstructured grids.

Parameters
cellthe cell
sourcethe dataset that contains the cell.
Returns
the volume

◆ allMomentsOrigRes()

static std::vector<vtkMomentsTensor> vtkMomentsHelper::allMomentsOrigRes ( int  dimension,
int  order,
int  fieldRank,
double  radius,
double  center[3],
vtkDataSet dataset,
std::string  nameOfPointData 
)
static

This function computes the moments at a given location and radius the moments are the projections of the function to the monomial basis they are evaluated using a numerical integration over the original unstructured dataset.

Parameters
dimension2D or 3D
orderthe maximal order up to which the moments are computed
fieldRank0 for scalar, 1 for vector and 2 for matrix
radiusthe integration radius at which the moments are computed
centerlocation where the moments are computed
datasetthe dataset of which the moments are computed
nameOfPointDatathe name of the array in the point data of which the momens are computed.
Returns
the moments

◆ allMomentsOrigResImageData()

static std::vector<vtkMomentsTensor> vtkMomentsHelper::allMomentsOrigResImageData ( int  dimension,
int  order,
int  fieldRank,
double  radius,
int  ptID,
vtkImageData dataset,
std::string  nameOfPointData 
)
static

This function computes the moments at a given location and radius the moments are the projections of the function to the monomial basis they are evaluated using a numerical integration over the original dataset if it is structured data.

Parameters
dimension2D or 3D
orderthe maximal order up to which the moments are computed
fieldRank0 for scalar, 1 for vector and 2 for matrix
radiusthe integration radius at which the moments are computed
ptIDpoint id of the location where the moments are computed
datasetthe dataset of which the moments are computed
nameOfPointDatathe name of the array in the point data of which the momens are computed.
Returns
the moments

◆ translationFactorAnalytic()

static double vtkMomentsHelper::translationFactorAnalytic ( double  radius,
int  dimension,
size_t  p,
size_t  q,
size_t  r 
)
static

This function computes the factor that needs to be removed for the translational normalization it corresponds to the moment of the function identical to one for the lowest orders, we know the analytic solution.

Parameters
radiusthe integration radius at which the moments are computed
dimension2D or 3D
pexponent in the basis function x^p*y^q*z^r
qexponent in the basis function x^p*y^q*z^r
rexponent in the basis function x^p*y^q*z^r
Returns
the translation factor

◆ translationFactor()

static double vtkMomentsHelper::translationFactor ( double  radius,
size_t  p,
size_t  q,
size_t  r,
vtkImageData stencil 
)
static

This function computes the factor that needs to be removed for the translational normalization it corresponds to the moment of the function identical to one if we do not know the analytic solution, we evaluate it numerically.

Parameters
centerlocation where the moments are computed
radiusthe integration radius at which the moments are computed
pexponent in the basis function x^p*y^q*z^r
qexponent in the basis function x^p*y^q*z^r
rexponent in the basis function x^p*y^q*z^r
stencilcontains the locations at which the dataset is evaluated for the integration
Returns
the translation factor

◆ BuildStencil()

static void vtkMomentsHelper::BuildStencil ( vtkImageData stencil,
double  radius,
int  numberOfIntegrationSteps,
int  dimension,
vtkDataSet source,
std::string  nameOfPointData 
)
static

This function generates the stencil, which contains the locations at which the dataset is evaluated for the integration.

Parameters
stencilcontains the locations at which the dataset is evaluated for the integration
radiusthe integration radius at which the moments are computed
numberOfIntegrationStepshow fine the discrete integration done in each dimension
dimension2D or 3D
sourcethe dataset
nameOfPointDatathe name of the array in the point data of which the momens are computed.
Returns
the moments

◆ CenterStencil()

static bool vtkMomentsHelper::CenterStencil ( double center,
vtkDataSet source,
vtkImageData stencil,
int  numberOfIntegrationSteps,
std::string  nameOfPointData 
)
static

This function moves the stencil to the current location, where the integration is supposed o be performed.

Parameters
centerthe location
sourcethe dataset
stencilcontains the locations at which the dataset is evaluated for the integration
numberOfIntegrationStepshow fine the discrete integration done in each dimension
Returns
the moments

◆ getTensorIndicesFromFieldIndex()

static std::vector<size_t> vtkMomentsHelper::getTensorIndicesFromFieldIndex ( size_t  index,
int  dimension,
int  order,
int  fieldRank 
)
static

Inverse function to getFieldIndexFromTensorIndices The output contains a vector with the tensor indices that describe the basis function that belongs to the given output array.

they are sorted by increasing order and then by the index as returned by vtkMomentsTensor.getIndices(i)

Parameters
indexthe index of this output field pointdata array
dimension2D or 3D
orderthe maximal order up to which the moments are computed
fieldRank0 for scalar, 1 for vector and 2 for matrix
Returns
the moments

◆ getTensorIndicesFromFieldIndexAsString()

static std::string vtkMomentsHelper::getTensorIndicesFromFieldIndexAsString ( size_t  index,
int  dimension,
int  order,
int  fieldRank 
)
static

The output contains the tensor indices that describe the basis function that belongs to the given output array as string.

Convenience function. they are sorted by increasing order and then by the index as returned by vtkMomentsTensor.getIndices(i)

Parameters
indexthe index of this output field pointdata array
dimension2D or 3D
orderthe maximal order up to which the moments are computed
fieldRank0 for scalar, 1 for vector and 2 for matrix
Returns
the moments

◆ getFieldIndexFromTensorIndices()

static size_t vtkMomentsHelper::getFieldIndexFromTensorIndices ( size_t  radiusIndex,
std::vector< size_t >  indices,
int  dimension,
int  fieldRank,
int  numberOfBasisFunctions 
)
static

Inverse function to getTensorIndicesFromFieldIndex given a vector with tensor indices and a radius, this function returns the index in the output of this algorithm that corresponds tothis basis function.

Parameters
radiusIndexindex of this radius in the radii vector
indicesthe given tensor indices
dimension2D or 3D
fieldRank0 for scalar, 1 for vector and 2 for matrix
numberOfBasisFunctionsnumber of basis functions in the source equals \sum_{i=0}^order dimension^o

◆ isCloseToEdge()

static bool vtkMomentsHelper::isCloseToEdge ( int  dimension,
int  ptId,
double  radius,
vtkImageData field 
)
static

the function returns true if the point lies within radius of the boundary of the dataset

Parameters
ptIdID of the point in question
fieldthe field that contains the point

◆ isEdge()

static bool vtkMomentsHelper::isEdge ( int  dimension,
int  ptId,
vtkImageData field 
)
static

the function returns true if the point lies on the boundary of the dataset

Parameters
ptIdID of the point in question
fieldthe field that contains the point

◆ getArrayIndex()

static vtkIdType vtkMomentsHelper::getArrayIndex ( std::vector< int coord,
std::vector< int dimensions 
)
static

Calculates the index of the coordinate in the 1D array that is treated as dimensions[0] x dimensions[1] x dimensions[2] matrix.

◆ getCoord()

static std::vector<int> vtkMomentsHelper::getCoord ( vtkIdType  index,
std::vector< int dimensions 
)
static

Calculates the coordinate as if we are in a dimensions[0] x dimensions[1] x dimensions[2] matrix based on the index of the 1D array.

◆ translateToOrigin()

static vtkImageData* vtkMomentsHelper::translateToOrigin ( vtkImageData data)
static

Translates the data to the origin (0, 0, 0)

◆ padField()

static vtkImageData* vtkMomentsHelper::padField ( vtkImageData field,
vtkImageData kernel,
int  dimension,
std::string  nameOfPointData 
)
static

Pad the field to a square where the size is max(field->Dimensions) + max(kernel->Dimensions)

◆ padKernel()

static vtkImageData* vtkMomentsHelper::padKernel ( vtkImageData kernel,
vtkImageData paddedField 
)
static

Pad the kernel to the same size as paddedField The center of the kernel is the origin of the final output and the rest is wrapped accordingly.


The documentation for this class was generated from the following file: