Logo ROOT   6.10/09
Reference Guide
PiecewiseInterpolation.h
Go to the documentation of this file.
1 /*****************************************************************************
2  * Project: RooFit *
3  * Package: RooFitCore *
4  * File: $Id: PiecewiseInterpolation.h,v 1.3 2007/05/11 09:11:30 verkerke Exp $
5  * Authors: *
6  * WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu *
7  * DK, David Kirkby, UC Irvine, dkirkby@uci.edu *
8  * *
9  * Copyright (c) 2000-2005, Regents of the University of California *
10  * and Stanford University. All rights reserved. *
11  * *
12  * Redistribution and use in source and binary forms, *
13  * with or without modification, are permitted according to the terms *
14  * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
15  *****************************************************************************/
16 #ifndef ROO_PIECEWISEINTERPOLATION
17 #define ROO_PIECEWISEINTERPOLATION
18 
19 #include "RooAbsReal.h"
20 #include "RooRealProxy.h"
21 #include "RooListProxy.h"
22 
23 #include "RooObjCacheManager.h"
24 
25 class RooRealVar;
26 class RooArgList ;
27 
29 public:
30 
32  PiecewiseInterpolation(const char *name, const char *title, const RooAbsReal& nominal, const RooArgList& lowSet, const RooArgList& highSet, const RooArgList& paramSet, Bool_t takeOwnerShip=kFALSE) ;
33  virtual ~PiecewiseInterpolation() ;
34 
35  PiecewiseInterpolation(const PiecewiseInterpolation& other, const char* name = 0);
36  virtual TObject* clone(const char* newname) const { return new PiecewiseInterpolation(*this, newname); }
37 
38  // virtual Double_t defaultErrorLevel() const ;
39 
40  // void printMetaArgs(std::ostream& os) const ;
41 
42  const RooArgList& lowList() const { return _lowSet ; }
43  const RooArgList& highList() const { return _highSet ; }
44  const RooArgList& paramList() const { return _paramSet ; }
45 
46  //virtual Bool_t forceAnalyticalInt(const RooAbsArg&) const { return kTRUE ; }
48 
49  Int_t getAnalyticalIntegralWN(RooArgSet& allVars, RooArgSet& analVars, const RooArgSet* normSet,const char* rangeName=0) const ;
50  Double_t analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName=0) const ;
51 
52  void setPositiveDefinite(bool flag=true){_positiveDefinite=flag;}
53 
54  void setInterpCode(RooAbsReal& param, int code);
55  void setAllInterpCodes(int code);
56  void printAllInterpCodes();
57 
58  virtual std::list<Double_t>* binBoundaries(RooAbsRealLValue& /*obs*/, Double_t /*xlo*/, Double_t /*xhi*/) const ;
59  virtual std::list<Double_t>* plotSamplingHint(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const ;
60  virtual Bool_t isBinnedDistribution(const RooArgSet& obs) const ;
61 
62 protected:
63 
64  class CacheElem : public RooAbsCacheElement {
65  public:
66  CacheElem() {} ;
67  virtual ~CacheElem() {} ;
69  RooArgList ret(_funcIntList) ;
70  ret.add(_lowIntList);
71  ret.add(_highIntList);
72  return ret ;
73  }
77  // will want std::vector<RooRealVar*> for low and high also
78  } ;
79  mutable RooObjCacheManager _normIntMgr ; // The integration cache manager
80 
81  RooRealProxy _nominal; // The nominal value
82  RooArgList _ownedList ; // List of owned components
83  RooListProxy _lowSet ; // Low-side variation
84  RooListProxy _highSet ; // High-side varaition
85  RooListProxy _paramSet ; // interpolation parameters
86  RooListProxy _normSet ; // interpolation parameters
87  Bool_t _positiveDefinite; // protect against negative and 0 bins.
88 
89  std::vector<int> _interpCode;
90 
91  Double_t evaluate() const;
92 
93  ClassDef(PiecewiseInterpolation,3) // Sum of RooAbsReal objects
94 };
95 
96 #endif
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &analVars, const RooArgSet *normSet, const char *rangeName=0) const
Advertise that all integrals can be handled internally.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
Double_t evaluate() const
Calculate and return current value of self.
const RooArgList & highList() const
virtual std::list< Double_t > * plotSamplingHint(RooAbsRealLValue &obs, Double_t xlo, Double_t xhi) const
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
virtual std::list< Double_t > * binBoundaries(RooAbsRealLValue &, Double_t, Double_t) const
WVE note: assumes nominal and alternates have identical structure, must add explicit check...
#define ClassDef(name, id)
Definition: Rtypes.h:297
std::vector< int > _interpCode
virtual ~PiecewiseInterpolation()
Destructor.
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
void setInterpCode(RooAbsReal &param, int code)
virtual RooArgList containedArgs(Action)
Double_t analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
Implement analytical integrations by doing appropriate weighting from component integrals functions t...
RooAbsCacheElement is the abstract base class for objects to be stored in RooAbsCache cache manager o...
RooObjCacheManager _normIntMgr
RooListProxy is the concrete proxy for RooArgList objects.
Definition: RooListProxy.h:25
const Bool_t kFALSE
Definition: RtypesCore.h:92
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
virtual Bool_t isBinnedDistribution(const RooArgSet &obs) const
WVE note: assumes nominal and alternates have identical structure, must add explicit check...
Class RooObjCacheManager is an implementation of class RooCacheManager<RooAbsCacheElement> and specia...
const RooArgList & paramList() const
Mother of all ROOT objects.
Definition: TObject.h:37
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
void setPositiveDefinite(bool flag=true)
RooRealProxy is the concrete proxy for RooAbsReal objects A RooRealProxy is the general mechanism to ...
Definition: RooRealProxy.h:23
Bool_t setBinIntegrator(RooArgSet &allVars)
virtual TObject * clone(const char *newname) const
const RooArgList & lowList() const