Main MRPT website > C++ reference
MRPT logo

CRejectionSamplingCapable.h

Go to the documentation of this file.
00001 /* +---------------------------------------------------------------------------+
00002    |          The Mobile Robot Programming Toolkit (MRPT) C++ library          |
00003    |                                                                           |
00004    |                   http://mrpt.sourceforge.net/                            |
00005    |                                                                           |
00006    |   Copyright (C) 2005-2011  University of Malaga                           |
00007    |                                                                           |
00008    |    This software was written by the Machine Perception and Intelligent    |
00009    |      Robotics Lab, University of Malaga (Spain).                          |
00010    |    Contact: Jose-Luis Blanco  <jlblanco@ctima.uma.es>                     |
00011    |                                                                           |
00012    |  This file is part of the MRPT project.                                   |
00013    |                                                                           |
00014    |     MRPT is free software: you can redistribute it and/or modify          |
00015    |     it under the terms of the GNU General Public License as published by  |
00016    |     the Free Software Foundation, either version 3 of the License, or     |
00017    |     (at your option) any later version.                                   |
00018    |                                                                           |
00019    |   MRPT is distributed in the hope that it will be useful,                 |
00020    |     but WITHOUT ANY WARRANTY; without even the implied warranty of        |
00021    |     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the         |
00022    |     GNU General Public License for more details.                          |
00023    |                                                                           |
00024    |     You should have received a copy of the GNU General Public License     |
00025    |     along with MRPT.  If not, see <http://www.gnu.org/licenses/>.         |
00026    |                                                                           |
00027    +---------------------------------------------------------------------------+ */
00028 #ifndef CRejectionSamplingCapable_H
00029 #define CRejectionSamplingCapable_H
00030 
00031 #include <mrpt/utils/utils_defs.h>
00032 #include <mrpt/bayes/CProbabilityParticle.h>
00033 #include <mrpt/random.h>
00034 
00035 namespace mrpt
00036 {
00037 namespace bayes
00038 {
00039         /** A base class for implementing rejection sampling in a generic state space.
00040          *   See the main method CRejectionSamplingCapable::rejectionSampling
00041          *  To use this class, create your own class as a child of this one and implement the desired
00042          *   virtual methods, and add any required internal data.
00043          */
00044         template <class TStateSpace>
00045         class CRejectionSamplingCapable
00046         {
00047         public:
00048                 typedef CProbabilityParticle<TStateSpace> TParticle;
00049 
00050         /** Virtual destructor
00051           */
00052         virtual ~CRejectionSamplingCapable()
00053                 {
00054                 }
00055 
00056                 /** Generates a set of N independent samples via rejection sampling.
00057                   * \param desiredSamples The number of desired samples to generate
00058                   * \param outSamples The output samples.
00059                   * \param timeoutTrials The maximum number of rejection trials for each generated sample (i.e. the maximum number of iterations). This can be used to set a limit to the time complexity of the algorithm for difficult probability densities.
00060                   *  All will have equal importance weights (a property of rejection sampling), although those samples
00061                   *   generated at timeout will have a different importance weights.
00062                   */
00063                 void rejectionSampling(
00064                         size_t                                                  desiredSamples,
00065                         std::vector<TParticle>                  &outSamples,
00066                         size_t                                                  timeoutTrials = 1000)
00067                 {
00068                         MRPT_START;
00069 
00070                         TStateSpace                                                     x;
00071                         typename std::vector<TParticle>::iterator       it;
00072 
00073                         // Set output size:
00074                         if ( outSamples.size() != desiredSamples )
00075                         {
00076                                 // Free old memory:
00077                                 for (it = outSamples.begin();it!=outSamples.end();it++)
00078                                         delete (it->d);
00079                                 outSamples.clear();
00080 
00081                                 // Reserve new memory:
00082                                 outSamples.resize( desiredSamples );
00083                                 for (it = outSamples.begin();it!=outSamples.end();it++)
00084                                         it->d = new TStateSpace;
00085                         }
00086 
00087                         // Rejection sampling loop:
00088                         double  acceptanceProb;
00089                         for (it = outSamples.begin();it!=outSamples.end();it++)
00090                         {
00091                                 size_t  timeoutCount = 0;
00092                                 double          bestLik = -1e250;
00093                                 TStateSpace     bestVal;
00094                                 do
00095                                 {
00096                                         RS_drawFromProposal( *it->d );
00097                                         acceptanceProb = RS_observationLikelihood( *it->d );
00098                                         ASSERT_(acceptanceProb>=0 && acceptanceProb<=1);
00099                                         if (acceptanceProb>bestLik)
00100                                         {
00101                                                 bestLik = acceptanceProb;
00102                                                 bestVal = *it->d;
00103                                         }
00104                                 } while (       acceptanceProb < mrpt::random::randomGenerator.drawUniform(0.0,0.999) &&
00105                                                         (++timeoutCount)<timeoutTrials );
00106 
00107                                 // Save weights:
00108                                 if (timeoutCount>=timeoutTrials)
00109                                 {
00110                                         it->log_w = log(bestLik);
00111                                         *it->d    = bestVal;
00112                                 }
00113                                 else
00114                                 {
00115                                         it->log_w = 0; // log(1.0);
00116                                 }
00117                         } // end for it
00118 
00119                         MRPT_END;
00120                 }
00121 
00122         protected:
00123                 /** Generates one sample, drawing from some proposal distribution.
00124                   */
00125                 virtual void RS_drawFromProposal( TStateSpace &outSample ) = 0;
00126 
00127                 /** Returns the NORMALIZED observation likelihood (linear, not exponential!!!) at a given point of the state space (values in the range [0,1]).
00128                   */
00129                 virtual double RS_observationLikelihood( const TStateSpace &x) = 0;
00130 
00131         }; // End of class def.
00132 
00133 } // End of namespace
00134 } // End of namespace
00135 
00136 #endif



Page generated by Doxygen 1.7.3 for MRPT 0.9.4 SVN: at Sat Mar 26 06:40:17 UTC 2011