core/vnl/algo/vnl_gaussian_kernel_1d.cxx
Go to the documentation of this file.
00001 // This is core/vnl/algo/vnl_gaussian_kernel_1d.cxx
00002 #ifdef VCL_NEEDS_PRAGMA_INTERFACE
00003 #pragma implementation
00004 #endif
00005 //:
00006 // \file
00007 // \author Andrew W. Fitzgibbon, Oxford RRG
00008 // \date   07 Aug 1997
00009 //
00010 //-----------------------------------------------------------------------------
00011 
00012 #include "vnl_gaussian_kernel_1d.h"
00013 #include <vcl_cmath.h>
00014 #include <vnl/vnl_math.h>
00015 
00016 // G(x) = 1/(sigma * sqrt(2*pi)) * exp(-0.5 * (x/sigma)^2)
00017 // x(g) = sigma * sqrt(-2 * log(g * sigma * sqrt(2*pi) ) )
00018 
00019 // Compute the x value at which a Gaussian becomes lower than cutoff.
00020 static inline
00021 double compute_width(double sigma, double cutoff)
00022 {
00023   return sigma * vcl_sqrt(-2 * vcl_log(cutoff * sigma * vcl_sqrt(2*vnl_math::pi)));
00024 }
00025 
00026 //: Construct a sampled 1D gaussian of standard deviation sigma.
00027 // The vector is normalized so that its sum is 0.5.
00028 vnl_gaussian_kernel_1d::vnl_gaussian_kernel_1d(double sigma, double cutoff):
00029   vec_((int)vcl_ceil(compute_width(sigma, cutoff)))
00030 {
00031   int wid = vec_.size();
00032   inscale_ = 0.5/(sigma * sigma);
00033   double area = 0;
00034   for (int i = 0; i < wid; ++i) {
00035     double v = G(i);
00036     area += v;
00037     vec_[i] = v;
00038   }
00039   vec_ *= (0.5/area);
00040 }
00041 
00042 double vnl_gaussian_kernel_1d::G(double x) const
00043 {
00044   return vcl_exp(-x*x * inscale_);
00045 }