LeastSquares.h
Go to the documentation of this file.
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2006-2009 Benoit Jacob <jacob.benoit.1@gmail.com>
5 //
6 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9 
10 #ifndef EIGEN2_LEASTSQUARES_H
11 #define EIGEN2_LEASTSQUARES_H
12 
13 namespace Eigen {
14 
84 template<typename VectorType>
85 void linearRegression(int numPoints,
86  VectorType **points,
87  VectorType *result,
88  int funcOfOthers )
89 {
90  typedef typename VectorType::Scalar Scalar;
92  const int size = points[0]->size();
93  result->resize(size);
94  HyperplaneType h(size);
95  fitHyperplane(numPoints, points, &h);
96  for(int i = 0; i < funcOfOthers; i++)
97  result->coeffRef(i) = - h.coeffs()[i] / h.coeffs()[funcOfOthers];
98  for(int i = funcOfOthers; i < size; i++)
99  result->coeffRef(i) = - h.coeffs()[i+1] / h.coeffs()[funcOfOthers];
100 }
101 
129 template<typename VectorType, typename HyperplaneType>
130 void fitHyperplane(int numPoints,
131  VectorType **points,
132  HyperplaneType *result,
133  typename NumTraits<typename VectorType::Scalar>::Real* soundness = 0)
134 {
135  typedef typename VectorType::Scalar Scalar;
138  ei_assert(numPoints >= 1);
139  int size = points[0]->size();
140  ei_assert(size+1 == result->coeffs().size());
141 
142  // compute the mean of the data
143  VectorType mean = VectorType::Zero(size);
144  for(int i = 0; i < numPoints; ++i)
145  mean += *(points[i]);
146  mean /= numPoints;
147 
148  // compute the covariance matrix
149  CovMatrixType covMat = CovMatrixType::Zero(size, size);
150  VectorType remean = VectorType::Zero(size);
151  for(int i = 0; i < numPoints; ++i)
152  {
153  VectorType diff = (*(points[i]) - mean).conjugate();
154  covMat += diff * diff.adjoint();
155  }
156 
157  // now we just have to pick the eigen vector with smallest eigen value
159  result->normal() = eig.eigenvectors().col(0);
160  if (soundness)
161  *soundness = eig.eigenvalues().coeff(0)/eig.eigenvalues().coeff(1);
162 
163  // let's compute the constant coefficient such that the
164  // plane pass trough the mean point:
165  result->offset() = - (result->normal().cwise()* mean).sum();
166 }
167 
168 } // end namespace Eigen
169 
170 #endif // EIGEN2_LEASTSQUARES_H
#define ei_assert
Computes eigenvalues and eigenvectors of selfadjoint matrices.
iterative scaling algorithm to equilibrate rows and column norms in matrices
Definition: matrix.hpp:471
void fitHyperplane(int numPoints, VectorType **points, HyperplaneType *result, typename NumTraits< typename VectorType::Scalar >::Real *soundness=0)
Definition: LeastSquares.h:130
void linearRegression(int numPoints, VectorType **points, VectorType *result, int funcOfOthers)
Definition: LeastSquares.h:85
const MatrixType & eigenvectors() const
Returns the eigenvectors of given matrix.
The matrix class, also used for vectors and row-vectors.
Definition: Matrix.h:127
#define EIGEN_STATIC_ASSERT_VECTOR_ONLY(TYPE)
Definition: StaticAssert.h:126
const RealVectorType & eigenvalues() const
Returns the eigenvalues of given matrix.


acado
Author(s): Milan Vukov, Rien Quirynen
autogenerated on Mon Jun 10 2019 12:34:48