CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Classes | Public Member Functions | Protected Member Functions | Protected Attributes
VerticalInterpPdf Class Reference

#include <VerticalInterpPdf.h>

Inheritance diagram for VerticalInterpPdf:

Classes

class  CacheElem
 

Public Member Functions

Double_t analyticalIntegralWN (Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
 
virtual Bool_t checkObservables (const RooArgSet *nset) const
 
virtual TObject * clone (const char *newname) const
 
const RooArgList & coefList () const
 
Double_t evaluate () const
 
virtual Bool_t forceAnalyticalInt (const RooAbsArg &) const
 
const RooArgList & funcList () const
 
Int_t getAnalyticalIntegralWN (RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet, const char *rangeName=0) const
 
 VerticalInterpPdf ()
 
 VerticalInterpPdf (const char *name, const char *title, const RooArgList &funcList, const RooArgList &coefList, Double_t quadraticRegion=0., Int_t quadraticAlgo=0)
 
 VerticalInterpPdf (const VerticalInterpPdf &other, const char *name=0)
 
virtual ~VerticalInterpPdf ()
 

Protected Member Functions

Double_t interpolate (Double_t coeff, Double_t central, RooAbsReal *fUp, RooAbsReal *fDown) const
 Iterator over coefficient list. More...
 

Protected Attributes

TIterator * _coefIter
 Iterator over FUNC list. More...
 
RooListProxy _coefList
 
TIterator * _funcIter
 
RooListProxy _funcList
 
RooObjCacheManager _normIntMgr
 
Int_t _quadraticAlgo
 
Double_t _quadraticRegion
 

Detailed Description

Vertical interpolation between multiple histograms (or pdfs). Based on RooRealSumPdf

Definition at line 12 of file VerticalInterpPdf.h.

Constructor & Destructor Documentation

VerticalInterpPdf::VerticalInterpPdf ( )

Definition at line 22 of file VerticalInterpPdf.cc.

Referenced by clone().

23 {
24  // Default constructor
25  // coverity[UNINIT_CTOR]
26  _funcIter = _funcList.createIterator() ;
27  _coefIter = _coefList.createIterator() ;
28  _quadraticRegion = 0;
29 }
RooListProxy _funcList
RooListProxy _coefList
TIterator * _coefIter
Iterator over FUNC list.
VerticalInterpPdf::VerticalInterpPdf ( const char *  name,
const char *  title,
const RooArgList &  funcList,
const RooArgList &  coefList,
Double_t  quadraticRegion = 0.,
Int_t  quadraticAlgo = 0 
)

Definition at line 33 of file VerticalInterpPdf.cc.

References _coefIter, _coefList, _funcIter, _funcList, and _quadraticAlgo.

33  :
34  RooAbsPdf(name,title),
35  _normIntMgr(this,10),
36  _funcList("!funcList","List of functions",this),
37  _coefList("!coefList","List of coefficients",this),
38  _quadraticRegion(quadraticRegion),
39  _quadraticAlgo(quadraticAlgo)
40 {
41 
42  if (inFuncList.getSize()!=2*inCoefList.getSize()+1) {
43  coutE(InputArguments) << "VerticalInterpPdf::VerticalInterpPdf(" << GetName()
44  << ") number of pdfs and coefficients inconsistent, must have Nfunc=1+2*Ncoef" << endl ;
45  assert(0);
46  }
47 
48  TIterator* funcIter = inFuncList.createIterator() ;
49  RooAbsArg* func;
50  while((func = (RooAbsArg*)funcIter->Next())) {
51  if (!dynamic_cast<RooAbsReal*>(func)) {
52  coutE(InputArguments) << "ERROR: VerticalInterpPdf::VerticalInterpPdf(" << GetName() << ") function " << func->GetName() << " is not of type RooAbsReal" << endl;
53  assert(0);
54  }
55  _funcList.add(*func) ;
56  }
57  delete funcIter;
58 
59  TIterator* coefIter = inCoefList.createIterator() ;
60  RooAbsArg* coef;
61  while((coef = (RooAbsArg*)coefIter->Next())) {
62  if (!dynamic_cast<RooAbsReal*>(coef)) {
63  coutE(InputArguments) << "ERROR: VerticalInterpPdf::VerticalInterpPdf(" << GetName() << ") coefficient " << coef->GetName() << " is not of type RooAbsReal" << endl;
64  assert(0);
65  }
66  _coefList.add(*coef) ;
67  }
68  delete coefIter;
69 
70  _funcIter = _funcList.createIterator() ;
71  _coefIter = _coefList.createIterator() ;
72 
73  if (_quadraticAlgo == -1) {
74  // multiplicative morphing: no way to do analytical integrals.
75  _forceNumInt = kTRUE;
76  } else if (_quadraticAlgo >= 100) {
77  _quadraticAlgo -= 100;
78  _forceNumInt = kTRUE;
79  }
80 }
RooListProxy _funcList
RooListProxy _coefList
RooObjCacheManager _normIntMgr
TIterator * _coefIter
Iterator over FUNC list.
VerticalInterpPdf::VerticalInterpPdf ( const VerticalInterpPdf other,
const char *  name = 0 
)

Definition at line 86 of file VerticalInterpPdf.cc.

References _coefIter, _coefList, _funcIter, and _funcList.

86  :
87  RooAbsPdf(other,name),
88  _normIntMgr(other._normIntMgr,this),
89  _funcList("!funcList",this,other._funcList),
90  _coefList("!coefList",this,other._coefList),
93 {
94  // Copy constructor
95 
96  _funcIter = _funcList.createIterator() ;
97  _coefIter = _coefList.createIterator() ;
98 }
RooListProxy _funcList
RooListProxy _coefList
RooObjCacheManager _normIntMgr
TIterator * _coefIter
Iterator over FUNC list.
VerticalInterpPdf::~VerticalInterpPdf ( )
virtual

Definition at line 103 of file VerticalInterpPdf.cc.

References _coefIter, and _funcIter.

104 {
105  // Destructor
106  delete _funcIter ;
107  delete _coefIter ;
108 }
TIterator * _coefIter
Iterator over FUNC list.

Member Function Documentation

Double_t VerticalInterpPdf::analyticalIntegralWN ( Int_t  code,
const RooArgSet *  normSet,
const char *  rangeName = 0 
) const

Definition at line 238 of file VerticalInterpPdf.cc.

References _coefIter, VerticalInterpPdf::CacheElem::_funcIntList, VerticalInterpPdf::CacheElem::_funcNormList, _normIntMgr, interpolate(), and relativeConstraints::value.

239 {
240  // Implement analytical integrations by deferring integration of component
241  // functions to integrators of components
242 
243  // Handle trivial passthrough scenario
244  if (code==0) return getVal(normSet2) ;
245 
246  RooAbsReal *coef;
247  Double_t value = 0;
248 
249  // WVE needs adaptation for rangeName feature
250  CacheElem* cache = (CacheElem*) _normIntMgr.getObjByIndex(code-1) ;
251 
252  TIterator* funcIntIter = cache->_funcIntList.createIterator() ;
253  RooAbsReal *funcInt = (RooAbsReal *) funcIntIter->Next();
254  Double_t central = funcInt->getVal();
255  value += central;
256 
257  _coefIter->Reset() ;
258  while((coef=(RooAbsReal*)_coefIter->Next())) {
259  Double_t coefVal = coef->getVal(normSet2) ;
260  RooAbsReal * funcIntUp = (RooAbsReal*)funcIntIter->Next() ;
261  RooAbsReal * funcIntDn = (RooAbsReal*)funcIntIter->Next() ;
262  value += interpolate(coefVal, central, funcIntUp, funcIntDn);
263  }
264 
265  delete funcIntIter ;
266 
267  Double_t normVal(1) ;
268  if (normSet2) {
269  normVal = 0 ;
270 
271  TIterator* funcNormIter = cache->_funcNormList.createIterator() ;
272 
273  RooAbsReal* funcNorm = (RooAbsReal*) funcNormIter->Next();
274  central = funcNorm->getVal(normSet2) ;
275 
276  _coefIter->Reset() ;
277  while((coef=(RooAbsReal*)_coefIter->Next())) {
278  RooAbsReal *funcNormUp = (RooAbsReal*)funcNormIter->Next() ;
279  RooAbsReal *funcNormDn = (RooAbsReal*)funcNormIter->Next() ;
280  Double_t coefVal = coef->getVal(normSet2) ;
281  normVal += interpolate(coefVal, central, funcNormUp, funcNormDn);
282  }
283 
284  delete funcNormIter ;
285  }
286 
287  return ( value > 0 ? value : 1E-9 ) / normVal;
288 }
Double_t interpolate(Double_t coeff, Double_t central, RooAbsReal *fUp, RooAbsReal *fDown) const
Iterator over coefficient list.
RooObjCacheManager _normIntMgr
TIterator * _coefIter
Iterator over FUNC list.
Bool_t VerticalInterpPdf::checkObservables ( const RooArgSet *  nset) const
virtual

Definition at line 150 of file VerticalInterpPdf.cc.

References _coefIter, _funcIter, and _quadraticAlgo.

151 {
152  // Check if FUNC is valid for given normalization set.
153  // Coeffient and FUNC must be non-overlapping, but func-coefficient
154  // pairs may overlap each other
155  //
156  // In the present implementation, coefficients may not be observables or derive
157  // from observables
158  if (_quadraticAlgo == -1) return false; // multiplicative morphing. we don't care.
159 
160 
161  _coefIter->Reset() ;
162  RooAbsReal* coef ;
163  while((coef=(RooAbsReal*)_coefIter->Next())) {
164  if (coef->dependsOn(*nset)) {
165  coutE(InputArguments) << "RooRealPdf::checkObservables(" << GetName() << "): ERROR coefficient " << coef->GetName()
166  << " depends on one or more of the following observables" ; nset->Print("1") ;
167  return true;
168  }
169  }
170 
171  _funcIter->Reset() ;
172  RooAbsReal* func ;
173  while((func = (RooAbsReal*)_funcIter->Next())) {
174  if (func->observableOverlaps(nset,*coef)) {
175  coutE(InputArguments) << "VerticalInterpPdf::checkObservables(" << GetName() << "): ERROR: coefficient " << coef->GetName()
176  << " and FUNC " << func->GetName() << " have one or more observables in common" << endl ;
177  return true;
178  }
179  }
180 
181  return false;
182 }
TIterator * _coefIter
Iterator over FUNC list.
virtual TObject* VerticalInterpPdf::clone ( const char *  newname) const
inlinevirtual

Definition at line 18 of file VerticalInterpPdf.h.

References VerticalInterpPdf().

18 { return new VerticalInterpPdf(*this,newname) ; }
const RooArgList& VerticalInterpPdf::coefList ( ) const
inline

Definition at line 29 of file VerticalInterpPdf.h.

References _coefList.

29 { return _coefList ; }
RooListProxy _coefList
Double_t VerticalInterpPdf::evaluate ( ) const

Definition at line 111 of file VerticalInterpPdf.cc.

References _coefIter, _funcIter, _quadraticAlgo, interpolate(), and relativeConstraints::value.

112 {
113  // Calculate the current value
114  Double_t value(0) ;
115 
116  // Do running sum of coef/func pairs, calculate lastCoef.
117  _funcIter->Reset() ;
118  _coefIter->Reset() ;
119  RooAbsReal* coef ;
120  RooAbsReal* func = (RooAbsReal*)_funcIter->Next();
121 
122  Double_t central = func->getVal();
123  value = central;
124 
125  if (_quadraticAlgo >= 0) {
126  // additive interpolation
127  while((coef=(RooAbsReal*)_coefIter->Next())) {
128  Double_t coefVal = coef->getVal() ;
129  RooAbsReal* funcUp = (RooAbsReal*)_funcIter->Next() ;
130  RooAbsReal* funcDn = (RooAbsReal*)_funcIter->Next() ;
131  value += interpolate(coefVal, central, funcUp, funcDn);
132  }
133  } else {
134  // multiplicative interpolation
135  while((coef=(RooAbsReal*)_coefIter->Next())) {
136  Double_t coefVal = coef->getVal() ;
137  RooAbsReal* funcUp = (RooAbsReal*)_funcIter->Next() ;
138  RooAbsReal* funcDn = (RooAbsReal*)_funcIter->Next() ;
139  value *= interpolate(coefVal, central, funcUp, funcDn);
140  }
141  }
142 
143  return value > 0 ? value : 1E-9 ;
144 }
Double_t interpolate(Double_t coeff, Double_t central, RooAbsReal *fUp, RooAbsReal *fDown) const
Iterator over coefficient list.
TIterator * _coefIter
Iterator over FUNC list.
virtual Bool_t VerticalInterpPdf::forceAnalyticalInt ( const RooAbsArg &  ) const
inlinevirtual

Definition at line 24 of file VerticalInterpPdf.h.

24 { return kTRUE ; }
const RooArgList& VerticalInterpPdf::funcList ( ) const
inline

Definition at line 28 of file VerticalInterpPdf.h.

References _funcList.

28 { return _funcList ; }
RooListProxy _funcList
Int_t VerticalInterpPdf::getAnalyticalIntegralWN ( RooArgSet &  allVars,
RooArgSet &  numVars,
const RooArgSet *  normSet,
const char *  rangeName = 0 
) const

Definition at line 188 of file VerticalInterpPdf.cc.

References VerticalInterpPdf::CacheElem::_funcIntList, _funcIter, VerticalInterpPdf::CacheElem::_funcNormList, and _normIntMgr.

190 {
191  // Advertise that all integrals can be handled internally.
192 
193  // Handle trivial no-integration scenario
194  if (allVars.getSize()==0) return 0 ;
195  if (_forceNumInt) return 0 ;
196 
197  // Select subset of allVars that are actual dependents
198  analVars.add(allVars) ;
199  RooArgSet* normSet = normSet2 ? getObservables(normSet2) : 0 ;
200 
201 
202  // Check if this configuration was created before
203  Int_t sterileIdx(-1) ;
204  CacheElem* cache = (CacheElem*) _normIntMgr.getObj(normSet,&analVars,&sterileIdx,(const char *)0);
205  if (cache) {
206  return _normIntMgr.lastIndex()+1 ;
207  }
208 
209  // Create new cache element
210  cache = new CacheElem ;
211 
212  // Make list of function projection and normalization integrals
213  _funcIter->Reset() ;
214  RooAbsReal *func ;
215  while((func=(RooAbsReal*)_funcIter->Next())) {
216  RooAbsReal* funcInt = func->createIntegral(analVars) ;
217  cache->_funcIntList.addOwned(*funcInt) ;
218  if (normSet && normSet->getSize()>0) {
219  RooAbsReal* funcNorm = func->createIntegral(*normSet) ;
220  cache->_funcNormList.addOwned(*funcNorm) ;
221  }
222  }
223 
224  // Store cache element
225  Int_t code = _normIntMgr.setObj(normSet,&analVars,(RooAbsCacheElement*)cache,0) ;
226 
227  if (normSet) {
228  delete normSet ;
229  }
230 
231  return code+1 ;
232 }
RooObjCacheManager _normIntMgr
Double_t VerticalInterpPdf::interpolate ( Double_t  coeff,
Double_t  central,
RooAbsReal *  fUp,
RooAbsReal *  fDown 
) const
protected

Iterator over coefficient list.

Definition at line 290 of file VerticalInterpPdf.cc.

References _quadraticAlgo, _quadraticRegion, and funct::pow().

Referenced by analyticalIntegralWN(), and evaluate().

291 {
292  if (_quadraticAlgo == -1) {
293  Double_t kappa = (coeff > 0 ? fUp->getVal()/central : central/fDn->getVal());
294  return pow(kappa,fabs(coeff));
295  }
296 
297  if (fabs(coeff) >= _quadraticRegion) {
298  return coeff * (coeff > 0 ? fUp->getVal() - central : central - fDn->getVal());
299  } else {
300  // quadratic interpolation coefficients between the three
301  if (_quadraticAlgo != 1) {
302  // quadratic interpolation null in zero and continuos at boundaries, but not differentiable at boundaries
303  // conditions:
304  // c_up (+_quadraticRegion) = +_quadraticRegion
305  // c_cen(+_quadraticRegion) = -_quadraticRegion
306  // c_dn (+_quadraticRegion) = 0
307  // c_up (-_quadraticRegion) = 0
308  // c_cen(-_quadraticRegion) = -_quadraticRegion
309  // c_dn (-_quadraticRegion) = +_quadraticRegion
310  // c_up(0) = c_dn(0) = c_cen(0) = 0
311  Double_t c_up = + coeff * (_quadraticRegion + coeff) / (2 * _quadraticRegion);
312  Double_t c_dn = - coeff * (_quadraticRegion - coeff) / (2 * _quadraticRegion);
313  Double_t c_cen = - coeff * coeff / _quadraticRegion;
314  return c_up * fUp->getVal() + c_dn * fDn->getVal() + c_cen * central;
315  } else { //if (_quadraticAlgo == 1) {
316  // quadratic interpolation that is everywhere differentiable, but it's not null in zero
317  // conditions on the function
318  // c_up (+_quadraticRegion) = +_quadraticRegion
319  // c_cen(+_quadraticRegion) = -_quadraticRegion
320  // c_dn (+_quadraticRegion) = 0
321  // c_up (-_quadraticRegion) = 0
322  // c_cen(-_quadraticRegion) = -_quadraticRegion
323  // c_dn (-_quadraticRegion) = +_quadraticRegion
324  // conditions on the derivatives
325  // c_up '(+_quadraticRegion) = +1
326  // c_cen'(+_quadraticRegion) = -1
327  // c_dn '(+_quadraticRegion) = 0
328  // c_up '(-_quadraticRegion) = 0
329  // c_cen'(-_quadraticRegion) = +1
330  // c_dn '(-_quadraticRegion) = -1
331  Double_t c_up = (_quadraticRegion + coeff) * (_quadraticRegion + coeff) / (4 * _quadraticRegion);
332  Double_t c_dn = (_quadraticRegion - coeff) * (_quadraticRegion - coeff) / (4 * _quadraticRegion);
333  Double_t c_cen = - c_up - c_dn;
334  return c_up * fUp->getVal() + c_dn * fDn->getVal() + c_cen * central;
335  }
336  }
337 
338 }
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40

Member Data Documentation

TIterator* VerticalInterpPdf::_coefIter
protected
RooListProxy VerticalInterpPdf::_coefList
protected

Definition at line 44 of file VerticalInterpPdf.h.

Referenced by coefList(), and VerticalInterpPdf().

TIterator* VerticalInterpPdf::_funcIter
protected
RooListProxy VerticalInterpPdf::_funcList
protected

Definition at line 43 of file VerticalInterpPdf.h.

Referenced by funcList(), and VerticalInterpPdf().

RooObjCacheManager VerticalInterpPdf::_normIntMgr
mutableprotected

Definition at line 41 of file VerticalInterpPdf.h.

Referenced by analyticalIntegralWN(), and getAnalyticalIntegralWN().

Int_t VerticalInterpPdf::_quadraticAlgo
protected

Definition at line 46 of file VerticalInterpPdf.h.

Referenced by checkObservables(), evaluate(), interpolate(), and VerticalInterpPdf().

Double_t VerticalInterpPdf::_quadraticRegion
protected

Definition at line 45 of file VerticalInterpPdf.h.

Referenced by interpolate().