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 itkANTSNeighborhoodCorrelationImageToImageMetricv4_h
19 #define itkANTSNeighborhoodCorrelationImageToImageMetricv4_h
20 
21 #include "itkImageToImageMetricv4.h"
22 #include "itkANTSNeighborhoodCorrelationImageToImageMetricv4GetValueAndDerivativeThreader.h"
23 
24 namespace itk {
25 
26 /** \class ANTSNeighborhoodCorrelationImageToImageMetricv4
27  *
28  * \brief Computes normalized cross correlation using a small neighborhood
29  * for each voxel between two images, with speed optimizations for dense
30  * registration.
31  *
32  * Please cite this reference for more details:
33  *
34  * Brian B. Avants, Nicholas J. Tustison, Gang Song, Philip A. Cook,
35  * Arno Klein, James C. Gee, A reproducible evaluation of ANTs similarity metric
36  * performance in brain image registration, NeuroImage, Volume 54, Issue 3,
37  * 1 February 2011, Pages 2033-2044, ISSN 1053-8119,
38  * DOI: 10.1016/j.neuroimage.2010.09.025.
39  *
40  * Around each voxel, the neighborhood is defined as a N-Dimensional
41  * rectangle centered at the voxel. The size of the rectangle is 2*radius+1.
42  * The normalized correlation between neighborhoods of fixed image and moving
43  * image are averaged over the whole image as the final metric.
44  * \note A radius less than 2 can be unstable. 2 is the default.
45  *
46  * This class uses a specific fast implementation that is described in the
47  * above paper. There are two particular speed-ups:
48  *
49  * 1) It is assumed that the derivative is only affected by changes in the
50  * transform at the center of the window. This is obviously not true but speeds
51  * the evaluation up considerably and works well in practice. This assumption
52  * is the main differentiation of this approach from a more generic one.
53  *
54  * 2) The evaluation uses on-the-fly queues with multi-threading and a sliding
55  * neighborhood window. This is described in the above paper and specifically
56  * optimized for dense registration.
57  *
58  *  Example of usage:
59  *
60  *  using MetricType = itk::ANTSNeighborhoodCorrelationImageToImageMetricv4
61  *    <ImageType, ImageType>;
62  *  using MetricTypePointer = MetricType::Pointer;
63  *  MetricTypePointer metric = MetricType::New();
64  *
65  *  // set all parameters
66  *  Size<Dimension> neighborhoodRadius;
67  *  neighborhoodRadius.Fill(2);
68  *  metric->SetRadius(neighborhood_radius);
69  *  metric->SetFixedImage(fixedImage);
70  *  metric->SetMovingImage(movingImage);
71  *  metric->SetFixedTransform(transformFix);
72  *  metric->SetMovingTransform(transformMov);
73  *
74  *  // initialization after parameters are set.
75  *  metric->Initialize();
76  *
77  *  // getting derivative and metric value
78  *  metric->GetValueAndDerivative(valueReturn, derivativeReturn);
79  *
80  *
81  * This class is templated over the type of the two input objects.
82  * This is the base class for a hierarchy of similarity metrics that may, in
83  * derived classes, operate on meshes, images, etc.  This class computes a
84  * value that measures the similarity between the two objects.
85  *
86  * \note Sparse sampling is not supported by this metric. An exception will be
87  * thrown if m_UseSampledPointSet is set. Support for sparse sampling
88  * will require a parallel implementation of the neighborhood scanning, which
89  * currently caches information as the neighborhood window moves.
90  *
91  * \ingroup ITKMetricsv4
92  */
93 template<typename TFixedImage, typename TMovingImage, typename TVirtualImage = TFixedImage,
94           typename TInternalComputationValueType = double,
95           typename TMetricTraits = DefaultImageToImageMetricTraitsv4<TFixedImage,TMovingImage,TVirtualImage,TInternalComputationValueType>
96           >
97 class ITK_TEMPLATE_EXPORT ANTSNeighborhoodCorrelationImageToImageMetricv4 :
98   public ImageToImageMetricv4<TFixedImage, TMovingImage, TVirtualImage, TInternalComputationValueType, TMetricTraits>
99 {
100 public:
101   ITK_DISALLOW_COPY_AND_ASSIGN(ANTSNeighborhoodCorrelationImageToImageMetricv4);
102 
103   /** Standard class type aliases. */
104   using Self = ANTSNeighborhoodCorrelationImageToImageMetricv4;
105   using Superclass = ImageToImageMetricv4<TFixedImage, TMovingImage, TVirtualImage,
106                              TInternalComputationValueType,TMetricTraits>;
107   using Pointer = SmartPointer<Self>;
108   using ConstPointer = SmartPointer<const Self>;
109 
110   /** Method for creation through the object factory. */
111   itkNewMacro(Self);
112 
113   /** Run-time type information (and related methods). */
114   itkTypeMacro(Self, Superclass);
115 
116   /** superclass types */
117   using MeasureType = typename Superclass::MeasureType;
118   using DerivativeType = typename Superclass::DerivativeType;
119   using DerivativeValueType = typename Superclass::DerivativeValueType;
120   using VirtualPointType = typename Superclass::VirtualPointType;
121   using FixedImagePointType = typename Superclass::FixedImagePointType;
122   using FixedImagePixelType = typename Superclass::FixedImagePixelType;
123   using FixedTransformType = typename Superclass::FixedTransformType;
124   using FixedImageGradientType = typename Superclass::FixedImageGradientType;
125   using FixedImageJacobianType = typename FixedTransformType::JacobianType;
126 
127   using MovingImagePointType = typename Superclass::MovingImagePointType;
128   using MovingImagePixelType = typename Superclass::MovingImagePixelType;
129   using MovingImageGradientType = typename Superclass::MovingImageGradientType;
130   using MovingTransformType = typename Superclass::MovingTransformType;
131   using MovingImageJacobianType = typename MovingTransformType::JacobianType;
132   using JacobianType = typename Superclass::JacobianType;
133 
134   using VirtualImageGradientType = typename Superclass::VirtualImageGradientType;
135 
136   using FixedImageType = typename Superclass::FixedImageType;
137   using MovingImageType = typename Superclass::MovingImageType;
138   using VirtualImageType = typename Superclass::VirtualImageType;
139   using FixedOutputPointType = typename Superclass::FixedOutputPointType;
140   using MovingOutputPointType = typename Superclass::MovingOutputPointType;
141 
142   using FixedTransformJacobianType = typename Superclass::FixedTransformType::JacobianType;
143   using MovingTransformJacobianType = typename Superclass::MovingTransformType::JacobianType;
144 
145   using NumberOfParametersType = typename Superclass::NumberOfParametersType;
146   using ImageDimensionType = typename Superclass::ImageDimensionType;
147 
148   using ImageRegionType = typename VirtualImageType::RegionType;
149   using RadiusType = typename VirtualImageType::SizeType;
150   using IndexType = typename VirtualImageType::IndexType;
151 
152   /* Image dimension accessors */
153   static constexpr ImageDimensionType FixedImageDimension = FixedImageType::ImageDimension;
154 
155   static constexpr ImageDimensionType MovingImageDimension = MovingImageType::ImageDimension;
156 
157   static constexpr ImageDimensionType VirtualImageDimension = VirtualImageType::ImageDimension;
158 
159   // Set the radius of the neighborhood window centered at each pixel.
160   // See the note above about using a radius less than 2.
161   itkSetMacro(Radius, RadiusType);
162 
163   // Get the Radius of the neighborhood window centered at each pixel
164   itkGetMacro(Radius, RadiusType);
165   itkGetConstMacro(Radius, RadiusType);
166 
167   void Initialize() override;
168 
169 protected:
170   ANTSNeighborhoodCorrelationImageToImageMetricv4();
171   ~ANTSNeighborhoodCorrelationImageToImageMetricv4() override = default;
172 
173   friend class ANTSNeighborhoodCorrelationImageToImageMetricv4GetValueAndDerivativeThreader< ThreadedImageRegionPartitioner< VirtualImageDimension >, Superclass, Self >;
174   using ANTSNeighborhoodCorrelationImageToImageMetricv4DenseGetValueAndDerivativeThreaderType =
175       ANTSNeighborhoodCorrelationImageToImageMetricv4GetValueAndDerivativeThreader< ThreadedImageRegionPartitioner< VirtualImageDimension >, Superclass, Self >;
176 
177   friend class ANTSNeighborhoodCorrelationImageToImageMetricv4GetValueAndDerivativeThreader< ThreadedIndexedContainerPartitioner, Superclass, Self >;
178   using ANTSNeighborhoodCorrelationImageToImageMetricv4SparseGetValueAndDerivativeThreaderType =
179       ANTSNeighborhoodCorrelationImageToImageMetricv4GetValueAndDerivativeThreader< ThreadedIndexedContainerPartitioner, Superclass, Self >;
180 
181   void PrintSelf(std::ostream & os, Indent indent) const override;
182 
183 private:
184   // Radius of the neighborhood window centered at each pixel
185   RadiusType m_Radius;
186 };
187 
188 } // end namespace itk
189 
190 #ifndef ITK_MANUAL_INSTANTIATION
191 #include "itkANTSNeighborhoodCorrelationImageToImageMetricv4.hxx"
192 #endif
193 
194 #endif
195