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