ITK  5.0.0
Insight Segmentation and Registration Toolkit
itkBSplineScatteredDataPointSetToImageFilter.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright Insight Software Consortium
4  *
5  * Licensed under the Apache License, Version 2.0 (the "License");
6  * you may not use this file except in compliance with the License.
7  * You may obtain a copy of the License at
8  *
9  * http://www.apache.org/licenses/LICENSE-2.0.txt
10  *
11  * Unless required by applicable law or agreed to in writing, software
12  * distributed under the License is distributed on an "AS IS" BASIS,
13  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14  * See the License for the specific language governing permissions and
15  * limitations under the License.
16  *
17  *=========================================================================*/
18 #ifndef itkBSplineScatteredDataPointSetToImageFilter_h
19 #define itkBSplineScatteredDataPointSetToImageFilter_h
20 
22 
25 #include "itkVectorContainer.h"
26 
27 #include "vnl/vnl_matrix.h"
28 
29 namespace itk
30 {
128 template< typename TInputPointSet, typename TOutputImage >
130  public PointSetToImageFilter< TInputPointSet, TOutputImage >
131 {
132 public:
133  ITK_DISALLOW_COPY_AND_ASSIGN(BSplineScatteredDataPointSetToImageFilter);
134 
140 
142  itkNewMacro( Self );
143 
146 
148  static constexpr unsigned int ImageDimension = TOutputImage::ImageDimension;
149 
150  using ImageType = TOutputImage;
151  using PointSetType = TInputPointSet;
152 
154  using PixelType = typename ImageType::PixelType;
157  using SizeType = typename ImageType::SizeType;
158  using IndexType = typename ImageType::IndexType;
159 
162  using PointSetPointer = typename PointSetType::Pointer;
163  using PointDataType = typename PointSetType::PixelType;
164  using PointDataContainerType = typename PointSetType::PointDataContainer;
165 
167  using RealType = float;
169 
177 
184 
185 
190  void SetSplineOrder( unsigned int );
191 
195  void SetSplineOrder( const ArrayType & );
196 
200  itkGetConstReferenceMacro( SplineOrder, ArrayType );
201 
206  itkSetMacro( NumberOfControlPoints, ArrayType );
207  itkGetConstReferenceMacro( NumberOfControlPoints, ArrayType );
209 
214  itkGetConstReferenceMacro( CurrentNumberOfControlPoints, ArrayType );
215 
220  void SetNumberOfLevels( unsigned int );
221 
226  void SetNumberOfLevels( const ArrayType & );
227 
232  itkGetConstReferenceMacro( NumberOfLevels, ArrayType );
233 
240  itkSetMacro( BSplineEpsilon, RealType );
241  itkGetConstMacro( BSplineEpsilon, RealType );
243 
260  itkSetMacro( CloseDimension, ArrayType );
261  itkGetConstReferenceMacro( CloseDimension, ArrayType );
263 
266  void SetPointWeights( WeightsContainerType *weights );
267 
271  itkSetMacro( GenerateOutputImage, bool );
272  itkGetConstReferenceMacro( GenerateOutputImage, bool );
273  itkBooleanMacro( GenerateOutputImage );
275 
278  {
279  return static_cast<PointDataImageType *>( this->ProcessObject::GetOutput( 1 ) );
280  }
281 
282 protected:
284  ~BSplineScatteredDataPointSetToImageFilter() override = default;
285 
286  void PrintSelf(std::ostream & os, Indent indent) const override;
287 
288  void ThreadedGenerateData( const RegionType &, ThreadIdType ) override;
289 
290  void DynamicThreadedGenerateData( const RegionType & ) override
291  {
292  itkExceptionMacro("This class requires threadId so it must use classic multi-threading model");
293  }
294 
295  void BeforeThreadedGenerateData() override;
296 
297  void AfterThreadedGenerateData() override;
298 
299  unsigned int SplitRequestedRegion( unsigned int, unsigned int, RegionType & ) override;
300 
301  void GenerateData() override;
302 
303 private:
304 
307  void RefineControlPointLattice();
308 
310  void UpdatePointSet();
311 
314  void GenerateOutputImage();
315 
317  void ThreadedGenerateDataForFitting( const RegionType &, ThreadIdType );
318 
320  void ThreadedGenerateDataForReconstruction( const RegionType &, ThreadIdType );
321 
324  void CollapsePhiLattice( PointDataImageType *, PointDataImageType *,
325  const RealType, const unsigned int );
326 
329  void SetPhiLatticeParametricDomainParameters();
330 
333  IndexType NumberToIndex( const unsigned int, const SizeType );
334 
335  bool m_DoMultilevel{ false };
336  bool m_GenerateOutputImage{ true };
337  bool m_UsePointWeights{ false };
338  unsigned int m_MaximumNumberOfLevels{ 1 };
339  unsigned int m_CurrentLevel{ 0 };
345 
347 
350 
351  vnl_matrix<RealType> m_RefinedLatticeCoefficients[ImageDimension];
352 
353  typename PointDataContainerType::Pointer m_InputPointData;
354  typename PointDataContainerType::Pointer m_OutputPointData;
355 
356  typename KernelType::Pointer m_Kernel[ImageDimension];
357 
362 
363  std::vector<RealImagePointer> m_OmegaLatticePerThread;
364  std::vector<PointDataImagePointer> m_DeltaLatticePerThread;
365 
366  RealType m_BSplineEpsilon{ static_cast< RealType >( 1e-3 ) };
367  bool m_IsFittingComplete{ false };
368 };
369 } // end namespace itk
370 
371 #ifndef ITK_MANUAL_INSTANTIATION
372 #include "itkBSplineScatteredDataPointSetToImageFilter.hxx"
373 #endif
374 
375 #endif
The base class for all process objects (source, filters, mappers) in the Insight data processing pipe...
Base class for all process objects that output image data.
BSpline kernel used for density estimation and nonparameteric regression.
typename TOutputImage::PointType PointType
BSpline kernel used for density estimation and nonparameteric regression.
typename OutputImageType::RegionType OutputImageRegionType
Image filter which provides a B-spline output approximation.
typename TOutputImage::SizeType SizeType
unsigned int ThreadIdType
Definition: itkIntTypes.h:99
Base class for filters that take a PointSet as input and produce an image as output. By default, if the user does not specify the size of the output image, the maximum size of the point-set&#39;s bounding box is used.
static constexpr double e
The base of the natural logarithm or Euler&#39;s number
Definition: itkMath.h:53
Define a front-end to the STL &quot;vector&quot; container that conforms to the IndexedContainerInterface.
Control indentation during Print() invocation.
Definition: itkIndent.h:49
Templated n-dimensional image class.
Definition: itkImage.h:75
DataObject * GetOutput(const DataObjectIdentifierType &key)