ITK  4.6.0
Insight Segmentation and Registration Toolkit
itkKalmanLinearEstimator.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 __itkKalmanLinearEstimator_h
19 #define __itkKalmanLinearEstimator_h
20 
21 #include "itkMacro.h"
22 
23 #include "vnl/vnl_vector_fixed.h"
24 #include "vnl/vnl_matrix_fixed.h"
25 
26 namespace itk
27 {
41 template< typename T, unsigned int VEstimatorDimension >
43 {
44 public:
47  itkStaticConstMacro(Dimension, unsigned int,
48  VEstimatorDimension);
49 
52  typedef vnl_vector_fixed< T, VEstimatorDimension > VectorType;
53 
56  typedef vnl_matrix_fixed< T, VEstimatorDimension, VEstimatorDimension > MatrixType;
57 
62  typedef T ValueType;
63 
68  void UpdateWithNewMeasure(const ValueType & newMeasure,
69  const VectorType & newPredictor);
70 
74  void ClearEstimation(void)
75  { m_Estimator = VectorType( T(0) ); }
76 
79  void ClearVariance(void)
80  {
81  m_Variance.set_identity();
82  }
83 
90  void SetVariance(const ValueType & var = 1.0)
91  {
92  m_Variance.set_identity();
93  m_Variance *= var;
94  }
96 
102  void SetVariance(const MatrixType & m)
103  { m_Variance = m; }
104 
107  const VectorType & GetEstimator(void) const
108  { return m_Estimator; }
109 
112  const MatrixType & GetVariance(void) const
113  { return m_Variance; }
114 
115 private:
120  void UpdateVariance(const VectorType &);
121 
125 
133 };
134 
135 template< typename T, unsigned int VEstimatorDimension >
136 void
139  const VectorType & newPredictor)
140 {
141  ValueType measurePrediction = dot_product(newPredictor, m_Estimator);
142 
143  ValueType errorMeasurePrediction = newMeasure - measurePrediction;
144 
145  VectorType Corrector = m_Variance * newPredictor;
146 
147  for ( unsigned int j = 0; j < VEstimatorDimension; j++ )
148  {
149  m_Estimator(j) += Corrector(j) * errorMeasurePrediction;
150  }
151 
152  UpdateVariance(newPredictor);
153 }
154 
155 template< typename T, unsigned int VEstimatorDimension >
156 void
158 ::UpdateVariance(const VectorType & newPredictor)
159 {
160  VectorType aux = m_Variance * newPredictor;
161 
162  ValueType denominator = 1.0 / ( 1.0 + dot_product(aux, newPredictor) );
163 
164  for ( unsigned int col = 0; col < VEstimatorDimension; col++ )
165  {
166  for ( unsigned int row = 0; row < VEstimatorDimension; row++ )
167  {
168  m_Variance(col, row) -= aux(col) * aux(row) * denominator;
169  }
170  }
171 }
172 } // end namespace itk
173 
174 #endif
const VectorType & GetEstimator(void) const
vnl_matrix_fixed< T, VEstimatorDimension, VEstimatorDimension > MatrixType
void UpdateWithNewMeasure(const ValueType &newMeasure, const VectorType &newPredictor)
const MatrixType & GetVariance(void) const
Implement a linear recursive estimator.
vnl_vector_fixed< T, VEstimatorDimension > VectorType
void SetVariance(const MatrixType &m)
void UpdateVariance(const VectorType &)
void SetVariance(const ValueType &var=1.0)
static const unsigned int Dimension