59 _mirrorLeft(false), _mirrorRight(false),
60 _asymLeft(false), _asymRight(false)
70 Mirror mirror,
double rho) :
72 _x(
"x",
"observable",this,
x),
77 _mirrorLeft(mirror==MirrorLeft || mirror==MirrorBoth || mirror==MirrorLeftAsymRight),
78 _mirrorRight(mirror==MirrorRight || mirror==MirrorBoth || mirror==MirrorAsymLeftRight),
79 _asymLeft(mirror==MirrorAsymLeft || mirror==MirrorAsymLeftRight || mirror==MirrorAsymBoth),
80 _asymRight(mirror==MirrorAsymRight || mirror==MirrorLeftAsymRight || mirror==MirrorAsymBoth),
99 Mirror mirror,
double rho) :
101 _x(
"x",
"Observable",this,xpdf),
106 _mirrorLeft(mirror==MirrorLeft || mirror==MirrorBoth || mirror==MirrorLeftAsymRight),
107 _mirrorRight(mirror==MirrorRight || mirror==MirrorBoth || mirror==MirrorAsymLeftRight),
108 _asymLeft(mirror==MirrorAsymLeft || mirror==MirrorAsymLeftRight || mirror==MirrorAsymBoth),
109 _asymRight(mirror==MirrorAsymRight || mirror==MirrorLeftAsymRight || mirror==MirrorAsymBoth),
126 RooAbsPdf(other,
name), _x(
"x",this,other._x), _nEvents(other._nEvents),
127 _dataPts(0), _dataWgts(0), _weights(0), _sumWgt(0),
128 _mirrorLeft( other._mirrorLeft ), _mirrorRight( other._mirrorRight ),
129 _asymLeft(other._asymLeft), _asymRight(other._asymRight),
172 inline bool operator()(
const struct Data&
a,
const struct Data&
b)
const
173 {
return a.x <
b.x; }
181 std::vector<Data> tmp;
183 double x0 = 0.,
x1 = 0.,
x2 = 0.;
187 for (
Int_t i = 0; i <
data.numEntries(); ++i) {
189 const double x = real.
getVal();
190 const double w =
data.weight();
209 std::sort(tmp.begin(), tmp.end(), cmp());
215 for (
unsigned i = 0; i < tmp.size(); ++i) {
221 std::vector<Data> tmp2;
243 const double xlo = std::min(
_hi,
245 const double xhi = std::max(
_lo,
247 if (xlo >= xhi)
continue;
255 for (
Int_t k = binlo; k <= binhi; ++k, chi += chi2incr) {
261 const double xlo = std::min(
_hi,
263 const double xhi = std::max(
_lo,
265 if (xlo >= xhi)
continue;
273 for (
Int_t k = binlo; k <= binhi; ++k, chi += chi2incr) {
280 const double xlo = std::min(
_hi,
282 const double xhi = std::max(
_lo,
284 if (xlo >= xhi)
continue;
292 for (
Int_t k = binlo; k <= binhi; ++k, chi += chi2incr) {
348 for (
Int_t i = imin + 2; i < imax; ++i)
363 }
else if (imin == imax) {
383 double max = -std::numeric_limits<double>::max();
400 for ( ; it < iend; ++it) {
401 const double r = (
x - *it) / sigmav;
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t hmin
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t b
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t r
Option_t Option_t TPoint TPoint const char x2
Option_t Option_t TPoint TPoint const char x1
TRObject operator()(const T1 &t1) const
bool contains(const RooAbsArg &var) const
Check if collection contains an argument with the same name as var.
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
virtual double getMax(const char *name=0) const
Get maximum of currently defined range.
virtual double getMin(const char *name=0) const
Get minimum of currently defined range.
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
bool matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
RooAbsArg * absArg() const
Return pointer to contained argument.
RooArgSet is a container object that can hold multiple RooAbsArg objects.
RooDataSet is a container class to hold unbinned data.
Class RooKeysPdf implements a one-dimensional kernel estimation p.d.f which model the distribution of...
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const override
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported.
RooKeysPdf()
coverity[UNINIT_CTOR]
double _lookupTable[_nPoints+1]
double maxVal(Int_t code) const override
Return maximum value for set of observables identified by code assigned in getMaxVal.
double g(double x, double sigma) const
void LoadDataSet(RooDataSet &data)
double analyticalIntegral(Int_t code, const char *rangeName=0) const override
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
double evaluate() const override
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
static const double _nSigma
Int_t getMaxVal(const RooArgSet &vars) const override
Advertise capability to determine maximum value of function for given set of observables.
RooRealVar represents a variable that can be changed from the outside.
double min(const char *rname=0) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.
double max(const char *rname=0) const
Query upper limit of range. This requires the payload to be RooAbsRealLValue or derived.
const T & arg() const
Return reference to object held in proxy.
const char * GetName() const override
Returns name of object.
RVec< PromoteType< T > > floor(const RVec< T > &v)
RVec< PromoteTypes< T0, T1 > > pow(const T0 &x, const RVec< T1 > &v)
RVec< PromoteType< T > > log(const RVec< T > &v)
RVec< PromoteType< T > > exp(const RVec< T > &v)
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
static uint64_t sum(uint64_t i)