starter/math/robust_estimators.h

   1 /*
   2  Copyright 2005, 2006 Computer Vision Lab,
   3  Ecole Polytechnique Federale de Lausanne (EPFL), Switzerland.
   4  Modified by Damian Stewart <damian@frey.co.nz> 2009-2010;
   5  modifications Copyright 2009, 2010 Damian Stewart <damian@frey.co.nz>.
   6
   7  Distributed under the terms of the GNU General Public License v3.
   8
   9  This file is part of The Artvertiser.
  10
  11  The Artvertiser is free software: you can redistribute it and/or modify
  12  it under the terms of the GNU General Public License as published by
  13  the Free Software Foundation, either version 3 of the License, or
  14  (at your option) any later version.
  15
  16  The Artvertiser is distributed in the hope that it will be useful,
  17  but WITHOUT ANY WARRANTY; without even the implied warranty of
  18  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  19  GNU General Public License for more details.
  20
  21  You should have received a copy of the GNU Lesser General Public License
  22  along with The Artvertiser.  If not, see <http://www.gnu.org/licenses/>.
  23  */
  24
  25 #ifndef ROBUST_ESTIMATORS_H
  26 #define ROBUST_ESTIMATORS_H
  27
  28 #include <math.h>
  29
  30 //! \ingroup starter
  31 //@{
  32
  33 double rho_tukey_without_sqrt(const double x2, const double c);
  34 double rho_tukey_without_sqrt_derivative(const double x2, const double c);
  35
  36 double rho_huber_without_sqrt(const double x2, const double c);
  37 double rho_huber_without_sqrt_derivative(const double x2, const double c);
  38
  39
  40 inline double rho_tukey_without_sqrt(const double x2, const double c)
  41 {
  42   double c2 = c * c;
  43
  44   if (x2 > c2)
  45     return c2 / 6;
  46   else
  47   {
  48     double t = 1. - x2 / c2;
  49
  50     return c2 / 6 * (1. - t * t * t);
  51   }
  52 }
  53
  54 inline double rho_tukey_without_sqrt_derivative(const double x2, const double c)
  55 {
  56   double c2 = c * c;
  57
  58   if (x2 > c2)
  59     return 0;
  60   else
  61   {
  62     double t = 1. - x2 / c2;
  63
  64     return t * t / 2.;
  65   }
  66 }
  67
  68 inline double rho_huber_without_sqrt(const double x2, const double c)
  69 {
  70   double c2 = c * c;
  71
  72   if (x2 < c2)
  73     return x2;
  74   else
  75     return 2 * c * sqrt(x2) - c2;
  76 }
  77
  78 inline double rho_huber_without_sqrt_derivative(const double x2, const double c)
  79 {
  80   double c2 = c * c;
  81
  82   if (x2 < c2)
  83     return 1;
  84   else
  85     return c / sqrt(x2);
  86 }
  87
  88 //@}
  89 #endif // ROBUST_ESTIMATORS_H