contrib/mul/vpdfl/vpdfl_gaussian_kernel_pdf_sampler.cxx
Go to the documentation of this file.
00001 // This is mul/vpdfl/vpdfl_gaussian_kernel_pdf_sampler.cxx
00002 #ifdef VCL_NEEDS_PRAGMA_INTERFACE
00003 #pragma implementation
00004 #endif
00005 //:
00006 // \file
00007 // \author Tim Cootes
00008 // \brief Sampler class for gaussian kernel PDF.
00009 
00010 #include "vpdfl_gaussian_kernel_pdf_sampler.h"
00011 
00012 #include <vcl_cassert.h>
00013 #include <vpdfl/vpdfl_gaussian_kernel_pdf.h>
00014 
00015 //=======================================================================
00016 // Dflt ctor
00017 //=======================================================================
00018 
00019 vpdfl_gaussian_kernel_pdf_sampler::vpdfl_gaussian_kernel_pdf_sampler():
00020  rng_(9667566ul)
00021 {
00022 }
00023 
00024 //=======================================================================
00025 // Destructor
00026 //=======================================================================
00027 
00028 vpdfl_gaussian_kernel_pdf_sampler::~vpdfl_gaussian_kernel_pdf_sampler()
00029 {
00030 }
00031 
00032 
00033 //=======================================================================
00034 // Method: is_a
00035 //=======================================================================
00036 
00037 vcl_string vpdfl_gaussian_kernel_pdf_sampler::is_a() const
00038 {
00039   return vcl_string("vpdfl_gaussian_kernel_pdf_sampler");
00040 }
00041 
00042 //=======================================================================
00043 // Method: is_class
00044 //=======================================================================
00045 
00046 bool vpdfl_gaussian_kernel_pdf_sampler::is_class(vcl_string const& s) const
00047 {
00048   return vpdfl_sampler_base::is_class(s) || s==vpdfl_gaussian_kernel_pdf_sampler::is_a();
00049 }
00050 
00051 //=======================================================================
00052 // Method: clone
00053 //=======================================================================
00054 
00055 vpdfl_sampler_base* vpdfl_gaussian_kernel_pdf_sampler::clone() const
00056 {
00057   return new vpdfl_gaussian_kernel_pdf_sampler(*this);
00058 }
00059 //=======================================================================
00060 
00061 void vpdfl_gaussian_kernel_pdf_sampler::reseed(unsigned long seed)
00062 {
00063   rng_.reseed(seed);
00064 }
00065 //=======================================================================
00066 
00067 //: Set model for which this is an instance
00068 // Error check that it is an axis gaussian.
00069 void vpdfl_gaussian_kernel_pdf_sampler::set_model(const vpdfl_pdf_base& model)
00070 {
00071   assert(model.is_class("vpdfl_gaussian_kernel_pdf"));
00072   // cannot use dynamic_cast<> without rtti - PVr
00073   // rtti currently turned off
00074   vpdfl_sampler_base::set_model(model);
00075 }
00076 
00077 void vpdfl_gaussian_kernel_pdf_sampler::sample_component(vnl_vector<double>& x,
00078                                                          int j)
00079 {
00080   const vpdfl_gaussian_kernel_pdf & kpdf = static_cast<const vpdfl_gaussian_kernel_pdf &>(model());
00081 
00082   int n_dims = kpdf.n_dims();
00083   x.set_size(n_dims);
00084 
00085   const double* m = kpdf.centre()[j].data_block();
00086   double w = kpdf.width()[j];
00087 
00088   double* x_data = x.data_block();
00089   for (int i=0;i<n_dims;++i)
00090     x_data[i] = m[i] + w*rng_.normal();
00091 }
00092 
00093 
00094 //=======================================================================
00095 
00096 void vpdfl_gaussian_kernel_pdf_sampler::sample(vnl_vector<double>& x)
00097 {
00098   const vpdfl_gaussian_kernel_pdf & kpdf = static_cast<const vpdfl_gaussian_kernel_pdf &>(model());
00099   int n = kpdf.centre().size();
00100 
00101   // Select component
00102   int j = 0;
00103   if (n>1) j=rng_.lrand32(0,n-1);
00104 
00105   sample_component(x,j);
00106 }
00107 
00108 //: Fill x with samples possibly chosen so as to represent the distribution
00109 //  Sample sequentially from each component.
00110 void vpdfl_gaussian_kernel_pdf_sampler::regular_samples(
00111                        vcl_vector<vnl_vector<double> >& x)
00112 {
00113   const vpdfl_gaussian_kernel_pdf & kpdf = static_cast<const vpdfl_gaussian_kernel_pdf &>(model());
00114   int n_k = kpdf.centre().size();
00115 
00116   int n_samples = x.size();
00117 
00118   for (int i=0;i<n_samples;++i)
00119   {
00120     // Ensure sample at centre.  (This biases towards the centres somewhat)
00121 //     if (i<n_k)
00122 //       x[i]=kpdf.centre()[i];
00123 //     else
00124       sample_component(x[i],i%n_k);
00125   }
00126 }
00127 
00128 //=======================================================================
00129 
00130 
00131 //: Return a reference to the pdf model
00132 // This is properly cast.
00133 const vpdfl_gaussian_kernel_pdf& vpdfl_gaussian_kernel_pdf_sampler::gaussian_kernel_pdf() const
00134 {
00135   return static_cast<const vpdfl_gaussian_kernel_pdf&>( model());
00136 }
00137