ITK
6.0.0
Insight Toolkit
|
#include <itkBSplineScatteredDataPointSetToImageFilter.h>
Image filter which provides a B-spline output approximation.
Given an n-D image with scattered data, this filter finds a fast approximation to that irregularly spaced data using uniform B-splines. The traditional method of inverting the observation matrix to find a least-squares fit is made obsolete. Therefore, memory issues are not a concern and inverting large matrices is not applicable. In addition, this allows fitting to be multi-threaded. This class generalizes from Lee's original paper to encompass n-D data in m-D parametric space and any feasible B-spline order as well as the option of specifying a confidence value for each point.
In addition to specifying the input point set, one must specify the number of control points. The specified number of control points must be greater than m_SplineOrder. If one wishes to use the multilevel component of this algorithm, one must also specify the number of levels in the hierarchy. If this is desired, the number of control points becomes the number of control points for the coarsest level. The algorithm then increases the number of control points at each level so that the B-spline n-D grid is refined to twice the previous level.
There are two parts to fitting scattered data: the parameterization assignment problem and the fitting problem given a parameterization. This filter only addresses the second problem in that the user must provide a parametric value for each scattered datum. Different parametric assignment schemes result in different B-spline object outputs.
This filter is general in that it accepts n-D scattered data in m-D parametric dimensions. Input to this filter is an m-D point set with a Vector data type of n dimensions. This means that the parametric values are stored in the points container of the point set whereas the scattered data are stored in the points data container of the point set.
Typical B-spline objects include curves, which have a parametric dimension of 1 and a data dimension of 2 or 3 (depending on the space in which the curve resides) and deformation fields which commonly have parametric and data dimensions of 2 or 3 (again depending on the space of the field). As an example, a curve through a set of 2D points has data dimension 2 and parametric dimension 1. The univariate curve could be represented as: <x(u),y(u)> Another example is a 3D deformation of 3D points, which has parametric dimension 3 and data dimension 3 and can be represented as: <dx(u,v,w), dy(u,v,w), dz(u,v,w)>. However, as mentioned before, the code is general such that, if the user wanted, she could model a time varying 3-D displacement field which resides in 4-D space as <dx(u, v, w, t), dy(u, v, w, t), dz(u, v, w, t)>.
The output is an image defining the sampled B-spline parametric domain where each pixel houses the sampled B-spline object value. For a curve fit to 3-D points, the output is a 1-D image where each voxel contains a vector with the approximated (x,y,z) location. The continuous, finite, rectilinear domain (as well as the sampling rate) is specified via the combination of the SetSpacing() and SetSize() functions. For a 2-D deformation on 2-D points, the output is a 2-D image where each voxel contains the approximated (dx, dy) vector.
The parameterization must be specified using SetPoint, where the actual coordinates of the point are set via SetPointData. For example, to compute a spline through the (ordered) 2D points (5,6) and (7,8), you should use:
This code was contributed in the Insight Journal paper: "N-D C^k B-Spline Scattered Data Approximation" by Nicholas J. Tustison, James C. Gee https://doi.org/10.54294/0d55to
Definition at line 132 of file itkBSplineScatteredDataPointSetToImageFilter.h.
Static Public Member Functions | |
static Pointer | New () |
Static Public Member Functions inherited from itk::PointSetToImageFilter< TInputPointSet, TOutputImage > | |
static Pointer | New () |
Static Public Member Functions inherited from itk::Object | |
static bool | GetGlobalWarningDisplay () |
static void | GlobalWarningDisplayOff () |
static void | GlobalWarningDisplayOn () |
static Pointer | New () |
static void | SetGlobalWarningDisplay (bool val) |
Static Public Member Functions inherited from itk::LightObject | |
static void | BreakOnError () |
static Pointer | New () |
Static Public Attributes | |
static constexpr unsigned int | ImageDimension = TOutputImage::ImageDimension |
Static Public Attributes inherited from itk::PointSetToImageFilter< TInputPointSet, TOutputImage > | |
static constexpr unsigned int | InputPointSetDimension = InputPointSetType::PointDimension |
static constexpr unsigned int | OutputImageDimension = TOutputImage::ImageDimension |
Static Public Attributes inherited from itk::ImageSource< TOutputImage > | |
static constexpr unsigned int | OutputImageDimension = TOutputImage::ImageDimension |
Private Member Functions | |
void | CollapsePhiLattice (PointDataImageType *, PointDataImageType *, const RealType, const unsigned int) |
void | GenerateOutputImage () |
IndexType | NumberToIndex (const unsigned int, const SizeType) |
void | RefineControlPointLattice () |
void | SetPhiLatticeParametricDomainParameters () |
void | ThreadedGenerateDataForFitting (const RegionType &, ThreadIdType) |
void | ThreadedGenerateDataForReconstruction (const RegionType &, ThreadIdType) |
void | ThreadedGenerateDataForUpdatingResidualValues (const RegionType &, ThreadIdType) |
Additional Inherited Members | |
Static Protected Member Functions inherited from itk::ImageSource< TOutputImage > | |
static const ImageRegionSplitterBase * | GetGlobalDefaultSplitter () |
static ITK_THREAD_RETURN_FUNCTION_CALL_CONVENTION | ThreaderCallback (void *arg) |
Static Protected Member Functions inherited from itk::ProcessObject | |
template<typename TSourceObject > | |
static void | MakeRequiredOutputs (TSourceObject &sourceObject, const DataObjectPointerArraySizeType numberOfRequiredOutputs) |
static constexpr float | progressFixedToFloat (uint32_t fixed) |
static uint32_t | progressFloatToFixed (float f) |
Protected Attributes inherited from itk::PointSetToImageFilter< TInputPointSet, TOutputImage > | |
DirectionType | m_Direction {} |
ValueType | m_InsideValue {} |
PointType | m_Origin {} |
ValueType | m_OutsideValue {} |
SizeType | m_Size {} |
SpacingType | m_Spacing {} |
Protected Attributes inherited from itk::ImageSource< TOutputImage > | |
bool | m_DynamicMultiThreading { true } |
Protected Attributes inherited from itk::ProcessObject | |
TimeStamp | m_OutputInformationMTime {} |
bool | m_Updating {} |
Protected Attributes inherited from itk::LightObject | |
std::atomic< int > | m_ReferenceCount {} |
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::ArrayType = FixedArray<unsigned int, Self::ImageDimension> |
Definition at line 179 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::ConstPointer = SmartPointer<const Self> |
Definition at line 142 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::ImageType = TOutputImage |
Definition at line 153 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::IndexType = typename ImageType::IndexType |
Definition at line 161 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::KernelOrder0Type = BSplineKernelFunction<0> |
Definition at line 184 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::KernelOrder1Type = BSplineKernelFunction<1> |
Definition at line 185 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::KernelOrder2Type = BSplineKernelFunction<2> |
Definition at line 186 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::KernelOrder3Type = BSplineKernelFunction<3> |
Definition at line 187 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::KernelType = CoxDeBoorBSplineKernelFunction<3> |
Interpolation kernel type (default spline order = 3).
Definition at line 183 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::OutputImageRegionType = RegionType |
Definition at line 159 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::PixelType = typename ImageType::PixelType |
Image type alias support
Definition at line 157 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::PointDataContainerPointer = typename PointDataContainerType::Pointer |
Definition at line 168 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::PointDataContainerType = typename PointSetType::PointDataContainer |
Definition at line 167 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::PointDataImagePointer = typename PointDataImageType::Pointer |
Definition at line 178 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::PointDataImageType = Image<PointDataType, Self::ImageDimension> |
Image types.
Definition at line 175 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::PointDataType = typename PointSetType::PixelType |
Definition at line 166 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::Pointer = SmartPointer<Self> |
Definition at line 141 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::PointSetPointer = typename PointSetType::Pointer |
Definition at line 165 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::PointSetType = TInputPointSet |
Definition at line 154 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::PointType = typename PointSetType::PointType |
PointSet type alias support
Definition at line 164 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::RealArrayType = FixedArray<RealType, Self::ImageDimension> |
Definition at line 180 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::RealImagePointer = typename RealImageType::Pointer |
Definition at line 177 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::RealImageType = Image<RealType, Self::ImageDimension> |
Definition at line 176 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::RealType = float |
Other type alias.
Definition at line 171 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::RegionType = typename ImageType::RegionType |
Definition at line 158 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::Self = BSplineScatteredDataPointSetToImageFilter |
Standard class type aliases.
Definition at line 139 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::SizeType = typename ImageType::SizeType |
Definition at line 160 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::Superclass = PointSetToImageFilter<TInputPointSet, TOutputImage> |
Definition at line 140 of file itkBSplineScatteredDataPointSetToImageFilter.h.
using itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::WeightsContainerType = VectorContainer<unsigned int, RealType> |
Definition at line 172 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
protected |
|
overrideprotecteddefault |
|
overrideprotectedvirtual |
If an imaging filter needs to perform processing after all processing threads have completed, the filter can can provide an implementation for AfterThreadedGenerateData(). The execution flow in the default GenerateData() method will be: 1) Allocate the output buffer 2) Call BeforeThreadedGenerateData() 3) Spawn threads, calling ThreadedGenerateData() in each thread. 4) Call AfterThreadedGenerateData() Note that this flow of control is only available if a filter provides a ThreadedGenerateData() method and NOT a GenerateData() method.
Reimplemented from itk::ImageSource< TOutputImage >.
|
overrideprotectedvirtual |
If an imaging filter needs to perform processing after the buffer has been allocated but before threads are spawned, the filter can can provide an implementation for BeforeThreadedGenerateData(). The execution flow in the default GenerateData() method will be: 1) Allocate the output buffer 2) Call BeforeThreadedGenerateData() 3) Spawn threads, calling ThreadedGenerateData() in each thread. 4) Call AfterThreadedGenerateData() Note that this flow of control is only available if a filter provides a ThreadedGenerateData() method and NOT a GenerateData() method.
Reimplemented from itk::ImageSource< TOutputImage >.
|
private |
Sub-function used by GenerateOutputImageFast() to generate the sampled B-spline object quickly.
|
inlineoverrideprotected |
Definition at line 303 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
overrideprotectedvirtual |
A version of GenerateData() specific for image processing filters. This implementation will split the processing across multiple threads. The buffer is allocated by this method. Then the BeforeThreadedGenerateData() method is called (if provided). Then, a series of threads are spawned each calling DynamicThreadedGenerateData(). After all the threads have completed processing, the AfterThreadedGenerateData() method is called (if provided). If an image processing filter cannot be threaded, the filter should provide an implementation of GenerateData(). That implementation is responsible for allocating the output buffer. If a filter can be threaded, it should NOT provide a GenerateData() method but should provide a DynamicThreadedGenerateData() instead.
Reimplemented from itk::PointSetToImageFilter< TInputPointSet, TOutputImage >.
|
private |
This function is not used as it requires an evaluation of all (SplineOrder+1)^ImageDimensions B-spline weights for each evaluation.
|
virtual |
Set/Get whether or not the sampled output B-spline object is constructed. The result of the fitting process is an n-D grid of control points which describe the continuous B-spline object.
|
virtual |
Set/Get the epsilon used for B-splines. The B-spline parametric domain in 1-D is defined on the half-closed interval [a,b). Extension to n-D is defined similarly. This presents some difficulty for defining the the image domain to be co-extensive with the parametric domain. We use the B-spline epsilon to push the edge of the image boundary inside the B-spline parametric domain.
|
virtual |
Set/Get the array to define the periodicity of the dimensions in the parametric space is to be. This array of 0/1 values defines whether a particular dimension of the parametric space is to be considered periodic or not. For example, if you are using interpolating along a 1D closed curve, the array type will have size 1, and you should set the first element of this array to the value "1". In the case that you were interpolating in a planar surface with cylindrical topology, the array type will have two components, and you should set to "1" the component that goes around the cylinder, and set to "0" the component that goes from the top of the cylinder to the bottom. This will indicate the periodicity of that parameter to the filter. Internally, in order to make periodic the domain of the parameter, the filter will reuse some of the points at the beginning of the domain as if they were also located at the end of the domain. The number of points to be reused will depend on the spline order. As a user, you don't need to replicate the points, the filter will do this for you.
|
virtual |
Get the number of current control points for each parametric dimension at the current fitting level. The B-spline mesh size is equal to the number of control points minus the spline order. Default = 4 in each dimension.
|
virtual |
Set/Get whether or not the sampled output B-spline object is constructed. The result of the fitting process is an n-D grid of control points which describe the continuous B-spline object.
|
overridevirtual |
Reimplemented from itk::PointSetToImageFilter< TInputPointSet, TOutputImage >.
|
virtual |
Set/Get the number of control points for each parametric dimension at the initial fitting level. The B-spline mesh size is equal to the number of control points minus the spline order. Default = 4 in each dimension.
|
virtual |
Get the number of fitting levels for all parametric dimensions. Starting with the mesh size implied by setting the number of control points, the mesh size is doubled at each fitting level. Default = 1 in all parametric dimensions.
|
inline |
Get the control point lattice produced by the fitting process.
Definition at line 287 of file itkBSplineScatteredDataPointSetToImageFilter.h.
References itk::ProcessObject::GetOutput().
|
virtual |
Get the spline order for all parametric dimensions. The spline order determines the continuity between B-spline elements and the degree of polynomial used to construct the B-spline elements. Default = 3.
|
static |
Method for creation through the object factory.
|
private |
Convert number to index given a size of image. Used to index the local control point neighborhoods.
|
overrideprotectedvirtual |
Methods invoked by Print() to print information about the object including superclasses. Typically not called by the user (use Print() instead) but used in the hierarchical print process to combine the output of several classes.
Reimplemented from itk::PointSetToImageFilter< TInputPointSet, TOutputImage >.
|
private |
Function used to propagate the fitting solution at one fitting level to the next level with the mesh resolution doubled.
|
virtual |
Set/Get the epsilon used for B-splines. The B-spline parametric domain in 1-D is defined on the half-closed interval [a,b). Extension to n-D is defined similarly. This presents some difficulty for defining the the image domain to be co-extensive with the parametric domain. We use the B-spline epsilon to push the edge of the image boundary inside the B-spline parametric domain.
|
virtual |
Set/Get the array to define the periodicity of the dimensions in the parametric space is to be. This array of 0/1 values defines whether a particular dimension of the parametric space is to be considered periodic or not. For example, if you are using interpolating along a 1D closed curve, the array type will have size 1, and you should set the first element of this array to the value "1". In the case that you were interpolating in a planar surface with cylindrical topology, the array type will have two components, and you should set to "1" the component that goes around the cylinder, and set to "0" the component that goes from the top of the cylinder to the bottom. This will indicate the periodicity of that parameter to the filter. Internally, in order to make periodic the domain of the parameter, the filter will reuse some of the points at the beginning of the domain as if they were also located at the end of the domain. The number of points to be reused will depend on the spline order. As a user, you don't need to replicate the points, the filter will do this for you.
|
virtual |
Set/Get whether or not the sampled output B-spline object is constructed. The result of the fitting process is an n-D grid of control points which describe the continuous B-spline object.
|
virtual |
Set/Get the number of control points for each parametric dimension at the initial fitting level. The B-spline mesh size is equal to the number of control points minus the spline order. Default = 4 in each dimension.
void itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::SetNumberOfLevels | ( | const ArrayType & | ) |
Set the number of fitting levels in each parametric dimension separately. Starting with the mesh size implied by setting the number of control points, the mesh size is doubled at each fitting level. Default = 1 in all parametric dimensions.
void itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::SetNumberOfLevels | ( | unsigned int | ) |
Set the number of fitting levels assuming the number of fitting levels is the same for each parametric dimension. Starting with the mesh size implied by setting the number of control points, the mesh size is doubled at each fitting level. Default = 1 in all parametric dimensions.
|
private |
Set the grid parametric domain parameters such as the origin, size, spacing, and direction.
void itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::SetPointWeights | ( | WeightsContainerType * | weights | ) |
A weighted fitting is possible where each input point is assigned a relative weighting.
void itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::SetSplineOrder | ( | const ArrayType & | ) |
Set the spline order for each parametric dimension separately. The spline order determines the continuity between B-spline elements and the degree of polynomial used to construct the B-spline elements. Default = 3.
void itk::BSplineScatteredDataPointSetToImageFilter< TInputPointSet, TOutputImage >::SetSplineOrder | ( | unsigned int | ) |
Set the spline order assuming it is the same in all parametric dimensions. The spline order determines the continuity between B-spline elements and the degree of polynomial used to construct the B-spline elements. Default = 3.
|
overrideprotected |
|
overrideprotected |
|
private |
Function used to generate the sampled B-spline object quickly.
|
private |
Function used to generate the sampled B-spline object quickly.
|
private |
Update the residuals for multi-level fitting.
|
staticconstexpr |
Extract dimension from the output image.
Definition at line 151 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 388 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 365 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 362 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 364 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 386 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 358 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 390 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 359 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 389 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 378 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 380 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 381 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 382 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 383 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 361 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 363 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 367 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 385 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 371 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 369 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 372 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 374 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 376 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 366 of file itkBSplineScatteredDataPointSetToImageFilter.h.
|
private |
Definition at line 360 of file itkBSplineScatteredDataPointSetToImageFilter.h.