Docs | All | Real | Category | PDF | DataSet | Plot | Container | Misc | Aux | User

RooFit Toolkit for Data Modeling
#include "RooNumConvPdf.hh"

RooNumConvPdf


class description - source file - inheritance tree (.pdf)

class RooNumConvPdf : public RooAbsPdf

Inheritance Chart:
TObject
<-
TNamed
RooPrintable
<-
RooAbsArg
<-
RooAbsReal
<-
RooAbsPdf
<-
RooNumConvPdf

    protected:
RooNumConvolution& conv() const virtual RooAbsGenContext* genContext(const RooArgSet& vars, const RooDataSet* prototype = 0, const RooArgSet* auxProto = 0, Bool_t verbose = kFALSE) const void initialize() const public:
RooNumConvPdf(const char* name, const char* title, RooRealVar& convVar, RooAbsPdf& pdf, RooAbsPdf& resmodel) RooNumConvPdf(const RooNumConvPdf& other, const char* name = "0") virtual ~RooNumConvPdf() static TClass* Class() void clearConvolutionWindow() virtual TObject* clone(const char* newname) const RooNumIntConfig& convIntConfig() virtual Double_t evaluate() const virtual TClass* IsA() const RooAbsReal& model() const RooAbsReal& pdf() const const TH2* profileData() const void setCallProfiling(Bool_t flag, Int_t nbinX = 40, Int_t nbinCall = 40, Int_t nCallHigh = 1000) void setCallWarning(Int_t threshold = 2000) void setConvolutionWindow(RooAbsReal& centerParam, RooAbsReal& widthParam, Double_t widthScaleFactor = 1) virtual void ShowMembers(TMemberInspector& insp, char* parent) virtual void Streamer(TBuffer& b) void StreamerNVirtual(TBuffer& b) RooRealVar& var() const

Data Members


    protected:
Bool_t _init RooNumConvolution* _conv RooRealProxy _origVar Original convolution variable RooRealProxy _origPdf Original input PDF RooRealProxy _origModel Original resolution model

Class Description

 Numeric 1-dimensional convolution operator PDF. This class can convolve any PDF
 with any other PDF

 This class should not be used blindly as numeric convolution is computing
 intensive and prone to stability fitting problems. If an analytic convolution
 can be calculated, you should use that or implement it if not available.
 RooNumConvPdf implements reasonable defaults that should convolve most
 functions reasonably well, but results strongly depend on the shape of your
 input PDFS so always check your result.

 The default integration engine for the numeric convolution is the
 adaptive Gauss-Kronrod method, which empirically seems the most robust
 for this task. You can override the convolution integration settings via
 the RooNumIntConfig object reference returned by the convIntConfig() member
 function

 By default the numeric convolution is integrated from -infinity to
 +infinity through a x -> 1/x coordinate transformation of the
 tails. For convolution with a very small bandwidth it may be
 advantageous (for both CPU consumption and stability) if the
 integration domain is limited to a finite range. The function
 setConvolutionWindow(mean,width,scale) allows to set a sliding
 window around the x value to be calculated taking a RooAbsReal
 expression for an offset and a width to be taken around the x
 value. These input expression can be RooFormulaVars or other
 function objects although the 3d 'scale' argument 'scale'
 multiplies the width RooAbsReal expression given in the 2nd
 argument, allowing for an appropriate window definition for most
 cases without need for a RooFormulaVar object: e.g. a Gaussian
 resolution PDF do setConvolutionWindow(gaussMean,gaussSigma,5)
 Note that for a 'wide' Gaussian the -inf to +inf integration
 may converge more quickly than that over a finite range!

 The default numeric precision is 1e-7, i.e. the global default for
 numeric integration but you should experiment with this value to
 see if it is sufficient for example by studying the number of function
 calls that MINUIT needs to fit your function as function of the
 convolution precision.

RooNumConvPdf(const char *name, const char *title, RooRealVar& convVar, RooAbsPdf& pdf, RooAbsPdf& resmodel) : RooAbsPdf(name,title), _init(kFALSE), _conv(0), _origVar("origVar","Original Convolution variable",this,convVar), _origPdf("origPdf","Original Input PDF",this,pdf), _origModel("origModel","Original Resolution model",this,resmodel)
 Constructor of convolution operator PDF

 convVar  :  convolution variable (on which both pdf and resmodel should depend)
 pdf      :  input 'physics' pdf
 resmodel :  input 'resultion' pdf

 output is pdf(x) (X) resmodel(x) = Int [ pdf(x') resmodel (x-x') ] dx'


RooNumConvPdf(const RooNumConvPdf& other, const char* name) : RooAbsPdf(other,name), _init(kFALSE), _origVar("origVar",this,other._origVar), _origPdf("origPdf",this,other._origPdf), _origModel("origModel",this,other._origModel)
 Copy constructor

~RooNumConvPdf()
 Destructor

Double_t evaluate() const

void initialize() const
 Save pointer to any prototype convolution object (only present if this object is made through
 a copy constructor)

RooAbsGenContext* genContext(const RooArgSet &vars, const RooDataSet *prototype, const RooArgSet* auxProto, Bool_t verbose) const



Inline Functions


                  TObject* clone(const char* newname) const
          RooNumIntConfig& convIntConfig()
                      void clearConvolutionWindow()
                      void setConvolutionWindow(RooAbsReal& centerParam, RooAbsReal& widthParam, Double_t widthScaleFactor = 1)
                      void setCallWarning(Int_t threshold = 2000)
                      void setCallProfiling(Bool_t flag, Int_t nbinX = 40, Int_t nbinCall = 40, Int_t nCallHigh = 1000)
                const TH2* profileData() const
               RooRealVar& var() const
               RooAbsReal& pdf() const
               RooAbsReal& model() const
        RooNumConvolution& conv() const
                   TClass* Class()
                   TClass* IsA() const
                      void ShowMembers(TMemberInspector& insp, char* parent)
                      void Streamer(TBuffer& b)
                      void StreamerNVirtual(TBuffer& b)
Last CVS Update: v 1.8 2005/06/20 15:44:55 wverkerke Top
Copyright © 2000-2005 University of California, Stanford University

Page maintained by Wouter Verkerke and David Kirkby

SourceForge.net Logo