Logo ROOT  
Reference Guide
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 "TBranch.h"
31#include "TMath.h"
32
33//#include <math.h>
34
35using namespace std;
36
38
39
40////////////////////////////////////////////////////////////////////////////////
41/// Constructor.
42/// \param[in] name
43/// \param[in] title
44/// \param[in] xx
45/// \param[in] yy
46/// \param[in] data
47/// \param[in] options
48/// \param[in] widthScaleFactor
49
50Roo2DKeysPdf::Roo2DKeysPdf(const char *name, const char *title,
51 RooAbsReal& xx, RooAbsReal & yy, RooDataSet& data, TString options, double widthScaleFactor):
52 RooAbsPdf(name,title),
53 x("x", "x dimension",this, xx),
54 y("y", "y dimension",this, yy)
55{
56 setWidthScaleFactor(widthScaleFactor);
57 loadDataSet(data, options);
58}
59
60
61////////////////////////////////////////////////////////////////////////////////
62/// Copy constructor.
63/// \param[in] other
64/// \param[in] name
65
66Roo2DKeysPdf::Roo2DKeysPdf(const Roo2DKeysPdf & other, const char* name) :
67 RooAbsPdf(other,name),
68 x("x", this, other.x),
69 y("y", this, other.y)
70{
71 if(_verbosedebug) { cout << "Roo2DKeysPdf::Roo2DKeysPdf copy ctor" << endl; }
72
73 _xMean = other._xMean;
74 _xSigma = other._xSigma;
75 _yMean = other._yMean;
76 _ySigma = other._ySigma;
77 _n = other._n;
78
82
83 _2pi = other._2pi;
84 _sqrt2pi = other._sqrt2pi;
85 _nEvents = other._nEvents;
86 _n16 = other._n16;
87 _debug = other._debug;
90
91 _lox = other._lox;
92 _hix = other._hix;
93 _loy = other._loy;
94 _hiy = other._hiy;
95 _xoffset = other._xoffset;
96 _yoffset = other._yoffset;
97
98 _x = new double[_nEvents];
99 _y = new double[_nEvents];
100 _hx = new double[_nEvents];
101 _hy = new double[_nEvents];
102
103 //copy the data and bandwidths
104 for(Int_t iEvt = 0; iEvt< _nEvents; iEvt++)
105 {
106 _x[iEvt] = other._x[iEvt];
107 _y[iEvt] = other._y[iEvt];
108 _hx[iEvt] = other._hx[iEvt];
109 _hy[iEvt] = other._hy[iEvt];
110 }
111}
112
113
114////////////////////////////////////////////////////////////////////////////////
115/// Destructor.
116
118 if(_verbosedebug) { cout << "Roo2DKeysPdf::Roo2KeysPdf dtor" << endl; }
119 delete[] _x;
120 delete[] _hx;
121 delete[] _y;
122 delete[] _hy;
123}
124
125
126////////////////////////////////////////////////////////////////////////////////
127/// Loads a new data set into the class instance.
128/// Returns 1 in case of error, 0 otherwise.
129/// \param[in] data
130/// \param[in] options
131
133{
134 if(_verbosedebug) { cout << "Roo2DKeysPdf::loadDataSet" << endl; }
135
136 setOptions(options);
137
138 if(_verbosedebug) { cout << "Roo2DKeysPdf::loadDataSet(RooDataSet& data, TString options)" << endl; }
139
140 _2pi = 2.0*TMath::Pi() ; //use pi from math.h
141 _sqrt2pi = sqrt(_2pi);
142 _nEvents = (Int_t)data.numEntries();
143 if(_nEvents == 0)
144 {
145 cout << "ERROR: Roo2DKeysPdf::loadDataSet The input data set is empty. Unable to begin generating the PDF" << endl;
146 return 1;
147 }
148 _n16 = TMath::Power(_nEvents, -0.166666666); // = (4/[n(dim(R) + 2)])^1/(dim(R)+4); dim(R) = 2
149
150 _lox = x.min();
151 _hix = x.max();
152 _loy = y.min();
153 _hiy = y.max();
154
155 _x = new double[_nEvents];
156 _y = new double[_nEvents];
157 _hx = new double[_nEvents];
158 _hy = new double[_nEvents];
159
160 double x0 = 0.0;
161 double x1 = 0.0;
162 double x_2 = 0.0;
163 double y0 = 0.0;
164 double y1 = 0.0;
165 double y_2 = 0.0;
166
167 //check that the data contain the variable we are interested in
168 Int_t bad = 0;
169 const RooAbsReal & xx = x.arg();
170 const RooAbsReal & yy = y.arg();
171 if(! (RooRealVar*)( (RooArgSet*)data.get(0) )->find( xx.GetName() ) )
172 {
173 cout << "Roo2DKeysPdf::Roo2DKeysPdf invalid RooAbsReal name: "<<xx.GetName()<<" not in the data set" <<endl;
174 bad = 1;
175 }
176 if(! (RooRealVar*)( (RooArgSet*)data.get(0) )->find( yy.GetName() ) )
177 {
178 cout << "Roo2DKeysPdf::Roo2DKeysPdf invalid RooAbsReal name: "<<yy.GetName()<<" not in the data set" << endl;
179 bad = 1;
180 }
181 if(bad)
182 {
183 cout << "Roo2DKeysPdf::Roo2DKeysPdf Unable to initialize object; incompatible RooDataSet doesn't contain"<<endl;
184 cout << " all of the RooAbsReal arguments"<<endl;
185 return 1;
186 }
187
188 //copy the data into local arrays
189 const RooArgSet * values = data.get();
190 const RooRealVar* X = ((RooRealVar*)(values->find(xx.GetName())) ) ;
191 const RooRealVar* Y = ((RooRealVar*)(values->find(yy.GetName())) ) ;
192
193 for (Int_t j=0;j<_nEvents;++j)
194 {
195 data.get(j) ;
196
197 _x[j] = X->getVal() ;
198 _y[j] = Y->getVal() ;
199
200 x0+=1; x1+=_x[j]; x_2+=_x[j]*_x[j];
201 y0+=1; y1+=_y[j]; y_2+=_y[j]*_y[j];
202 }
203
204 //==========================================//
205 //calculate the mean and sigma for the data //
206 //==========================================//
207 if(_nEvents == 0)
208 {
209 cout << "Roo2DKeysPdf::Roo2DKeysPdf Empty data set was used; can't generate a PDF"<<endl;
210 }
211
212 _xMean = x1/x0;
214
215 _yMean = y1/y0;
217
219
220 //calculate the PDF
222}
223
224
225////////////////////////////////////////////////////////////////////////////////
226
228{
229 if(_verbosedebug) { cout << "Roo2DKeysPdf::setOptions" << endl; }
230
231 options.ToLower();
232 if( options.Contains("a") ) _BandWidthType = 0;
233 else _BandWidthType = 1;
234 if( options.Contains("n") ) _BandWidthType = 1;
235 else _BandWidthType = 0;
236 if( options.Contains("m") ) _MirrorAtBoundary = 1;
237 else _MirrorAtBoundary = 0;
238 if( options.Contains("d") ) _debug = 1;
239 else _debug = 0;
240 if( options.Contains("v") ) { _debug = 1; _verbosedebug = 1; }
241 else _verbosedebug = 0;
242 if( options.Contains("vv") ) { _vverbosedebug = 1; }
243 else _vverbosedebug = 0;
244
245 if( _debug )
246 {
247 cout << "Roo2DKeysPdf::setOptions(TString options) options = "<< options << endl;
248 cout << "\t_BandWidthType = " << _BandWidthType << endl;
249 cout << "\t_MirrorAtBoundary = " << _MirrorAtBoundary << endl;
250 cout << "\t_debug = " << _debug << endl;
251 cout << "\t_verbosedebug = " << _verbosedebug << endl;
252 cout << "\t_vverbosedebug = " << _vverbosedebug << endl;
253 }
254}
255
256
257////////////////////////////////////////////////////////////////////////////////
258
260{
261 cout << "Roo2DKeysPdf::getOptions(void)" << endl;
262 cout << "\t_BandWidthType = " << _BandWidthType << endl;
263 cout << "\t_MirrorAtBoundary = " << _MirrorAtBoundary << endl;
264 cout << "\t_debug = " << _debug << endl;
265 cout << "\t_verbosedebug = " << _verbosedebug << endl;
266 cout << "\t_vverbosedebug = " << _vverbosedebug << endl;
267}
268
269
270////////////////////////////////////////////////////////////////////////////////
271/// Calculates the kernel bandwidth for x & y and the probability look up table _p[i][j]
272/// \param[in] kernel
273
275{
276 if(_verbosedebug) { cout << "Roo2DKeysPdf::calculateBandWidth(Int_t kernel)" << endl; }
277 if(kernel != -999)
278 {
279 _BandWidthType = kernel;
280 }
281
282 double h = 0.0;
283
284 double sigSum = _xSigma*_xSigma + _ySigma*_ySigma;
285 double sqrtSum = sqrt( sigSum );
286 double sigProd = _ySigma*_xSigma;
287 if(sigProd != 0.0) h = _n16*sqrt( sigSum/sigProd );
288 if(sqrtSum == 0)
289 {
290 cout << "Roo2DKeysPdf::calculateBandWidth The sqr(variance sum) == 0.0. " << " Your dataset represents a delta function."<<endl;
291 return 1;
292 }
293
294 double hXSigma = h * _xSigma;
295 double hYSigma = h * _ySigma;
296 double xhmin = hXSigma * sqrt(2.)/10; //smallest anticipated bandwidth
297 double yhmin = hYSigma * sqrt(2.)/10;
298
299 //////////////////////////////////////
300 //calculate bandwidths from the data//
301 //////////////////////////////////////
302 if(_BandWidthType == 1) //calculate a trivial bandwidth
303 {
304 cout << "Roo2DKeysPdf::calculateBandWidth Using a normal bandwidth (same for a given dimension) based on"<<endl;
305 cout << " h_j = n^{-1/6}*sigma_j for the j^th dimension and n events * "<<_widthScaleFactor<<endl;
306 double hxGaussian = _n16 * _xSigma * _widthScaleFactor;
307 double hyGaussian = _n16 * _ySigma * _widthScaleFactor;
308 for(Int_t j=0;j<_nEvents;++j)
309 {
310 _hx[j] = hxGaussian;
311 _hy[j] = hyGaussian;
312 if(_hx[j]<xhmin) _hx[j] = xhmin;
313 if(_hy[j]<yhmin) _hy[j] = yhmin;
314 }
315 }
316 else //use an adaptive bandwidth to reduce the dependence on global data distribution
317 {
318 cout << "Roo2DKeysPdf::calculateBandWidth Using an adaptive bandwidth (in general different for all events) [default]"<<endl;
319 cout << " scaled by a factor of "<<_widthScaleFactor<<endl;
320 double xnorm = h * TMath::Power(_xSigma/sqrtSum, 1.5) * _widthScaleFactor;
321 double ynorm = h * TMath::Power(_ySigma/sqrtSum, 1.5) * _widthScaleFactor;
322 for(Int_t j=0;j<_nEvents;++j)
323 {
324 double f_ti = TMath::Power( g(_x[j], _x, hXSigma, _y[j], _y, hYSigma), -0.25 ) ;
325 _hx[j] = xnorm * f_ti;
326 _hy[j] = ynorm * f_ti;
327 if(_hx[j]<xhmin) _hx[j] = xhmin;
328 if(_hy[j]<yhmin) _hy[j] = yhmin;
329 }
330 }
331
332 return 0;
333}
334
335
336////////////////////////////////////////////////////////////////////////////////
337/// Evaluates the kernel estimation for x,y, interpolating between the points if necessary
338///
339/// Uses the caching intrinsic in RFC to bypass the grid and remove
340/// the grid and extrapolation approximation in the kernel estimation method
341/// implementation.
342
344{
345 if(_vverbosedebug) { cout << "Roo2DKeysPdf::evaluate()" << endl; }
346 return evaluateFull(x,y);
347}
348
349
350////////////////////////////////////////////////////////////////////////////////
351/// Evaluates the sum of the product of the 2D kernels
352/// for use in calculating the fixed kernel estimate, f,
353/// given the bandwidths _hx[j] and _hy[j].
354///
355/// _n is calculated once in the constructor.
356/// \param[in] thisX
357/// \param[in] thisY
358
359double Roo2DKeysPdf::evaluateFull(double thisX, double thisY) const
360{
361 if( _vverbosedebug ) { cout << "Roo2DKeysPdf::evaluateFull()" << endl; }
362
363 double f=0.0;
364
365 double rx2, ry2, zx, 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
377 zx += highBoundaryCorrection(thisX, _hx[j], x.max(), _x[j])
378 + lowBoundaryCorrection(thisX, _hx[j], x.min(), _x[j]);
379 zy += highBoundaryCorrection(thisY, _hy[j], y.max(), _y[j])
380 + lowBoundaryCorrection(thisY, _hy[j], y.min(), _y[j]);
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) { cout << "Roo2DKeysPdf::highBoundaryCorrection" << 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) { cout << "Roo2DKeysPdf::lowBoundaryCorrection" << 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) cout << "The Bandwidth Type selected is Trivial" << endl;
469 else cout << "The Bandwidth Type selected is Adaptive" << 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 cout << "Roo2DKeysPdf::getMean unknown axis "<<axis<<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 cout << "Roo2DKeysPdf::getSigma unknown axis "<<axis<<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";
510 TString nName = name;
511 nName += "_Ntuple";
512 writeHistToFile( outputFile, histName);
513 writeNTupleToFile( outputFile, nName);
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 TFile * file = 0;
526 cout << "Roo2DKeysPdf::writeHistToFile This member function is temporarily disabled" <<endl;
527 //make sure that any existing file is not over written
528 file = new TFile(outputFile, "UPDATE");
529 if (!file)
530 {
531 cout << "Roo2DKeysPdf::writeHistToFile unable to open file "<< outputFile <<endl;
532 return;
533 }
534
535
536 const RooAbsReal & xx = x.arg();
537 const RooAbsReal & yy = y.arg();
538 RooArgSet values( RooArgList( xx, yy ));
539 RooRealVar * xArg = ((RooRealVar*)(values.find(xx.GetName())) ) ;
540 RooRealVar * yArg = ((RooRealVar*)(values.find(yy.GetName())) ) ;
541
542 TH2F * hist = (TH2F*)xArg->createHistogram("hist", *yArg);
543 hist = (TH2F*)this->fillHistogram(hist, RooArgList(*xArg, *yArg) );
544 hist->SetName(histName);
545
546 file->Write();
547 file->Close();
548}
549
550
551////////////////////////////////////////////////////////////////////////////////
552/// Saves the data and calculated bandwidths to a file,
553/// as a record of what produced the PDF and to give a reduced
554/// data set in order to facilitate re-calculation in the future.
555/// \param[in] outputFile Name of the file where to store the data
556/// \param[in] name Name of the tree which will contain the data
557
558void Roo2DKeysPdf::writeNTupleToFile(char * outputFile, const char * name) const
559{
560 TFile * file = 0;
561
562 //make sure that any existing file is not over written
563 file = new TFile(outputFile, "UPDATE");
564 if (!file)
565 {
566 cout << "Roo2DKeysPdf::writeNTupleToFile unable to open file "<< outputFile <<endl;
567 return;
568 }
569 RooAbsReal & xArg = (RooAbsReal&)x.arg();
570 RooAbsReal & yArg = (RooAbsReal&)y.arg();
571
572 double theX, theY, hx/*, hy*/;
573 TString label = name;
574 label += " the source data for 2D Keys PDF";
575 TTree * _theTree = new TTree(name, label);
576 if(!_theTree) { cout << "Unable to get a TTree for output" << endl; return; }
577 _theTree->SetAutoSave(1000000000); // autosave when 1 Gbyte written
578
579 //name the TBranches the same as the RooAbsReal's
580 const char * xname = xArg.GetName();
581 const char * yname = yArg.GetName();
582 if (!strcmp(xname,"")) xname = "x";
583 if (!strcmp(yname,"")) yname = "y";
584
585 _theTree->Branch(xname, &theX, " x/D");
586 _theTree->Branch(yname, &theY, " y/D");
587 _theTree->Branch("hx", &hx, " hx/D");
588 _theTree->Branch("hy", &hx, " hy/D");
589
590 for(Int_t iEvt = 0; iEvt < _nEvents; iEvt++)
591 {
592 theX = _x[iEvt];
593 theY = _y[iEvt];
594 hx = _hx[iEvt];
595 hx = _hy[iEvt];
596 _theTree->Fill();
597 }
598 file->Write();
599 file->Close();
600}
601
602
603////////////////////////////////////////////////////////////////////////////////
604/// Prints out _p[_nPoints][_nPoints] indicating the domain limits.
605/// \param[out] out Output stream where to print
606
607void Roo2DKeysPdf::PrintInfo(ostream & out) const
608{
609 out << "Roo2DKeysPDF instance domain information:"<<endl;
610 out << "\tX_min = " << _lox <<endl;
611 out << "\tX_max = " << _hix <<endl;
612 out << "\tY_min = " << _loy <<endl;
613 out << "\tY_max = " << _hiy <<endl;
614
615 out << "Data information:" << endl;
616 out << "\t<x> = " << _xMean <<endl;
617 out << "\tsigma(x) = " << _xSigma <<endl;
618 out << "\t<y> = " << _yMean <<endl;
619 out << "\tsigma(y) = " << _ySigma <<endl;
620
621 out << "END of info for Roo2DKeys pdf instance"<< endl;
622}
#define d(i)
Definition: RSha256.hxx:102
#define f(i)
Definition: RSha256.hxx:104
#define h(i)
Definition: RSha256.hxx:106
int Int_t
Definition: RtypesCore.h:45
#define ClassImp(name)
Definition: Rtypes.h:375
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.
Definition: Roo2DKeysPdf.h:26
RooRealProxy y
Definition: Roo2DKeysPdf.h:79
double getMean(const char *axis) const
double _yoffset
Definition: Roo2DKeysPdf.h:111
Roo2DKeysPdf(const char *name, const char *title, RooAbsReal &xx, RooAbsReal &yy, RooDataSet &data, TString options="a", double widthScaleFactor=1.0)
Constructor.
double _widthScaleFactor
Definition: Roo2DKeysPdf.h:112
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
Definition: Roo2DKeysPdf.h:115
void writeToFile(char *outputFile, const char *name) const
Int_t _verbosedebug
Definition: Roo2DKeysPdf.h:118
double _ySigma
Definition: Roo2DKeysPdf.h:103
double * _hy
Definition: Roo2DKeysPdf.h:98
double _sqrt2pi
Definition: Roo2DKeysPdf.h:106
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
Definition: Roo2DKeysPdf.h:116
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
Definition: Roo2DKeysPdf.h:119
double * _x
Definition: Roo2DKeysPdf.h:95
double * _y
Definition: Roo2DKeysPdf.h:97
double * _hx
Definition: Roo2DKeysPdf.h:96
double getSigma(const char *axis) const
double _xoffset
Definition: Roo2DKeysPdf.h:110
~Roo2DKeysPdf() override
Destructor.
void setWidthScaleFactor(double widthScaleFactor)
Definition: Roo2DKeysPdf.h:124
double _xSigma
Definition: Roo2DKeysPdf.h:101
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
Definition: Roo2DKeysPdf.h:78
RooAbsArg * find(const char *name) const
Find object with given name in list.
TH1 * createHistogram(const char *name, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:64
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:94
TH1 * fillHistogram(TH1 *hist, const RooArgList &plotVars, double scaleFactor=1, const RooArgSet *projectedVars=0, bool scaling=true, const RooArgSet *condObs=0, 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:57
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:55
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
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.
A ROOT file is a suite of consecutive data records (TKey instances) with a well defined format.
Definition: TFile.h:54
void SetName(const char *name) override
Change the name of this histogram.
Definition: TH1.cxx:8835
2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:257
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
Basic string class.
Definition: TString.h:136
void ToLower()
Change string to lower-case.
Definition: TString.cxx:1150
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:624
A TTree represents a columnar dataset.
Definition: TTree.h:79
virtual Int_t Fill()
Fill all branches.
Definition: TTree.cxx:4571
virtual void SetAutoSave(Long64_t autos=-300000000)
In case of a program crash, it will be possible to recover the data in the tree up to the last AutoSa...
Definition: TTree.cxx:8299
TBranch * Branch(const char *name, T *obj, Int_t bufsize=32000, Int_t splitlevel=99)
Add a new branch, and infer the data type from the type of obj being passed.
Definition: TTree.h:350
RVec< PromoteType< T > > exp(const RVec< T > &v)
Definition: RVec.hxx:1744
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
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Returns x raised to the power y.
Definition: TMath.h:718
constexpr Double_t Pi()
Definition: TMath.h:37
Definition: file.py:1