39 #include "Utilities/func/MathUtil.h" 59 calc.
SetTh23(TMath::DegToRad() * 43);
65 bool savePredToFile=
false,
66 bool loadSamplesFromFile=
true,
67 bool saveSamplesToFile=
false)
69 auto TFText = [](
bool val) {
return val ?
"TRUE" :
"FALSE"; };
71 std::cout <<
"Loading predictions from file: " << TFText(loadPredFromFile) <<
std::endl;
73 std::cout <<
"Loading MCMC samples from file: " << TFText(loadSamplesFromFile) <<
std::endl;
74 std::cout <<
"Saving MCMC samples to file: " << TFText(saveSamplesToFile) <<
std::endl;
85 std::unique_ptr<ana::PredictionExtrap>
pred;
86 if (!loadPredFromFile)
89 pred = std::make_unique<ana::PredictionNoExtrap>(
103 if (savePredToFile && !loadPredFromFile)
106 pred->SaveTo(&outf,
"pred");
124 std::unique_ptr<ana::MCMCSamples> samples;
125 if (!loadSamplesFromFile)
129 config.num_warmup = 500;
130 config.num_samples = 2000;
134 fitter.SetStanConfig(config);
136 fitter.Fit(&calc, systs);
137 samples = fitter.ExtractSamples();
139 if (saveSamplesToFile)
142 samples->SaveTo(&outf,
"samples");
153 &fitDmSq32Scaled_UniformPrior, 30, 2.2, 2.8);
158 marker.SetMarkerColor(kMagenta);
159 marker.SetMarkerSize(3);
161 c.SaveAs(
"/nova/ana/users/jwolcott/scratch/test_stanfit_statsonly_surface_contour.png");
163 auto marg_dm2 = samples->MarginalizeTo(&fitDmSq32Scaled_UniformPrior);
165 auto h_marg_dm2 = marg_dm2.ToTH1(bin_marg_dm2);
167 c.SaveAs(
"/nova/ana/users/jwolcott/scratch/test_stanfit_statsonly_dm2_marg.png");
176 c.SaveAs(
"/nova/ana/users/jwolcott/scratch/test_stanfit_statsonly_bestfitpred.png");
virtual void SetL(double L)=0
void SetDmsq32(const T &dmsq32) override
virtual void SetDmsq21(const T &dmsq21)=0
Simple record of shifts applied to systematic parameters.
TH2 * QuantileSurface(Quantile quantile) const
Configuration parameters for the Stan MCMC fitter, bundled up together for easy storage and parameter...
Helper struct for the cache. Might not need this.
Divide bin contents by bin widths.
virtual void SetTh13(const T &th13)=0
const std::string SAVED_SAMPLES_FILE
double GetBestFitY() const
const std::string SAVED_PRED_FILE
T sqr(T x)
More efficient square function than pow(x,2)
virtual void SetDmsq32(const T &dmsq32)=0
Representation of a spectrum in any variable, with associated POT.
void DrawBestFit(Color_t color, Int_t marker=kFullCircle) const
Draw the best fit point.
static std::unique_ptr< MCMCSamples > LoadFrom(TDirectory *dir, const std::string &name)
Load from file.
For nominal spectra and reweighting systs (xsec/flux)
Spectrum FakeData(double pot) const
Synonymous with AsimovData(). Retained for compatibility.
void SetValue(osc::IOscCalcAdjustable *osc, double val) const override
void Draw() const
Draw the surface itself.
void DrawContour(TH2 *fc, Style_t style, Color_t color, double minchi=-1)
const ConstrainedFitVarWithPrior fitDmSq32Scaled_UniformPrior
void ResetCalculator(osc::IOscCalcAdjustable &calc)
const ConstrainedFitVarWithPrior fitDmSq32Scaled_UniformPrior & kFitDmSq32Scaled
const HistAxis kNumuCCOptimisedAxis("Reconstructed Neutrino Energy (GeV)", kNumuCCEOptimisedBinning, kCCE)
HistAxis that implements optimised numuCCE from L. Vinton. See docdb 16332. This was close to 'custC'...
const SystShifts kNoShift
const ConstrainedFitVarWithPrior fitSsqTh23_UniformPriorSsqTh23
virtual void SetRho(double rho)=0
void SetTh23(const T &th23) override
const ConstrainedFitVarWithPrior fitSsqTh23_UniformPriorSsqTh23 & kFitSinSqTheta23
virtual void SetTh23(const T &th23)=0
stan::math::var PriorUniformInFitVar(const stan::math::var &, const osc::IOscCalcAdjustableStan *)
Prior uniform in the variable that is being fitted.
Version of FitVarWithPrior for use with constrained FitVar_StanSupports.
TH1 * DataMCComparison(const Spectrum &data, const Spectrum &mc, EBinType bintype)
void test_stanfit_statsonly(bool loadPredFromFile=true, bool savePredToFile=false, bool loadSamplesFromFile=true, bool saveSamplesToFile=false)
static Binning Simple(int n, double lo, double hi, const std::vector< std::string > &labels={})
Fitter type that bolts the Stan fitting tools onto CAFAna.
virtual void SetTh12(const T &th12)=0
virtual void SetdCP(const T &dCP)=0
Compare a single data spectrum to the MC + cosmics expectation.
double GetBestFitX() const