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

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

Inheritance diagram for ana::FrequentistSurface:

Public Member Functions
	FrequentistSurface (const IExperiment expt, osc::IOscCalcAdjustable calc, const FitAxis &xax, const FitAxis &yax, 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)

	FrequentistSurface (const IExperiment expt, osc::IOscCalcAdjustable 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 std::string &name) 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< FrequentistSurface >	LoadFrom (TDirectory *dir, const std::string &name)

static std::unique_ptr< FrequentistSurface >	LoadFromMulti (std::vector< TFile * > files, std::string label)

static std::unique_ptr< FrequentistSurface >	LoadFromMulti (std::string wildcard, std::string label)

Protected Member Functions
	FrequentistSurface ()

void	CreateHistograms (const FitAxis &xax, const FitAxis &yax, 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::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)

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)

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)

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
bool	fParallel

MinuitFitter::FitOpts	fFitOpts

std::vector< TH2 * >	fProfHists

double	fBestLikelihood

double	fBestFitX

double	fBestFitY

TH2F *	fHist

bool	fLogX

bool	fLogY

std::vector< double >	fSeedValues

std::vector< int >	fBinMask

Detailed Description

Log-likelihood scan across two parameters.

Definition at line 27 of file FrequentistSurface.h.

Constructor & Destructor Documentation

ana::FrequentistSurface::FrequentistSurface	(	const IExperiment *	expt,
		osc::IOscCalcAdjustable *	calc,
		const FitAxis &	xax,
		const FitAxis &	yax,
		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

expt	The experiment object to draw $\chi^2$ values from
calc	Values for oscillation parameters to be held fixed
xax	Definition of the x-axis
yax	Definition of the y-axis
profVars	Oscillation parameters to profile over
profSysts	Systematic parameters to profile over
seedPts	Try all combinations of these params as seeds
systSeedPts	Try all of these systematic combinations as seeds
parallel	Use all the cores on this machine? Be careful...

Definition at line 32 of file FrequentistSurface.cxx.

References CreateHistograms(), FillSurface(), FindMinimum(), ana::ISurface::fLogX, ana::ISurface::fLogY, ana::ISurface::fSeedValues, ana::FitAxis::islog, registry_explorer::v, and ana::FitAxis::var.

     : fParallel(parallel), fFitOpts(opts)
   {
     fLogX = xax.islog;
     fLogY = yax.islog;
 
     CreateHistograms(xax, yax, profVars, profSysts);
 
     for(const IFitVar* v: profVars) fSeedValues.push_back(v->GetValue( calc));
 
     FillSurface(expt, calc, xax.var, yax.var, profVars, profSysts, seedPts, systSeedPts);
 
     FindMinimum(expt, calc, xax.var, yax.var, profVars, profSysts, seedPts, systSeedPts);
   }

ana::FrequentistSurface::FrequentistSurface	(	const IExperiment *	expt,
		osc::IOscCalcAdjustable *	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`
	)

inline

Parameters

expt	The experiment object to draw $\chi^2$ values from
calc	Values for oscillation parameters to be held fixed
xvar	Oscillation parameter to place on the x axis
nbinsx	Number of bins along x axis
xmin	Minimum value of x axis
xmax	Maximum value of x axis
nbinsy	Number of bins along y axis
ymin	Minimum value of y axis
ymax	Maximum value of y axis
profVars	Oscillation parameters to profile over
profSysts	Systematic parameters to profile over
seedPts	Try all combinations of these params as seeds
systSeedPts	Try all of these systematic combinations as seeds
parallel	Use all the cores on this machine? Be careful...

Definition at line 63 of file FrequentistSurface.h.

References FrequentistSurface(), ana::MinuitFitter::kNormal, plot_validation_datamc::opts, and ~FrequentistSurface().

                                                                    {},
                        const std::vector<const ISyst*>& profSysts = {},
                        const SeedList& seedPts = SeedList(),
                        const std::vector<SystShifts>& systSeedPts = {},
                        bool parallel = false,
                        MinuitFitter::FitOpts opts = MinuitFitter::kNormal)
       : FrequentistSurface(expt, calc,
                            FitAxis(xvar, nbinsx, xmin, xmax, false),
                            FitAxis(yvar, nbinsy, ymin, ymax, false),
                            profVars, profSysts, seedPts, systSeedPts, parallel, opts)
     {
     }

ana::FrequentistSurface::~FrequentistSurface ( )

virtual

Definition at line 57 of file FrequentistSurface.cxx.

References CreateHistograms().

Referenced by FrequentistSurface().

58 {

59 }

ana::FrequentistSurface::FrequentistSurface ( )

inlineprotected

Definition at line 91 of file FrequentistSurface.h.

References CreateHistograms(), FillSurface(), FillSurfacePoint(), FindMinimum(), ProgressBarTitle(), string, submit_syst::x, and submit_syst::y.

Referenced by FrequentistSurface().

91 {};

Member Function Documentation

double ana::ISurface::BestLikelihood ( ) const

inlineinherited

Definition at line 29 of file ISurface.h.

References ana::ISurface::fBestLikelihood.

29 {return fBestLikelihood;}

ana::ISurface::fBestLikelihood

double fBestLikelihood

Definition: ISurface.h:65

double ana::ISurface::BinCenterX ( int bin ) const

protectedinherited

Definition at line 265 of file ISurface.cxx.

References visualisationForPaperMasterPlot::ax, ana::ISurface::fHist, and ana::ISurface::fLogX.

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

   {
     const TAxis* ax = fHist->GetXaxis();
     return fLogX ? ax->GetBinCenterLog(bin) : ax->GetBinCenter(bin);
   }

double ana::ISurface::BinCenterY ( int bin ) const

protectedinherited

Definition at line 272 of file ISurface.cxx.

References visualisationForPaperMasterPlot::ax, ana::ISurface::fHist, and ana::ISurface::fLogY.

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

   {
     const TAxis* ax = fHist->GetYaxis();
     return fLogY ? ax->GetBinCenterLog(bin) : ax->GetBinCenter(bin);
   }

void ana::ISurface::CheckMask ( const std::string & func ) const

protectedinherited

Definition at line 279 of file ISurface.cxx.

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

Referenced by ana::ISurface::DrawBestFit(), ana::ISurface::DrawContour(), and ana::ISurface::GetGraphs().

   {
     if(!fBinMask.empty()) {
       std::cout << "Cannot call " << func << "() on a partial surface!" << std::endl;
       abort();
     }
   }

void ana::FrequentistSurface::CreateHistograms	(	const FitAxis &	xax,
		const FitAxis &	yax,
		const std::vector< const IFitVar * > &	profVars,
		const std::vector< const ISyst * > &	profSysts
	)

protected

Definition at line 63 of file FrequentistSurface.cxx.

References ana::ExpandedHistogram(), ana::ISurface::fHist, fProfHists, MECModelEnuComparisons::i, ana::FitAxis::islog, ana::IFitVar::LatexName(), ana::FitAxis::max, ana::FitAxis::min, ana::FitAxis::nbins, ProgressBarTitle(), string, plotROC::title, and ana::FitAxis::var.

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

   {
     fHist = ExpandedHistogram(";"+xax.var->LatexName()+";"+yax.var->LatexName(),
                               xax.nbins, xax.min, xax.max, xax.islog,
                               yax.nbins, yax.min, yax.max, yax.islog);
 
     for(unsigned int i = 0; i < profVars.size()+profSysts.size(); ++i){
       std::string title;
       if(i < profVars.size())
         title = profVars[i]->LatexName();
       else
         title = profSysts[i-profVars.size()]->LatexName();
 
       fProfHists.push_back(ExpandedHistogram(title+";"+xax.var->LatexName()+";"+yax.var->LatexName(),
                                              xax.nbins, xax.min, xax.max, xax.islog,
                                              yax.nbins, yax.min, yax.max, yax.islog));
     }
 
   }

void ana::ISurface::Draw ( ) const

inherited

Draw the surface itself.

Definition at line 21 of file ISurface.cxx.

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

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

   {
     // Can be useful to draw a partial surface for debugging
     //    CheckMask();
 
     EnsureAxes();
 
     fHist->Draw("colz same");
 
     // colz obliterated them
     gPad->RedrawAxis();
 
     gPad->Update();
   }

void ana::ISurface::DrawBestFit	(	Color_t	color,
		Int_t	marker = `kFullCircle`
	)		const

inherited

Draw the best fit point.

Definition at line 37 of file ISurface.cxx.

References ana::ISurface::CheckMask(), ana::ISurface::EnsureAxes(), ana::ISurface::fBestFitX, and ana::ISurface::fBestFitY.

Referenced by cc(), demo5(), demo_CPT(), drawSensitivity(), ana::ISurface::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().

   {
     CheckMask("DrawBestFit");
     EnsureAxes();
 
     TMarker *mark = new TMarker(fBestFitX, fBestFitY, marker);
     mark->SetMarkerSize(1.5);
     mark->SetMarkerColor(color);
     mark->Draw();
     gPad->Update();
   }

void ana::ISurface::DrawContour	(	TH2 *	fc,
		Style_t	style,
		Color_t	color,
		double	minchi = `-1`
	)

inherited

Parameters

fc	Surface to compare against for this significance level
style	Line style for TAttLine, solid, dotted, dashed etc
color	Line 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 ana::ISurface::CheckMask(), ana::ISurface::EnsureAxes(), MECModelEnuComparisons::g, ana::ISurface::GetGraphs(), and APDGainPoints::gs.

Referenced by cc(), demo5(), Draw2DSurface(), DrawContours(), drawSensitivity(), ana::ISurface::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().

   {
     CheckMask("DrawContour");
     EnsureAxes();
 
     std::vector<TGraph *> gs = GetGraphs(fc, minchi);
 
     for (TGraph *g: gs)
     {
       g->SetLineWidth(3);//2);
       g->SetLineStyle(style);
       g->SetLineColor(color);
       g->Draw("l");
     }
 
     gPad->Update();
   }

void ana::ISurface::EnsureAxes ( ) const

protectedinherited

Definition at line 70 of file ISurface.cxx.

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

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

   {
     // Could have a file temporarily open
     DontAddDirectory guard;
 
     // If this pad has already been drawn in, already has axes
     if (gPad && !gPad->GetListOfPrimitives()->IsEmpty()) return;
 
     const TAxis *ax = fHist->GetXaxis();
     const TAxis *ay = fHist->GetYaxis();
     const double Nx = ax->GetNbins();
     const double Ny = ay->GetNbins();
 
     // Axes with limits where the user originally requested them, which we
     // adjusted to be the centres of the first and last bins.
     TH2 *axes = new TH2C(UniqueName().c_str(),
                          TString::Format(";%s;%s",
                                          ax->GetTitle(), ay->GetTitle()),
                          Nx - 1, BinCenterX(1), BinCenterX(Nx),
                          Ny - 1, BinCenterY(1), BinCenterY(Ny));
     axes->Draw();
 
     if(fHist){
       // "colz same" will reuse axis's min and max, so set them helpfully here
       axes->SetMinimum(fHist->GetMinimum());
       axes->SetMaximum(fHist->GetMaximum());
     }
 
     axes->SetTitle(fHist->GetTitle());
     axes->GetXaxis()->SetLabelSize(ax->GetLabelSize());
     axes->GetYaxis()->SetLabelSize(ay->GetLabelSize());
     axes->GetXaxis()->SetLabelOffset(ax->GetLabelOffset());
     axes->GetYaxis()->SetLabelOffset(ay->GetLabelOffset());
     axes->GetXaxis()->CenterTitle();
     axes->GetYaxis()->CenterTitle();
 
     if(fLogX) gPad->SetLogx();
     if(fLogY) gPad->SetLogy();
 
     gPad->Update();
   }

void ana::FrequentistSurface::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
	)

protectedvirtual

Reimplemented in ana::SurfaceKrige.

Definition at line 108 of file FrequentistSurface.cxx.

References ana::ThreadPool::AddTask(), bin, ana::ISurface::BinCenterX(), ana::ISurface::BinCenterY(), om::cerr, ana::IExperiment::ChiSq(), ana::Progress::Done(), e, allTimeWatchdog::endl, ana::ISurface::fBinMask, ana::ISurface::fHist, FillSurfacePoint(), ana::ThreadPool::Finish(), genie::utils::style::Format(), fParallel, site_stats_from_log::job, ana::JobNumber(), ana::ThreadPool::NThreads(), ana::NumJobs(), ana::IFitVar::Penalty(), time_estimates::pool, cacheDefinitionData::prog, ProgressBarTitle(), ana::RunningOnGrid(), ana::Progress::SetProgress(), ana::IFitVar::ShortName(), ana::ThreadPool::ShowProgress(), fillBadChanDBTables::step, update_sam_good_runs_metadata::stride, string, submit_syst::x, and submit_syst::y.

Referenced by FrequentistSurface().

   {
     // Nothing created during surface filling belongs in a
     // directory. Unfortunately the local guards in Spectrum etc are racey when
     // run in parallel. But this should cover the whole lot safely.
     DontAddDirectory guard;
 
     const std::string progTitle = ProgressBarTitle(xvar, yvar, profVars, profSysts);
 
     Progress *prog = 0;
     // Difficult to keep a progress bar properly up to date when threaded
     if (!fParallel) prog = new Progress(progTitle);
     ThreadPool *pool = 0;
 
     if(fParallel){
       // A hack/workaround needed for parallel running:
       //
       // Give all the constituents of the Prediction a chance to do lazy
       // initialization, before they race themselves trying to do it in
       // parallel.
       expt->ChiSq(calc);
 
       pool = new ThreadPool;
       pool->ShowProgress(progTitle + TString::Format(" using %d threads", pool->NThreads()).Data());
     }
 
     const int Nx = fHist->GetNbinsX();
     const int Ny = fHist->GetNbinsY();
 
     // Fill bins in "random" order so that the progress bar is accurate
     // straight away instead of possibly being misled by whatever atypical
     // points we start with. This step is a prime which guarantees we get every
     // cell.
     int step = 7919;
     // Very unlikely (Nx or Ny is a multiple of step), but just to be safe.
     if ((Nx * Ny) % step == 0) step = 1;
 
     int bin = 0;
     int neval = 0;
 
     // Allow the surface to be parallelised across multiple jobs by splitting up
     // the full surface into patches, and only running bins that fall inside a
     // certain patch
     int first = 0, last = Nx * Ny;
     if (RunningOnGrid() && NumJobs()>1) {
       double stride = double(Nx*Ny) / double(NumJobs());
       int job = JobNumber();
       first = job * stride;
       last = (job+1) * stride;
     }
 
     do{
       if (neval < first) {
         ++neval;
         bin = (bin + step) % (Nx * Ny);
         continue;
       } else if (neval >= last) {
         break;
       }
 
       const int x = bin % Nx + 1;
       const int y = bin / Nx + 1;
 
       const double xv = BinCenterX(x);
       const double yv = BinCenterY(y);
 
       if (xvar->Penalty(xv, calc) > 1e-10)
       {
         std::cerr << "Warning! " << xvar->ShortName() << " = " << xv
                   << " has penalty of " << xvar->Penalty(xv, calc)
                   << " that could have been applied in surface. "
                   << "This should never happen." << std::endl;
       }
       if (yvar->Penalty(yv, calc) > 1e-10)
       {
         std::cerr << "Warning! " << yvar->ShortName() << " = " << yv
                   << " has penalty of " << yvar->Penalty(yv, calc)
                   << " that could have been applied in surface. "
                   << "This should never happen." << std::endl;
       }
 
       ThreadPool::func_t task = [=](){
         FillSurfacePoint(expt, calc,
                          xvar, xv, yvar, yv,
                          profVars, profSysts, seedPts, systSeedPts);
       };
 
       ++neval;
       if(fParallel){
         pool->AddTask(task);
       }
       else{
         task(); // Just do it straight away
         prog->SetProgress(double(neval-first) / double(last-first));
       }
 
       if (RunningOnGrid() && NumJobs()>1) fBinMask.push_back(bin);
 
       bin = (bin + step) % (Nx * Ny);
     } while (bin != 0);
 
 
     if(fParallel){
       pool->Finish();
       delete pool;
     }
     else{
       prog->Done();
       delete prog;
     }
   }

double ana::FrequentistSurface::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
	)

protected

Definition at line 227 of file FrequentistSurface.cxx.

References calc, ana::IExperiment::ChiSq(), osc::_IOscCalcAdjustable< 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().

   {
     if(fParallel){
       // Need to take our own copy so that we don't get overwritten by someone
       // else's changes.
       calc = calc->Copy();
     }
 
     xvar->SetValue(calc, x);
     yvar->SetValue(calc, y);
 
     //Make sure that the profiled values of fitvars do not persist between steps.
     for(int i = 0; i < (int)fSeedValues.size(); ++i) profVars[i]->SetValue( calc, fSeedValues[i] );
 
     double chi;
     if(profVars.empty() && profSysts.empty()){
       chi = expt->ChiSq(calc);
     }
     else{
       MinuitFitter fitter(expt, profVars, profSysts);
       fitter.SetFitOpts(fFitOpts);
       SystShifts bestSysts;
       chi = fitter.Fit(calc, bestSysts, seedPts, systSeedPts, MinuitFitter::kQuiet)->EvalMetricVal();
 
       for(unsigned int i = 0; i < profVars.size(); ++i){
         fProfHists[i]->Fill(x, y, profVars[i]->GetValue(calc));
       }
       for(unsigned int j = 0; j < profSysts.size(); ++j){
         fProfHists[j+profVars.size()]->Fill(x, y, bestSysts.GetShift(profSysts[j]));
       }
     }
 
     fHist->Fill(x, y, chi);
 
     if(fParallel) delete calc;
 
     return chi;
   }

void ana::ISurface::FillSurfObj	(	ISurface &	surf,
		TDirectory *	dir
	)

staticprotectedinherited

Definition at line 229 of file ISurface.cxx.

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

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

   {
 
     const TVectorD v = *(TVectorD *) dir->Get("minValues");
     const TVectorD s = *(TVectorD *) dir->Get("seeds");
     TVectorD* m = (TVectorD *) dir->Get("mask");
 
     surf.fHist = (TH2F *) dir->Get("hist");
     surf.fBestLikelihood = v[0];
     surf.fBestFitX = v[1];
     surf.fBestFitY = v[2];
 
     for (int idx = 0; idx < s.GetNrows(); ++idx)
       surf.fSeedValues.push_back(s[idx]);
 
     if (m) {
       for (int idx = 0; idx < m->GetNrows(); ++idx)
         surf.fBinMask.push_back((*m)[idx]);
     }
 
   }

void ana::FrequentistSurface::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
	)

protected

Definition at line 275 of file FrequentistSurface.cxx.

References bin, ana::ISurface::BinCenterX(), ana::ISurface::BinCenterY(), confusionMatrixTree::count, ana::ISurface::fBestFitX, ana::ISurface::fBestFitY, ana::ISurface::fBestLikelihood, ana::ISurface::fBinMask, fFitOpts, ana::ISurface::fHist, ana::IFitter::Fit(), ana::ISurface::fSeedValues, ana::IFitVar::GetValue(), MECModelEnuComparisons::i, makeTrainCVSamples::int, ana::SystShifts::Nominal(), ana::NumJobs(), ana::RunningOnGrid(), ana::MinuitFitter::SetFitOpts(), ana::IFitVar::SetValue(), submit_syst::x, and submit_syst::y.

Referenced by FrequentistSurface().

   {
     // Location of the best minimum found from filled surface
     double minchi = 1e10;
     int minx = fHist->GetNbinsX()/2;
     int miny = fHist->GetNbinsY()/2;
     for(int x = 1; x <= fHist->GetNbinsX(); ++x){
       for(int y = 1; y <= fHist->GetNbinsY(); ++y){
         int bin = ((y-1)*fHist->GetNbinsX())+(x-1);
         const double chi = fHist->GetBinContent(x, y);
         if(RunningOnGrid() && NumJobs()>1
           && !std::count(fBinMask.begin(), fBinMask.end(), bin)) continue;
         if (chi < minchi){
           minchi = chi;
           minx = x;
           miny = y;
         }
       }
     }
 
     std::vector<const IFitVar*> allVars = {xvar, yvar};
     allVars.insert(allVars.end(), profVars.begin(), profVars.end());
     MinuitFitter fit(expt, allVars, profSysts);
     fit.SetFitOpts(fFitOpts);
     // Seed from best grid point
     xvar->SetValue(calc, BinCenterX(minx));
     yvar->SetValue(calc, BinCenterY(miny));
     for(int i = 0; i < (int)fSeedValues.size(); ++i) profVars[i]->SetValue( calc, fSeedValues[i] );
     SystShifts systSeed = SystShifts::Nominal();
     fBestLikelihood = fit.Fit(calc, systSeed, seedPts)->EvalMetricVal();
     fBestFitX = xvar->GetValue(calc);
     fBestFitY = yvar->GetValue(calc);
 
     for(int x = 0; x < fHist->GetNbinsX()+2; ++x){
       for(int y = 0; y < fHist->GetNbinsY()+2; ++y){
         fHist->SetBinContent(x, y, fHist->GetBinContent(x, y)-fBestLikelihood);
       }
     }
 
     fHist->SetMinimum(0);
   }

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;}

ana::ISurface::fBestFitX

double fBestFitX

Definition: ISurface.h:66

double ana::ISurface::GetBestFitY ( ) const

inlineinherited

Definition at line 31 of file ISurface.h.

References bin, ana::ISurface::BinCenterX(), ana::ISurface::BinCenterY(), plot_validation_datamc_2018::color, dir, ana::ISurface::Draw(), ana::ISurface::DrawBestFit(), ana::ISurface::DrawContour(), ana::ISurface::fBestFitY, ana::ISurface::FillSurfObj(), ana::ISurface::Flat, ana::ISurface::GetGraphs(), ana::ISurface::SaveToHelper(), ana::ISurface::SetTitle(), submit_syst::str, demo5::surf, and ana::ISurface::ToTH2().

Referenced by test_stanfit_statsonly().

31 {return fBestFitY;}

ana::ISurface::fBestFitY

double fBestFitY

Definition: ISurface.h:67

std::vector< TGraph * > ana::ISurface::GetGraphs	(	TH2 *	fc,
		double	minchi = `-1`
	)

inherited

For expert use, custom painting of contours.

Definition at line 113 of file ISurface.cxx.

References ana::ISurface::CheckMask(), 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().

   {
     CheckMask("GetGraphs");
 
     std::vector<TGraph*> ret;
 
     if(minchi < 0) minchi = fBestLikelihood;
     std::unique_ptr<TH2> surf((TH2*)fHist->Clone(UniqueName().c_str()));
     surf->Add(fc, -1);
 
     TVirtualPad* bak = gPad;
 
     const bool wasbatch = gROOT->IsBatch();
     gROOT->SetBatch(); // User doesn't want to see our temporary canvas
     TCanvas tmp;
 
     gStyle->SetOptStat(0);
 
     const double level = minchi-fBestLikelihood;
     surf->SetContour(1, &level);
     surf->Draw("cont list");
 
     tmp.Update();
     tmp.Paint();
 
     gROOT->SetBatch(wasbatch);
     gPad = bak;
 
     // The graphs we need (contained inside TLists, contained inside
     // TObjArrays) are in the list of specials. But we need to be careful about
     // types, because other stuff can get in here too (TDatabasePDG for
     // example).
     TCollection* specs = gROOT->GetListOfSpecials();
 
     TIter nextSpec(specs);
     while(TObject* spec = nextSpec()){
       if(!spec->InheritsFrom(TObjArray::Class())) continue;
       TObjArray* conts = (TObjArray*)spec;
 
       if(conts->IsEmpty()) continue;
 
       if(!conts->At(0)->InheritsFrom(TList::Class())) continue;
       TList* cont = (TList*)conts->At(0);
 
       TIter nextObj(cont);
       // Contour could be split into multiple pieces
       while(TObject* obj = nextObj()){
         if(!obj->InheritsFrom(TGraph::Class())) continue;
 
         ret.push_back((TGraph*)obj->Clone(UniqueName().c_str()));
       } // end for obj
     } // end for spec
 
     return ret;
   }

std::vector<TH2*> ana::FrequentistSurface::GetProfiledHists ( )

inline

Maps of the values taken on by the profiled parameters.

Definition at line 83 of file FrequentistSurface.h.

References dir, compareCafs::files, fProfHists, label, LoadFrom(), LoadFromMulti(), SaveTo(), string, and wildcard().

Referenced by run_joint_fit_2020_contours().

83 {return fProfHists;}

ana::FrequentistSurface::fProfHists

std::vector< TH2 * > fProfHists

Definition: FrequentistSurface.h:131

std::unique_ptr< FrequentistSurface > ana::FrequentistSurface::LoadFrom	(	TDirectory *	dir,
		const std::string &	name
	)

static

Definition at line 350 of file FrequentistSurface.cxx.

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

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

   {
     dir = dir->GetDirectory(name.c_str()); // switch to subdir
     assert(dir);
 
     DontAddDirectory guard;
 
     TObjString *tag = (TObjString *) dir->Get("type");
     assert(tag);
     assert(tag->GetString() == "FrequentistSurface" || tag->GetString() == "Surface");
     delete tag;
 
     std::unique_ptr<FrequentistSurface> surf(new FrequentistSurface);
     ISurface::FillSurfObj(*surf, dir);
 
     for (std::size_t idx = 0; idx < surf->fSeedValues.size(); ++idx)
     {
       // Search for old "marg" name here too for backwards compatibility
       TH2 *h = (TH2 *) dir->Get(TString::Format("profHists/hist%lu", idx));
       if (h)
         surf->fProfHists.push_back(h);
       else
         surf->fProfHists.push_back((TH2 *) dir->Get(TString::Format("margHists/hist%lu", idx)));
     }
 
     const TObjString* logx = (TObjString*)dir->Get("logx");
     const TObjString* logy = (TObjString*)dir->Get("logy");
     // Tolerate missing log tags for backwards compatibility
     surf->fLogX = (logx && logx->GetString() == "yes");
     surf->fLogY = (logy && logy->GetString() == "yes");
 
     delete dir;
 
     return surf;
   }

std::unique_ptr< FrequentistSurface > ana::FrequentistSurface::LoadFromMulti	(	std::vector< TFile * >	files,
		std::string	label
	)

static

Definition at line 387 of file FrequentistSurface.cxx.

References ana::assert(), bin, confusionMatrixTree::count, om::cout, allTimeWatchdog::endl, MakeMiniprodValidationCuts::f, ana::ISurface::fBestLikelihood, ana::ISurface::fBinMask, fProfHists, make_syst_table_plots::h, MECModelEnuComparisons::i, calib::j, LoadFrom(), nbins, runNovaSAM::ret, ana::ISurface::ToTH2(), submit_syst::x, and submit_syst::y.

Referenced by GetProfiledHists(), LoadFromMulti(), MergeSurface(), MergeSurfaces(), and test_surf_stride().

   {
     std::vector<std::unique_ptr<FrequentistSurface>> surfs;
     for (TFile* f : files) {
       surfs.push_back(FrequentistSurface::LoadFrom(f, label));
     }
 
     int Nx = surfs[0]->fHist->GetNbinsX();
     int Ny = surfs[0]->fHist->GetNbinsY();
     size_t nbins = Nx * Ny;
     std::vector<int> binMask;
 
     // Loop over the surfaces to check all bins are present
     for (size_t i = 0; i < surfs.size(); ++i) {
       for (int bin : surfs[i]->fBinMask) {
         if (std::count(binMask.begin(), binMask.end(), bin))
           assert(false && "Repeated bin found in surfaces being merged.");
         binMask.push_back(bin);
       }
     }
     if (binMask.size() != nbins) {
       std::cout << "Missing bins found in surfaces being merged. "
         << "Are you sure you included all files for this surface?"
         << std::endl;
       assert(false);
     }
 
     DontAddDirectory guard;
 
     // Create return surface and initialise with first in list
     std::unique_ptr<FrequentistSurface> ret(new FrequentistSurface);
     ret->fHist = new TH2F(*(TH2F*)surfs[0]->ToTH2(0));
     ret->fBestLikelihood = surfs[0]->fBestLikelihood;
     ret->fBestFitX = surfs[0]->fBestFitX;
     ret->fBestFitY = surfs[0]->fBestFitY;
     ret->fLogX = surfs[0]->fLogX;
     ret->fLogY = surfs[0]->fLogY;
     ret->fSeedValues = surfs[0]->fSeedValues;
     for (TH2* h : surfs[0]->fProfHists)
       ret->fProfHists.push_back((TH2F*)h->Clone());
 
     // Loop over other surfaces and add them in
     for (size_t i = 1; i < surfs.size(); ++i) {
       ret->fHist->Add(surfs[i]->ToTH2(0));
       if (surfs[i]->fBestLikelihood < ret->fBestLikelihood) {
         ret->fBestLikelihood = surfs[i]->fBestLikelihood;
         ret->fBestFitX = surfs[i]->fBestFitX;
         ret->fBestFitY = surfs[i]->fBestFitY;
       }
       for (size_t j = 0; j < ret->fProfHists.size(); ++j)
         ret->fProfHists[j]->Add(surfs[i]->fProfHists[j]);
     }
 
     // Scale hist by global minimum
     for(int x = 0; x < ret->fHist->GetNbinsX()+2; ++x){
       for(int y = 0; y < ret->fHist->GetNbinsY()+2; ++y){
         ret->fHist->SetBinContent(x, y, ret->fHist->GetBinContent(x, y)-ret->fBestLikelihood);
       }
     }
 
     ret->fHist->SetMinimum(0);
 
     return ret;
   }

std::unique_ptr< FrequentistSurface > ana::FrequentistSurface::LoadFromMulti	(	std::string	wildcard,
		std::string	label
	)

static

Definition at line 454 of file FrequentistSurface.cxx.

References MakeMiniprodValidationCuts::f, plotROC::fileName, compareCafs::files, LoadFromMulti(), runNovaSAM::ret, string, and ana::Wildcard().

   {
     std::vector<std::string> fileNames = Wildcard(wildcard);
     std::vector<TFile*> files;
     for (std::string fileName : fileNames) {
       files.push_back(TFile::Open(fileName.c_str(), "read"));
     }
     std::unique_ptr<FrequentistSurface> ret = LoadFromMulti(files, label);
     for(TFile* f: files) delete f;
     return ret;
   }

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 87 of file FrequentistSurface.cxx.

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

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

   {
     std::string title = "Filling surface for "+yvar->ShortName()+" vs "+xvar->ShortName();
 
     if(!profVars.empty() || !profSysts.empty()){
       title += " (profiling ";
 
       for(const IFitVar* v: profVars) title += v->ShortName() + ", ";
       for(const ISyst* s: profSysts)  title += s->ShortName() + ", ";
 
       // Always have one superfluous ", " at the end
       title.resize(title.size()-2);
       title += ")";
     }
 
     return title;
   }

void ana::FrequentistSurface::SaveTo	(	TDirectory *	dir,
		const std::string &	name
	)		const

Definition at line 324 of file FrequentistSurface.cxx.

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

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

   {
     TDirectory *tmp = gDirectory;
 
     dir = dir->mkdir(name.c_str()); // switch to subdir
     dir->cd();
 
     TObjString("FrequentistSurface").Write("type");
 
     ISurface::SaveToHelper(dir);
 
     TDirectory *profDir = dir->mkdir("profHists");
     int idx = 0;
     for (auto it: FrequentistSurface::fProfHists)
     {
       profDir->cd();
       it->Write(TString::Format("hist%d", idx++));
     }
 
     dir->Write();
     delete dir;
 
     tmp->cd();
   }

void ana::ISurface::SaveToHelper ( TDirectory * dir ) const

protectedinherited

dir should already be the appropriate sub-directory

Definition at line 198 of file ISurface.cxx.

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

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

   {
     TDirectory *oldDir = gDirectory;
 
     dir->cd();
 
     TVectorD v(3);
     v[0] = fBestLikelihood;
     v[1] = fBestFitX;
     v[2] = fBestFitY;
     v.Write("minValues");
 
     fHist->Write("hist");
 
     TVectorD s(fSeedValues.size(), &fSeedValues[0]);
     s.Write("seeds");
 
     if (!fBinMask.empty()) {
       std::vector<double> tmp(fBinMask.begin(), fBinMask.end());
       TVectorD m(tmp.size(), &tmp[0]);
       m.Write("mask");
     }
 
     TObjString(fLogX ? "yes" : "no").Write("logx");
     TObjString(fLogY ? "yes" : "no").Write("logy");
 
     if (oldDir)
       oldDir->cd();
   }

void ana::ISurface::SetTitle ( const char * str )

inherited

Definition at line 192 of file ISurface.cxx.

References ana::ISurface::fHist.

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

   {
     fHist->SetTitle(str);
   }

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::ISurface::GetBestFitY(), joint_fit_2017_contours(), joint_fit_2018_contours(), joint_fit_2019_contours(), joint_fit_future_contour_univ(), LoadFromMulti(), MakeCAFSensitivities_for_FNEX(), Plot2DSlice(), Plot2DSlices(), plot_3flavor_withsysts(), run_joint_fit_2020_contours(), sensitivity2018(), sensitivity2020(), and test_stanfit_withsysts().

   {
     // Could have a file temporarily open
     DontAddDirectory guard;
 
     TH2 *ret = new TH2F(*fHist);
 
     if (minchi >= 0)
     {
       for (int x = 0; x < ret->GetNbinsX() + 2; ++x)
       {
         for (int y = 0; y < ret->GetNbinsY() + 2; ++y)
         {
           ret->SetBinContent(x, y, ret->GetBinContent(x, y) + fBestLikelihood - minchi);
         }
       }
     }
 
     return ret;
   }

Member Data Documentation

double ana::ISurface::fBestFitX

protectedinherited

Definition at line 66 of file ISurface.h.

Referenced by ana::BayesianSurface::BuildHist(), ana::ISurface::DrawBestFit(), ana::ISurface::FillSurfObj(), FindMinimum(), ana::ISurface::GetBestFitX(), and ana::ISurface::SaveToHelper().

double ana::ISurface::fBestFitY

protectedinherited

Definition at line 67 of file ISurface.h.

Referenced by ana::BayesianSurface::BuildHist(), ana::ISurface::DrawBestFit(), ana::ISurface::FillSurfObj(), FindMinimum(), ana::ISurface::GetBestFitY(), and ana::ISurface::SaveToHelper().

double ana::ISurface::fBestLikelihood

protectedinherited

Definition at line 65 of file ISurface.h.

Referenced by ana::ISurface::BestLikelihood(), ana::BayesianSurface::BuildHist(), ana::ISurface::FillSurfObj(), FindMinimum(), ana::ISurface::GetGraphs(), LoadFromMulti(), ana::ISurface::SaveToHelper(), and ana::ISurface::ToTH2().

std::vector<int> ana::ISurface::fBinMask

protectedinherited

Definition at line 71 of file ISurface.h.

Referenced by ana::ISurface::CheckMask(), FillSurface(), ana::ISurface::FillSurfObj(), FindMinimum(), LoadFromMulti(), and ana::ISurface::SaveToHelper().

MinuitFitter::FitOpts ana::FrequentistSurface::fFitOpts

protected

Definition at line 128 of file FrequentistSurface.h.

Referenced by FillSurfacePoint(), and FindMinimum().

TH2F* ana::ISurface::fHist

protectedinherited

Definition at line 68 of file ISurface.h.

Referenced by ana::ISurface::BinCenterX(), ana::ISurface::BinCenterY(), ana::BayesianSurface::BuildHist(), CreateHistograms(), ana::ISurface::Draw(), ana::ISurface::EnsureAxes(), FillSurface(), FillSurfacePoint(), ana::ISurface::FillSurfObj(), FindMinimum(), ana::Flat(), ana::ISurface::GetGraphs(), ana::BayesianSurface::QuantileSurface(), ana::ISurface::SaveToHelper(), ana::ISurface::SetTitle(), and ana::ISurface::ToTH2().

bool ana::ISurface::fLogX

protectedinherited

Definition at line 69 of file ISurface.h.

Referenced by ana::BayesianSurface::BayesianSurface(), ana::ISurface::BinCenterX(), ana::BayesianSurface::BuildHist(), ana::ISurface::EnsureAxes(), FrequentistSurface(), and ana::ISurface::SaveToHelper().

bool ana::ISurface::fLogY

protectedinherited

Definition at line 69 of file ISurface.h.

Referenced by ana::BayesianSurface::BayesianSurface(), ana::ISurface::BinCenterY(), ana::BayesianSurface::BuildHist(), ana::ISurface::EnsureAxes(), FrequentistSurface(), and ana::ISurface::SaveToHelper().

bool ana::FrequentistSurface::fParallel

protected

Definition at line 126 of file FrequentistSurface.h.

Referenced by FillSurface(), and FillSurfacePoint().

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

protected

Definition at line 131 of file FrequentistSurface.h.

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

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

protectedinherited

Definition at line 70 of file ISurface.h.

Referenced by FillSurfacePoint(), ana::ISurface::FillSurfObj(), FindMinimum(), FrequentistSurface(), and ana::ISurface::SaveToHelper().

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

/cvmfs/nova-development.opensciencegrid.org/novasoft/releases/N21-01-25/CAFAna/Fit/FrequentistSurface.h
/cvmfs/nova-development.opensciencegrid.org/novasoft/releases/N21-01-25/CAFAna/Fit/FrequentistSurface.cxx

Public Member Functions

Static Public Member Functions

Protected Member Functions

Static Protected Member Functions

Protected Attributes

Detailed Description

Constructor & Destructor Documentation

Member Function Documentation

Member Data Documentation