Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooBarlowBeestonLL.cxx
Go to the documentation of this file.
1// @(#)root/roostats:$Id$
2// Author: Kyle Cranmer, George Lewis
3/*************************************************************************
4 * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
5 * All rights reserved. *
6 * *
7 * For the licensing terms see $ROOTSYS/LICENSE. *
8 * For the list of contributors see $ROOTSYS/README/CREDITS. *
9 *************************************************************************/
10
11//////////////////////////////////////////////////////////////////////////////
12/** \class RooStats::HistFactory::RooBarlowBeestonLL
13 * \ingroup HistFactory
14//
15// Class RooBarlowBeestonLL implements the profile likelihood estimator for
16// a given likelihood and set of parameters of interest. The value return by
17// RooBarlowBeestonLL is the input likelihood nll minimized w.r.t all nuisance parameters
18// (which are all parameters except for those listed in the constructor) minus
19// the -log(L) of the best fit. Note that this function is slow to evaluate
20// as a MIGRAD minimization step is executed for each function evaluation
21*/
22
23#include <stdexcept>
24#include <cmath>
25#include <iostream>
26
28#include "RooAbsReal.h"
29#include "RooMsgService.h"
30#include "RooRealVar.h"
31
33#include "RooCategory.h"
34#include "RooSimultaneous.h"
35#include "RooArgList.h"
36
39
40#include <TMath.h>
41
42
43
44
45////////////////////////////////////////////////////////////////////////////////
46
48 RooAbsReal& nllIn /*, const RooArgSet& observables*/) :
49 RooAbsReal(name,title),
50 _nll("input","-log(L) function",this,nllIn)
51 // _obs("paramOfInterest","Parameters of interest",this),
52 // _par("nuisanceParam","Nuisance parameters",this,false,false),
53{
54 // Constructor of profile likelihood given input likelihood nll w.r.t
55 // the given set of variables. The input log likelihood is minimized w.r.t
56 // to all other variables of the likelihood at each evaluation and the
57 // value of the global log likelihood minimum is always subtracted.
58
59 // Determine actual parameters and observables
60 /*
61 std::unique_ptr<RooArgSet> actualObs{nllIn.getObservables(observables)};
62 std::unique_ptr<RooArgSet> actualPars{nllIn.getParameters(observables)};
63
64 _obs.add(*actualObs) ;
65 _par.add(*actualPars) ;
66 */
67}
68
69
70
71////////////////////////////////////////////////////////////////////////////////
72
75 _nll("nll",this,other._nll),
76 // _obs("obs",this,other._obs),
77 // _par("par",this,other._par),
78 _paramFixed(other._paramFixed)
79{
80 // Copy constructor
81
82 // _paramAbsMin.addClone(other._paramAbsMin) ;
83 // _obsAbsMin.addClone(other._obsAbsMin) ;
84
85}
86
87
88////////////////////////////////////////////////////////////////////////////////
89
91 for (auto const *var : static_range_cast<RooRealVar *>(*bin_center)) {
92 RooRealVar* target = static_cast<RooRealVar*>(observables->find(var->GetName())) ;
93 target->setVal(var->getVal()) ;
94 }
95}
96
97
98////////////////////////////////////////////////////////////////////////////////
99
101 bool verbose=false;
102
103 if(!_data) {
104 std::cout << "Error: Must initialize data before initializing cache" << std::endl;
105 throw std::runtime_error("Uninitialized Data");
106 }
107 if(!_pdf) {
108 std::cout << "Error: Must initialize model pdf before initializing cache" << std::endl;
109 throw std::runtime_error("Uninitialized model pdf");
110 }
111
112 // Get the data bin values for all channels and bins
113 std::map< std::string, std::vector<double> > ChannelBinDataMap;
115
116 // Get a list of constraint terms
118 RooArgList constraints;
119 RooArgSet* obsSet = std::unique_ptr<RooArgSet>{_pdf->getObservables(*_data)}.release();
120 FactorizeHistFactoryPdf(*obsSet, *_pdf, obsTerms, constraints);
121
122 if( obsTerms.empty() ) {
123 std::cout << "Error: Found no observable terms with pdf: " << _pdf->GetName()
124 << " using dataset: " << _data->GetName() << std::endl;
125 return;
126 }
127 if( constraints.empty() ) {
128 std::cout << "Error: Found no constraint terms with pdf: " << _pdf->GetName()
129 << " using dataset: " << _data->GetName() << std::endl;
130 return;
131 }
132
133 // Loop over the channels
134 auto simPdf = static_cast<RooSimultaneous*>(_pdf);
135 auto channelCat = static_cast<RooCategory const*>(&simPdf->indexCat());
136 for (const auto& nameIdx : *channelCat) {
137
138 // Warning: channel cat name is not necessarily the same name
139 // as the pdf's (for example, if someone does edits)
140 RooAbsPdf* channelPdf = simPdf->getPdf(nameIdx.first.c_str());
141 std::string channel_name = channelPdf->GetName();
142
143 // First, we check if this channel uses Stat Uncertainties:
144 RooArgList* gammas = nullptr;
145 ParamHistFunc* param_func=nullptr;
147 if( ! hasStatUncert ) {
148 if(verbose) {
149 std::cout << "Channel: " << channel_name
150 << " doesn't have statistical uncertainties"
151 << std::endl;
152 }
153 continue;
154 }
155 else {
156 if(verbose) std::cout << "Found ParamHistFunc: " << param_func->GetName() << std::endl;
157 }
158
159 // Now, loop over the bins in this channel
160 // To Do: Check that the index convention
161 // still works for 2-d (ie matches the
162 // convention in ParamHistFunc, etc)
163 int num_bins = param_func->numBins();
164
165 // Initialize the vector to the number of bins, allowing
166 // us to skip gamma's if they are constant
167
168 std::vector<BarlowCache> temp_cache( num_bins );
170
171 for( Int_t bin_index = 0; bin_index < num_bins; ++bin_index ) {
172
173 // Create a cache object
174 BarlowCache cache;
175
176 // Get the gamma for this bin, and skip if it's constant
177 RooRealVar* gamma_stat = dynamic_cast<RooRealVar*>(&(param_func->getParameter(bin_index)));
178 if(!gamma_stat) throw std::runtime_error("ParamHistFunc contains non-RooRealVar, not supported in RooBarlowBeestonLL");
179 if( gamma_stat->isConstant() ) {
180 if(verbose) std::cout << "Ignoring constant gamma: " << gamma_stat->GetName() << std::endl;
181 continue;
182 }
183 else {
184 cache.hasStatUncert=true;
186 cache.gamma = gamma_stat;
187 _statUncertParams.insert( gamma_stat->GetName() );
188 }
189
190 // Store a snapshot of the bin center
191 RooArgSet* bin_center = (RooArgSet*) param_func->get( bin_index )->snapshot();
192 cache.bin_center = bin_center;
193 cache.observables = obsSet;
194
195 cache.binVolume = param_func->binVolume();
196
197 // Delete this part, simply a comment
198 RooArgList obs_list( *(cache.bin_center) );
199
200 // Get the gamma's constraint term
201 RooAbsReal* pois_mean = nullptr;
202 RooRealVar* tau = nullptr;
204 if( !tau || !pois_mean ) {
205 std::cout << "Failed to find pois mean or tau parameter for " << gamma_stat->GetName() << std::endl;
206 }
207 else {
208 if (verbose) {
209 std::cout << "Found pois mean and tau for parameter: " << gamma_stat->GetName() << " tau: " << tau->GetName()
210 << " " << tau->getVal() << " pois_mean: " << pois_mean->GetName() << " " << pois_mean->getVal()
211 << std::endl;
212 }
213 }
214
215 cache.tau = tau;
216 cache.nom_pois_mean = pois_mean;
217
218 // Get the RooRealSumPdf
220 if( sum_pdf == nullptr ) {
221 std::cout << "Failed to find RooRealSumPdf in channel " << channel_name
222 << ", therefor skipping this channel for analytic uncertainty minimization"
223 << std::endl;
225 break;
226 }
227 cache.sumPdf = sum_pdf;
228
229 // And set the data value for this bin
231 std::cout << "Error: channel with name: " << channel_name
232 << " not found in BinDataMap" << std::endl;
233 throw std::runtime_error("BinDataMap");
234 }
235 double nData = ChannelBinDataMap[channel_name].at(bin_index);
236 cache.nData = nData;
237
238 temp_cache.at(bin_index) = cache;
239 //_barlowCache[channel_name].at(bin_index) = cache;
240
241 } // End: Loop over bins
242
244 std::cout << "Adding channel: " << channel_name
245 << " to the barlow cache" << std::endl;
246 _barlowCache[channel_name] = temp_cache;
247 }
248
249
250 } // End: Loop over channels
251
252
253
254 // Successfully initialized the cache
255 // Printing some info
256 /*
257 std::map< std::string, std::vector< BarlowCache > >::iterator iter_cache;
258 for( iter_cache = _barlowCache.begin(); iter_cache != _barlowCache.end(); ++iter_cache ) {
259
260 std::string channel_name = (*iter_cache).first;
261 std::vector< BarlowCache >& channel_cache = (*iter_cache).second;
262
263
264 for( unsigned int i = 0; i < channel_cache.size(); ++i ) {
265 BarlowCache& bin_cache = channel_cache.at(i);
266
267
268 RooRealVar* gamma = bin_cache.gamma;
269 RooRealVar* tau = bin_cache.tau;
270 RooAbsReal* pois_mean = bin_cache.nom_pois_mean;
271 RooAbsPdf* sum_pdf = (RooAbsPdf*) bin_cache.sumPdf;
272 double binVolume = bin_cache.binVolume;
273
274
275 if( !bin_cache.hasStatUncert ) {
276 std::cout << "Barlow Cache for Channel: " << channel_name
277 << " Bin: " << i
278 << " NO STAT UNCERT"
279 << std::endl;
280 }
281 else {
282 std::cout << "Barlow Cache for Channel: " << channel_name
283 << " Bin: " << i
284 << " gamma: " << gamma->GetName()
285 << " tau: " << tau->GetName()
286 << " pois_mean: " << pois_mean->GetName()
287 << " sum_pdf: " << sum_pdf->GetName()
288 << " binVolume: " << binVolume
289 << std::endl;
290 }
291
292 }
293 }
294 */
295
296}
297
298
299////////////////////////////////////////////////////////////////////////////////
300
303 bool stripDisconnected) const {
305
307 toRemove.reserve( _statUncertParams.size());
308
309 for (auto const& arg : outputSet) {
310
311 // If there is a gamma in the name,
312 // strip it from the list of dependencies
313
314 if( _statUncertParams.find(arg->GetName()) != _statUncertParams.end() ) {
315 toRemove.add( *arg );
316 }
317 }
318
319 for( auto& arg : toRemove) outputSet.remove( *arg, true );
320
321 return errorInBaseCall || false;
322
323}
324
325
326/*
327////////////////////////////////////////////////////////////////////////////////
328
329const RooArgSet& RooStats::HistFactory::RooBarlowBeestonLL::bestFitParams() const
330{
331 validateAbsMin() ;
332 return _paramAbsMin ;
333}
334
335
336////////////////////////////////////////////////////////////////////////////////
337
338const RooArgSet& RooStats::HistFactory::RooBarlowBeestonLL::bestFitObs() const
339{
340 validateAbsMin() ;
341 return _obsAbsMin ;
342}
343*/
344
345
346
347////////////////////////////////////////////////////////////////////////////////
348/// Optimized implementation of createProfile for profile likelihoods.
349/// Return profile of original function in terms of stated parameters
350/// of interest rather than profiling recursively.
351
352/*
353RooAbsReal* RooStats::HistFactory::RooBarlowBeestonLL::createProfile(const RooArgSet& paramsOfInterest)
354{
355 return nll().createProfile(paramsOfInterest) ;
356}
357*/
358
359
360/*
361void RooStats::HistFactory::RooBarlowBeestonLL::FactorizePdf(const RooArgSet &observables, RooAbsPdf &pdf, RooArgList &obsTerms, RooArgList &constraints) const {
362 // utility function to factorize constraint terms from a pdf
363 // (from G. Petrucciani)
364 const std::type_info & id = typeid(pdf);
365 if (id == typeid(RooProdPdf)) {
366 RooProdPdf *prod = dynamic_cast<RooProdPdf *>(&pdf);
367 RooArgList list(prod->pdfList());
368 for (int i = 0, n = list.size(); i < n; ++i) {
369 RooAbsPdf *pdfi = (RooAbsPdf *) list.at(i);
370 FactorizePdf(observables, *pdfi, obsTerms, constraints);
371 }
372 } else if (id == typeid(RooSimultaneous) ) { //|| id == typeid(RooSimultaneousOpt)) {
373 RooSimultaneous *sim = dynamic_cast<RooSimultaneous *>(&pdf);
374 std::unique_ptr<RooAbsCategoryLValue> cat{(RooAbsCategoryLValue *) sim->indexCat().Clone()};
375 for (int ic = 0, nc = cat->numBins((const char *)0); ic < nc; ++ic) {
376 cat->setBin(ic);
377 FactorizePdf(observables, *sim->getPdf(cat->getLabel()), obsTerms, constraints);
378 }
379 } else if (pdf.dependsOn(observables)) {
380 if (!obsTerms.contains(pdf)) obsTerms.add(pdf);
381 } else {
382 if (!constraints.contains(pdf)) constraints.add(pdf);
383 }
384}
385*/
386
387
388
389////////////////////////////////////////////////////////////////////////////////
390
392{
393 // Loop over the channels (keys to the map)
394 //clock_t time_before_setVal, time_after_setVal;
395 //time_before_setVal=clock();
396 std::map< std::string, std::vector< BarlowCache > >::iterator iter_cache;
397 for( iter_cache = _barlowCache.begin(); iter_cache != _barlowCache.end(); ++iter_cache ) {
398
399 std::string channel_name = (*iter_cache).first;
400 std::vector< BarlowCache >& channel_cache = (*iter_cache).second;
401
402 /* Slower way to find the channel vector:
403 // Get the vector of bin uncertainty caches for this channel
404 if( _barlowCache.find( channel_name ) == _barlowCache.end() ) {
405 std::cout << "Error: channel: " << channel_name
406 << " not found in barlow Cache" << std::endl;
407 throw std::runtime_error("Channel not in barlow cache");
408 }
409
410 std::vector< BarlowCache >& channel_cache = _barlowCache[ channel_name ];
411 */
412
413 // Loop over the bins in the cache
414 // Set all gamma's to 0
415 for( unsigned int i = 0; i < channel_cache.size(); ++i ) {
417 if( !bin_cache.hasStatUncert ) continue;
418 RooRealVar* gamma = bin_cache.gamma;
419 gamma->setVal(0.0);
420 }
421 std::vector< double > nu_b_vec( channel_cache.size() );
422 for( unsigned int i = 0; i < channel_cache.size(); ++i ) {
424 if( !bin_cache.hasStatUncert ) continue;
425
426 RooAbsPdf* sum_pdf = static_cast<RooAbsPdf*>(bin_cache.sumPdf);
427 RooArgSet* obsSet = bin_cache.observables;
428 double binVolume = bin_cache.binVolume;
429
430 bin_cache.SetBinCenter();
431 double nu_b = sum_pdf->getVal(*obsSet)*sum_pdf->expectedEvents(*obsSet)*binVolume;
432 nu_b_vec.at(i) = nu_b;
433 }
434
435 // Loop over the bins in the cache
436 // Set all gamma's to 1
437 for( unsigned int i = 0; i < channel_cache.size(); ++i ) {
439 if( !bin_cache.hasStatUncert ) continue;
440 RooRealVar* gamma = bin_cache.gamma;
441 gamma->setVal(1.0);
442 }
443 std::vector< double > nu_b_stat_vec( channel_cache.size() );
444 for( unsigned int i = 0; i < channel_cache.size(); ++i ) {
446 if( !bin_cache.hasStatUncert ) continue;
447
448 RooAbsPdf* sum_pdf = static_cast<RooAbsPdf*>(bin_cache.sumPdf);
449 RooArgSet* obsSet = bin_cache.observables;
450 double binVolume = bin_cache.binVolume;
451
452 bin_cache.SetBinCenter();
453 double nu_b_stat = sum_pdf->getVal(*obsSet)*sum_pdf->expectedEvents(*obsSet)*binVolume - nu_b_vec.at(i);
454 nu_b_stat_vec.at(i) = nu_b_stat;
455 }
456 //time_after_setVal=clock();
457
458 // Done with the first loops.
459 // Now evaluating the function
460
461 //clock_t time_before_eval, time_after_eval;
462
463 // Loop over the bins in the cache
464 //time_before_eval=clock();
465 for( unsigned int i = 0; i < channel_cache.size(); ++i ) {
466
468
469 if( !bin_cache.hasStatUncert ) {
470 //std::cout << "Bin: " << i << " of " << channel_cache.size()
471 // << " in channel: " << channel_name
472 // << " doesn't have stat uncertainties" << std::endl;
473 continue;
474 }
475
476 // Set the observable to the bin center
477 bin_cache.SetBinCenter();
478
479 // Get the cached objects
480 RooRealVar* gamma = bin_cache.gamma;
481 RooRealVar* tau = bin_cache.tau;
482 RooAbsReal* pois_mean = bin_cache.nom_pois_mean;
483 //RooAbsPdf* sum_pdf = (RooAbsPdf*) bin_cache.sumPdf;
484 //RooArgSet* obsSet = bin_cache.observables;
485 //double binVolume = bin_cache.binVolume;
486
487 // Get the values necessary for
488 // the analytic minimization
489 double nu_b = nu_b_vec.at(i);
490 double nu_b_stat = nu_b_stat_vec.at(i);
491
492 double tau_val = tau->getVal();
493 double nData = bin_cache.nData;
494 double m_val = pois_mean->getVal();
495
496 // Initialize the minimized value of gamma
497 double gamma_hat_hat = 1.0;
498
499 // Check that the quadratic term is > 0
500 if(nu_b_stat > 0.00000001) {
501
503 double B = nu_b*tau_val + nu_b*nu_b_stat - nData*nu_b_stat - m_val*nu_b_stat;
504 double C = -1*m_val*nu_b;
505
506 double discrim = B*B-4*A*C;
507
508 if( discrim < 0 ) {
509 std::cout << "Warning: Discriminant (B*B - 4AC) < 0" << std::endl;
510 std::cout << "Warning: Taking B*B - 4*A*C == 0" << std::endl;
511 discrim=0;
512 //throw std::runtime_error("BarlowBeestonLL::evaluate() : B*B - 4AC < 0");
513 }
514 if( A <= 0 ) {
515 std::cout << "Warning: A <= 0" << std::endl;
516 throw std::runtime_error("BarlowBeestonLL::evaluate() : A < 0");
517 }
518
519 gamma_hat_hat = ( -1*B + std::sqrt(discrim) ) / (2*A);
520 }
521
522 // If the quadratic term is 0, we simply
523 // use a linear equation
524 else {
525 gamma_hat_hat = m_val/tau_val;
526 }
527
528 // Check for NAN
530 std::cout << "ERROR: gamma hat hat is NAN" << std::endl;
531 throw std::runtime_error("BarlowBeestonLL::evaluate() : gamma hat hat is NAN");
532 }
533
534 if( gamma_hat_hat <= 0 ) {
535 std::cout << "WARNING: gamma hat hat <= 0. Setting to 0" << std::endl;
536 gamma_hat_hat = 0;
537 }
538
539 /*
540 std::cout << "n: " << bin_cache.nData << " "
541 << "nu_stat: " << nu_b_stat << " "
542 << "nu: " << nu_b << " "
543 << "tau: " << tau->getVal() << " "
544 << "m: " << pois_mean->getVal() << " "
545 << "A: " << A << " "
546 << "B: " << B << " "
547 << "C: " << C << " "
548 << "gamma hat hat: " << gamma_hat_hat
549 << std::endl;
550 */
551
552 gamma->setVal( gamma_hat_hat );
553
554 }
555
556 //time_after_eval=clock();
557
558 //float time_setVal = ((float) time_after_setVal - (float) time_before_setVal) / ((float) CLOCKS_PER_SEC);
559 //float time_eval = ((float) time_after_eval - (float) time_before_eval) / ((float) CLOCKS_PER_SEC);
560
561 /*
562 std::cout << "Barlow timing for channel: " << channel_name
563 << " SetVal: " << time_setVal
564 << " Eval: " << time_eval
565 << std::endl;
566 */
567 }
568
569
570 return _nll;
571
572}
573
574
575
576/*
577////////////////////////////////////////////////////////////////////////////////
578/// Check that parameters and likelihood value for 'best fit' are still valid. If not,
579/// because the best fit has never been calculated, or because constant parameters have
580/// changed value or parameters have changed const/float status, the minimum is recalculated
581
582void RooStats::HistFactory::RooBarlowBeestonLL::validateAbsMin() const
583{
584 // Check if constant status of any of the parameters have changed
585 if (_absMinValid) {
586 _piter->Reset() ;
587 RooAbsArg* par ;
588 while((par=(RooAbsArg*)_piter->Next())) {
589 if (_paramFixed[par->GetName()] != par->isConstant()) {
590 cxcoutI(Minimization) << "RooStats::HistFactory::RooBarlowBeestonLL::evaluate(" << GetName() << ") constant status of parameter " << par->GetName() << " has changed from "
591 << (_paramFixed[par->GetName()]?"fixed":"floating") << " to " << (par->isConstant()?"fixed":"floating")
592 << ", recalculating absolute minimum" << std::endl ;
593 _absMinValid = false ;
594 break ;
595 }
596 }
597 }
598
599
600 // If we don't have the absolute minimum w.r.t all observables, calculate that first
601 if (!_absMinValid) {
602
603 cxcoutI(Minimization) << "RooStats::HistFactory::RooBarlowBeestonLL::evaluate(" << GetName() << ") determining minimum likelihood for current configurations w.r.t all observable" << std::endl ;
604
605
606 // Save current values of non-marginalized parameters
607 std::unique_ptr<RooArgSet> obsStart{(RooArgSet*) _obs.snapshot(false)};
608
609 // Start from previous global minimum
610 if (_paramAbsMin.size()>0) {
611 const_cast<RooSetProxy&>(_par).assignValueOnly(_paramAbsMin) ;
612 }
613 if (_obsAbsMin.size()>0) {
614 const_cast<RooSetProxy&>(_obs).assignValueOnly(_obsAbsMin) ;
615 }
616
617 // Find minimum with all observables floating
618 const_cast<RooSetProxy&>(_obs).setAttribAll("Constant",false) ;
619 _minuit->migrad() ;
620
621 // Save value and remember
622 _absMin = _nll ;
623 _absMinValid = true ;
624
625 // Save parameter values at abs minimum as well
626 _paramAbsMin.removeAll() ;
627
628 // Only store non-constant parameters here!
629 std::unique_ptr<RooArgSet> tmp{_par.selectByAttrib("Constant",false)};
630 _paramAbsMin.addClone(*tmp) ;
631
632 _obsAbsMin.addClone(_obs) ;
633
634 // Save constant status of all parameters
635 _piter->Reset() ;
636 RooAbsArg* par ;
637 while((par=(RooAbsArg*)_piter->Next())) {
638 _paramFixed[par->GetName()] = par->isConstant() ;
639 }
640
641 if (dologI(Minimization)) {
642 cxcoutI(Minimization) << "RooStats::HistFactory::RooBarlowBeestonLL::evaluate(" << GetName() << ") minimum found at (" ;
643
644 RooAbsReal* arg ;
645 bool first=true ;
646 _oiter->Reset() ;
647 while ((arg=(RooAbsReal*)_oiter->Next())) {
648 ccxcoutI(Minimization) << (first?"":", ") << arg->GetName() << "=" << arg->getVal() ;
649 first=false ;
650 }
651 ccxcoutI(Minimization) << ")" << std::endl ;
652 }
653
654 // Restore original parameter values
655 const_cast<RooSetProxy&>(_obs) = *obsStart ;
656
657 }
658}
659*/
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 char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t target
char name[80]
Definition TGX11.cxx:110
A class which maps the current values of a RooRealVar (or a set of RooRealVars) to one of a number of...
const_iterator begin() const
const_iterator end() const
RooFit::OwningPtr< RooArgSet > getParameters(const RooAbsData *data, bool stripDisconnected=true) const
Create a list of leaf nodes in the arg tree starting with ourself as top node that don't match any of...
Abstract interface for all probability density functions.
Definition RooAbsPdf.h:40
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:59
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:103
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
Object to represent discrete states.
Definition RooCategory.h:28
Variable that can be changed from the outside.
Definition RooRealVar.h:37
Facilitates simultaneous fitting of multiple PDFs to subsets of a given dataset.
Class RooBarlowBeestonLL implements the profile likelihood estimator for a given likelihood and set o...
double evaluate() const override
Optimized implementation of createProfile for profile likelihoods.
bool getParameters(const RooArgSet *depList, RooArgSet &outputSet, bool stripDisconnected=true) const override
Fills a list with leaf nodes in the arg tree starting with ourself as top node that don't match any o...
RooBarlowBeestonLL(const char *name, const char *title, RooAbsReal &nll)
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
bool getStatUncertaintyFromChannel(RooAbsPdf *channel, ParamHistFunc *&paramfunc, RooArgList *gammaList)
void FactorizeHistFactoryPdf(const RooArgSet &, RooAbsPdf &, RooArgList &, RooArgList &)
void getDataValuesForObservables(std::map< std::string, std::vector< double > > &ChannelBinDataMap, RooAbsData *data, RooAbsPdf *simPdf)
RooAbsPdf * getSumPdfFromChannel(RooAbsPdf *channel)
int getStatUncertaintyConstraintTerm(RooArgList *constraints, RooRealVar *gamma_stat, RooAbsReal *&pois_mean, RooRealVar *&tau)
Bool_t IsNaN(Double_t x)
Definition TMath.h:896