Public Member Functions | Static Public Member Functions | Protected Member Functions | Static Protected Member Functions | Protected Attributes | Friends | List of all members
ana::FrequentistSurface Class Reference

Log-likelihood scan across two parameters. More...

#include "/cvmfs/nova-development.opensciencegrid.org/novasoft/releases/N20-08-04/CAFAna/Fit/FrequentistSurface.h"

Inheritance diagram for ana::FrequentistSurface:
ana::ISurface ana::SurfaceKrige

Public Member Functions

 FrequentistSurface (const IExperiment *expt, osc::IOscCalculatorAdjustable *calc, const IFitVar *xvar, int nbinsx, double xmin, double xmax, const IFitVar *yvar, int nbinsy, double ymin, double ymax, const std::vector< const IFitVar * > &profVars={}, const std::vector< const ISyst * > &profSysts={}, const SeedList &seedPts=SeedList(), const std::vector< SystShifts > &systSeedPts={}, bool parallel=false, MinuitFitter::FitOpts opts=MinuitFitter::kNormal)
 
virtual ~FrequentistSurface ()
 
std::vector< TH2 * > GetProfiledHists ()
 Maps of the values taken on by the profiled parameters. More...
 
void SaveTo (TDirectory *dir) const
 
double BestLikelihood () const
 
double GetBestFitX () const
 
double GetBestFitY () const
 
void Draw () const
 Draw the surface itself. More...
 
void DrawBestFit (Color_t color, Int_t marker=kFullCircle) const
 Draw the best fit point. More...
 
void DrawContour (TH2 *fc, Style_t style, Color_t color, double minchi=-1)
 
std::vector< TGraph * > GetGraphs (TH2 *fc, double minchi=-1)
 For expert use, custom painting of contours. More...
 
TH2 * ToTH2 (double minchi=-1) const
 
void SetTitle (const char *str)
 

Static Public Member Functions

static std::unique_ptr< FrequentistSurfaceLoadFrom (TDirectory *dir)
 

Protected Member Functions

 FrequentistSurface ()
 
void CreateHistograms (const IFitVar *xvar, int nbinsx, double xmin, double xmax, const IFitVar *yvar, int nbinsy, double ymin, double ymax, const std::vector< const IFitVar * > &profVars, const std::vector< const ISyst * > &profSysts)
 
std::string ProgressBarTitle (const IFitVar *xvar, const IFitVar *yvar, const std::vector< const IFitVar * > &profVars, const std::vector< const ISyst * > &profSysts) const
 
virtual void FillSurface (const IExperiment *expt, osc::IOscCalculatorAdjustable *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)
 
double FillSurfacePoint (const IExperiment *expt, osc::IOscCalculatorAdjustable *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)
 
void FindMinimum (const IExperiment *expt, osc::IOscCalculatorAdjustable *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)
 
void EnsureAxes () const
 

Static Protected Member Functions

static void FillSurfObj (ISurface &surf, TDirectory *dir)
 

Protected Attributes

bool fParallel
 
MinuitFitter::FitOpts fFitOpts
 
std::vector< TH2 * > fProfHists
 
double fBestLikelihood
 
double fBestFitX
 
double fBestFitY
 
TH2F * fHist
 
std::vector< double > fSeedValues
 

Friends

class NumuSurface
 
class NueSurface
 

Detailed Description

Log-likelihood scan across two parameters.

Definition at line 26 of file FrequentistSurface.h.

Constructor & Destructor Documentation

ana::FrequentistSurface::FrequentistSurface ( const IExperiment expt,
osc::IOscCalculatorAdjustable calc,
const IFitVar xvar,
int  nbinsx,
double  xmin,
double  xmax,
const IFitVar yvar,
int  nbinsy,
double  ymin,
double  ymax,
const std::vector< const IFitVar * > &  profVars = {},
const std::vector< const ISyst * > &  profSysts = {},
const SeedList seedPts = SeedList(),
const std::vector< SystShifts > &  systSeedPts = {},
bool  parallel = false,
MinuitFitter::FitOpts  opts = MinuitFitter::kNormal 
)
Parameters
exptThe experiment object to draw $ \chi^2 $ values from
calcValues for oscillation parameters to be held fixed
xvarOscillation parameter to place on the x axis
nbinsxNumber of bins along x axis
xminMinimum value of x axis
xmaxMaximum value of x axis
nbinsyNumber of bins along y axis
yminMinimum value of y axis
ymaxMaximum value of y axis
profVarsOscillation parameters to profile over
profSystsSystematic parameters to profile over
seedPtsTry all combinations of these params as seeds
systSeedPtsTry all of these systematic combinations as seeds
parallelUse all the cores on this machine? Be careful...

Definition at line 31 of file FrequentistSurface.cxx.

References CreateHistograms(), FillSurface(), FindMinimum(), ana::ISurface::fSeedValues, and registry_explorer::v.

41  : fParallel(parallel), fFitOpts(opts)
42  {
44  yvar, nbinsy, ymin, ymax,
45  profVars, profSysts);
46 
47  for(const IFitVar* v: profVars) fSeedValues.push_back(v->GetValue( calc));
48 
49  FillSurface(expt, calc, xvar, yvar, profVars, profSysts, seedPts, systSeedPts);
50 
51  FindMinimum(expt, calc, xvar, yvar, profVars, profSysts, seedPts, systSeedPts);
52  }
void CreateHistograms(const IFitVar *xvar, int nbinsx, double xmin, double xmax, const IFitVar *yvar, int nbinsy, double ymin, double ymax, const std::vector< const IFitVar * > &profVars, const std::vector< const ISyst * > &profSysts)
virtual void FillSurface(const IExperiment *expt, osc::IOscCalculatorAdjustable *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::map< std::string, double > xmax
MinuitFitter::FitOpts fFitOpts
Double_t ymax
Definition: plot.C:25
expt
Definition: demo5.py:34
std::vector< double > fSeedValues
Definition: ISurface.h:63
Int_t nbinsx
Definition: plot.C:23
void FindMinimum(const IExperiment *expt, osc::IOscCalculatorAdjustable *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)
Double_t ymin
Definition: plot.C:24
ana::FrequentistSurface::~FrequentistSurface ( )
virtual

Definition at line 55 of file FrequentistSurface.cxx.

References CreateHistograms().

56  {
57  }
ana::FrequentistSurface::FrequentistSurface ( )
inlineprotected

Member Function Documentation

double ana::ISurface::BestLikelihood ( ) const
inlineinherited
void ana::FrequentistSurface::CreateHistograms ( const IFitVar xvar,
int  nbinsx,
double  xmin,
double  xmax,
const IFitVar yvar,
int  nbinsy,
double  ymin,
double  ymax,
const std::vector< const IFitVar * > &  profVars,
const std::vector< const ISyst * > &  profSysts 
)
protected

Definition at line 61 of file FrequentistSurface.cxx.

References ana::ExpandedHistogram(), ana::ISurface::fHist, fProfHists, MECModelEnuComparisons::i, ana::IFitVar::LatexName(), nbinsx, ProgressBarTitle(), plotROC::title, xmax, make_mec_shifts_plots::xmin, ymax, and ymin.

Referenced by FrequentistSurface(), and ~FrequentistSurface().

65  {
66  fHist = ExpandedHistogram(";"+xvar->LatexName()+";"+yvar->LatexName(),
67  nbinsx, xmin, xmax,
68  nbinsy, ymin, ymax);
69 
70  for(unsigned int i = 0; i < profVars.size()+profSysts.size(); ++i){
72  if(i < profVars.size())
73  title = profVars[i]->LatexName();
74  else
75  title = profSysts[i-profVars.size()]->LatexName();
76 
77  fProfHists.push_back(ExpandedHistogram(title+";"+xvar->LatexName()+";"+yvar->LatexName(),
78  nbinsx, xmin, xmax,
79  nbinsy, ymin, ymax));
80  }
81 
82  }
std::map< std::string, double > xmax
std::vector< TH2 * > fProfHists
Double_t ymax
Definition: plot.C:25
TH2F * fHist
Definition: ISurface.h:62
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
Int_t nbinsx
Definition: plot.C:23
Double_t ymin
Definition: plot.C:24
TH2F * ExpandedHistogram(const std::string &title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax)
Definition: Utilities.cxx:198
void ana::ISurface::Draw ( ) const
inherited

Draw the surface itself.

Definition at line 20 of file ISurface.cxx.

References ana::ISurface::EnsureAxes(), and ana::ISurface::fHist.

Referenced by cc(), demo5(), ana::NueSurface::Draw(), ana::NumuSurface::Draw(), ana::ISurface::GetBestFitY(), modularextrap_demo_nue(), modularextrap_demo_numu(), Plot2DSlices(), test_ana(), test_external_constraints(), test_stanfit_statsonly(), and test_stanfit_withsysts().

21  {
22  EnsureAxes();
23 
24  fHist->Draw("colz same");
25 
26  // colz obliterated them
27  gPad->RedrawAxis();
28 
29  gPad->Update();
30  }
TH2F * fHist
Definition: ISurface.h:62
void EnsureAxes() const
Definition: ISurface.cxx:64
void ana::ISurface::DrawBestFit ( Color_t  color,
Int_t  marker = kFullCircle 
) const
inherited
void ana::ISurface::DrawContour ( TH2 *  fc,
Style_t  style,
Color_t  color,
double  minchi = -1 
)
inherited
Parameters
fcSurface to compare against for this significance level
styleLine style for TAttLine, solid, dotted, dashed etc
colorLine color for TAttLine
minchi$\chi^2$ of best fit to compare against. Default: best fit from this surface.

Definition at line 45 of file ISurface.cxx.

References ana::ISurface::EnsureAxes(), MECModelEnuComparisons::g, ana::ISurface::GetGraphs(), and APDGainPoints::gs.

Referenced by cc(), demo5(), Draw2DSurface(), ana::NueSurface::DrawContour(), ana::NumuSurface::DrawContour(), drawSensitivity(), ana::ISurface::GetBestFitY(), getSensitivity(), make_contours(), make_nom_expt(), modularextrap_demo_nue(), modularextrap_demo_numu(), Plot2DSlice(), Plot2DSlices(), plotContProf(), plots(), starPlot(), syst_test(), template_basic(), template_GENIE_systs(), template_nonGENIE_systs(), test_ana(), test_external_constraints(), test_stanfit_statsonly(), and test_stanfit_withsysts().

47  {
48  EnsureAxes();
49 
50  std::vector<TGraph *> gs = GetGraphs(fc, minchi);
51 
52  for (TGraph *g: gs)
53  {
54  g->SetLineWidth(3);//2);
55  g->SetLineStyle(style);
56  g->SetLineColor(color);
57  g->Draw("l");
58  }
59 
60  gPad->Update();
61  }
std::vector< TGraph * > GetGraphs(TH2 *fc, double minchi=-1)
For expert use, custom painting of contours.
Definition: ISurface.cxx:115
void EnsureAxes() const
Definition: ISurface.cxx:64
void ana::ISurface::EnsureAxes ( ) const
protectedinherited

Definition at line 64 of file ISurface.cxx.

References visualisationForPaperMasterPlot::ax, file_size_ana::axes, ana::ISurface::fHist, genie::utils::style::Format(), and ana::UniqueName().

Referenced by ana::ISurface::Draw(), ana::ISurface::DrawBestFit(), ana::ISurface::DrawContour(), and ana::NueSurface::DrawValley().

65  {
66  // Could have a file temporarily open
67  DontAddDirectory guard;
68 
69  // If this pad has already been drawn in, already has axes
70  if (gPad && !gPad->GetListOfPrimitives()->IsEmpty()) return;
71 
72  // Old, hackier solution
73  /*
74  std::cout << gPad->GetListOfPrimitives()->GetEntries() << std::endl;
75  // Which pads have we already drawn axes in? Never draw axes in them
76  // again. Unfortunately UniqueID() never seems to be set. If that's the
77  // case, set it to a random value and hope...
78  static std::set<UInt_t> already;
79  if(already.count(gPad->GetUniqueID())) return;
80  if(gPad->GetUniqueID() == 0) gPad->SetUniqueID(rand());
81  already.insert(gPad->GetUniqueID());
82  */
83 
84  const TAxis *ax = fHist->GetXaxis();
85  const TAxis *ay = fHist->GetYaxis();
86  const double Nx = ax->GetNbins();
87  const double Ny = ay->GetNbins();
88 
89  // Axes with limits where the user originally requested them, which we
90  // adjusted to be the centres of the first and last bins.
91  TH2 *axes = new TH2C(UniqueName().c_str(),
92  TString::Format(";%s;%s",
93  ax->GetTitle(), ay->GetTitle()),
94  Nx - 1, ax->GetBinCenter(1), ax->GetBinCenter(Nx),
95  Ny - 1, ay->GetBinCenter(1), ay->GetBinCenter(Ny));
96  axes->Draw();
97 
98  if(fHist){
99  // "colz same" will reuse axis's min and max, so set them helpfully here
100  axes->SetMinimum(fHist->GetMinimum());
101  axes->SetMaximum(fHist->GetMaximum());
102  }
103 
104  axes->SetTitle(fHist->GetTitle());
105  axes->GetXaxis()->SetLabelSize(ax->GetLabelSize());
106  axes->GetYaxis()->SetLabelSize(ay->GetLabelSize());
107  axes->GetXaxis()->SetLabelOffset(ax->GetLabelOffset());
108  axes->GetYaxis()->SetLabelOffset(ay->GetLabelOffset());
109  axes->GetXaxis()->CenterTitle();
110  axes->GetYaxis()->CenterTitle();
111  gPad->Update();
112  }
TH2F * fHist
Definition: ISurface.h:62
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
Definition: Style.cxx:154
std::string UniqueName()
Return a different string each time, for creating histograms.
Definition: Utilities.cxx:30
void ana::FrequentistSurface::FillSurface ( const IExperiment expt,
osc::IOscCalculatorAdjustable 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 
)
protectedvirtual

Reimplemented in ana::SurfaceKrige.

Definition at line 107 of file FrequentistSurface.cxx.

References ana::ThreadPool::AddMemberTask(), bin, om::cerr, ana::Progress::Done(), e, allTimeWatchdog::endl, ana::ISurface::fHist, FillSurfacePoint(), ana::ThreadPool::Finish(), fParallel, Get, ana::IFitVar::Penalty(), time_estimates::pool, cacheDefinitionData::prog, ProgressBarTitle(), ana::Progress::SetProgress(), ana::IFitVar::ShortName(), ana::ThreadPool::ShowProgress(), fillBadChanDBTables::step, submit_syst::x, and submit_syst::y.

Referenced by FrequentistSurface().

114  {
115  if (fParallel && !(profVars.empty() && profSysts.empty()))
116  {
117  // Minuit calls this at some point, and it creates a static. Which
118  // doesn't like happening on multiple threads at once. So do it upfront
119  // while we're still single-threaded.
121  }
122 
123  // Nothing created during surface filling belongs in a
124  // directory. Unfortunately the local guards in Spectrum etc are racey when
125  // run in parallel. But this should cover the whole lot safely.
126  DontAddDirectory guard;
127 
128  const std::string progTitle = ProgressBarTitle(xvar, yvar, profVars, profSysts);
129 
130  Progress *prog = 0;
131  // Difficult to keep a progress bar properly up to date when threaded
132  if (!fParallel) prog = new Progress(progTitle);
133  ThreadPool *pool = 0;
134  if (fParallel)
135  {
136  pool = new ThreadPool;
137  pool->ShowProgress(progTitle);
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  do
155  {
156  const int x = bin % Nx + 1;
157  const int y = bin / Nx + 1;
158 
159  const double xv = fHist->GetXaxis()->GetBinCenter(x);
160  const double yv = fHist->GetYaxis()->GetBinCenter(y);
161 
162  if (xvar->Penalty(xv, calc) > 1e-10)
163  {
164  std::cerr << "Warning! " << xvar->ShortName() << " = " << xv
165  << " has penalty of " << xvar->Penalty(xv, calc)
166  << " that could have been applied in surface. "
167  << "This should never happen." << std::endl;
168  }
169  if (yvar->Penalty(yv, calc) > 1e-10)
170  {
171  std::cerr << "Warning! " << yvar->ShortName() << " = " << yv
172  << " has penalty of " << yvar->Penalty(yv, calc)
173  << " that could have been applied in surface. "
174  << "This should never happen." << std::endl;
175  }
176 
177  if (fParallel)
178  {
179  pool->AddMemberTask(this, &FrequentistSurface::FillSurfacePoint,
180  expt, calc,
181  xvar, xv, yvar, yv,
182  profVars, profSysts, seedPts, systSeedPts);
183  } else
184  {
185  FillSurfacePoint(expt, calc,
186  xvar, xv, yvar, yv,
187  profVars, profSysts, seedPts, systSeedPts);
188  ++neval;
189  prog->SetProgress(neval / double(Nx * Ny));
190  }
191 
192  bin = (bin + step) % (Nx * Ny);
193  } while (bin != 0);
194 
195 
196  if (fParallel)
197  {
198  pool->Finish();
199  delete pool;
200  } else
201  {
202  prog->Done();
203  delete prog;
204  }
205  }
double FillSurfacePoint(const IExperiment *expt, osc::IOscCalculatorAdjustable *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)
OStream cerr
Definition: OStream.cxx:7
cout<< "Opened file "<< fin<< " ixs= "<< ixs<< endl;if(ixs==0) hhh=(TH1F *) fff-> Get("h1")
Definition: AddMC.C:8
std::string ProgressBarTitle(const IFitVar *xvar, const IFitVar *yvar, const std::vector< const IFitVar * > &profVars, const std::vector< const ISyst * > &profSysts) const
expt
Definition: demo5.py:34
TH2F * fHist
Definition: ISurface.h:62
float bin[41]
Definition: plottest35.C:14
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
Float_t e
Definition: plot.C:35
double ana::FrequentistSurface::FillSurfacePoint ( const IExperiment expt,
osc::IOscCalculatorAdjustable 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 
)
protected

Definition at line 245 of file FrequentistSurface.cxx.

References calc, ana::IExperiment::ChiSq(), osc::_IOscCalculatorAdjustable< T >::Copy(), fFitOpts, ana::ISurface::fHist, ana::IFitter::Fit(), fParallel, fProfHists, ana::ISurface::fSeedValues, ana::SystShifts::GetShift(), MECModelEnuComparisons::i, makeTrainCVSamples::int, calib::j, ana::IFitter::kQuiet, ana::MinuitFitter::SetFitOpts(), and ana::IFitVar::SetValue().

Referenced by FillSurface(), and FrequentistSurface().

253  {
254  osc::IOscCalculatorAdjustable* calcNoHash = 0; // specific to parallel mode
255 
256  if(fParallel){
257  // Need to take our own copy so that we don't get overwritten by someone
258  // else's changes.
259  calc = calc->Copy();
260  }
261 
262  xvar->SetValue(calc, x);
263  yvar->SetValue(calc, y);
264 
265  if (fParallel){
266  calcNoHash = new OscCalcNoHash(calc);
267  }
268 
269  //Make sure that the profiled values of fitvars do not persist between steps.
270  for(int i = 0; i < (int)fSeedValues.size(); ++i) profVars[i]->SetValue( calc, fSeedValues[i] );
271 
272  double chi;
273  if(profVars.empty() && profSysts.empty()){
274  chi = expt->ChiSq(fParallel ? calcNoHash : calc);
275  }
276  else{
277  MinuitFitter fitter(expt, profVars, profSysts);
278  fitter.SetFitOpts(fFitOpts);
279  SystShifts bestSysts;
280  chi = fitter.Fit(calc, bestSysts, seedPts, systSeedPts, MinuitFitter::kQuiet);
281 
282  for(unsigned int i = 0; i < profVars.size(); ++i){
283  fProfHists[i]->Fill(x, y, profVars[i]->GetValue(calc));
284  }
285  for(unsigned int j = 0; j < profSysts.size(); ++j){
286  fProfHists[j+profVars.size()]->Fill(x, y, bestSysts.GetShift(profSysts[j]));
287  }
288  }
289 
290  fHist->Fill(x, y, chi);
291 
292  if(fParallel){
293  delete calc;
294  delete calcNoHash;
295  }
296  return chi;
297  }
MinuitFitter::FitOpts fFitOpts
std::vector< TH2 * > fProfHists
expt
Definition: demo5.py:34
TH2F * fHist
Definition: ISurface.h:62
const double j
Definition: BetheBloch.cxx:29
osc::OscCalculatorDumb calc
std::vector< double > fSeedValues
Definition: ISurface.h:63
virtual _IOscCalculatorAdjustable< T > * Copy() const =0
void ana::ISurface::FillSurfObj ( ISurface surf,
TDirectory *  dir 
)
staticprotectedinherited

Definition at line 219 of file ISurface.cxx.

References ana::ISurface::fBestFitX, ana::ISurface::fBestFitY, ana::ISurface::fBestLikelihood, ana::ISurface::fHist, ana::ISurface::fSeedValues, compare_h5_caf::idx, and registry_explorer::v.

Referenced by ana::ISurface::GetBestFitY(), LoadFrom(), and ana::BayesianSurface::LoadFrom().

220  {
221 
222  const TVectorD v = *(TVectorD *) dir->Get("minValues");
223  const TVectorD s = *(TVectorD *) dir->Get("seeds");
224 
225  surf.fHist = (TH2F *) dir->Get("hist");
226  surf.fBestLikelihood = v[0];
227  surf.fBestFitX = v[1];
228  surf.fBestFitY = v[2];
229 
230  for (int idx = 0; idx < s.GetNrows(); ++idx)
231  surf.fSeedValues.push_back(s[idx]);
232 
233  }
const XML_Char * s
Definition: expat.h:262
TDirectory * dir
Definition: macro.C:5
surf
Definition: demo5.py:35
void ana::FrequentistSurface::FindMinimum ( const IExperiment expt,
osc::IOscCalculatorAdjustable 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 
)
protected

Definition at line 301 of file FrequentistSurface.cxx.

References ana::ISurface::fBestFitX, ana::ISurface::fBestFitY, ana::ISurface::fBestLikelihood, fFitOpts, ana::ISurface::fHist, ana::IFitter::Fit(), ana::ISurface::fSeedValues, ana::IFitVar::GetValue(), MECModelEnuComparisons::i, makeTrainCVSamples::int, ana::SystShifts::Nominal(), ana::MinuitFitter::SetFitOpts(), ana::IFitVar::SetValue(), submit_syst::x, and submit_syst::y.

Referenced by FrequentistSurface().

308  {
309  // Location of the best minimum found from filled surface
310  double minchi = 1e10;
311  int minx = fHist->GetNbinsX()/2;
312  int miny = fHist->GetNbinsY()/2;
313  for(int x = 1; x <= fHist->GetNbinsX(); ++x){
314  for(int y = 1; y <= fHist->GetNbinsY(); ++y){
315  const double chi = fHist->GetBinContent(x, y);
316  if(chi < minchi){
317  minchi = chi;
318  minx = x;
319  miny = y;
320  }
321  }
322  }
323 
324  std::vector<const IFitVar*> allVars = {xvar, yvar};
325  allVars.insert(allVars.end(), profVars.begin(), profVars.end());
326  MinuitFitter fit(expt, allVars, profSysts);
327  fit.SetFitOpts(fFitOpts);
328  // Seed from best grid point
329  xvar->SetValue(calc, fHist->GetXaxis()->GetBinCenter(minx));
330  yvar->SetValue(calc, fHist->GetYaxis()->GetBinCenter(miny));
331  for(int i = 0; i < (int)fSeedValues.size(); ++i) profVars[i]->SetValue( calc, fSeedValues[i] );
332  SystShifts systSeed = SystShifts::Nominal();
333  fBestLikelihood = fit.Fit(calc, systSeed, seedPts);
334  fBestFitX = xvar->GetValue(calc);
335  fBestFitY = yvar->GetValue(calc);
336 
337  for(int x = 0; x < fHist->GetNbinsX()+2; ++x){
338  for(int y = 0; y < fHist->GetNbinsY()+2; ++y){
339  fHist->SetBinContent(x, y, fHist->GetBinContent(x, y)-fBestLikelihood);
340  }
341  }
342 
343  fHist->SetMinimum(0);
344  }
MinuitFitter::FitOpts fFitOpts
static SystShifts Nominal()
Definition: SystShifts.h:34
expt
Definition: demo5.py:34
TH2F * fHist
Definition: ISurface.h:62
double fBestFitY
Definition: ISurface.h:61
std::vector< double > fSeedValues
Definition: ISurface.h:63
double fBestFitX
Definition: ISurface.h:60
double fBestLikelihood
Definition: ISurface.h:59
double ana::ISurface::GetBestFitX ( ) const
inlineinherited

Definition at line 30 of file ISurface.h.

References ana::ISurface::fBestFitX.

Referenced by test_stanfit_statsonly().

30 {return fBestFitX;}
double fBestFitX
Definition: ISurface.h:60
double ana::ISurface::GetBestFitY ( ) const
inlineinherited
std::vector< TGraph * > ana::ISurface::GetGraphs ( TH2 *  fc,
double  minchi = -1 
)
inherited

For expert use, custom painting of contours.

Definition at line 115 of file ISurface.cxx.

References ana::ISurface::fBestLikelihood, ana::ISurface::fHist, runNovaSAM::ret, PandAna.Demos.tute_pid_validation::specs, demo5::surf, tmp, and ana::UniqueName().

Referenced by ana::ISurface::DrawContour(), and ana::ISurface::GetBestFitY().

116  {
117  std::vector<TGraph*> ret;
118 
119  if(minchi < 0) minchi = fBestLikelihood;
120  std::unique_ptr<TH2> surf((TH2*)fHist->Clone(UniqueName().c_str()));
121  surf->Add(fc, -1);
122 
123  TVirtualPad* bak = gPad;
124 
125  const bool wasbatch = gROOT->IsBatch();
126  gROOT->SetBatch(); // User doesn't want to see our temporary canvas
127  TCanvas tmp;
128 
129  gStyle->SetOptStat(0);
130 
131  const double level = minchi-fBestLikelihood;
132  surf->SetContour(1, &level);
133  surf->Draw("cont list");
134 
135  tmp.Update();
136  tmp.Paint();
137 
138  gROOT->SetBatch(wasbatch);
139  gPad = bak;
140 
141  // The graphs we need (contained inside TLists, contained inside
142  // TObjArrays) are in the list of specials. But we need to be careful about
143  // types, because other stuff can get in here too (TDatabasePDG for
144  // example).
145  TCollection* specs = gROOT->GetListOfSpecials();
146 
147  TIter nextSpec(specs);
148  while(TObject* spec = nextSpec()){
149  if(!spec->InheritsFrom(TObjArray::Class())) continue;
150  TObjArray* conts = (TObjArray*)spec;
151 
152  if(conts->IsEmpty()) continue;
153 
154  if(!conts->At(0)->InheritsFrom(TList::Class())) continue;
155  TList* cont = (TList*)conts->At(0);
156 
157  TIter nextObj(cont);
158  // Contour could be split into multiple pieces
159  while(TObject* obj = nextObj()){
160  if(!obj->InheritsFrom(TGraph::Class())) continue;
161 
162  ret.push_back((TGraph*)obj->Clone(UniqueName().c_str()));
163  } // end for obj
164  } // end for spec
165 
166  return ret;
167  }
Float_t tmp
Definition: plot.C:36
TH2F * fHist
Definition: ISurface.h:62
double fBestLikelihood
Definition: ISurface.h:59
std::string UniqueName()
Return a different string each time, for creating histograms.
Definition: Utilities.cxx:30
surf
Definition: demo5.py:35
std::vector<TH2*> ana::FrequentistSurface::GetProfiledHists ( )
inline

Maps of the values taken on by the profiled parameters.

Definition at line 60 of file FrequentistSurface.h.

References dir, fProfHists, LoadFrom(), and SaveTo().

Referenced by run_joint_fit_2020_contours().

60 {return fProfHists;}
std::vector< TH2 * > fProfHists
std::unique_ptr< FrequentistSurface > ana::FrequentistSurface::LoadFrom ( TDirectory *  dir)
static

Definition at line 367 of file FrequentistSurface.cxx.

References ana::assert(), ana::ISurface::FillSurfObj(), genie::utils::style::Format(), make_syst_table_plots::h, compare_h5_caf::idx, demo5::surf, and getGoodRuns4SAM::tag.

Referenced by compare_fits(), demoFitContours(), DrawExtrapSurface(), GetProfiledHists(), joint_fit_2018_contours(), joint_fit_2019_contours(), median_contours(), plot_joint_fit_2020_contours(), Plotting_CompareMultipleContours(), sensitivity2018(), and sensitivity2020().

368  {
369  DontAddDirectory guard;
370 
371  TObjString *tag = (TObjString *) dir->Get("type");
372  assert(tag);
373  assert(tag->GetString() == "FrequentistSurface" || tag->GetString() == "Surface");
374 
375  std::unique_ptr<FrequentistSurface> surf(new FrequentistSurface);
377 
378  for (std::size_t idx = 0; idx < surf->fSeedValues.size(); ++idx)
379  {
380  // Search for old "marg" name here too for backwards compatibility
381  TH2 *h = (TH2 *) dir->Get(TString::Format("profHists/hist%lu", idx));
382  if (h)
383  surf->fProfHists.push_back(h);
384  else
385  surf->fProfHists.push_back((TH2 *) dir->Get(TString::Format("margHists/hist%lu", idx)));
386  }
387 
388  return surf;
389  }
TDirectory * dir
Definition: macro.C:5
assert(nhit_max >=nhit_nbins)
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
Definition: Style.cxx:154
static void FillSurfObj(ISurface &surf, TDirectory *dir)
Definition: ISurface.cxx:219
surf
Definition: demo5.py:35
std::string ana::FrequentistSurface::ProgressBarTitle ( const IFitVar xvar,
const IFitVar yvar,
const std::vector< const IFitVar * > &  profVars,
const std::vector< const ISyst * > &  profSysts 
) const
protected

Definition at line 86 of file FrequentistSurface.cxx.

References ana::IFitVar::ShortName(), plotROC::title, and registry_explorer::v.

Referenced by CreateHistograms(), FillSurface(), and FrequentistSurface().

89  {
90  std::string title = "Filling surface for "+yvar->ShortName()+" vs "+xvar->ShortName();
91 
92  if(!profVars.empty() || !profSysts.empty()){
93  title += " (profiling ";
94 
95  for(const IFitVar* v: profVars) title += v->ShortName() + ", ";
96  for(const ISyst* s: profSysts) title += s->ShortName() + ", ";
97 
98  // Always have one superfluous ", " at the end
99  title.resize(title.size()-2);
100  title += ")";
101  }
102 
103  return title;
104  }
const XML_Char * s
Definition: expat.h:262
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
void ana::FrequentistSurface::SaveTo ( TDirectory *  dir) const

Definition at line 347 of file FrequentistSurface.cxx.

References genie::utils::style::Format(), fProfHists, compare_h5_caf::idx, it, ana::ISurface::SaveTo(), and tmp.

Referenced by getContProf(), getContProf_Sensitivity(), GetProfiledHists(), joint_fit_future_contour_univ(), make_surfprof(), make_surfprof_sensitivity(), MakeExtrapSurface(), run_joint_fit_2020_contours(), sensitivity2018(), and sensitivity2020().

348  {
349  TDirectory *tmp = gDirectory;
350  dir->cd();
351  TObjString("FrequentistSurface").Write("type");
352 
354 
355  TDirectory *profDir = dir->mkdir("profHists");
356  int idx = 0;
358  {
359  profDir->cd();
360  it->Write(TString::Format("hist%d", idx++));
361  }
362 
363  tmp->cd();
364  }
set< int >::iterator it
Float_t tmp
Definition: plot.C:36
std::vector< TH2 * > fProfHists
void SaveTo(TDirectory *dir) const
Definition: ISurface.cxx:198
TDirectory * dir
Definition: macro.C:5
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
Definition: Style.cxx:154
void ana::ISurface::SetTitle ( const char *  str)
inherited
TH2 * ana::ISurface::ToTH2 ( double  minchi = -1) const
inherited

Definition at line 170 of file ISurface.cxx.

References ana::ISurface::fBestLikelihood, ana::ISurface::fHist, runNovaSAM::ret, submit_syst::x, and submit_syst::y.

Referenced by CAF_makeCAFSensitivities_for_FNEX(), demoFitContours(), ana::Flat(), ana::ISurface::GetBestFitY(), joint_fit_2017_contours(), joint_fit_2018_contours(), joint_fit_2019_contours(), joint_fit_future_contour_univ(), MakeCAFSensitivities_for_FNEX(), Plot2DSlice(), Plot2DSlices(), run_joint_fit_2020_contours(), sensitivity2018(), sensitivity2020(), and test_stanfit_withsysts().

171  {
172  // Could have a file temporarily open
173  DontAddDirectory guard;
174 
175  TH2 *ret = new TH2F(*fHist);
176 
177  if (minchi >= 0)
178  {
179  for (int x = 0; x < ret->GetNbinsX() + 2; ++x)
180  {
181  for (int y = 0; y < ret->GetNbinsY() + 2; ++y)
182  {
183  ret->SetBinContent(x, y, ret->GetBinContent(x, y) + fBestLikelihood - minchi);
184  }
185  }
186  }
187 
188  return ret;
189  }
TH2F * fHist
Definition: ISurface.h:62
double fBestLikelihood
Definition: ISurface.h:59

Friends And Related Function Documentation

friend class NueSurface
friend

Definition at line 30 of file FrequentistSurface.h.

friend class NumuSurface
friend

Definition at line 29 of file FrequentistSurface.h.

Member Data Documentation

double ana::ISurface::fBestFitX
protectedinherited
double ana::ISurface::fBestFitY
protectedinherited
double ana::ISurface::fBestLikelihood
protectedinherited
MinuitFitter::FitOpts ana::FrequentistSurface::fFitOpts
protected

Definition at line 104 of file FrequentistSurface.h.

Referenced by FillSurfacePoint(), and FindMinimum().

TH2F* ana::ISurface::fHist
protectedinherited
bool ana::FrequentistSurface::fParallel
protected

Definition at line 102 of file FrequentistSurface.h.

Referenced by FillSurface(), and FillSurfacePoint().

std::vector<TH2*> ana::FrequentistSurface::fProfHists
protected

Definition at line 107 of file FrequentistSurface.h.

Referenced by CreateHistograms(), FillSurfacePoint(), GetProfiledHists(), and SaveTo().

std::vector<double> ana::ISurface::fSeedValues
protectedinherited

The documentation for this class was generated from the following files: