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

#include "/cvmfs/nova-development.opensciencegrid.org/novasoft/releases/N21-01-15/CAFAna/Fit/ISurface.h"

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

Public Member Functions

 ISurface ()
 
virtual ~ISurface ()
 
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)
 

Protected Member Functions

double BinCenterX (int bin) const
 
double BinCenterY (int bin) const
 
void SaveToHelper (TDirectory *dir) const
 dir should already be the appropriate sub-directory More...
 
void EnsureAxes () const
 
void CheckMask (const std::string &func) const
 

Static Protected Member Functions

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

Protected Attributes

double fBestLikelihood
 
double fBestFitX
 
double fBestFitY
 
TH2F * fHist
 
bool fLogX
 
bool fLogY
 
std::vector< double > fSeedValues
 
std::vector< intfBinMask
 

Friends

TH2 * Flat (double, const ISurface &)
 Helper function for the gaussian approximation surfaces. More...
 

Detailed Description

Definition at line 17 of file ISurface.h.

Constructor & Destructor Documentation

ana::ISurface::ISurface ( )
inline

Definition at line 20 of file ISurface.h.

21  : fBestLikelihood(std::numeric_limits<double>::signaling_NaN()),
22  fBestFitX(std::numeric_limits<double>::signaling_NaN()),
23  fBestFitY(std::numeric_limits<double>::signaling_NaN()),
24  fHist(nullptr), fLogX(false), fLogY(false)
25  {}
TH2F * fHist
Definition: ISurface.h:68
bool fLogX
Definition: ISurface.h:69
double fBestFitY
Definition: ISurface.h:67
bool fLogY
Definition: ISurface.h:69
double fBestFitX
Definition: ISurface.h:66
double fBestLikelihood
Definition: ISurface.h:65
virtual ana::ISurface::~ISurface ( )
inlinevirtual

Definition at line 27 of file ISurface.h.

27 {}

Member Function Documentation

double ana::ISurface::BestLikelihood ( ) const
inline

Definition at line 29 of file ISurface.h.

References fBestLikelihood.

29 {return fBestLikelihood;}
double fBestLikelihood
Definition: ISurface.h:65
double ana::ISurface::BinCenterX ( int  bin) const
protected

Definition at line 265 of file ISurface.cxx.

References visualisationForPaperMasterPlot::ax, fHist, and fLogX.

Referenced by ana::BayesianSurface::BuildHist(), EnsureAxes(), ana::FrequentistSurface::FillSurface(), ana::FrequentistSurface::FindMinimum(), and GetBestFitY().

266  {
267  const TAxis* ax = fHist->GetXaxis();
268  return fLogX ? ax->GetBinCenterLog(bin) : ax->GetBinCenter(bin);
269  }
TH2F * fHist
Definition: ISurface.h:68
bool fLogX
Definition: ISurface.h:69
float bin[41]
Definition: plottest35.C:14
double ana::ISurface::BinCenterY ( int  bin) const
protected

Definition at line 272 of file ISurface.cxx.

References visualisationForPaperMasterPlot::ax, fHist, and fLogY.

Referenced by ana::BayesianSurface::BuildHist(), EnsureAxes(), ana::FrequentistSurface::FillSurface(), ana::FrequentistSurface::FindMinimum(), and GetBestFitY().

273  {
274  const TAxis* ax = fHist->GetYaxis();
275  return fLogY ? ax->GetBinCenterLog(bin) : ax->GetBinCenter(bin);
276  }
TH2F * fHist
Definition: ISurface.h:68
float bin[41]
Definition: plottest35.C:14
bool fLogY
Definition: ISurface.h:69
void ana::ISurface::CheckMask ( const std::string func) const
protected

Definition at line 279 of file ISurface.cxx.

References om::cout, allTimeWatchdog::endl, and fBinMask.

Referenced by DrawBestFit(), DrawContour(), and GetGraphs().

280  {
281  if(!fBinMask.empty()) {
282  std::cout << "Cannot call " << func << "() on a partial surface!" << std::endl;
283  abort();
284  }
285  }
double func(double x, double y)
OStream cout
Definition: OStream.cxx:6
std::vector< int > fBinMask
Definition: ISurface.h:71
void ana::ISurface::Draw ( ) const

Draw the surface itself.

Definition at line 21 of file ISurface.cxx.

References EnsureAxes(), and fHist.

Referenced by cc(), demo5(), GetBestFitY(), modularextrap_demo_nue(), modularextrap_demo_numu(), Plot2DSlices(), plot_3flavor_withsysts(), test_ana(), test_stanfit_statsonly(), test_stanfit_withsysts(), and test_surf_stride().

22  {
23  // Can be useful to draw a partial surface for debugging
24  // CheckMask();
25 
26  EnsureAxes();
27 
28  fHist->Draw("colz same");
29 
30  // colz obliterated them
31  gPad->RedrawAxis();
32 
33  gPad->Update();
34  }
TH2F * fHist
Definition: ISurface.h:68
void EnsureAxes() const
Definition: ISurface.cxx:70
void ana::ISurface::DrawBestFit ( Color_t  color,
Int_t  marker = kFullCircle 
) const

Draw the best fit point.

Definition at line 37 of file ISurface.cxx.

References CheckMask(), EnsureAxes(), fBestFitX, and fBestFitY.

Referenced by cc(), demo5(), demo_CPT(), drawSensitivity(), GetBestFitY(), getSensitivity(), make_contours(), modularextrap_demo_nue(), modularextrap_demo_numu(), Plot2DSlices(), plot_3flavor_withsysts(), starPlot(), syst_test(), template_basic(), template_GENIE_systs(), template_nonGENIE_systs(), test_ana(), test_stanfit_statsonly(), and test_stanfit_withsysts().

38  {
39  CheckMask("DrawBestFit");
40  EnsureAxes();
41 
42  TMarker *mark = new TMarker(fBestFitX, fBestFitY, marker);
43  mark->SetMarkerSize(1.5);
44  mark->SetMarkerColor(color);
45  mark->Draw();
46  gPad->Update();
47  }
void CheckMask(const std::string &func) const
Definition: ISurface.cxx:279
double fBestFitY
Definition: ISurface.h:67
void EnsureAxes() const
Definition: ISurface.cxx:70
double fBestFitX
Definition: ISurface.h:66
void ana::ISurface::DrawContour ( TH2 *  fc,
Style_t  style,
Color_t  color,
double  minchi = -1 
)
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 50 of file ISurface.cxx.

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

Referenced by cc(), demo5(), Draw2DSurface(), DrawContours(), drawSensitivity(), GetBestFitY(), getSensitivity(), make_contours(), make_nom_expt(), modularextrap_demo_nue(), modularextrap_demo_numu(), Plot2DSlice(), Plot2DSlices(), plot_3flavor_withsysts(), plotContProf(), plots(), starPlot(), syst_test(), template_basic(), template_GENIE_systs(), template_nonGENIE_systs(), test_ana(), test_stanfit_statsonly(), and test_stanfit_withsysts().

52  {
53  CheckMask("DrawContour");
54  EnsureAxes();
55 
56  std::vector<TGraph *> gs = GetGraphs(fc, minchi);
57 
58  for (TGraph *g: gs)
59  {
60  g->SetLineWidth(3);//2);
61  g->SetLineStyle(style);
62  g->SetLineColor(color);
63  g->Draw("l");
64  }
65 
66  gPad->Update();
67  }
void CheckMask(const std::string &func) const
Definition: ISurface.cxx:279
std::vector< TGraph * > GetGraphs(TH2 *fc, double minchi=-1)
For expert use, custom painting of contours.
Definition: ISurface.cxx:113
void EnsureAxes() const
Definition: ISurface.cxx:70
void ana::ISurface::EnsureAxes ( ) const
protected

Definition at line 70 of file ISurface.cxx.

References visualisationForPaperMasterPlot::ax, file_size_ana::axes, BinCenterX(), BinCenterY(), fHist, fLogX, fLogY, genie::utils::style::Format(), and ana::UniqueName().

Referenced by Draw(), DrawBestFit(), and DrawContour().

71  {
72  // Could have a file temporarily open
73  DontAddDirectory guard;
74 
75  // If this pad has already been drawn in, already has axes
76  if (gPad && !gPad->GetListOfPrimitives()->IsEmpty()) return;
77 
78  const TAxis *ax = fHist->GetXaxis();
79  const TAxis *ay = fHist->GetYaxis();
80  const double Nx = ax->GetNbins();
81  const double Ny = ay->GetNbins();
82 
83  // Axes with limits where the user originally requested them, which we
84  // adjusted to be the centres of the first and last bins.
85  TH2 *axes = new TH2C(UniqueName().c_str(),
86  TString::Format(";%s;%s",
87  ax->GetTitle(), ay->GetTitle()),
88  Nx - 1, BinCenterX(1), BinCenterX(Nx),
89  Ny - 1, BinCenterY(1), BinCenterY(Ny));
90  axes->Draw();
91 
92  if(fHist){
93  // "colz same" will reuse axis's min and max, so set them helpfully here
94  axes->SetMinimum(fHist->GetMinimum());
95  axes->SetMaximum(fHist->GetMaximum());
96  }
97 
98  axes->SetTitle(fHist->GetTitle());
99  axes->GetXaxis()->SetLabelSize(ax->GetLabelSize());
100  axes->GetYaxis()->SetLabelSize(ay->GetLabelSize());
101  axes->GetXaxis()->SetLabelOffset(ax->GetLabelOffset());
102  axes->GetYaxis()->SetLabelOffset(ay->GetLabelOffset());
103  axes->GetXaxis()->CenterTitle();
104  axes->GetYaxis()->CenterTitle();
105 
106  if(fLogX) gPad->SetLogx();
107  if(fLogY) gPad->SetLogy();
108 
109  gPad->Update();
110  }
double BinCenterX(int bin) const
Definition: ISurface.cxx:265
double BinCenterY(int bin) const
Definition: ISurface.cxx:272
TH2F * fHist
Definition: ISurface.h:68
bool fLogX
Definition: ISurface.h:69
bool fLogY
Definition: ISurface.h:69
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:29
void ana::ISurface::FillSurfObj ( ISurface surf,
TDirectory *  dir 
)
staticprotected

Definition at line 229 of file ISurface.cxx.

References fBestFitX, fBestFitY, fBestLikelihood, fBinMask, fHist, fSeedValues, compare_h5_caf::idx, m, and registry_explorer::v.

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

230  {
231 
232  const TVectorD v = *(TVectorD *) dir->Get("minValues");
233  const TVectorD s = *(TVectorD *) dir->Get("seeds");
234  TVectorD* m = (TVectorD *) dir->Get("mask");
235 
236  surf.fHist = (TH2F *) dir->Get("hist");
237  surf.fBestLikelihood = v[0];
238  surf.fBestFitX = v[1];
239  surf.fBestFitY = v[2];
240 
241  for (int idx = 0; idx < s.GetNrows(); ++idx)
242  surf.fSeedValues.push_back(s[idx]);
243 
244  if (m) {
245  for (int idx = 0; idx < m->GetNrows(); ++idx)
246  surf.fBinMask.push_back((*m)[idx]);
247  }
248 
249  }
const XML_Char * s
Definition: expat.h:262
TDirectory * dir
Definition: macro.C:5
surf
Definition: demo5.py:35
double ana::ISurface::GetBestFitX ( ) const
inline

Definition at line 30 of file ISurface.h.

References fBestFitX.

Referenced by test_stanfit_statsonly().

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

For expert use, custom painting of contours.

Definition at line 113 of file ISurface.cxx.

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

Referenced by DrawContour(), and GetBestFitY().

114  {
115  CheckMask("GetGraphs");
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  }
void CheckMask(const std::string &func) const
Definition: ISurface.cxx:279
Float_t tmp
Definition: plot.C:36
TH2F * fHist
Definition: ISurface.h:68
double fBestLikelihood
Definition: ISurface.h:65
std::string UniqueName()
Return a different string each time, for creating histograms.
Definition: Utilities.cxx:29
surf
Definition: demo5.py:35
void ana::ISurface::SaveToHelper ( TDirectory *  dir) const
protected

dir should already be the appropriate sub-directory

Definition at line 198 of file ISurface.cxx.

References fBestFitX, fBestFitY, fBestLikelihood, fBinMask, fHist, fLogX, fLogY, fSeedValues, m, tmp, and registry_explorer::v.

Referenced by GetBestFitY(), ana::BayesianSurface::SaveTo(), and ana::FrequentistSurface::SaveTo().

199  {
200  TDirectory *oldDir = gDirectory;
201 
202  dir->cd();
203 
204  TVectorD v(3);
205  v[0] = fBestLikelihood;
206  v[1] = fBestFitX;
207  v[2] = fBestFitY;
208  v.Write("minValues");
209 
210  fHist->Write("hist");
211 
212  TVectorD s(fSeedValues.size(), &fSeedValues[0]);
213  s.Write("seeds");
214 
215  if (!fBinMask.empty()) {
216  std::vector<double> tmp(fBinMask.begin(), fBinMask.end());
217  TVectorD m(tmp.size(), &tmp[0]);
218  m.Write("mask");
219  }
220 
221  TObjString(fLogX ? "yes" : "no").Write("logx");
222  TObjString(fLogY ? "yes" : "no").Write("logy");
223 
224  if (oldDir)
225  oldDir->cd();
226  }
Float_t tmp
Definition: plot.C:36
const XML_Char * s
Definition: expat.h:262
TH2F * fHist
Definition: ISurface.h:68
bool fLogX
Definition: ISurface.h:69
double fBestFitY
Definition: ISurface.h:67
bool fLogY
Definition: ISurface.h:69
std::vector< double > fSeedValues
Definition: ISurface.h:70
TDirectory * dir
Definition: macro.C:5
double fBestFitX
Definition: ISurface.h:66
std::vector< int > fBinMask
Definition: ISurface.h:71
double fBestLikelihood
Definition: ISurface.h:65
void ana::ISurface::SetTitle ( const char *  str)

Definition at line 192 of file ISurface.cxx.

References fHist.

Referenced by GetBestFitY(), and plots().

193  {
194  fHist->SetTitle(str);
195  }
TH2F * fHist
Definition: ISurface.h:68
TH2 * ana::ISurface::ToTH2 ( double  minchi = -1) const

Definition at line 170 of file ISurface.cxx.

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

Referenced by CAF_makeCAFSensitivities_for_FNEX(), demoFitContours(), GetBestFitY(), joint_fit_2017_contours(), joint_fit_2018_contours(), joint_fit_2019_contours(), joint_fit_future_contour_univ(), ana::FrequentistSurface::LoadFromMulti(), MakeCAFSensitivities_for_FNEX(), Plot2DSlice(), Plot2DSlices(), plot_3flavor_withsysts(), 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:68
double fBestLikelihood
Definition: ISurface.h:65

Friends And Related Function Documentation

TH2* Flat ( double  ,
const ISurface  
)
friend

Helper function for the gaussian approximation surfaces.

Definition at line 253 of file ISurface.cxx.

Referenced by ana::Gaussian2Sigma1D(), ana::Gaussian2Sigma2D(), ana::Gaussian3Sigma1D(), ana::Gaussian3Sigma2D(), ana::Gaussian68Percent1D(), ana::Gaussian68Percent2D(), ana::Gaussian90Percent1D(), ana::Gaussian90Percent2D(), ana::Gaussian95Percent1D(), ana::Gaussian95Percent2D(), ana::Gaussian99Percent1D(), ana::Gaussian99Percent2D(), GetBestFitY(), and ana::BayesianSurface::QuantileSurface().

254  {
255  TH2* h = new TH2F(*s.fHist);
256 
257  for(int x = 0; x < h->GetNbinsX()+2; ++x)
258  for(int y = 0; y < h->GetNbinsY()+2; ++y)
259  h->SetBinContent(x, y, level);
260 
261  return h;
262  }
const XML_Char * s
Definition: expat.h:262

Member Data Documentation

double ana::ISurface::fBestFitX
protected
double ana::ISurface::fBestFitY
protected
double ana::ISurface::fBestLikelihood
protected
std::vector<int> ana::ISurface::fBinMask
protected
TH2F* ana::ISurface::fHist
protected
bool ana::ISurface::fLogX
protected
bool ana::ISurface::fLogY
protected
std::vector<double> ana::ISurface::fSeedValues
protected

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