Logo ROOT  
Reference Guide
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
Loading...
Searching...
No Matches
Roo2DKeysPdf.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitModels *
4 * @(#)root/roofit:$Id$
5 * Authors: *
6 * AB, Adrian Bevan, Liverpool University, bevan@slac.stanford.edu *
7 * *
8 * Copyright (c) 2000-2005, Regents of the University of California, *
9 * Liverpool University, *
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
17/** \class Roo2DKeysPdf
18 \ingroup Roofit
19
20Two-dimensional kernel estimation PDF.
21
22<b>This function has been superseded by the more general RooNDKeysPdf.</b>
23*/
24
25#include "Roo2DKeysPdf.h"
26#include "RooRealVar.h"
27#include "TTree.h"
28#include "TH2.h"
29#include "TFile.h"
30#include "TMath.h"
31
32using std::ostream;
33
35
36
37////////////////////////////////////////////////////////////////////////////////
38/// Constructor.
39/// \param[in] name
40/// \param[in] title
41/// \param[in] xx
42/// \param[in] yy
43/// \param[in] data
44/// \param[in] options
45/// \param[in] widthScaleFactor
46
47Roo2DKeysPdf::Roo2DKeysPdf(const char *name, const char *title,
49 RooAbsPdf(name,title),
50 x("x", "x dimension",this, xx),
51 y("y", "y dimension",this, yy)
52{
54 loadDataSet(data, options);
55}
56
57
58////////////////////////////////////////////////////////////////////////////////
59/// Copy constructor.
60/// \param[in] other
61/// \param[in] name
62
65 x("x", this, other.x),
66 y("y", this, other.y)
67{
68 if(_verbosedebug) { std::cout << "Roo2DKeysPdf::Roo2DKeysPdf copy ctor" << std::endl; }
69
70 _xMean = other._xMean;
71 _xSigma = other._xSigma;
72 _yMean = other._yMean;
73 _ySigma = other._ySigma;
74 _n = other._n;
75
76 _BandWidthType = other._BandWidthType;
77 _MirrorAtBoundary = other._MirrorAtBoundary;
78 _widthScaleFactor = other._widthScaleFactor;
79
80 _2pi = other._2pi;
81 _sqrt2pi = other._sqrt2pi;
82 _nEvents = other._nEvents;
83 _n16 = other._n16;
84 _debug = other._debug;
85 _verbosedebug = other._verbosedebug;
86 _vverbosedebug = other._vverbosedebug;
87
88 _lox = other._lox;
89 _hix = other._hix;
90 _loy = other._loy;
91 _hiy = other._hiy;
92 _xoffset = other._xoffset;
93 _yoffset = other._yoffset;
94
95 _x = new double[_nEvents];
96 _y = new double[_nEvents];
97 _hx = new double[_nEvents];
98 _hy = new double[_nEvents];
99
100 //copy the data and bandwidths
101 for(Int_t iEvt = 0; iEvt< _nEvents; iEvt++)
102 {
103 _x[iEvt] = other._x[iEvt];
104 _y[iEvt] = other._y[iEvt];
105 _hx[iEvt] = other._hx[iEvt];
106 _hy[iEvt] = other._hy[iEvt];
107 }
108}
109
110
111////////////////////////////////////////////////////////////////////////////////
112/// Destructor.
113
115 if(_verbosedebug) { std::cout << "Roo2DKeysPdf::Roo2KeysPdf dtor" << std::endl; }
116 delete[] _x;
117 delete[] _hx;
118 delete[] _y;
119 delete[] _hy;
120}
121
122
123////////////////////////////////////////////////////////////////////////////////
124/// Loads a new data set into the class instance.
125/// Returns 1 in case of error, 0 otherwise.
126/// \param[in] data
127/// \param[in] options
128
130{
131 if(_verbosedebug) { std::cout << "Roo2DKeysPdf::loadDataSet" << std::endl; }
132
133 setOptions(options);
134
135 if(_verbosedebug) { std::cout << "Roo2DKeysPdf::loadDataSet(RooDataSet& data, TString options)" << std::endl; }
136
137 _2pi = 2.0*TMath::Pi() ; //use pi from math.h
138 _sqrt2pi = sqrt(_2pi);
139 _nEvents = (Int_t)data.numEntries();
140 if(_nEvents == 0)
141 {
142 std::cout << "ERROR: Roo2DKeysPdf::loadDataSet The input data set is empty. Unable to begin generating the PDF" << std::endl;
143 return 1;
144 }
145 _n16 = std::pow(_nEvents, -1./6.); // = (4/[n(dim(R) + 2)])^1/(dim(R)+4); dim(R) = 2
146
147 _lox = x.min();
148 _hix = x.max();
149 _loy = y.min();
150 _hiy = y.max();
151
152 _x = new double[_nEvents];
153 _y = new double[_nEvents];
154 _hx = new double[_nEvents];
155 _hy = new double[_nEvents];
156
157 double x0 = 0.0;
158 double x1 = 0.0;
159 double x_2 = 0.0;
160 double y0 = 0.0;
161 double y1 = 0.0;
162 double y_2 = 0.0;
163
164 //check that the data contain the variable we are interested in
165 Int_t bad = 0;
166 const RooAbsReal & xx = x.arg();
167 const RooAbsReal & yy = y.arg();
168 if(! static_cast<RooRealVar*>((const_cast<RooArgSet *>(data.get(0)))->find( xx.GetName() )) )
169 {
170 std::cout << "Roo2DKeysPdf::Roo2DKeysPdf invalid RooAbsReal name: "<<xx.GetName()<<" not in the data set" << std::endl;
171 bad = 1;
172 }
173 if(! static_cast<RooRealVar*>((const_cast<RooArgSet *>(data.get(0)))->find( yy.GetName() )) )
174 {
175 std::cout << "Roo2DKeysPdf::Roo2DKeysPdf invalid RooAbsReal name: "<<yy.GetName()<<" not in the data set" << std::endl;
176 bad = 1;
177 }
178 if(bad)
179 {
180 std::cout << "Roo2DKeysPdf::Roo2DKeysPdf Unable to initialize object; incompatible RooDataSet doesn't contain"<< std::endl;
181 std::cout << " all of the RooAbsReal arguments"<< std::endl;
182 return 1;
183 }
184
185 //copy the data into local arrays
186 const RooArgSet * values = data.get();
187 auto X = static_cast<RooRealVar const*>(values->find(xx.GetName()));
188 auto Y = static_cast<RooRealVar const*>(values->find(yy.GetName()));
189
190 for (Int_t j=0;j<_nEvents;++j)
191 {
192 data.get(j) ;
193
194 _x[j] = X->getVal() ;
195 _y[j] = Y->getVal() ;
196
197 x0+=1; x1+=_x[j]; x_2+=_x[j]*_x[j];
198 y0+=1; y1+=_y[j]; y_2+=_y[j]*_y[j];
199 }
200
201 //==========================================//
202 //calculate the mean and sigma for the data //
203 //==========================================//
204 if(_nEvents == 0)
205 {
206 std::cout << "Roo2DKeysPdf::Roo2DKeysPdf Empty data set was used; can't generate a PDF"<< std::endl;
207 }
208
209 _xMean = x1/x0;
211
212 _yMean = y1/y0;
214
216
217 //calculate the PDF
219}
220
221
222////////////////////////////////////////////////////////////////////////////////
223
225{
226 if(_verbosedebug) { std::cout << "Roo2DKeysPdf::setOptions" << std::endl; }
227
228 options.ToLower();
229 if( options.Contains("a") ) _BandWidthType = 0;
230 else _BandWidthType = 1;
231 if( options.Contains("n") ) _BandWidthType = 1;
232 else _BandWidthType = 0;
233 if( options.Contains("m") ) _MirrorAtBoundary = 1;
234 else _MirrorAtBoundary = 0;
235 if( options.Contains("d") ) _debug = 1;
236 else _debug = 0;
237 if( options.Contains("v") ) { _debug = 1; _verbosedebug = 1; }
238 else _verbosedebug = 0;
239 if( options.Contains("vv") ) { _vverbosedebug = 1; }
240 else _vverbosedebug = 0;
241
242 if( _debug )
243 {
244 std::cout << "Roo2DKeysPdf::setOptions(TString options) options = "<< options << std::endl;
245 std::cout << "\t_BandWidthType = " << _BandWidthType << std::endl;
246 std::cout << "\t_MirrorAtBoundary = " << _MirrorAtBoundary << std::endl;
247 std::cout << "\t_debug = " << _debug << std::endl;
248 std::cout << "\t_verbosedebug = " << _verbosedebug << std::endl;
249 std::cout << "\t_vverbosedebug = " << _vverbosedebug << std::endl;
250 }
251}
252
253
254////////////////////////////////////////////////////////////////////////////////
255
257{
258 std::cout << "Roo2DKeysPdf::getOptions(void)" << std::endl;
259 std::cout << "\t_BandWidthType = " << _BandWidthType << std::endl;
260 std::cout << "\t_MirrorAtBoundary = " << _MirrorAtBoundary << std::endl;
261 std::cout << "\t_debug = " << _debug << std::endl;
262 std::cout << "\t_verbosedebug = " << _verbosedebug << std::endl;
263 std::cout << "\t_vverbosedebug = " << _vverbosedebug << std::endl;
264}
265
266
267////////////////////////////////////////////////////////////////////////////////
268/// Calculates the kernel bandwidth for x & y and the probability look up table _p[i][j]
269/// \param[in] kernel
270
272{
273 if(_verbosedebug) { std::cout << "Roo2DKeysPdf::calculateBandWidth(Int_t kernel)" << std::endl; }
274 if(kernel != -999)
275 {
277 }
278
279 double h = 0.0;
280
282 double sqrtSum = sqrt( sigSum );
283 double sigProd = _ySigma*_xSigma;
284 if(sigProd != 0.0) h = _n16*sqrt( sigSum/sigProd );
285 if(sqrtSum == 0)
286 {
287 std::cout << "Roo2DKeysPdf::calculateBandWidth The sqr(variance sum) == 0.0. " << " Your dataset represents a delta function."<< std::endl;
288 return 1;
289 }
290
291 double hXSigma = h * _xSigma;
292 double hYSigma = h * _ySigma;
293 double xhmin = hXSigma * sqrt(2.)/10; //smallest anticipated bandwidth
294 double yhmin = hYSigma * sqrt(2.)/10;
295
296 //////////////////////////////////////
297 //calculate bandwidths from the data//
298 //////////////////////////////////////
299 if(_BandWidthType == 1) //calculate a trivial bandwidth
300 {
301 std::cout << "Roo2DKeysPdf::calculateBandWidth Using a normal bandwidth (same for a given dimension) based on"<< std::endl;
302 std::cout << " h_j = n^{-1/6}*sigma_j for the j^th dimension and n events * "<<_widthScaleFactor<< std::endl;
305 for(Int_t j=0;j<_nEvents;++j)
306 {
307 _hx[j] = hxGaussian;
308 _hy[j] = hyGaussian;
309 if(_hx[j]<xhmin) _hx[j] = xhmin;
310 if(_hy[j]<yhmin) _hy[j] = yhmin;
311 }
312 }
313 else //use an adaptive bandwidth to reduce the dependence on global data distribution
314 {
315 std::cout << "Roo2DKeysPdf::calculateBandWidth Using an adaptive bandwidth (in general different for all events) [default]"<< std::endl;
316 std::cout << " scaled by a factor of "<<_widthScaleFactor<< std::endl;
317 double xnorm = h * std::pow(_xSigma/sqrtSum, 1.5) * _widthScaleFactor;
318 double ynorm = h * std::pow(_ySigma/sqrtSum, 1.5) * _widthScaleFactor;
319 for(Int_t j=0;j<_nEvents;++j)
320 {
321 double f_ti = std::pow( g(_x[j], _x, hXSigma, _y[j], _y, hYSigma), -0.25 ) ;
322 _hx[j] = xnorm * f_ti;
323 _hy[j] = ynorm * f_ti;
324 if(_hx[j]<xhmin) _hx[j] = xhmin;
325 if(_hy[j]<yhmin) _hy[j] = yhmin;
326 }
327 }
328
329 return 0;
330}
331
332
333////////////////////////////////////////////////////////////////////////////////
334/// Evaluates the kernel estimation for x,y, interpolating between the points if necessary
335///
336/// Uses the caching intrinsic in RFC to bypass the grid and remove
337/// the grid and extrapolation approximation in the kernel estimation method
338/// implementation.
339
341{
342 if(_vverbosedebug) { std::cout << "Roo2DKeysPdf::evaluate()" << std::endl; }
343 return evaluateFull(x,y);
344}
345
346
347////////////////////////////////////////////////////////////////////////////////
348/// Evaluates the sum of the product of the 2D kernels
349/// for use in calculating the fixed kernel estimate, f,
350/// given the bandwidths _hx[j] and _hy[j].
351///
352/// _n is calculated once in the constructor.
353/// \param[in] thisX
354/// \param[in] thisY
355
356double Roo2DKeysPdf::evaluateFull(double thisX, double thisY) const
357{
358 if( _vverbosedebug ) { std::cout << "Roo2DKeysPdf::evaluateFull()" << std::endl; }
359
360 double f=0.0;
361
362 double rx2;
363 double ry2;
364 double zx;
365 double zy;
367 {
368 for (Int_t j = 0; j < _nEvents; ++j)
369 {
370 rx2 = 0.0; ry2 = 0.0; zx = 0.0; zy = 0.0;
371 if(_hx[j] != 0.0) rx2 = (thisX - _x[j])/_hx[j];
372 if(_hy[j] != 0.0) ry2 = (thisY - _y[j])/_hy[j];
373
374 if(_hx[j] != 0.0) zx = exp(-0.5*rx2*rx2)/_hx[j];
375 if(_hy[j] != 0.0) zy = exp(-0.5*ry2*ry2)/_hy[j];
376
381 f += zy * zx;
382 // f += _n * zy * zx; // ooops this is a normalisation factor :(
383 }
384 }
385 else
386 {
387 for (Int_t j = 0; j < _nEvents; ++j)
388 {
389 rx2 = 0.0; ry2 = 0.0; zx = 0.0; zy = 0.0;
390 if(_hx[j] != 0.0) rx2 = (thisX - _x[j])/_hx[j];
391 if(_hy[j] != 0.0) ry2 = (thisY - _y[j])/_hy[j];
392
393 if(_hx[j] != 0.0) zx = exp(-0.5*rx2*rx2)/_hx[j];
394 if(_hy[j] != 0.0) zy = exp(-0.5*ry2*ry2)/_hy[j];
395 f += zy * zx;
396 // f += _n * zy * zx; // ooops this is a normalisation factor :(
397 }
398 }
399 return f;
400}
401
402
403////////////////////////////////////////////////////////////////////////////////
404/// Apply the mirror at boundary correction to a dimension given the space position to evaluate
405/// at (thisVar), the bandwidth at this position (thisH), the boundary (high/low) and the
406/// value of the data kernel that this correction is being applied to tVar (i.e. the _x[ix] etc.).
407/// \param[in] thisVar
408/// \param[in] thisH
409/// \param[in] high
410/// \param[in] tVar
411
412double Roo2DKeysPdf::highBoundaryCorrection(double thisVar, double thisH, double high, double tVar) const
413{
414 if(_vverbosedebug) { std::cout << "Roo2DKeysPdf::highBoundaryCorrection" << std::endl; }
415
416 if(thisH == 0.0) return 0.0;
417 double correction = (thisVar + tVar - 2.0* high )/thisH;
418 return exp(-0.5*correction*correction)/thisH;
419}
420
421
422////////////////////////////////////////////////////////////////////////////////
423
424double Roo2DKeysPdf::lowBoundaryCorrection(double thisVar, double thisH, double low, double tVar) const
425{
426 if(_vverbosedebug) { std::cout << "Roo2DKeysPdf::lowBoundaryCorrection" << std::endl; }
427
428 if(thisH == 0.0) return 0.0;
429 double correction = (thisVar + tVar - 2.0* low )/thisH;
430 return exp(-0.5*correction*correction)/thisH;
431}
432
433
434////////////////////////////////////////////////////////////////////////////////
435/// Calculates f(t_i) for the bandwidths.
436/// \f$ g = 1/(N_{evt} * \sigma_j * \sqrt{2\pi})*\sum_{all evts}{\prod d K[ \exp({-(xd - ti)/\sigma_{j}d^2}) ]}\f$
437/// \param[in] varMean1
438/// \param[in] _var1
439/// \param[in] sigma1
440/// \param[in] varMean2
441/// \param[in] _var2
442/// \param[in] sigma2
443
444double Roo2DKeysPdf::g(double varMean1, double * _var1, double sigma1, double varMean2, double * _var2, double sigma2) const
445{
446 if((_nEvents == 0.0) || (sigma1 == 0.0) || (sigma2 == 0)) return 0.0;
447
448 double c1 = -1.0/(2.0*sigma1*sigma1);
449 double c2 = -1.0/(2.0*sigma2*sigma2);
450 double d = 4.0*c1*c2 /(_sqrt2pi*_nEvents);
451 double z = 0.0;
452
453 for (Int_t i = 0; i < _nEvents; ++i)
454 {
455 double r1 = _var1[i] - varMean1;
456 double r2 = _var2[i] - varMean2;
457 z += exp( c1 * r1*r1 ) * exp( c2 * r2*r2 );
458 }
459 z = z*d;
460 return z;
461}
462
463
464////////////////////////////////////////////////////////////////////////////////
465
467{
468 if(_BandWidthType == 1) std::cout << "The Bandwidth Type selected is Trivial" << std::endl;
469 else std::cout << "The Bandwidth Type selected is Adaptive" << std::endl;
470
471 return _BandWidthType;
472}
473
474
475////////////////////////////////////////////////////////////////////////////////
476
477double Roo2DKeysPdf::getMean(const char * axis) const
478{
479 if(!strcmp(axis,x.GetName()) || !strcmp(axis,"x") || !strcmp(axis,"X")) return _xMean;
480 else if(!strcmp(axis,y.GetName()) || !strcmp(axis,"y") || !strcmp(axis,"Y")) return _yMean;
481 else
482 {
483 std::cout << "Roo2DKeysPdf::getMean unknown axis "<<axis<< std::endl;
484 }
485 return 0.0;
486}
487
488
489////////////////////////////////////////////////////////////////////////////////
490
491double Roo2DKeysPdf::getSigma(const char * axis) const
492{
493 if(!strcmp(axis,x.GetName()) || !strcmp(axis,"x") || !strcmp(axis,"X")) return _xSigma;
494 else if(!strcmp(axis,y.GetName()) || !strcmp(axis,"y") || !strcmp(axis,"Y")) return _ySigma;
495 else
496 {
497 std::cout << "Roo2DKeysPdf::getSigma unknown axis "<<axis<< std::endl;
498 }
499 return 0.0;
500}
501
502
503
504////////////////////////////////////////////////////////////////////////////////
505
506void Roo2DKeysPdf::writeToFile(char * outputFile, const char * name) const
507{
508 TString histName = name;
509 histName += "_hist";
511 nName += "_Ntuple";
512 writeHistToFile( outputFile, histName);
514}
515
516
517////////////////////////////////////////////////////////////////////////////////
518/// Plots the PDF as a histogram and saves it to a file, so that it can be loaded in
519/// as a Roo2DHist PDF in the future to save on calculation time.
520/// \param[in] outputFile Name of the file where to store the PDF
521/// \param[in] histName PDF histogram name
522
523void Roo2DKeysPdf::writeHistToFile(char * outputFile, const char * histName) const
524{
525 std::cout << "Roo2DKeysPdf::writeHistToFile This member function is temporarily disabled" << std::endl;
526 //make sure that any existing file is not over written
527 std::unique_ptr<TFile> file{TFile::Open(outputFile, "UPDATE")};
528 if (!file)
529 {
530 std::cout << "Roo2DKeysPdf::writeHistToFile unable to open file "<< outputFile << std::endl;
531 return;
532 }
533
534
535 const RooAbsReal & xx = x.arg();
536 const RooAbsReal & yy = y.arg();
537 RooArgSet values( RooArgList( xx, yy ));
538 RooRealVar * xArg = (static_cast<RooRealVar*>(values.find(xx.GetName())) ) ;
539 RooRealVar * yArg = (static_cast<RooRealVar*>(values.find(yy.GetName())) ) ;
540
541 TH2F * hist = (TH2F*)xArg->createHistogram("hist", *yArg);
542 hist = static_cast<TH2F*>(this->fillHistogram(hist, RooArgList(*xArg, *yArg) ));
543 hist->SetName(histName);
544
545 file->Write();
546 file->Close();
547}
548
549
550////////////////////////////////////////////////////////////////////////////////
551/// Saves the data and calculated bandwidths to a file,
552/// as a record of what produced the PDF and to give a reduced
553/// data set in order to facilitate re-calculation in the future.
554/// \param[in] outputFile Name of the file where to store the data
555/// \param[in] name Name of the tree which will contain the data
556
557void Roo2DKeysPdf::writeNTupleToFile(char * outputFile, const char * name) const
558{
559 TFile * file = nullptr;
560
561 //make sure that any existing file is not over written
562 file = new TFile(outputFile, "UPDATE");
563 if (!file)
564 {
565 std::cout << "Roo2DKeysPdf::writeNTupleToFile unable to open file "<< outputFile << std::endl;
566 return;
567 }
570
571 double theX;
572 double theY;
573 double hx;
574 TString label = name;
575 label += " the source data for 2D Keys PDF";
576 TTree * _theTree = new TTree(name, label);
577 if(!_theTree) { std::cout << "Unable to get a TTree for output" << std::endl; return; }
578 _theTree->SetAutoSave(1000000000); // autosave when 1 Gbyte written
579
580 //name the TBranches the same as the RooAbsReal's
581 const char * xname = xArg.GetName();
582 const char * yname = yArg.GetName();
583 if (!strcmp(xname,"")) xname = "x";
584 if (!strcmp(yname,"")) yname = "y";
585
586 _theTree->Branch(xname, &theX, " x/D");
587 _theTree->Branch(yname, &theY, " y/D");
588 _theTree->Branch("hx", &hx, " hx/D");
589 _theTree->Branch("hy", &hx, " hy/D");
590
591 for(Int_t iEvt = 0; iEvt < _nEvents; iEvt++)
592 {
593 theX = _x[iEvt];
594 theY = _y[iEvt];
595 hx = _hx[iEvt];
596 hx = _hy[iEvt];
597 _theTree->Fill();
598 }
599 file->Write();
600 file->Close();
601}
602
603
604////////////////////////////////////////////////////////////////////////////////
605/// Prints out _p[_nPoints][_nPoints] indicating the domain limits.
606/// \param[out] out Output stream where to print
607
608void Roo2DKeysPdf::PrintInfo(ostream & out) const
609{
610 out << "Roo2DKeysPDF instance domain information:"<< std::endl;
611 out << "\tX_min = " << _lox << std::endl;
612 out << "\tX_max = " << _hix << std::endl;
613 out << "\tY_min = " << _loy << std::endl;
614 out << "\tY_max = " << _hiy << std::endl;
615
616 out << "Data information:" << std::endl;
617 out << "\t<x> = " << _xMean << std::endl;
618 out << "\tsigma(x) = " << _xSigma << std::endl;
619 out << "\t<y> = " << _yMean << std::endl;
620 out << "\tsigma(y) = " << _ySigma << std::endl;
621
622 out << "END of info for Roo2DKeys pdf instance"<< std::endl;
623}
#define d(i)
Definition RSha256.hxx:102
#define f(i)
Definition RSha256.hxx:104
#define g(i)
Definition RSha256.hxx:105
#define h(i)
Definition RSha256.hxx:106
int Int_t
Definition RtypesCore.h:45
#define ClassImp(name)
Definition Rtypes.h:374
#define X(type, name)
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
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 x1
Option_t Option_t TPoint TPoint const char y1
char name[80]
Definition TGX11.cxx:110
Two-dimensional kernel estimation PDF.
RooRealProxy y
double getMean(const char *axis) const
Roo2DKeysPdf(const char *name, const char *title, RooAbsReal &xx, RooAbsReal &yy, RooDataSet &data, TString options="a", double widthScaleFactor=1.0)
Constructor.
double _widthScaleFactor
void writeNTupleToFile(char *outputFile, const char *name) const
Saves the data and calculated bandwidths to a file, as a record of what produced the PDF and to give ...
Int_t _BandWidthType
void writeToFile(char *outputFile, const char *name) const
Int_t _verbosedebug
double * _hy
double lowBoundaryCorrection(double thisVar, double thisH, double low, double tVar) const
double evaluate() const override
Evaluates the kernel estimation for x,y, interpolating between the points if necessary.
Int_t getBandWidthType() const
Int_t _MirrorAtBoundary
void PrintInfo(std::ostream &) const
Prints out _p[_nPoints][_nPoints] indicating the domain limits.
Int_t loadDataSet(RooDataSet &data, TString options)
Loads a new data set into the class instance.
double g(double var1, double *_var1, double sigma1, double var2, double *_var2, double sigma2) const
Calculates f(t_i) for the bandwidths.
Int_t _vverbosedebug
double * _x
double * _y
double * _hx
double getSigma(const char *axis) const
~Roo2DKeysPdf() override
Destructor.
void setWidthScaleFactor(double widthScaleFactor)
double highBoundaryCorrection(double thisVar, double thisH, double high, double tVar) const
Apply the mirror at boundary correction to a dimension given the space position to evaluate at (thisV...
void writeHistToFile(char *outputFile, const char *histName) const
Plots the PDF as a histogram and saves it to a file, so that it can be loaded in as a Roo2DHist PDF i...
Int_t calculateBandWidth(Int_t kernel=-999)
Calculates the kernel bandwidth for x & y and the probability look up table _p[i][j].
void setOptions(TString options)
void getOptions(void) const
double evaluateFull(double thisX, double thisY) const
Evaluates the sum of the product of the 2D kernels for use in calculating the fixed kernel estimate,...
RooRealProxy x
RooAbsArg * find(const char *name) const
Find object with given name in list.
Abstract interface for all probability density functions.
Definition RooAbsPdf.h:42
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:59
TH1 * fillHistogram(TH1 *hist, const RooArgList &plotVars, double scaleFactor=1, const RooArgSet *projectedVars=nullptr, bool scaling=true, const RooArgSet *condObs=nullptr, bool setError=true) const
Fill the ROOT histogram 'hist' with values sampled from this function at the bin centers.
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition RooArgList.h:22
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:24
Container class to hold unbinned data.
Definition RooDataSet.h:34
Variable that can be changed from the outside.
Definition RooRealVar.h:37
double max(const char *rname=nullptr) 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.
double min(const char *rname=nullptr) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.
A ROOT file is an on-disk file, usually with extension .root, that stores objects in a file-system-li...
Definition TFile.h:53
static TFile * Open(const char *name, Option_t *option="", const char *ftitle="", Int_t compress=ROOT::RCompressionSetting::EDefaults::kUseCompiledDefault, Int_t netopt=0)
Create / open a file.
Definition TFile.cxx:4086
Int_t Write(const char *name=nullptr, Int_t opt=0, Int_t bufsiz=0) override
Write memory objects to this file.
Definition TFile.cxx:2433
void Close(Option_t *option="") override
Close a file.
Definition TFile.cxx:949
void SetName(const char *name) override
Change the name of this histogram.
Definition TH1.cxx:8962
2-D histogram with a float per channel (see TH1 documentation)
Definition TH2.h:303
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
Basic string class.
Definition TString.h:139
void ToLower()
Change string to lower-case.
Definition TString.cxx:1182
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:632
A TTree represents a columnar dataset.
Definition TTree.h:79
Double_t y[n]
Definition legend1.C:17
return c1
Definition legend1.C:41
Double_t x[n]
Definition legend1.C:17
return c2
Definition legend2.C:14
constexpr Double_t Pi()
Definition TMath.h:37