itkStatisticalShapePointPenalty.h

Go to the documentation of this file.

/*=========================================================================
 *
 *  Copyright UMC Utrecht and contributors
 *
 *  Licensed under the Apache License, Version 2.0 (the "License");
 *  you may not use this file except in compliance with the License.
 *  You may obtain a copy of the License at
 *
 *        http://www.apache.org/licenses/LICENSE-2.0.txt
 *
 *  Unless required by applicable law or agreed to in writing, software
 *  distributed under the License is distributed on an "AS IS" BASIS,
 *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 *  See the License for the specific language governing permissions and
 *  limitations under the License.
 *
 *=========================================================================*/
#ifndef itkStatisticalShapePointPenalty_h
#define itkStatisticalShapePointPenalty_h
 
#include "itkSingleValuedPointSetToPointSetMetric.h"
 
#include "itkPoint.h"
#include "itkPointSet.h"
#include "itkImage.h"
#include "itkArray.h"
#include <itkVariableSizeMatrix.h>
 
#include <vnl/vnl_matrix.h>
#include <vnl/vnl_math.h>
#include <vnl/vnl_vector.h>
#include <vnl/algo/vnl_real_eigensystem.h>
#include <vnl/algo/vnl_symmetric_eigensystem.h>
// #include <vnl/algo/vnl_svd.h>
#include <vnl/algo/vnl_svd_economy.h>
 
#include <string>
 
namespace itk
{
 
template <typename TFixedPointSet, typename TMovingPointSet>
class ITK_TEMPLATE_EXPORT StatisticalShapePointPenalty
  : public SingleValuedPointSetToPointSetMetric<TFixedPointSet, TMovingPointSet>
{
public:
  ITK_DISALLOW_COPY_AND_MOVE(StatisticalShapePointPenalty);

  using Self = StatisticalShapePointPenalty;
  using Superclass = SingleValuedPointSetToPointSetMetric<TFixedPointSet, TMovingPointSet>;
  using Pointer = SmartPointer<Self>;
  using ConstPointer = SmartPointer<const Self>;

  itkNewMacro(Self);

  itkOverrideGetNameOfClassMacro(StatisticalShapePointPenalty);

  using typename Superclass::TransformType;
  using typename Superclass::TransformPointer;
  using typename Superclass::ParametersType;
  using typename Superclass::TransformJacobianType;
  using typename Superclass::NonZeroJacobianIndicesType;
 
  using typename Superclass::MeasureType;
  using typename Superclass::DerivativeType;
  using typename Superclass::DerivativeValueType;
  using typename Superclass::FixedPointSetType;
  using typename Superclass::MovingPointSetType;
  using typename Superclass::FixedPointSetConstPointer;
  using typename Superclass::MovingPointSetConstPointer;
 
  using typename Superclass::PointIterator;
 
  using typename Superclass::InputPointType;
  using typename Superclass::OutputPointType;
 
  using CoordinateType = typename OutputPointType::CoordinateType;
  using VnlVectorType = vnl_vector<CoordinateType>;
  using VnlMatrixType = vnl_matrix<CoordinateType>;
  // typedef itk::Array<VnlVectorType *> ProposalDerivativeType;
  using ProposalDerivativeType = typename std::vector<VnlVectorType *>;
  // typedef typename vnl_vector<VnlVectorType *> ProposalDerivativeType; //Cannot be linked
  using PCACovarianceType = vnl_svd_economy<CoordinateType>;

  void
  Initialize() override;

  MeasureType
  GetValue(const ParametersType & parameters) const override;

  void
  GetDerivative(const ParametersType & parameters, DerivativeType & Derivative) const override;

  void
  GetValueAndDerivative(const ParametersType & parameters,
                        MeasureType &          Value,
                        DerivativeType &       Derivative) const override;

  itkSetClampMacro(ShrinkageIntensity, MeasureType, 0.0, 1.0);
  itkGetConstMacro(ShrinkageIntensity, MeasureType);
 
  itkSetMacro(ShrinkageIntensityNeedsUpdate, bool);
  itkBooleanMacro(ShrinkageIntensityNeedsUpdate);

  itkSetClampMacro(BaseVariance, MeasureType, -1.0, NumericTraits<MeasureType>::max());
  itkGetConstMacro(BaseVariance, MeasureType);
 
  itkSetMacro(BaseVarianceNeedsUpdate, bool);
  itkBooleanMacro(BaseVarianceNeedsUpdate);
 
  itkSetClampMacro(CentroidXVariance, MeasureType, -1.0, NumericTraits<MeasureType>::max());
  itkGetConstMacro(CentroidXVariance, MeasureType);
 
  itkSetClampMacro(CentroidYVariance, MeasureType, -1.0, NumericTraits<MeasureType>::max());
  itkGetConstMacro(CentroidYVariance, MeasureType);
 
  itkSetClampMacro(CentroidZVariance, MeasureType, -1.0, NumericTraits<MeasureType>::max());
  itkGetConstMacro(CentroidZVariance, MeasureType);
 
  itkSetClampMacro(SizeVariance, MeasureType, -1.0, NumericTraits<MeasureType>::max());
  itkGetConstMacro(SizeVariance, MeasureType);
 
  itkSetMacro(VariancesNeedsUpdate, bool);
  itkBooleanMacro(VariancesNeedsUpdate);
 
  itkSetClampMacro(CutOffValue, MeasureType, 0.0, NumericTraits<MeasureType>::max());
  itkGetConstMacro(CutOffValue, MeasureType);
 
  itkSetClampMacro(CutOffSharpness,
                   MeasureType,
                   NumericTraits<MeasureType>::NonpositiveMin(),
                   NumericTraits<MeasureType>::max());
  itkGetConstMacro(CutOffSharpness, MeasureType);
 
  itkSetMacro(ShapeModelCalculation, int);
  itkGetConstReferenceMacro(ShapeModelCalculation, int);
 
  itkSetMacro(NormalizedShapeModel, bool);
  itkGetConstReferenceMacro(NormalizedShapeModel, bool);
  itkBooleanMacro(NormalizedShapeModel);
 
  itkSetConstObjectMacro(EigenVectors, vnl_matrix<double>);
  itkSetConstObjectMacro(EigenValues, vnl_vector<double>);
  itkSetConstObjectMacro(MeanVector, vnl_vector<double>);
 
  itkSetConstObjectMacro(CovarianceMatrix, vnl_matrix<double>);
 
protected:
  StatisticalShapePointPenalty();
  ~StatisticalShapePointPenalty() override;

  void
  PrintSelf(std::ostream & os, Indent indent) const override;
 
private:
  void
  FillProposalVector(const OutputPointType & fixedPoint, const unsigned int vertexindex) const;
 
  void
  FillProposalDerivative(const OutputPointType & fixedPoint, const unsigned int vertexindex) const;
 
  void
  UpdateCentroidAndAlignProposalVector(const unsigned int shapeLength) const;
 
  void
  UpdateCentroidAndAlignProposalDerivative(const unsigned int shapeLength) const;
 
  void
  UpdateL2(const unsigned int shapeLength) const;
 
  void
  NormalizeProposalVector(const unsigned int shapeLength) const;
 
  void
  UpdateL2AndNormalizeProposalDerivative(const unsigned int shapeLength) const;
 
  void
  CalculateValue(MeasureType &   value,
                 VnlVectorType & differenceVector,
                 VnlVectorType & centerrotated,
                 VnlVectorType & eigrot) const;
 
  void
  CalculateDerivative(DerivativeType &      derivative,
                      const MeasureType &   value,
                      const VnlVectorType & differenceVector,
                      const VnlVectorType & eigrot,
                      const unsigned int    shapeLength) const;
 
  void
  CalculateCutOffValue(MeasureType & value) const;
 
  void
  CalculateCutOffDerivative(typename DerivativeType::element_type & derivativeElement, const MeasureType & value) const;
 
  const VnlVectorType * m_MeanVector{};
  const VnlMatrixType * m_CovarianceMatrix{};
  const VnlMatrixType * m_EigenVectors{};
  const VnlVectorType * m_EigenValues{};
 
  VnlMatrixType * m_InverseCovarianceMatrix{};
 
  double m_CentroidXVariance{};
  double m_CentroidXStd{};
  double m_CentroidYVariance{};
  double m_CentroidYStd{};
  double m_CentroidZVariance{};
  double m_CentroidZStd{};
  double m_SizeVariance{};
  double m_SizeStd{};
 
  bool m_ShrinkageIntensityNeedsUpdate{};
  bool m_BaseVarianceNeedsUpdate{};
  bool m_VariancesNeedsUpdate{};
 
  VnlVectorType * m_EigenValuesRegularized{};
 
  mutable ProposalDerivativeType * m_ProposalDerivative{};
  unsigned int                     m_ProposalLength{};
  bool                             m_NormalizedShapeModel{};
  int                              m_ShapeModelCalculation{};
  double                           m_ShrinkageIntensity{};
  double                           m_BaseVariance{};
  double                           m_BaseStd{};
  mutable VnlVectorType            m_ProposalVector{};
  mutable VnlVectorType            m_MeanValues{};
 
  double m_CutOffValue{};
  double m_CutOffSharpness{};
};
 
} // end namespace itk
 
#ifndef ITK_MANUAL_INSTANTIATION
#  include "itkStatisticalShapePointPenalty.hxx"
#endif
 
#endif