library: libRooFit
#include "RooAbsAnaConvPdf.h"

RooAbsAnaConvPdf


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

class RooAbsAnaConvPdf : public RooAbsPdf

Inheritance Chart:
TObject
<-
TNamed
RooPrintable
<-
RooAbsArg
<-
RooAbsReal
<-
RooAbsPdf
<-
RooAbsAnaConvPdf
<-
RooBCPEffDecay
RooBCPGenDecay
RooBDecay
RooBMixDecay
RooDecay
RooNonCPEigenDecay
 
    This is an abstract class, constructors will not be documented.
    Look at the header to check for available constructors.


    protected:
Bool_t changeModel(const RooResolutionModel& newModel) const RooRealVar* convVar() const virtual Double_t evaluate() const virtual RooAbsGenContext* genContext(const RooArgSet& vars, const RooDataSet* prototype = 0, const RooArgSet* auxProto = 0, Bool_t verbose = kFALSE) const void makeCoefVarList() const RooArgSet* parseIntegrationRequest(const RooArgSet& intSet, Int_t& coefCode, RooArgSet* analVars = 0) const virtual void syncNormalizationPostHook(RooAbsReal* norm, const RooArgSet* nset) const virtual Bool_t syncNormalizationPreHook(RooAbsReal* norm, const RooArgSet* nset) const public:
virtual ~RooAbsAnaConvPdf() virtual Double_t analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName = "0") const static TClass* Class() virtual Double_t coefAnalyticalIntegral(Int_t coef, Int_t code, const char* rangeName = "0") const virtual Double_t coefficient(Int_t basisIndex) const Int_t declareBasis(const char* expression, const RooArgList& params) virtual Bool_t forceAnalyticalInt(const RooAbsArg& dep) const virtual Int_t getAnalyticalIntegralWN(RooArgSet& allVars, RooArgSet& analVars, const RooArgSet* normSet, const char* rangeName = "0") const virtual Int_t getCoefAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName = "0") const Double_t getCoefNorm(Int_t coefIdx, const RooArgSet& nset, const char* rangeName) const Double_t getCoefNorm(Int_t coefIdx, const RooArgSet* nset = 0, const char* rangeName = "0") const virtual TClass* IsA() const virtual Bool_t isDirectGenSafe(const RooAbsArg& arg) const virtual void printToStream(ostream& stream, RooPrintable::PrintOption opt = Standard, TString indent = ) const virtual void ShowMembers(TMemberInspector& insp, char* parent) virtual void Streamer(TBuffer& b) void StreamerNVirtual(TBuffer& b)

Data Members


    protected:
Bool_t _isCopy RooResolutionModel* _model Original resolution model RooRealVar* _convVar Convolution variable RooListProxy _convSet Set of (resModel (x) basisFunc) convolution objects RooArgList _basisList List of created basis functions RooArgSet* _convNormSet Subset of last normalization that applies to convolutions TIterator* _convSetIter ! Iterator over _convNormSet RooArgList _coefVarList RooNormListManager _coefNormMgr ! Coefficient normalization manager RooAICRegistry _codeReg

Class Description

  RooAbsAnaConvPdf is the base class of for PDFs that represents a
  physics model that can be analytically convolved with a resolution model

  To achieve factorization between the physics model and the resolution
  model, each physics model must be able to be written in the form
           _ _                 _              _
    Phys(x,a,b) = Sum_k coef_k(a) * basis_k(x,b)

  where basis_k are a limited number of functions in terms of the variable
  to be convoluted and coef_k are coefficients independent of the convolution
  variable.

  Classes derived from RooResolutionModel implement
         _ _                        _                  _
   R_k(x,b,c) = Int(dx') basis_k(x',b) * resModel(x-x',c)

  which RooAbsAnaConvPdf uses to construct the pdf for [ Phys (x) R ] :
          _ _ _                 _          _ _
    PDF(x,a,b,c) = Sum_k coef_k(a) * R_k(x,b,c)

  A minimal implementation of a RooAbsAnaConvPdf physics model consists of

  - A constructor that declares the required basis functions using the declareBasis() method.
    The declareBasis() function assigns a unique identifier code to each declare basis

  - An implementation of coefficient(Int_t code) returning the coefficient value for each
    declared basis function

  Optionally, analytical integrals can be provided for the coefficient functions. The
  interface for this is quite similar to that for integrals of regular PDFs. Two functions,

   Int_t getCoefAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars)
   Double_t coefAnalyticalIntegral(Int_t coef, Int_t code),

  advertise the coefficient integration capabilities and implement them respectively.
  Please see RooAbsPdf for additional details. Advertised analytical integrals must be
  valid for all coefficients.

~RooAbsAnaConvPdf()
 Destructor

Int_t declareBasis(const char* expression, const RooArgList& params)
 Declare a basis function for use in this physics model. The string expression
 must be a valid RooFormulVar expression representing the basis function, referring
 to the convolution variable as '@0', and any additional parameters (supplied in
 'params' as '@1','@2' etc.

 The return value is a unique identifier code, that will be passed to coefficient()
 to identify the basis function for which the coefficient is requested. If the
 resolution model used does not support the declared basis function, code -1 is
 returned.


Bool_t changeModel(const RooResolutionModel& newModel)
 Change the resolution model to given model

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

Bool_t isDirectGenSafe(const RooAbsArg& arg) const
 All direct generation of convolution arg if model is truth model

const RooRealVar* convVar() const
 Return a pointer to the convolution variable instance used in the resolution model

Double_t evaluate() const
 Calculate the current unnormalized value of the PDF

 PDF = sum_k coef_k * [ basis_k (x) ResModel ]


Int_t getAnalyticalIntegralWN(RooArgSet& allVars, RooArgSet& analVars, const RooArgSet* normSet2, const char* /*rangeName*/) const
 Handle trivial no-integration scenario

Double_t analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName) const
 Return analytical integral defined by given scenario code.

 For unnormalized integrals this is
                    _                _
   PDF = sum_k Int(dx) coef_k * Int(dy) [ basis_k (x) ResModel ].
       _
 where x is the set of coefficient dependents to be integrated
 and y the set of basis function dependents to be integrated.

 For normalized integrals this becomes

         sum_k Int(dx) coef_k * Int(dy) [ basis_k (x) ResModel ].
  PDF =  --------------------------------------------------------
         sum_k Int(dv) coef_k * Int(dw) [ basis_k (x) ResModel ].

 where x is the set of coefficient dependents to be integrated,
 y the set of basis function dependents to be integrated,
 v is the set of coefficient dependents over which is normalized and
 w is the set of basis function dependents over which is normalized.

 Set x must be contained in v and set y must be contained in w.


Int_t getCoefAnalyticalIntegral(RooArgSet& /*allVars*/, RooArgSet& /*analVars*/, const char* /*rangeName*/) const
 Default implementation of function advertising integration capabilities: no integrals
 are advertised.

Double_t coefAnalyticalIntegral(Int_t coef, Int_t code, const char* /*rangeName*/) const
 Default implementation of function implementing advertised integrals. Only
 the pass-through scenario (no integration) is implemented.

Bool_t forceAnalyticalInt(const RooAbsArg& /*dep*/) const
 This function forces RooRealIntegral to offer all integration dependents
 to RooAbsAnaConvPdf::getAnalyticalIntegralWN() for consideration for
 analytical integration, if RRI considers this to be unsafe (e.g. due
 to hidden Jacobian terms).

 RooAbsAnaConvPdf will not attempt to actually integrate all these dependents
 but feed them to the resolution models integration interface, which will
 make the final determination on how to integrate these dependents.

Double_t getCoefNorm(Int_t coefIdx, const RooArgSet* nset, const char* rangeName) const

void makeCoefVarList() const
 Build complete list of coefficient variables

Bool_t syncNormalizationPreHook(RooAbsReal* /*norm*/,const RooArgSet* nset) const
 Overload of hook function in RooAbsPdf::syncNormalization(). This functions serves
 two purposes:

   - Modify default normalization behaviour of RooAbsPdf: integration requests over
     unrelated variables are properly executed (introducing a trivial multiplication
     for each unrelated dependent). This is necessary if composite resolution models
     are used in which the components do not necessarily all have the same set
     of dependents.

   - Built the sub set of normalization dependents that is contained the basis function
     resolution model convolution (to be used in syncNormalizationPostHook().

void syncNormalizationPostHook(RooAbsReal* /*norm*/,const RooArgSet* /*nset*/) const
 Overload of hook function in RooAbsPdf::syncNormalization(). This function propagates
 the syncNormalization() call to all basis-function/resolution-model convolution component
 objects and fixes the physics models client-server links by adding each variable that
 serves any of the convolution objects normalizations. PDFs by default have all client-server
 links that control the unnormalized value (as returned by evaluate()), but convoluted PDFs
 have a non-trivial normalization term that may introduce dependencies on additional server
 that exclusively appear in the normalization.

void printToStream(ostream& os, PrintOption opt, TString indent) const
 Print info about this object to the specified stream. In addition to the info
 from RooAbsPdf::printToStream() we add:

   Verbose : detailed information on convolution integrals



Inline Functions


           Double_t getCoefNorm(Int_t coefIdx, const RooArgSet* nset = 0, const char* rangeName = "0") const
           Double_t coefficient(Int_t basisIndex) const
         RooArgSet* parseIntegrationRequest(const RooArgSet& intSet, Int_t& coefCode, RooArgSet* analVars = 0) const
            TClass* Class()
            TClass* IsA() const
               void ShowMembers(TMemberInspector& insp, char* parent)
               void Streamer(TBuffer& b)
               void StreamerNVirtual(TBuffer& b)


Last update: Tue Jun 28 18:08:46 2005
Copyright (c) 2000-2005, Regents of the University of California *


ROOT page - Class index - Class Hierarchy - Top of the page

This page has been automatically generated. If you have any comments or suggestions about the page layout send a mail to ROOT support, or contact the developers with any questions or problems regarding ROOT.