FrequentistSurface.cxx
Go to the documentation of this file.
2 
3 #include <iostream>
4 #include <functional>
5 #include <limits>
6 
7 #include "TCanvas.h"
8 #include "TGraph.h"
9 #include "TH2F.h"
10 #include "TObjArray.h"
11 #include "TPad.h"
12 #include "TROOT.h"
13 #include "TStyle.h"
14 #include "TKey.h"
15 #include "TVectorD.h"
16 #include "TObjString.h"
17 //#include "TCollection.h"
18 
19 #include "CAFAna/Core/Var.h"
23 #include "CAFAna/Core/IFitVar.h"
24 #include "CAFAna/Core/Progress.h"
25 #include "CAFAna/Core/ThreadPool.h"
26 #include "CAFAna/Core/Utilities.h"
27 
28 #include "OscLib/IOscCalc.h"
29 
30 namespace ana
31 {
32  //----------------------------------------------------------------------
35  const FitAxis& xax,
36  const FitAxis& yax,
37  const std::vector<const IFitVar*>& profVars,
38  const std::vector<const ISyst*>& profSysts,
39  const SeedList& seedPts,
40  const std::vector<SystShifts>& systSeedPts,
41  bool parallel,
43  : fParallel(parallel), fFitOpts(opts)
44  {
45  fLogX = xax.islog;
46  fLogY = yax.islog;
47 
48  CreateHistograms(xax, yax, profVars, profSysts);
49 
50  for(const IFitVar* v: profVars) fSeedValues.push_back(v->GetValue( calc));
51 
52  FillSurface(expt, calc, xax.var, yax.var, profVars, profSysts, seedPts, systSeedPts);
53 
54  FindMinimum(expt, calc, xax.var, yax.var, profVars, profSysts, seedPts, systSeedPts);
55  }
56 
57  //---------------------------------------------------------------------
59  {
60  }
61 
62  //---------------------------------------------------------------------
64  CreateHistograms(const FitAxis& xax, const FitAxis& yax,
65  const std::vector<const IFitVar*>& profVars,
66  const std::vector<const ISyst*>& profSysts)
67  {
68  fHist = ExpandedHistogram(";"+xax.var->LatexName()+";"+yax.var->LatexName(),
69  xax.nbins, xax.min, xax.max, xax.islog,
70  yax.nbins, yax.min, yax.max, yax.islog);
71 
72  for(unsigned int i = 0; i < profVars.size()+profSysts.size(); ++i){
74  if(i < profVars.size())
75  title = profVars[i]->LatexName();
76  else
77  title = profSysts[i-profVars.size()]->LatexName();
78 
79  fProfHists.push_back(ExpandedHistogram(title+";"+xax.var->LatexName()+";"+yax.var->LatexName(),
80  xax.nbins, xax.min, xax.max, xax.islog,
81  yax.nbins, yax.min, yax.max, yax.islog));
82  }
83 
84  }
85 
86  //---------------------------------------------------------------------
88  ProgressBarTitle(const IFitVar* xvar, const IFitVar* yvar,
89  const std::vector<const IFitVar*>& profVars,
90  const std::vector<const ISyst*>& profSysts) const
91  {
92  std::string title = "Filling surface for "+yvar->ShortName()+" vs "+xvar->ShortName();
93 
94  if(!profVars.empty() || !profSysts.empty()){
95  title += " (profiling ";
96 
97  for(const IFitVar* v: profVars) title += v->ShortName() + ", ";
98  for(const ISyst* s: profSysts) title += s->ShortName() + ", ";
99 
100  // Always have one superfluous ", " at the end
101  title.resize(title.size()-2);
102  title += ")";
103  }
104 
105  return title;
106  }
107 
108  //---------------------------------------------------------------------
111  const IFitVar *xvar, const IFitVar *yvar,
112  const std::vector<const IFitVar *> &profVars,
113  const std::vector<const ISyst *> &profSysts,
114  const SeedList& seedPts,
115  const std::vector<SystShifts> &systSeedPts)
116  {
117  // Nothing created during surface filling belongs in a
118  // directory. Unfortunately the local guards in Spectrum etc are racey when
119  // run in parallel. But this should cover the whole lot safely.
120  DontAddDirectory guard;
121 
122  const std::string progTitle = ProgressBarTitle(xvar, yvar, profVars, profSysts);
123 
124  Progress *prog = 0;
125  // Difficult to keep a progress bar properly up to date when threaded
126  if (!fParallel) prog = new Progress(progTitle);
127  ThreadPool *pool = 0;
128 
129  if(fParallel){
130  // A hack/workaround needed for parallel running:
131  //
132  // Give all the constituents of the Prediction a chance to do lazy
133  // initialization, before they race themselves trying to do it in
134  // parallel.
135  expt->ChiSq(calc);
136 
137  pool = new ThreadPool;
138  pool->ShowProgress(progTitle + TString::Format(" using %d threads", pool->NThreads()).Data());
139  }
140 
141  const int Nx = fHist->GetNbinsX();
142  const int Ny = fHist->GetNbinsY();
143 
144  // Fill bins in "random" order so that the progress bar is accurate
145  // straight away instead of possibly being misled by whatever atypical
146  // points we start with. This step is a prime which guarantees we get every
147  // cell.
148  int step = 7919;
149  // Very unlikely (Nx or Ny is a multiple of step), but just to be safe.
150  if ((Nx * Ny) % step == 0) step = 1;
151 
152  int bin = 0;
153  int neval = 0;
154 
155  // Allow the surface to be parallelised across multiple jobs by splitting up
156  // the full surface into patches, and only running bins that fall inside a
157  // certain patch
158  int first = 0, last = Nx * Ny;
159  if (RunningOnGrid() && NumJobs()>1) {
160  double stride = double(Nx*Ny) / double(NumJobs());
161  int job = JobNumber();
162  first = job * stride;
163  last = (job+1) * stride;
164  }
165 
166  do{
167  if (neval < first) {
168  ++neval;
169  bin = (bin + step) % (Nx * Ny);
170  continue;
171  } else if (neval >= last) {
172  break;
173  }
174 
175  const int x = bin % Nx + 1;
176  const int y = bin / Nx + 1;
177 
178  const double xv = BinCenterX(x);
179  const double yv = BinCenterY(y);
180 
181  if (xvar->Penalty(xv, calc) > 1e-10)
182  {
183  std::cerr << "Warning! " << xvar->ShortName() << " = " << xv
184  << " has penalty of " << xvar->Penalty(xv, calc)
185  << " that could have been applied in surface. "
186  << "This should never happen." << std::endl;
187  }
188  if (yvar->Penalty(yv, calc) > 1e-10)
189  {
190  std::cerr << "Warning! " << yvar->ShortName() << " = " << yv
191  << " has penalty of " << yvar->Penalty(yv, calc)
192  << " that could have been applied in surface. "
193  << "This should never happen." << std::endl;
194  }
195 
196  ThreadPool::func_t task = [=](){
197  FillSurfacePoint(expt, calc,
198  xvar, xv, yvar, yv,
199  profVars, profSysts, seedPts, systSeedPts);
200  };
201 
202  ++neval;
203  if(fParallel){
204  pool->AddTask(task);
205  }
206  else{
207  task(); // Just do it straight away
208  prog->SetProgress(double(neval-first) / double(last-first));
209  }
210 
211  if (RunningOnGrid() && NumJobs()>1) fBinMask.push_back(bin);
212 
213  bin = (bin + step) % (Nx * Ny);
214  } while (bin != 0);
215 
216 
217  if(fParallel){
218  pool->Finish();
219  delete pool;
220  }
221  else{
222  prog->Done();
223  delete prog;
224  }
225  }
226 
227  //----------------------------------------------------------------------
230  const IFitVar* xvar, double x,
231  const IFitVar* yvar, double y,
232  const std::vector<const IFitVar*>& profVars,
233  const std::vector<const ISyst*>& profSysts,
234  const SeedList& seedPts,
235  const std::vector<SystShifts>& systSeedPts)
236  {
237  if(fParallel){
238  // Need to take our own copy so that we don't get overwritten by someone
239  // else's changes.
240  calc = calc->Copy();
241  }
242 
243  xvar->SetValue(calc, x);
244  yvar->SetValue(calc, y);
245 
246  //Make sure that the profiled values of fitvars do not persist between steps.
247  for(int i = 0; i < (int)fSeedValues.size(); ++i) profVars[i]->SetValue( calc, fSeedValues[i] );
248 
249  expt->Reset();
250 
251  double chi;
252  if(profVars.empty() && profSysts.empty()){
253  chi = expt->ChiSq(calc);
254  }
255  else{
256  MinuitFitter fitter(expt, profVars, profSysts);
257  fitter.SetFitOpts(fFitOpts);
258  SystShifts bestSysts;
259  chi = fitter.Fit(calc, bestSysts, seedPts, systSeedPts, MinuitFitter::kQuiet)->EvalMetricVal();
260 
261  for(unsigned int i = 0; i < profVars.size(); ++i){
262  fProfHists[i]->Fill(x, y, profVars[i]->GetValue(calc));
263  }
264  for(unsigned int j = 0; j < profSysts.size(); ++j){
265  fProfHists[j+profVars.size()]->Fill(x, y, bestSysts.GetShift(profSysts[j]));
266  }
267  }
268 
269  fHist->Fill(x, y, chi);
270 
271  if(fParallel) delete calc;
272 
273  return chi;
274  }
275 
276 
277  //---------------------------------------------------------------------
280  const IFitVar* xvar, const IFitVar* yvar,
281  const std::vector<const IFitVar*>& profVars,
282  const std::vector<const ISyst*>& profSysts,
283  const SeedList& seedPts,
284  const std::vector<SystShifts>& systSeedPts)
285  {
286  // Location of the best minimum found from filled surface
287  double minchi = 1e10;
288  int minx = fHist->GetNbinsX()/2;
289  int miny = fHist->GetNbinsY()/2;
290  for(int x = 1; x <= fHist->GetNbinsX(); ++x){
291  for(int y = 1; y <= fHist->GetNbinsY(); ++y){
292  int bin = ((y-1)*fHist->GetNbinsX())+(x-1);
293  const double chi = fHist->GetBinContent(x, y);
294  if(RunningOnGrid() && NumJobs()>1
295  && !std::count(fBinMask.begin(), fBinMask.end(), bin)) continue;
296  if (chi < minchi){
297  minchi = chi;
298  minx = x;
299  miny = y;
300  }
301  }
302  }
303 
304  std::vector<const IFitVar*> allVars = {xvar, yvar};
305  allVars.insert(allVars.end(), profVars.begin(), profVars.end());
306  MinuitFitter fit(expt, allVars, profSysts);
307  fit.SetFitOpts(fFitOpts);
308  expt->Reset();
309  // Seed from best grid point
310  xvar->SetValue(calc, BinCenterX(minx));
311  yvar->SetValue(calc, BinCenterY(miny));
312  for(int i = 0; i < (int)fSeedValues.size(); ++i) profVars[i]->SetValue( calc, fSeedValues[i] );
313  SystShifts systSeed = SystShifts::Nominal();
314  fBestLikelihood = fit.Fit(calc, systSeed, seedPts)->EvalMetricVal();
315  fBestFitX = xvar->GetValue(calc);
316  fBestFitY = yvar->GetValue(calc);
317 
318  for(int x = 0; x < fHist->GetNbinsX()+2; ++x){
319  for(int y = 0; y < fHist->GetNbinsY()+2; ++y){
320  fHist->SetBinContent(x, y, fHist->GetBinContent(x, y)-fBestLikelihood);
321  }
322  }
323 
324  fHist->SetMinimum(0);
325  }
326 
327  //----------------------------------------------------------------------
328  void FrequentistSurface::SaveTo(TDirectory* dir, const std::string& name) const
329  {
330  TDirectory *tmp = gDirectory;
331 
332  dir = dir->mkdir(name.c_str()); // switch to subdir
333  dir->cd();
334 
335  TObjString("FrequentistSurface").Write("type");
336 
338 
339  TDirectory *profDir = dir->mkdir("profHists");
340  int idx = 0;
342  {
343  profDir->cd();
344  it->Write(TString::Format("hist%d", idx++));
345  }
346 
347  dir->Write();
348  delete dir;
349 
350  tmp->cd();
351  }
352 
353  //----------------------------------------------------------------------
354  std::unique_ptr<FrequentistSurface> FrequentistSurface::LoadFrom(TDirectory* dir, const std::string& name)
355  {
356  dir = dir->GetDirectory(name.c_str()); // switch to subdir
357  assert(dir);
358 
359  DontAddDirectory guard;
360 
361  TObjString *tag = (TObjString *) dir->Get("type");
362  assert(tag);
363  assert(tag->GetString() == "FrequentistSurface" || tag->GetString() == "Surface");
364  delete tag;
365 
366  std::unique_ptr<FrequentistSurface> surf(new FrequentistSurface);
367  ISurface::FillSurfObj(*surf, dir);
368 
369  for (std::size_t idx = 0; idx < surf->fSeedValues.size(); ++idx)
370  {
371  // Search for old "marg" name here too for backwards compatibility
372  TH2 *h = (TH2 *) dir->Get(TString::Format("profHists/hist%lu", idx));
373  if (h)
374  surf->fProfHists.push_back(h);
375  else
376  surf->fProfHists.push_back((TH2 *) dir->Get(TString::Format("margHists/hist%lu", idx)));
377  }
378 
379  const TObjString* logx = (TObjString*)dir->Get("logx");
380  const TObjString* logy = (TObjString*)dir->Get("logy");
381  // Tolerate missing log tags for backwards compatibility
382  surf->fLogX = (logx && logx->GetString() == "yes");
383  surf->fLogY = (logy && logy->GetString() == "yes");
384 
385  delete dir;
386 
387  return surf;
388  }
389 
390  //----------------------------------------------------------------------
391  std::unique_ptr<FrequentistSurface> FrequentistSurface::LoadFromMulti(
392  std::vector<TFile*> files, std::string label)
393  {
394  std::vector<std::unique_ptr<FrequentistSurface>> surfs;
395  for (TFile* f : files) {
396  surfs.push_back(FrequentistSurface::LoadFrom(f, label));
397  }
398 
399  int Nx = surfs[0]->fHist->GetNbinsX();
400  int Ny = surfs[0]->fHist->GetNbinsY();
401  size_t nbins = Nx * Ny;
402  std::vector<int> binMask;
403 
404  // Loop over the surfaces to check all bins are present
405  for (size_t i = 0; i < surfs.size(); ++i) {
406  for (int bin : surfs[i]->fBinMask) {
407  if (std::count(binMask.begin(), binMask.end(), bin))
408  assert(false && "Repeated bin found in surfaces being merged.");
409  binMask.push_back(bin);
410  }
411  }
412  if (binMask.size() != nbins) {
413  std::cout << "Missing bins found in surfaces being merged. "
414  << "Are you sure you included all files for this surface?"
415  << std::endl;
416  assert(false);
417  }
418 
419  DontAddDirectory guard;
420 
421  // Create return surface and initialise with first in list
422  std::unique_ptr<FrequentistSurface> ret(new FrequentistSurface);
423  ret->fHist = (TH2F*)surfs[0]->ToTH2(0)->Clone();
424  for (TH2* h : surfs[0]->fProfHists)
425  ret->fProfHists.push_back((TH2F*)h->Clone());
426 
427  ret->fLogX = surfs[0]->fLogX;
428  ret->fLogY = surfs[0]->fLogY;
429  ret->fSeedValues = surfs[0]->fSeedValues;
430 
431  // Find job that contains minimum bin (and therefore true best fit)
432  int best = -1;
433  double minChi = std::numeric_limits<double>::max();
434  for (size_t i = 0; i < surfs.size(); ++i) {
435  TH2F* hist = (TH2F*)surfs[i]->ToTH2(0);
436  for (size_t bin : surfs[i]->fBinMask) {
437  const int x = bin % Nx + 1;
438  const int y = bin / Nx + 1;
439  double val = hist->GetBinContent(x, y);
440  if (val < minChi) {
441  minChi = val;
442  best = i;
443  }
444  }
445  }
446 
447  ret->fBestLikelihood = surfs[best]->fBestLikelihood;
448  ret->fBestFitX = surfs[best]->fBestFitX;
449  ret->fBestFitY = surfs[best]->fBestFitY;
450 
451  // Loop over other surfaces and add them in
452  for (size_t i = 1; i < surfs.size(); ++i) {
453  ret->fHist->Add(surfs[i]->ToTH2(0));
454  for (size_t j = 0; j < ret->fProfHists.size(); ++j)
455  ret->fProfHists[j]->Add(surfs[i]->fProfHists[j]);
456  }
457 
458  // Scale hist by global minimum
459  for(int x = 0; x < ret->fHist->GetNbinsX()+2; ++x){
460  for(int y = 0; y < ret->fHist->GetNbinsY()+2; ++y){
461  ret->fHist->SetBinContent(x, y, ret->fHist->GetBinContent(x, y)-ret->fBestLikelihood);
462  }
463  }
464 
465  ret->fHist->SetMinimum(0);
466 
467  return ret;
468  }
469 
470  //----------------------------------------------------------------------
471  std::unique_ptr<FrequentistSurface> FrequentistSurface::LoadFromMulti(
473  {
474  std::vector<std::string> fileNames = Wildcard(wildcard);
475  std::vector<TFile*> files;
476  for (std::string fileName : fileNames) {
477  files.push_back(TFile::Open(fileName.c_str(), "read"));
478  }
479  std::unique_ptr<FrequentistSurface> ret = LoadFromMulti(files, label);
480  for(TFile* f: files) delete f;
481  return ret;
482  }
483 
484  // See eg the statistics section of the PDG
485  TH2* Gaussian68Percent2D(const FrequentistSurface& s){return Flat(2.30, s);}
486  TH2* Gaussian90Percent2D(const FrequentistSurface& s){return Flat(4.61, s);}
487  TH2* Gaussian95Percent2D(const FrequentistSurface& s){return Flat(5.99, s);}
488  TH2* Gaussian2Sigma2D (const FrequentistSurface& s){return Flat(6.18, s);}
489  TH2* Gaussian99Percent2D(const FrequentistSurface& s){return Flat(9.21, s);}
490  TH2* Gaussian3Sigma2D (const FrequentistSurface& s){return Flat(11.83, s);}
491 
492  TH2* Gaussian68Percent1D(const FrequentistSurface& s){return Flat(1.00, s);}
493  TH2* Gaussian90Percent1D(const FrequentistSurface& s){return Flat(2.71, s);}
494  TH2* Gaussian95Percent1D(const FrequentistSurface& s){return Flat(3.84, s);}
495  TH2* Gaussian2Sigma1D (const FrequentistSurface& s){return Flat(4.00, s);}
496  TH2* Gaussian99Percent1D(const FrequentistSurface& s){return Flat(6.63, s);}
497  TH2* Gaussian3Sigma1D (const FrequentistSurface& s){return Flat(9.00, s);}
498 
499 
500 } // namespace
const XML_Char * name
Definition: expat.h:151
void Finish()
Wait for all threads to complete before returning.
Definition: ThreadPool.cxx:50
const std::string & LatexName() const
Definition: IFitVar.h:38
virtual double ChiSq(osc::IOscCalcAdjustable *osc, const SystShifts &syst=SystShifts::Nominal()) const
Definition: IExperiment.h:18
TH2 * Gaussian90Percent1D(const FrequentistSurface &s)
Up-value surface for 90% confidence in 1D in gaussian approximation.
fileName
Definition: plotROC.py:78
size_t JobNumber()
What&#39;s the process number for a grid job?
Definition: Utilities.cxx:437
virtual _IOscCalcAdjustable< T > * Copy() const =0
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
set< int >::iterator it
MinuitFitter::FitOpts fFitOpts
double BinCenterX(int bin) const
Definition: ISurface.cxx:265
double FillSurfacePoint(const IExperiment *expt, osc::IOscCalcAdjustable *calc, const IFitVar *xvar, double x, const IFitVar *yvar, double y, const std::vector< const IFitVar * > &profVars, const std::vector< const ISyst * > &profSysts, const SeedList &seedPts, const std::vector< SystShifts > &systSeedPts)
bool RunningOnGrid()
Is this a grid (condor) job?
Definition: Utilities.cxx:358
std::function< void(void)> func_t
The type of the user&#39;s worker functions.
Definition: ThreadPool.h:33
Simple record of shifts applied to systematic parameters.
Definition: SystShifts.h:20
TH2 * Gaussian68Percent2D(const FrequentistSurface &s)
Up-value surface for 68% confidence in 2D in gaussian approximation.
const IFitVar * var
Definition: FitAxis.h:17
double BinCenterY(int bin) const
Definition: ISurface.cxx:272
OStream cerr
Definition: OStream.cxx:7
TH2 * Gaussian2Sigma2D(const FrequentistSurface &s)
Up-value surface for 2 sigma confidence in 2D in gaussian approximation.
Float_t tmp
Definition: plot.C:36
TH2 * Gaussian2Sigma1D(const FrequentistSurface &s)
Up-value surface for 2 sigma confidence in 1D in gaussian approximation.
static SystShifts Nominal()
Definition: SystShifts.h:34
osc::OscCalcDumb calc
std::string ProgressBarTitle(const IFitVar *xvar, const IFitVar *yvar, const std::vector< const IFitVar * > &profVars, const std::vector< const ISyst * > &profSysts) const
Encapsulate code to systematically shift a caf::SRProxy.
Definition: ISyst.h:14
int NThreads() const
Definition: ThreadPool.h:40
bool logy
TH2F * ExpandedHistogram(const std::string &title, int nbinsx, double xmin, double xmax, bool xlog, int nbinsy, double ymin, double ymax, bool ylog)
Internal helper for Surface and FCSurface.
Definition: Utilities.cxx:148
const char * label
std::vector< TH2 * > fProfHists
Log-likelihood scan across two parameters.
const XML_Char * s
Definition: expat.h:262
TH2 * Gaussian95Percent2D(const FrequentistSurface &s)
Up-value surface for 95% confidence in 2D in gaussian approximation.
const int nbins
Definition: cellShifts.C:15
expt
Definition: demo5.py:34
std::vector< std::string > Wildcard(const std::string &wildcardString)
Find files matching a UNIX glob, plus expand environment variables.
Definition: UtilsExt.cxx:268
TH2F * fHist
Definition: ISurface.h:68
double min
Definition: FitAxis.h:19
void FindMinimum(const IExperiment *expt, osc::IOscCalcAdjustable *calc, const IFitVar *xvar, const IFitVar *yvar, const std::vector< const IFitVar * > &profVars, const std::vector< const ISyst * > &profSysts, const SeedList &seedPts, const std::vector< SystShifts > &systSeedPts)
std::vector< std::string > wildcard(const std::string &wildcardString)
Definition: convert.C:9
virtual std::unique_ptr< IFitSummary > Fit(osc::IOscCalcAdjustable *seed, SystShifts &bestSysts=junkShifts, const SeedList &seedPts=SeedList(), const std::vector< SystShifts > &systSeedPts={}, Verbosity verb=kVerbose) const
Master fitting method. Depends on FitHelper and FitHelperSeeded.
Definition: IFitter.cxx:69
bool fLogX
Definition: ISurface.h:69
virtual double Penalty(double, osc::IOscCalcAdjustable *) const
Definition: IFitVar.h:35
TH2 * Gaussian3Sigma2D(const FrequentistSurface &s)
Up-value surface for 3 sigma confidence in 2D in gaussian approximation.
bool islog
Definition: FitAxis.h:21
int nbins
Definition: FitAxis.h:18
T GetShift(const ISyst *syst) const
TH2 * Gaussian99Percent2D(const FrequentistSurface &s)
Up-value surface for 99% confidence in 2D in gaussian approximation.
static std::unique_ptr< FrequentistSurface > LoadFromMulti(std::vector< TFile * > files, std::string label)
double max
Definition: FitAxis.h:20
double fBestFitY
Definition: ISurface.h:67
TH2 * Gaussian3Sigma1D(const FrequentistSurface &s)
Up-value surface for 3 sigma confidence in 1D in gaussian approximation.
const double j
Definition: BetheBloch.cxx:29
void SaveToHelper(TDirectory *dir) const
dir should already be the appropriate sub-directory
Definition: ISurface.cxx:198
float bin[41]
Definition: plottest35.C:14
void AddTask(const func_t &func)
Definition: ThreadPool.cxx:89
bool fLogY
Definition: ISurface.h:69
OStream cout
Definition: OStream.cxx:6
std::vector< double > fSeedValues
Definition: ISurface.h:70
bool logx
TH2 * ToTH2(double minchi=-1) const
Definition: ISurface.cxx:170
void SetProgress(double frac)
Update the progress fraction between zero and one.
Definition: Progress.cxx:41
const std::string & ShortName() const
Definition: IFitVar.h:37
TDirectory * dir
Definition: macro.C:5
double fBestFitX
Definition: ISurface.h:66
size_t NumJobs()
Definition: Utilities.cxx:448
virtual void SetValue(osc::IOscCalcAdjustable *osc, double val) const =0
TH2 * Gaussian90Percent2D(const FrequentistSurface &s)
Up-value surface for 90% confidence in 2D in gaussian approximation.
void SaveTo(TDirectory *dir, const std::string &name) const
Base class defining interface for experiments.
Definition: IExperiment.h:14
friend TH2 * Flat(double, const ISurface &)
Helper function for the gaussian approximation surfaces.
Definition: ISurface.cxx:253
std::vector< int > fBinMask
Definition: ISurface.h:71
A very simple thread pool for use by Surface.
Definition: ThreadPool.h:18
assert(nhit_max >=nhit_nbins)
A simple ascii-art progress bar.
Definition: Progress.h:9
void ShowProgress(const std::string &title)
Definition: ThreadPool.cxx:43
Interface definition for fittable variables.
Definition: IFitVar.h:17
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
Definition: Style.cxx:154
TH2 * Gaussian99Percent1D(const FrequentistSurface &s)
Up-value surface for 99% confidence in 1D in gaussian approximation.
virtual void Reset() const
Definition: IExperiment.h:32
static void FillSurfObj(ISurface &surf, TDirectory *dir)
Definition: ISurface.cxx:229
void CreateHistograms(const FitAxis &xax, const FitAxis &yax, const std::vector< const IFitVar * > &profVars, const std::vector< const ISyst * > &profSysts)
Prevent histograms being added to the current directory.
Definition: UtilsExt.h:46
Float_t e
Definition: plot.C:35
TH2 * Gaussian95Percent1D(const FrequentistSurface &s)
Up-value surface for 95% confidence in 1D in gaussian approximation.
T max(sqlite3 *const db, std::string const &table_name, std::string const &column_name)
Definition: statistics.h:66
void Done()
Call this when action is completed.
Definition: Progress.cxx:90
static std::unique_ptr< FrequentistSurface > LoadFrom(TDirectory *dir, const std::string &name)
double fBestLikelihood
Definition: ISurface.h:65
TH2 * Gaussian68Percent1D(const FrequentistSurface &s)
Up-value surface for 68% confidence in 1D in gaussian approximation.
Collect information describing the axis of a fit variable.
Definition: FitAxis.h:10
surf
Definition: demo5.py:35
void SetFitOpts(FitOpts opts)
virtual double GetValue(const osc::IOscCalcAdjustable *osc) const =0
virtual void FillSurface(const IExperiment *expt, osc::IOscCalcAdjustable *calc, const IFitVar *xvar, const IFitVar *yvar, const std::vector< const IFitVar * > &profVars, const std::vector< const ISyst * > &profSysts, const SeedList &seedPts, const std::vector< SystShifts > &systSeedPts)
Perform MINUIT fits in one or two dimensions.
Definition: MinuitFitter.h:17
enum BeamMode string