Logo ROOT   6.10/09
Reference Guide
RooIntegralMorph.h
Go to the documentation of this file.
1 /*****************************************************************************
2  * Project: RooFit *
3  * *
4  * Copyright (c) 2000-2007, Regents of the University of California *
5  * and Stanford University. All rights reserved. *
6  * *
7  * Redistribution and use in source and binary forms, *
8  * with or without modification, are permitted according to the terms *
9  * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
10  *****************************************************************************/
11 
12 #ifndef ROOLINEARMORPH
13 #define ROOLINEARMORPH
14 
15 #include "RooAbsCachedPdf.h"
16 #include "RooRealProxy.h"
17 #include "RooCategoryProxy.h"
18 #include "RooAbsReal.h"
19 #include "RooAbsCategory.h"
20 class RooBrentRootFinder ;
21 
22 class TH1D ;
23 
25 public:
26  RooIntegralMorph() : _cache(nullptr) {
27  // coverity[UNINIT_CTOR]
28  } ;
29  RooIntegralMorph(const char *name, const char *title,
30  RooAbsReal& _pdf1,
31  RooAbsReal& _pdf2,
32  RooAbsReal& _x,
34  RooIntegralMorph(const RooIntegralMorph& other, const char* name=0) ;
35  virtual TObject* clone(const char* newname) const { return new RooIntegralMorph(*this,newname); }
36  inline virtual ~RooIntegralMorph() { }
37 
39  // P.d.f is self normalized
40  return kTRUE ;
41  }
42  void setCacheAlpha(Bool_t flag) {
43  // Activate caching of p.d.f. shape for all values of alpha as well
44  _cacheMgr.sterilize() ; _cacheAlpha = flag ;
45  }
46  Bool_t cacheAlpha() const {
47  // If true caching of p.d.f for all alpha values is active
48  return _cacheAlpha ;
49  }
50 
51  virtual void preferredObservableScanOrder(const RooArgSet& obs, RooArgSet& orderedObs) const ;
52 
54  public:
56  ~MorphCacheElem() ;
57  void calculate(TIterator* iter) ;
59 
60  protected:
61 
62  void findRange() ;
64  Int_t binX(Double_t x) ;
65  void fillGap(Int_t ixlo, Int_t ixhi,Double_t splitPoint=0.5) ;
66  void interpolateGap(Int_t ixlo, Int_t ixhi) ;
67 
70  RooAbsPdf* _pdf1 ; // PDF1
71  RooAbsPdf* _pdf2 ; // PDF2
72  RooRealVar* _x ; // X
73  RooAbsReal* _alpha ; // ALPHA
74  RooAbsReal* _c1 ; // CDF of PDF 1
75  RooAbsReal* _c2 ; // CDF of PDF 2
76  RooAbsFunc* _cb1 ; // Binding of CDF1
77  RooAbsFunc* _cb2 ; // Binding of CDF2
78  RooBrentRootFinder* _rf1 ; // ROOT finder on CDF1
79  RooBrentRootFinder* _rf2 ; // ROOT finder of CDF2 ;
80 
81  std::vector<Double_t> _yatX ; //
82  std::vector<Double_t> _calcX; //
85 
87 
88  } ;
89 
90 protected:
91 
92  friend class MorphCacheElem ;
93  virtual PdfCacheElem* createCache(const RooArgSet* nset) const ;
94  virtual const char* inputBaseName() const ;
95  virtual RooArgSet* actualObservables(const RooArgSet& nset) const ;
96  virtual RooArgSet* actualParameters(const RooArgSet& nset) const ;
97  virtual void fillCacheObject(PdfCacheElem& cache) const ;
98 
99  RooRealProxy pdf1 ; // First input shape
100  RooRealProxy pdf2 ; // Second input shape
101  RooRealProxy x ; // Observable
102  RooRealProxy alpha ; // Interpolation parameter
103  Bool_t _cacheAlpha ; // If true, both (x,alpha) are cached
104  mutable MorphCacheElem* _cache ; // Current morph cache element in use
105 
106 
107  Double_t evaluate() const ;
108 
109 private:
110 
111  ClassDef(RooIntegralMorph,1) // Linear shape interpolation operator p.d.f
112 };
113 
114 #endif
void setCacheAlpha(Bool_t flag)
MorphCacheElem * _cache
virtual TObject * clone(const char *newname) const
std::vector< Double_t > _calcX
std::vector< Double_t > _yatX
virtual RooArgSet * actualObservables(const RooArgSet &nset) const
Observable to be cached for given choice of normalization.
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
MorphCacheElem(RooIntegralMorph &self, const RooArgSet *nset)
Construct of cache element, copy relevant input from RooIntegralMorph, create the cdfs from the input...
Iterator abstract base class.
Definition: TIterator.h:30
Bool_t cacheAlpha() const
virtual RooArgSet * actualParameters(const RooArgSet &nset) const
Parameters of the cache.
#define ClassDef(name, id)
Definition: Rtypes.h:297
Bool_t selfNormalized() const
Int_t binX(Double_t x)
Return the bin number enclosing the given x value.
Class RooIntegralMorph is an implementation of the histogram interpolation technique described by Ale...
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
Double_t calcX(Double_t y, Bool_t &ok)
Calculate the x value of the output p.d.f at the given cdf value y.
virtual const char * inputBaseName() const
Return base name component for cache components in this case a string encoding the names of both end ...
void findRange()
Determine which range of y values can be mapped to x values from the numeric inversion of the input c...
Implement the abstract 1-dimensional root finding interface using the Brent-Decker method...
tomato 1-D histogram with a double per channel (see TH1 documentation)}
Definition: TH1.h:594
virtual void fillCacheObject(PdfCacheElem &cache) const
Fill the cache with the interpolated shape.
const Bool_t kFALSE
Definition: RtypesCore.h:92
void fillGap(Int_t ixlo, Int_t ixhi, Double_t splitPoint=0.5)
Fill all empty histogram bins between bins ixlo and ixhi.
RooAbsCachedPdf is the abstract base class for p.d.f.s that need or want to cache their evaluate() ou...
double Double_t
Definition: RtypesCore.h:55
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
RooRealProxy pdf1
Double_t y[n]
Definition: legend1.C:17
RooRealProxy alpha
virtual void preferredObservableScanOrder(const RooArgSet &obs, RooArgSet &orderedObs) const
Indicate to the RooAbsCachedPdf base class that for the filling of the cache the traversal of the x s...
Mother of all ROOT objects.
Definition: TObject.h:37
virtual RooArgList containedArgs(Action)
Return all RooAbsArg components contained in this cache.
Double_t evaluate() const
Dummy.
void interpolateGap(Int_t ixlo, Int_t ixhi)
Fill empty histogram bins between ixlo and ixhi with values obtained from linear interpolation of ixl...
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
RooRealProxy is the concrete proxy for RooAbsReal objects A RooRealProxy is the general mechanism to ...
Definition: RooRealProxy.h:23
virtual ~RooIntegralMorph()
RooObjCacheManager _cacheMgr
virtual PdfCacheElem * createCache(const RooArgSet *nset) const
Create and return a derived MorphCacheElem.
const Bool_t kTRUE
Definition: RtypesCore.h:91
Abstract interface for evaluating a real-valued function of one real variable and performing numerica...
Definition: RooAbsFunc.h:23
void calculate(TIterator *iter)
Calculate shape of p.d.f for x,alpha values defined by dIter iterator over cache histogram.