Page MenuHomeHEPForge
Files F8308808

No OneTemporary
Actions

View Options

diff --git a/parameters/config.xml b/parameters/config.xml
index f093aee..01df49d 100644
--- a/parameters/config.xml
+++ b/parameters/config.xml
@@ -1,219 +1,226 @@
<nuisance>
<!-- # ###################################################### -->
<!-- # NUISANCE CONFIGURATION OPTIONS -->
<!-- # This file is read in by default at runtime -->
<!-- # If you want to override on a case by case bases use -q at runtime -->
<!-- # ###################################################### -->
<!-- # MAIN Configs -->
<!-- # ###################################################### -->
<!-- # Logger goes from -->
<!-- # 1 Quiet -->
<!-- # 2 Fitter -->
<!-- # 3 Samples -->
<!-- # 4 Reconfigure Loops -->
<!-- # 5 Every Event print out (SHOUT) -->
<!-- # -1 DEBUGGING -->
<config VERBOSITY='4'/>
<!-- # ERROR goes from -->
<!-- # 0 NONE -->
<!-- # 1 FATAL -->
<!-- # 2 WARN -->
<config ERROR='2'/>
<config TRACE='0'/>
<config cores='1' />
<config spline_test_throws='50' />
<config spline_cores='1' />
<config spline_chunks='20' />
<config spline_procchunk='-1' />
<config Electron_NThetaBins='4' />
<config Electron_NEnergyBins='4' />
<config Electron_ThetaWidth='1.0' />
<config Electron_EnergyWidth='0.10' />
<config RemoveFSIParticles='0' />
<config RemoveUndefParticles='0' />
<config RemoveNuclearParticles='0'/>
<config MINERvASaveExtraCCQE='0' />
<!-- # Input Configs -->
<!-- # ###################################################### -->
<!-- # Default Requirements file for the externalDataFitter Package -->
<!-- # MAX Events : -1 is no cut. Events will be scaled automatically to give good xsec predictions. -->
<config MAXEVENTS='-1'/>
<!-- Include empty stacks in the THStack -->
<config includeemptystackhists='0'/>
<!-- # Turn on/off event manager -->
<!-- # EventManager enables us to only loop number of events once for multiple projections of the same measurements -->
<!-- # e.g. MiniBooNE CC1pi+ Q2 and MiniBooNE CC1pi+ Tmu would ordinarily require 2 reconfigures, but with this enabled it requires only one -->
<config EventManager='1'/>
<!-- # Event Directories -->
<!-- # Can setup default directories and use @EVENT_DIR/path to link to it -->
<config EVENT_DIR='/data2/stowell/NIWG/'/>
<config NEUT_EVENT_DIR='/data2/stowell/NIWG/neut/fit_samples_neut5.3.3/'/>
<config GENIE_EVENT_DIR='/data2/stowell/NIWG/genie/fit_samples_R.2.10.0/'/>
<config NUWRO_EVENT_DIR='/data2/stowell/NIWG/nuwro/fit_samples/'/>
<config GIBUU_EVENT_DIR='/data/GIBUU/DIR/'/>
<config SaveNuWroExtra='0' />
<!-- # In PrepareGENIE the reconstructed splines can be saved into the file -->
<config save_genie_splines='1'/>
<!-- # In InputHandler the option to regenerate NuWro flux/xsec plots is available -->
<!-- # Going to move this to its own app soon -->
<!-- # DEVEL CONFIG OPTION, don't touch! -->
<config CacheSize='0'/>
<!-- # ReWeighting Configuration Options -->
<!-- # ###################################################### -->
<!-- # Convert Dials in output statements using dial conversion card -->
<config convert_dials='0'/>
<!-- # Vetos can be used to specify RW dials NOT to be loaded in -->
<!-- # Useful if one specific one has an issue -->
<config FitWeight_fNIWGRW_veto=''/>
<config FitWeight_fNuwroRW_veto=''/>
<config FitWeight_fNeutRW_veto=''/>
<config FitWeight_fGenieRW_veto=''/>
<!-- # Output Options -->
<!-- # ###################################################### -->
<!-- # Save Nominal prediction with all rw engines at default -->
<config savenominal='0'/>
<!-- # Save prefit with values at starting values -->
<config saveprefit='0'/>
<!-- # Here's the full list of drawing options -->
<!-- DATA/MC/EVT/FINE/RATIO/MODES/SHAPE/WGHT/WEIGHTS/FLUX/XSEC/MASK/COV/INCOV/DECOMP/CANVPDG/CANVMC/SETTINGS'/>
<!-- #config drawopts DATA/MC/EVT/FINE/RATIO/MODES/SHAPE/RESIDUAL/MATRIX/FLUX/MASK/MAP -->
<!-- #config drawopts DATA/MC -->
<config drawopts='DATA/MC/EVT/FINE/RATIO/MODES/SHAPE/WEIGHTS/FLUX/XSEC/MASK/COV/INCOV/DECOMP/CANVPDG/CANVMC/SETTINGS'/>
<config InterpolateSigmaQ0Histogram='1' />
<config InterpolateSigmaQ0HistogramRes='100' />
<config InterpolateSigmaQ0HistogramThrow='1' />
<config InterpolateSigmaQ0HistogramNTHROWS='100000' />
<!-- # Save the shape scaling applied with option SHAPE into the main MC hist -->
<config saveshapescaling='0'/>
<config CorrectGENIEMECNorm='1'/>
<!-- # Set style of 1D output histograms -->
<config linecolour='1'/>
<config linestyle='1'/>
<config linewidth='1'/>
<!-- # For GenericFlux -->
<config isLiteMode='0'/>
<!-- # Statistical Options -->
<!-- # ###################################################### -->
<!-- # Add MC Statistical error to likelihoods -->
<config addmcerror='0'/>
<!-- # NUISMIN Configurations -->
<!-- # ###################################################### -->
<config MAXCALLS='1000000'/>
<config MAXITERATIONS='1000000'/>
<config TOLERANCE='0.001'/>
<!-- # Number of events required in low stats routines -->
<config LOWSTATEVENTS='25000'/>
<!-- # Error band generator configs -->
<!-- # ###################################################### -->
<!-- # For -f ErrorBands creates error bands for given measurements -->
<!-- # How many throws do we want (The higher the better precision) -->
<config error_throws='250'/>
<!-- # Are we throwing uniform or according to Gaussian? -->
<!-- # Only use uniform if wanting to study the limits of a dial. -->
<config error_uniform='0'/>
<config WriteSeperateStacks='1'/>
<!-- # Other Individual Case Configs -->
<!-- # ###################################################### -->
<!-- # Covariance throw options for fake data studies with MiniBooNE data. -->
<config thrown_covariance='FULL'/>
<config throw_mc_stat='0.0'/>
<config throw_diag_syst='0'/>
<config throw_corr_syst='0'/>
<config throw_mc_stat='0'/>
<!-- # Apply a shift to the muon momentum before calculation of Q2 -->
<config muon_momentum_shift='0.0'/>
<config muon_momentum_throw='0'/>
<!-- # MINERvA Specific Configs -->
<config MINERvA_CCinc_XSec_2DEavq3_nu_hadron_cut='0'/>
<config MINERvA_CCinc_XSec_2DEavq3_nu_useq3true='0'/>
<config Modes_split_PN_NN='0'/>
<config SignalReconfigures='0'/>
<!-- # SciBooNE specific -->
<config SciBarDensity='1.04'/>
<config SciBarRecoDist='12.0'/>
<config PenetratingMuonEnergy='1.4'/>
<config FlatEfficiency='1.0'/>
<config NumRangeSteps='50'/>
<config UseProton='true'/>
<config UseZackEff='false'/>
<config MINERvADensity='1.04'/>
<config MINERvARecoDist='10.0'/>
<!-- Different way of reweighting in GENIE -->
<config GENIEWeightEngine_CCRESMode="kModeMaMv"/>
<!--
+<config GENIEWeightEngine_CCRESMode="kModeMa"/>
+-->
<config GENIEWeightEngine_CCQEMode="kModeMa"/>
+<!--
+<config GENIEWeightEngine_CCQEMode="kModeMa"/> Normal MAQE Llewelyn-Smith
+<config GENIEWeightEngine_CCQEMode="kModeNormAndMaShape"/> Taking shape into account (don't really use this...)
+<config GENIEWeightEngine_CCQEMode="kModeZExp"/> Using z-expansion
-->
+
<!-- CCQE/2p2h/1pi Gaussian enhancement method -->
<!-- Apply tilt-shift weights or normal Gaussian parameters-->
<config Gaussian_Enhancement="Normal" />
<!--
<config Gaussian_Enhancement="Tilt-Shift" />
-->
<config NToyThrows='100000' />
<!-- Use NOvA Weights or not -->
<config NOvA_Weights="false" />
<!-- Tune name, for GENIE v3+ -->
<config GENIETune="G18_02a_00_000" />
<config GENIEEventGeneratorList="Default" />
<!--
<config GENIEXSecModelCCRES="genie::ReinSehgalRESPXSec" />
<config GENIEXSecModelCOH="genie::ReinSehgalCOHPiPXSec" />
<config GENIEXSecModelCCQE="genie::LwlynSmithQELCCPXSec" />
-->
<!--
<config GENIEXSecModelCCQE="genie::NievesQELCCPXSec" />
<config GENIEXSecModelCCRES="genie::BergerSehgalRESPXSec2014" />
<config GENIEXSecModelCOH="genie::BergerSehgalCOHPiPXSec2015" />
-->
<config UseShapeCovar="0" />
<config CC0piNBINS="156" />
</nuisance>
diff --git a/src/Reweight/GENIEWeightEngine.cxx b/src/Reweight/GENIEWeightEngine.cxx
index 31afebb..ceeed8f 100644
--- a/src/Reweight/GENIEWeightEngine.cxx
+++ b/src/Reweight/GENIEWeightEngine.cxx
@@ -1,277 +1,348 @@
#include "GENIEWeightEngine.h"
#ifdef __GENIE_EMP_MECRW_ENABLED
#include "ReWeight/GReWeightXSecEmpiricalMEC.h"
#endif
GENIEWeightEngine::GENIEWeightEngine(std::string name) {
#ifdef __GENIE_ENABLED__
// Setup the NEUT Reweight engien
fCalcName = name;
LOG(FIT) << "Setting up GENIE RW : " << fCalcName << std::endl;
// Create RW Engine suppressing cout
StopTalking();
fGenieRW = new genie::rew::GReWeight();
// Get List of Vetos (Just for debugging)
- std::string rw_engine_list =
- FitPar::Config().GetParS("FitWeight_fGenieRW_veto");
+ std::string rw_engine_list = FitPar::Config().GetParS("FitWeight_fGenieRW_veto");
bool xsec_ncel = rw_engine_list.find("xsec_ncel") == std::string::npos;
bool xsec_ccqe = rw_engine_list.find("xsec_ccqe") == std::string::npos;
bool xsec_coh = rw_engine_list.find("xsec_coh") == std::string::npos;
bool xsec_nnres = rw_engine_list.find("xsec_nonresbkg") == std::string::npos;
bool xsec_nudis = rw_engine_list.find("nuclear_dis") == std::string::npos;
- bool xsec_resdec =
- rw_engine_list.find("hadro_res_decay") == std::string::npos;
+ bool xsec_resdec = rw_engine_list.find("hadro_res_decay") == std::string::npos;
bool xsec_fzone = rw_engine_list.find("hadro_intranuke") == std::string::npos;
bool xsec_intra = rw_engine_list.find("hadro_fzone") == std::string::npos;
bool xsec_agky = rw_engine_list.find("hadro_agky") == std::string::npos;
bool xsec_qevec = rw_engine_list.find("xsec_ccqe_vec") == std::string::npos;
bool xsec_dis = rw_engine_list.find("xsec_dis") == std::string::npos;
bool xsec_nc = rw_engine_list.find("xsec_nc") == std::string::npos;
bool xsec_ccres = rw_engine_list.find("xsec_ccres") == std::string::npos;
bool xsec_ncres = rw_engine_list.find("xsec_ncres") == std::string::npos;
bool xsec_nucqe = rw_engine_list.find("nuclear_qe") == std::string::npos;
- bool xsec_qeaxial =
- rw_engine_list.find("xsec_ccqe_axial") == std::string::npos;
+ bool xsec_qeaxial = rw_engine_list.find("xsec_ccqe_axial") == std::string::npos;
#ifdef __GENIE_EMP_MECRW_ENABLED
bool xsec_empMEC = rw_engine_list.find("xsec_empMEC") == std::string::npos;
#endif
#ifndef GENIE_PRE_R3
genie::RunOpt* grunopt = genie::RunOpt::Instance();
grunopt->EnableBareXSecPreCalc(true);
grunopt->SetEventGeneratorList(Config::GetParS("GENIEEventGeneratorList"));
grunopt->SetTuneName(Config::GetParS("GENIETune"));
grunopt->BuildTune();
std::string genv = std::string(getenv("GENIE")) + "/config/Messenger_laconic.xml";
genie::utils::app_init::MesgThresholds(genv);
#endif
// Now actually add the RW Calcs
if (xsec_ncel)
fGenieRW->AdoptWghtCalc("xsec_ncel", new genie::rew::GReWeightNuXSecNCEL);
if (xsec_ccqe) {
fGenieRW->AdoptWghtCalc("xsec_ccqe", new genie::rew::GReWeightNuXSecCCQE);
// (dynamic_cast<GReWeightNuXSecCCQE*> (fGenieRW->WghtCalc("xsec_ccqe")))
// ->SetXSecModel( FitPar::Config().GetParS("GENIEXSecModelCCQE") );
}
#ifdef __GENIE_EMP_MECRW_ENABLED
if (xsec_empMEC) {
fGenieRW->AdoptWghtCalc("xsec_empMEC",
new genie::rew::GReWeightXSecEmpiricalMEC);
}
#endif
if (xsec_coh) {
fGenieRW->AdoptWghtCalc("xsec_coh", new genie::rew::GReWeightNuXSecCOH());
// (dynamic_cast<GReWeightNuXSecCOH*> (fGenieRW->WghtCalc("xsec_coh")))
// ->SetXSecModel( FitPar::Config().GetParS("GENIEXSecModelCOH") );
}
if (xsec_nnres)
fGenieRW->AdoptWghtCalc("xsec_nonresbkg",
new genie::rew::GReWeightNonResonanceBkg);
if (xsec_nudis)
fGenieRW->AdoptWghtCalc("nuclear_dis", new genie::rew::GReWeightDISNuclMod);
if (xsec_resdec)
fGenieRW->AdoptWghtCalc("hadro_res_decay",
new genie::rew::GReWeightResonanceDecay);
if (xsec_fzone)
fGenieRW->AdoptWghtCalc("hadro_fzone", new genie::rew::GReWeightFZone);
if (xsec_intra)
fGenieRW->AdoptWghtCalc("hadro_intranuke", new genie::rew::GReWeightINuke);
if (xsec_agky)
fGenieRW->AdoptWghtCalc("hadro_agky", new genie::rew::GReWeightAGKY);
if (xsec_qevec)
fGenieRW->AdoptWghtCalc("xsec_ccqe_vec",
new genie::rew::GReWeightNuXSecCCQEvec);
#if __GENIE_VERSION__ >= 212
if (xsec_qeaxial)
fGenieRW->AdoptWghtCalc("xsec_ccqe_axial",
new genie::rew::GReWeightNuXSecCCQEaxial);
#endif
if (xsec_dis)
fGenieRW->AdoptWghtCalc("xsec_dis", new genie::rew::GReWeightNuXSecDIS);
if (xsec_nc)
fGenieRW->AdoptWghtCalc("xsec_nc", new genie::rew::GReWeightNuXSecNC);
if (xsec_ccres) {
#if __GENIE_VERSION__ < 213
fGenieRW->AdoptWghtCalc("xsec_ccres", new genie::rew::GReWeightNuXSecCCRES);
#else
fGenieRW->AdoptWghtCalc(
"xsec_ccres",
new genie::rew::GReWeightNuXSecCCRES(
FitPar::Config().GetParS("GENIEXSecModelCCRES"), "Default"));
#endif
}
if (xsec_ncres)
fGenieRW->AdoptWghtCalc("xsec_ncres", new genie::rew::GReWeightNuXSecNCRES);
if (xsec_nucqe)
fGenieRW->AdoptWghtCalc("nuclear_qe", new genie::rew::GReWeightFGM);
- if (xsec_ccqe) {
- GReWeightNuXSecCCQE *rwccqe =
- dynamic_cast<GReWeightNuXSecCCQE *>(fGenieRW->WghtCalc("xsec_ccqe"));
- rwccqe->SetMode(GReWeightNuXSecCCQE::kModeMa);
- }
+#if __GENIE_VERSION__ >= 212
+ // Set the CCQE reweighting style
+ GReWeightNuXSecCCQE * rwccqe = dynamic_cast<GReWeightNuXSecCCQE *> (fGenieRW->WghtCalc("xsec_ccqe"));
+
+ // For MaCCQE reweighting
+ std::string ccqetype = FitPar::Config().GetParS("GENIEWeightEngine_CCQEMode");
+ if (ccqetype == "kModeMa") {
+ LOG(FIT) << "Setting GENIE ReWeight CCQE to kModeMa" << std::endl;
+ rwccqe->SetMode(GReWeightNuXSecCCQE::kModeMa);
+ } else if (ccqetype == "kModeNormAndMaShape") {
+ LOG(FIT) << "Setting GENIE ReWeight CCQE to kModeNormAndMaShape" << std::endl;
+ rwccqe->SetMode(GReWeightNuXSecCCQE::kModeNormAndMaShape);
+ // For z-expansion reweighting
+ } else if (ccqetype == "kModeZExp") {
+ LOG(FIT) << "Setting GENIE ReWeight CCQE to kModeZExp" << std::endl;
+ rwccqe->SetMode(GReWeightNuXSecCCQE::kModeZExp);
+ } else {
+ ERR(FTL) << "Did not find specified GENIE ReWeight CCQE mode" << std::endl;
+ ERR(FTL) << "You provided: " << ccqetype << " in parameters/config.xml" << std::endl;
+ throw;
+ }
+
+ // Check the UserPhysicsOptions too!
+ AlgConfigPool *Pool = genie::AlgConfigPool::Instance();
+ Registry* full = Pool->GlobalParameterList();
+ std::string name_ax = full->GetAlg("AxialFormFactorModel").name;
+ std::string config_ax = full->GetAlg("AxialFormFactorModel").config;
+ if (name_ax == "genie::DipoleAxialFormFactorModel" && ccqetype == "kModeZExp") {
+ ERR(FTL) << "Trying to run Z Expansion reweighting with Llewelyn-Smith model." << std::endl;
+ ERR(FTL) << "Please check your " << std::getenv("GENIE") << "/config/UserPhysicsOptions.xml to match generated" << std::endl;
+ ERR(FTL) << "You're telling me " << name_ax << "/" << config_ax << std::endl;
+ ERR(FTL) << "Also check your " << std::getenv("NUISANCE") << "/parameters/config.xml GENIEWeightEngine_CCQEMode: " << ccqetype << std::endl;
+ throw;
+ }
+
+ if (name_ax == "genie::ZExpAxialFormFactorModel" && ccqetype != "kModeZExp") {
+ ERR(FTL) << "Trying to run Llewelyn-Smith reweighting with Z Expansion model." << std::endl;
+ ERR(FTL) << "Please change your " << std::getenv("GENIE") << "/config/UserPhysicsOptions.xml to match generated" << std::endl;
+ ERR(FTL) << "You're telling me " << name_ax << "/" << config_ax << std::endl;
+ ERR(FTL) << "Also check your " << std::getenv("NUISANCE") << "/parameters/config.xml GENIEWeightEngine_CCQEMode: " << ccqetype << std::endl;
+ throw;
+ }
+
+ std::string name_qelcc = full->GetAlg("XSecModel@genie::EventGenerator/QEL-CC").name;
+ std::string config_qelcc = full->GetAlg("XSecModel@genie::EventGenerator/QEL-CC").config;
+ if (config_qelcc == "Default" && ccqetype == "kModeZExp") {
+ ERR(FTL) << "Trying to run Z Expansion reweighting with Llewelyn-Smith model." << std::endl;
+ ERR(FTL) << "Please change your " << std::getenv("GENIE") << "/config/UserPhysicsOptions.xml to match generated" << std::endl;
+ ERR(FTL) << "You're telling me " << name_qelcc << "/" << config_qelcc << std::endl;
+ ERR(FTL) << "Also check your " << std::getenv("NUISANCE") << "/parameters/config.xml GENIEWeightEngine_CCQEMode: " << ccqetype << std::endl;
+ throw;
+ }
+
+ if (config_qelcc == "ZExp" && ccqetype != "kModeZExp") {
+ ERR(FTL) << "Trying to run Llewelyn-Smith reweighting with Z Expansion model." << std::endl;
+ ERR(FTL) << "Please change your " << std::getenv("GENIE") << "/config/UserPhysicsOptions.xml to match generated" << std::endl;
+ ERR(FTL) << "You're telling me " << name_qelcc << "/" << config_qelcc << std::endl;
+ ERR(FTL) << "Also check your " << std::getenv("NUISANCE") << "/parameters/config.xml GENIEWeightEngine_CCQEMode: " << ccqetype << std::endl;
+ throw;
+ }
+#endif
if (xsec_ccres) {
// Default to include shape and normalization changes for CCRES (can be
// changed downstream if desired)
GReWeightNuXSecCCRES *rwccres =
- dynamic_cast<GReWeightNuXSecCCRES *>(fGenieRW->WghtCalc("xsec_ccres"));
+ dynamic_cast<GReWeightNuXSecCCRES *>(fGenieRW->WghtCalc("xsec_ccres"));
std::string marestype =
- FitPar::Config().GetParS("GENIEWeightEngine_CCRESMode");
+ FitPar::Config().GetParS("GENIEWeightEngine_CCRESMode");
if (!marestype.compare("kModeNormAndMaMvShape")) {
rwccres->SetMode(GReWeightNuXSecCCRES::kModeNormAndMaMvShape);
} else if (!marestype.compare("kModeMaMv")) {
rwccres->SetMode(GReWeightNuXSecCCRES::kModeMaMv);
} else {
THROW("Unkown MARES Mode in GENIE Weight Engine : " << marestype);
}
}
if (xsec_ncres) {
// Default to include shape and normalization changes for NCRES (can be
// changed downstream if desired)
GReWeightNuXSecNCRES *rwncres =
- dynamic_cast<GReWeightNuXSecNCRES *>(fGenieRW->WghtCalc("xsec_ncres"));
+ dynamic_cast<GReWeightNuXSecNCRES *>(fGenieRW->WghtCalc("xsec_ncres"));
rwncres->SetMode(GReWeightNuXSecNCRES::kModeMaMv);
}
if (xsec_dis) {
// Default to include shape and normalization changes for DIS (can be
// changed downstream if desired)
GReWeightNuXSecDIS *rwdis =
- dynamic_cast<GReWeightNuXSecDIS *>(fGenieRW->WghtCalc("xsec_dis"));
+ dynamic_cast<GReWeightNuXSecDIS *>(fGenieRW->WghtCalc("xsec_dis"));
rwdis->SetMode(GReWeightNuXSecDIS::kModeABCV12u);
// Set Abs Twk Config
fIsAbsTwk = (FitPar::Config().GetParB("setabstwk"));
}
// allow cout again
StartTalking();
#else
ERR(FTL) << "GENIE RW NOT ENABLED" << std::endl;
#endif
};
void GENIEWeightEngine::IncludeDial(std::string name, double startval) {
#ifdef __GENIE_ENABLED__
// Get First enum
int nuisenum = Reweight::ConvDial(name, kGENIE);
+ // Check ZExp sillyness in GENIE
+ // If ZExpansion parameters are used we need to set a different mode in GENIE ReWeight...
+ // GENIE doesn't have a setter either...
+#if __GENIE_VERSION__ >= 212
+ std::string ccqetype = FitPar::Config().GetParS("GENIEWeightEngine_CCQEMode");
+ if (ccqetype != "kModeZExp" && (name == "ZExpA1CCQE" || name == "ZExpA2CCQE" || name == "ZExpA3CCQE" || name == "ZExpA4CCQE")) {
+ ERR(FTL) << "Found a Z-expansion parameter in GENIE although the GENIE ReWeighting engine is set to use Llewelyn-Smith and MaQE!" << std::endl;
+ ERR(FTL) << "Change your GENIE UserPhysicsOptions.xml in " << std::getenv("GENIE") << "/config/UserPhysicsOptions.xml to match requirements" << std::endl;
+ throw;
+ }
+
+ if ((ccqetype != "kModeMa" && ccqetype != "kModeMaNormAndMaShape") && (name == "MaCCQE")) {
+ ERR(FTL) << "Found MaCCQE parameter in GENIE although the GENIE ReWeighting engine is set to not use this!" << std::endl;
+ ERR(FTL) << "Change your GENIE UserPhysicsOptions.xml in " << std::getenv("GENIE") << "/config/UserPhysicsOptions.xml to match requirements" << std::endl;
+ throw;
+ }
+#endif
+
// Setup Maps
fEnumIndex[nuisenum]; // = std::vector<size_t>(0);
fNameIndex[name]; // = std::vector<size_t>(0);
// Split by commas
std::vector<std::string> allnames = GeneralUtils::ParseToStr(name, ",");
for (uint i = 0; i < allnames.size(); i++) {
std::string singlename = allnames[i];
// Get RW
genie::rew::GSyst_t rwsyst = GSyst::FromString(singlename);
// Fill Maps
int index = fValues.size();
fValues.push_back(0.0);
fGENIESysts.push_back(rwsyst);
// Initialize dial
- std::cout << "Registering " << singlename << " from " << name << std::endl;
+ LOG(FIT) << "Registering " << singlename << " from " << name << std::endl;
fGenieRW->Systematics().Init(fGENIESysts[index]);
// If Absolute
if (fIsAbsTwk) {
GSystUncertainty::Instance()->SetUncertainty(rwsyst, 1.0, 1.0);
}
// Setup index
fEnumIndex[nuisenum].push_back(index);
fNameIndex[name].push_back(index);
}
// Set Value if given
if (startval != -999.9) {
SetDialValue(nuisenum, startval);
}
#endif
};
void GENIEWeightEngine::SetDialValue(int nuisenum, double val) {
#ifdef __GENIE_ENABLED__
std::vector<size_t> indices = fEnumIndex[nuisenum];
for (uint i = 0; i < indices.size(); i++) {
fValues[indices[i]] = val;
fGenieRW->Systematics().Set(fGENIESysts[indices[i]], val);
}
#endif
}
void GENIEWeightEngine::SetDialValue(std::string name, double val) {
#ifdef __GENIE_ENABLED__
std::vector<size_t> indices = fNameIndex[name];
for (uint i = 0; i < indices.size(); i++) {
fValues[indices[i]] = val;
fGenieRW->Systematics().Set(fGENIESysts[indices[i]], val);
}
#endif
}
void GENIEWeightEngine::Reconfigure(bool silent) {
#ifdef __GENIE_ENABLED__
// Hush now...
if (silent)
StopTalking();
// Reconf
fGenieRW->Reconfigure();
fGenieRW->Print();
// Shout again
if (silent)
StartTalking();
#endif
}
double GENIEWeightEngine::CalcWeight(BaseFitEvt *evt) {
double rw_weight = 1.0;
#ifdef __GENIE_ENABLED__
// Make nom weight
if (!evt) {
THROW("evt not found : " << evt);
}
// Skip Non GENIE
if (evt->fType != kGENIE)
return 1.0;
if (!(evt->genie_event)) {
THROW("evt->genie_event not found!" << evt->genie_event);
}
if (!(evt->genie_event->event)) {
THROW("evt->genie_event->event GHepRecord not found!"
- << (evt->genie_event->event));
+ << (evt->genie_event->event));
}
if (!fGenieRW) {
THROW("GENIE RW Not Found!" << fGenieRW);
}
rw_weight = fGenieRW->CalcWeight(*(evt->genie_event->event));
// std::cout << "Returning GENIE Weight for electron scattering = " <<
// rw_weight << std::endl;
#endif
// Return rw_weight
return rw_weight;
}
diff --git a/src/Reweight/GENIEWeightEngine.h b/src/Reweight/GENIEWeightEngine.h
index e0a7f7e..85866ce 100644
--- a/src/Reweight/GENIEWeightEngine.h
+++ b/src/Reweight/GENIEWeightEngine.h
@@ -1,88 +1,89 @@
#ifndef WEIGHT_ENGINE_GENIE_H
#define WEIGHT_ENGINE_GENIE_H
#include "FitLogger.h"
#ifdef __GENIE_ENABLED__
#ifdef GENIE_PRE_R3
+#include "Algorithm/AlgConfigPool.h"
#include "EVGCore/EventRecord.h"
#include "GHEP/GHepRecord.h"
#include "ReWeight/GSyst.h"
#include "ReWeight/GSystUncertainty.h"
#include "Ntuple/NtpMCEventRecord.h"
#include "ReWeight/GReWeight.h"
#include "ReWeight/GReWeightAGKY.h"
#include "ReWeight/GReWeightDISNuclMod.h"
#include "ReWeight/GReWeightFGM.h"
#include "ReWeight/GReWeightFZone.h"
#include "ReWeight/GReWeightINuke.h"
#include "ReWeight/GReWeightNonResonanceBkg.h"
#include "ReWeight/GReWeightNuXSecCCQE.h"
#include "ReWeight/GReWeightNuXSecCCQEvec.h"
#include "ReWeight/GReWeightNuXSecCCRES.h"
#include "ReWeight/GReWeightNuXSecCOH.h"
#include "ReWeight/GReWeightNuXSecDIS.h"
#include "ReWeight/GReWeightNuXSecNC.h"
#include "ReWeight/GReWeightNuXSecNCEL.h"
#include "ReWeight/GReWeightNuXSecNCRES.h"
#include "ReWeight/GReWeightResonanceDecay.h"
#if __GENIE_VERSION__ >= 212
#include "ReWeight/GReWeightNuXSecCCQEaxial.h"
#endif
#else
#include "Framework/EventGen/EventRecord.h"
#include "Framework/GHEP/GHepRecord.h"
#include "Framework/GHEP/GHepParticle.h"
#include "Framework/Ntuple/NtpMCEventRecord.h"
#include "Framework/Utils/RunOpt.h"
#include "Framework/Utils/AppInit.h"
#include "RwFramework/GSyst.h"
#include "RwFramework/GSystUncertainty.h"
#include "RwFramework/GReWeight.h"
#include "RwCalculators/GReWeightAGKY.h"
#include "RwCalculators/GReWeightDISNuclMod.h"
#include "RwCalculators/GReWeightFGM.h"
#include "RwCalculators/GReWeightFZone.h"
#include "RwCalculators/GReWeightINuke.h"
#include "RwCalculators/GReWeightNonResonanceBkg.h"
#include "RwCalculators/GReWeightNuXSecCCQE.h"
#include "RwCalculators/GReWeightNuXSecCCQEvec.h"
#include "RwCalculators/GReWeightNuXSecCCRES.h"
#include "RwCalculators/GReWeightNuXSecCOH.h"
#include "RwCalculators/GReWeightNuXSecDIS.h"
#include "RwCalculators/GReWeightNuXSecNC.h"
#include "RwCalculators/GReWeightNuXSecNCEL.h"
#include "RwCalculators/GReWeightNuXSecNCRES.h"
#include "RwCalculators/GReWeightResonanceDecay.h"
#include "RwCalculators/GReWeightNuXSecCCQEaxial.h"
#endif
using namespace genie;
using namespace genie::rew;
#endif
#include "GeneratorUtils.h"
#include "WeightEngineBase.h"
#include "FitWeight.h"
class GENIEWeightEngine : public WeightEngineBase {
public:
GENIEWeightEngine(std::string name);
~GENIEWeightEngine() {};
void IncludeDial(std::string name, double startval);
void SetDialValue(int rwenum, double val);
void SetDialValue(std::string name, double val);
void Reconfigure(bool silent = false);
double CalcWeight(BaseFitEvt* evt);
inline bool NeedsEventReWeight() { return true; };
#ifdef __GENIE_ENABLED__
std::vector<genie::rew::GSyst_t> fGENIESysts;
genie::rew::GReWeight* fGenieRW; //!< Genie RW Object
#endif
};
#endif

File Metadata

Mime Type
text/x-diff
Expires
Sat, Dec 21, 1:12 PM (1 d, 1 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4017282
Default Alt Text
(25 KB)

Event Timeline