root/Math/DistSampler.h

0001 // @(#)root/mathcore:$Id$
0002 // Author: L. Moneta Fri Sep 22 15:06:47 2006
0003
0004 /**********************************************************************
0005  *                                                                    *
0006  * Copyright (c) 2006  LCG ROOT Math Team, CERN/PH-SFT                *
0007  *                                                                    *
0008  *                                                                    *
0009  **********************************************************************/
0010
0011 // Header file for class DistSampler
0012
0013 #ifndef ROOT_Math_DistSampler
0014 #define ROOT_Math_DistSampler
0015
0016 #include "Math/IFunctionfwd.h"
0017
0018 #include "Math/WrappedFunction.h"
0019
0020 #include <algorithm>
0021 #include <vector>
0022 #include <cassert>
0023
0024 class TRandom;
0025
0026 namespace ROOT {
0027
0028    namespace Fit {
0029
0030       class DataRange;
0031       class BinData;
0032       class UnBinData;
0033    }
0034
0035    namespace Math {
0036
0037       class DistSamplerOptions;
0038
0039 /**
0040    @defgroup Random Interface classes for Random number generation
0041
0042    Pseudo-random numbers generator classes and for generation of random number distributions.
0043    These classes implement several pseudo-random number generators and method for generation of random numbers
0044    according to arbitrary distributions
0045
0046    @ingroup MathCore
0047
0048 */
0049
0050 /**
0051    Interface class for generic sampling of a distribution,
0052    i.e. generating random numbers according to arbitrary distributions
0053
0054    @ingroup Random
0055 */
0056
0057
0058 class DistSampler {
0059
0060 public:
0061
0062    /// default constructor
0063    DistSampler() : fOwnFunc(false), fRange(nullptr), fFunc(nullptr) {}
0064
0065
0066    /// virtual destructor
0067    virtual ~DistSampler();
0068
0069
0070
0071    /// set the parent function distribution to use for sampling (generic case)
0072    template<class Function>
0073    void SetFunction(Function & func, unsigned int dim) {
0074       WrappedMultiFunction<Function &> wf(func, dim);
0075       fData.resize(dim);
0076       // need to clone to avoid temporary
0077       DoSetFunction(wf,true);
0078    }
0079
0080    /// set the parent function distribution to use for random sampling (one dim case)
0081    virtual void SetFunction(const ROOT::Math::IGenFunction & func)  {
0082       SetFunction<const ROOT::Math::IGenFunction>(func, 1);
0083    }
0084
0085
0086    /// set the parent function distribution to use for random sampling (multi-dim case)
0087    virtual void SetFunction(const ROOT::Math::IMultiGenFunction & func)  {
0088       DoSetFunction(func,false);
0089    }
0090
0091    /// return the dimension of the parent distribution (and the data)
0092    unsigned int NDim() const { return fData.size(); }
0093
0094
0095    /**
0096       Initialize the sampling generator with the given algorithm.
0097       Implemented by the derived classes who needs it
0098       (like UnuranSampler).
0099       If nothing is specified use default algorithm
0100       from DistSamplerOptions::SetDefaultAlgorithm
0101    */
0102    virtual bool Init(const char * =""/* algorithm */) { return true;}
0103
0104    /**
0105       Initialize the generators with the given DistSamplerOption object.
0106       The string will include the algorithm and in case additional options
0107       which can be interpreted by a re-implemented method in the derived class.
0108       The default implementation just calls the above method
0109       passing just the algorithm name
0110    */
0111    virtual bool Init(const DistSamplerOptions & opt );
0112
0113
0114    /**
0115        Set the random engine to be used.
0116        To be implemented by the derived classes who provides
0117        random sampling
0118    */
0119    virtual void SetRandom(TRandom *  ) {}
0120
0121    /**
0122        Set the random seed for the TRandom instances used by the sampler
0123        classes.
0124        To be implemented by the derived classes who provides random sampling
0125    */
0126    virtual void SetSeed(unsigned int /*seed*/ ) {}
0127
0128    /**
0129       Get the random engine used by the sampler.
0130       To be implemented by the derived classes who needs it
0131       Returns zero by default
0132     */
0133    virtual TRandom * GetRandom() { return nullptr; }
0134
0135    /// Set the range in a given dimension.
0136    void SetRange(double xmin, double xmax, int icoord = 0);
0137
0138    /// Set the range for all dimensions.
0139    void SetRange(const double * xmin, const double * xmax);
0140    /// Set the range for all dimensions (use std::vector)
0141    void SetRange(const std::vector<double> & xmin, const std::vector<double> & xmax){
0142       assert(xmin.size() >= NDim() && xmax.size() >= NDim());
0143       SetRange(xmin.data(),xmax.data());
0144    }
0145
0146    /// Set the range using the ROOT::Fit::DataRange class.
0147    void SetRange(const ROOT::Fit::DataRange & range);
0148
0149    /// Set the mode of the distribution (1D case).
0150    /// It could be useful or needed by some sampling methods.
0151    /// It is implemented by derived classes if needed (e.g. TUnuranSampler)
0152    virtual void SetMode(double  ) {}
0153
0154    /// Set the mode of the distribution (Multi-dim case).
0155    virtual void SetMode(const std::vector<double> &) {}
0156
0157    /// Set the normalization area of distribution.
0158    /// Implemented by derived classes if needed
0159    virtual void SetArea(double) {}
0160
0161    /// Use the log of the provided pdf.
0162    /// Implemented by the derived classes
0163    virtual void SetUseLogPdf(bool = true) {}
0164
0165    /// Set usage of Derivative of PDF.
0166    /// Can be implemented by derived class
0167    virtual void SetDPdf(const ROOT::Math::IGenFunction & ) {}
0168
0169    /// Set usage of Cumulative of PDF.
0170    /// Can be implemented by derived class
0171    virtual void SetCdf(const ROOT::Math::IGenFunction &) {}
0172
0173    /// Get the parent distribution function (must be called after setting the function).
0174    const ROOT::Math::IMultiGenFunction & ParentPdf() const {
0175       return *fFunc;
0176    }
0177
0178    /// Check if there is a parent distribution defined.
0179    bool HasParentPdf() const { return fFunc != nullptr; }
0180
0181    /**
0182       Sample one event in one dimension.
0183       Specialized implementation could be provided by the derived classes
0184    */
0185    virtual double Sample1D() {
0186       Sample(&fData[0]);
0187       return fData[0];
0188    }
0189
0190    /**
0191       Sample one event and return an array x with
0192       sample coordinates values.
0193     */
0194    const double *  Sample() {
0195       Sample(&fData[0]);
0196       return &fData.front();
0197    }
0198
0199    /**
0200       Sample one event in multi-dimension by filling the given array.
0201       Return false if the sampling failed.
0202       Abstract method to be re-implemented by the derived classes
0203    */
0204    virtual bool Sample(double * x) = 0;
0205
0206    /**
0207       Sample one bin given an estimate of the pdf in the bin.
0208       (this can be function value at the center or its integral in the bin
0209       divided by the bin width)
0210       By default do not do random sample, just return the function values
0211       Typically Poisson statistics will be used
0212     */
0213    virtual bool SampleBin(double prob, double & value, double * error = nullptr) {
0214       value = prob;
0215       if (error) *error = 0;
0216       return true;
0217    }
0218    /**
0219       Sample a set of bins given a vector of probabilities
0220       Typically multinomial statistics will be used and the sum of the probabilities
0221       will be equal to the total number of events to be generated
0222       For sampling the bins independently, SampleBin should be used
0223     */
0224    virtual bool SampleBins(unsigned int n, const double * prob, double * values, double * errors  = nullptr)  {
0225       std::copy(prob,prob+n, values); // default impl returns prob values (Asimov data)
0226       if (errors) std::fill(errors,errors+n,0);
0227       return true;
0228    }
0229
0230
0231    /**
0232       Generate a un-binned data set by filling the given data set object.
0233       If the data set object is not empty, the new generated data will be appended to the
0234       existing one.
0235    */
0236    virtual bool Generate(unsigned int nevt, ROOT::Fit::UnBinData & data);
0237    /**
0238       Generate a vector of events by filling the passed data vector.
0239       The flag eventRow indicates how the events are arranged in the multi-dim case.
0240       The can be arranged in rows or in columns.
0241       With eventRow=false events are the columns in data: {x1,x2,.....,xn},{y1,....yn}
0242       With eventRow=true  events are rows in data: {x1,y1},{x2,y2},.....{xn,yn}
0243    */
0244    virtual bool Generate(unsigned int nevt, double * data, bool eventRow = false);
0245
0246    /**
0247       Generate a binned data set.
0248       A range must have been set before (otherwise inf is returned)
0249       and the bins are equidistant in the previously defined range
0250       bin center values must be present in given data set
0251       If the sampler is implemented by a random one, the entries
0252       will be binned according to the Poisson distribution
0253       It is assumed the distribution is normalized, otherwise the nevt must be scaled
0254       accordingly. The expected value/bin nexp  = f(x_i) * binArea/ nevt
0255       Extend control if use a fixed (i.e. multinomial statistics) or floating total number of events
0256    */
0257    virtual bool Generate(unsigned int nevt, const  int * nbins, ROOT::Fit::BinData & data, bool extend = true, bool expErr = true);
0258    /**
0259       Same as before but passing the range in case of 1 dim data.
0260     */
0261    bool Generate(unsigned int nevt, int nbins, double xmin, double xmax, ROOT::Fit::BinData & data, bool extend = true, bool expErr = true ) {
0262       SetRange(xmin,xmax);
0263       int nbs[1]; nbs[0] = nbins;
0264       return Generate(nevt, nbs, data, extend, expErr);
0265    }
0266
0267
0268 protected:
0269
0270    // internal method to set the function
0271    virtual void DoSetFunction(const ROOT::Math::IMultiGenFunction & func, bool copy);
0272    // internal method to set the dimension
0273    virtual void DoSetDimension(unsigned int ndim);
0274    // check if generator have been initialized correctly and one can start generating
0275    bool IsInitialized() ;
0276    /// return the data range of the Pdf . Must be called after setting the function
0277    const ROOT::Fit::DataRange & PdfRange() const {
0278       assert(fRange);
0279       return *fRange;
0280    }
0281
0282 private:
0283
0284    // private methods
0285
0286    bool fOwnFunc;                               ///< flag to indicate if the function is owned
0287    mutable std::vector<double> fData;           ///<! internal array used to cached the sample data
0288    ROOT::Fit::DataRange    *   fRange;          ///< data range
0289    const ROOT::Math::IMultiGenFunction * fFunc; ///< internal function (ND)
0290
0291
0292 };
0293
0294    } // end namespace Math
0295
0296 } // end namespace ROOT
0297
0298
0299 #endif /* ROOT_Math_DistSampler */