Public Member Functions | Private Member Functions | List of all members
rwgt::NuReweight Class Reference

#include <NuReweight.h>

Inheritance diagram for rwgt::NuReweight:
rwgt::GENIEReweight

Public Member Functions

 NuReweight ()
 <constructor More...
 
 ~NuReweight ()
 
double CalcWeight (const simb::MCTruth &truth, const simb::GTruth &gtruth) const
 
- Public Member Functions inherited from rwgt::GENIEReweight
 GENIEReweight ()
 <constructor More...
 
 ~GENIEReweight ()
 Set the nominal values for the reweight parameters. More...
 
void AddReweightValue (ReweightLabel_t rLabel, double value)
 Change a reweight parameter. If it hasn't been added yet add it. More...
 
void ChangeParameterValue (ReweightLabel_t rLabel, double value)
 Configure the weight calculators. More...
 
double NominalParameterValue (ReweightLabel_t rLabel)
 Return the configured value of the given parameter. More...
 
double ReweightParameterValue (ReweightLabel_t rLabel)
 Add reweight parameters to the list. More...
 
genie::rew::GReWeight * WeightCalculator ()
 
void Configure ()
 Reconfigure the weight calculators. More...
 
void Reconfigure ()
 Simple Configuration functions for configuring a single weight calculator. More...
 
void ReweightNCEL (double ma, double eta)
 Simple Configurtion of the CCQE axial weight calculator. More...
 
void ReweightQEMA (double ma)
 Simple Configuration of the CCQE vector weight calculator. More...
 
void ReweightQEVec (double mv)
 
void ReweightQEZExp (double norm, double a1, double a2, double a3, double a4)
 Simple Configuration of the CC Resonance weight calculator. More...
 
void ReweightResGanged (double ma, double mv=0.0)
 Simple Configuration of the Coherant weight calculator. More...
 
void ReweightCCRes (double ma, double mv=0.0)
 Simple Configurtion of the NC Resonance weight calculator. More...
 
void ReweightNCRes (double ma, double mv=0.0)
 Simple Configuration of the NC and CC Resonance weight calculator with the axial mass parameter for NC/CC ganged together. More...
 
void ReweightCoh (double ma, double r0)
 Simple Configuration of the Non-Resonance Background weight calculator. More...
 
void ReweightNonResRvp1pi (double sigma)
 Simple Configuration of the Non-Resonance Background weight calculator. More...
 
void ReweightNonResRvbarp1pi (double sigma)
 Simple Configuration of the Non-Resonance Background weight calculator. Here it is being configured for v+p and vbar + n (2 pi) type interactions. More...
 
void ReweightNonResRvp2pi (double sigma)
 Simple Configuration of the Non-Resonance Background weight calculator. More...
 
void ReweightNonResRvbarp2pi (double sigma)
 Simple Configuration of the Resonance decay model weight calculator. More...
 
void ReweightResDecay (double gamma, double eta, double theta)
 Simple Configuration of the Total NC cross section. More...
 
void ReweightNC (double norm)
 Simple Configuration of the DIS FF model weight calculator. More...
 
void ReweightDIS (double aht, double bht, double cv1u, double cv2u)
 Simple Configuration of the DIS nuclear model. More...
 
void ReweightDISnucl (bool mode)
 Simple Configuration of the DIS AGKY hadronization model. More...
 
void ReweightAGKY (double xF, double pT)
 Simple Configuration of the Intranuke Nuclear model. More...
 
void ReweightFormZone (double sigma)
 Simple Configuration of the Fermigas model reweight calculator. More...
 
void ReweightFGM (double kF, double sf)
 End of Simple Reweight Configurations. More...
 
void ReweightIntraNuke (ReweightLabel_t name, double sigma)
 Simple Configuration of the Formation Zone reweight calculator. More...
 
void ReweightIntraNuke (int name, double sigma)
 
void MaQEshape ()
 
void MaQErate ()
 
void CCRESshape ()
 
void CCRESrate ()
 
void NCRESshape ()
 
void NCRESrate ()
 
void DIS_BYshape ()
 
void DIS_BYrate ()
 
void UseSigmaDef ()
 
void UseStandardDef ()
 
void SetNominalValues ()
 Return the nominal value for the given parameter. More...
 
double CalculateSigma (ReweightLabel_t label, double value)
 Calculate the weights. More...
 
double CalculateWeight (const genie::EventRecord &evr) const
 
void ConfigureNCEL ()
 Configure the MaQE weight calculator. More...
 
void ConfigureQEMA ()
 Configure the QE vector FF weight calculator. More...
 
void ConfigureQEVec ()
 Configure the CCRES calculator. More...
 
void ConfigureCCRes ()
 Configure the NCRES calculator. More...
 
void ConfigureNCRes ()
 Configure the ResBkg (kno) weight calculator. More...
 
void ConfigureResBkg ()
 Configure the ResDecay weight calculator. More...
 
void ConfgureResDecay ()
 Configure the NC weight calculator. More...
 
void ConfigureNC ()
 Configure the DIS (Bodek-Yang) weight calculator. More...
 
void ConfigureDIS ()
 Configure the Coherant model weight calculator. More...
 
void ConfigureCoh ()
 Configure the hadronization (AGKY) weight calculator. More...
 
void ConfigureAGKY ()
 Configure the DIS nuclear model weight calculator. More...
 
void ConfigureDISNucMod ()
 Configure the FG model weight calculator. More...
 
void ConfigureFGM ()
 Configure the Formation Zone weight calculator. More...
 
void ConfigureFZone ()
 Configure the intranuke weight calculator. More...
 
void ConfigureINuke ()
 configure the weight parameters being used More...
 
void ConfigureParameters ()
 

Private Member Functions

genie::EventRecord RetrieveGHEP (const simb::MCTruth &truth, const simb::GTruth &gtruth) const
 

Additional Inherited Members

- Protected Attributes inherited from rwgt::GENIEReweight
bool fReweightNCEL
 
bool fReweightQEMA
 
bool fReweightQEVec
 
bool fReweightCCRes
 
bool fReweightNCRes
 
bool fReweightResBkg
 
bool fReweightResDecay
 
bool fReweightNC
 
bool fReweightDIS
 
bool fReweightCoh
 
bool fReweightAGKY
 
bool fReweightDISNucMod
 
bool fReweightFGM
 
bool fReweightFZone
 
bool fReweightINuke
 
bool fReweightZexp
 
bool fReweightMEC
 
bool fMaQEshape
 
bool fMaCCResShape
 
bool fMaNCResShape
 
bool fDISshape
 
bool fUseSigmaDef
 
std::vector< int > fReWgtParameterName
 
std::vector< double > fReWgtParameterValue
 
std::map< int, double > fNominalParameters
 
genie::rew::GReWeight * fWcalc
 

Detailed Description

Definition at line 15 of file NuReweight.h.

Constructor & Destructor Documentation

rwgt::NuReweight::NuReweight ( )

<constructor

destructor

Definition at line 118 of file NuReweight.cxx.

rwgt::NuReweight::~NuReweight ( )

Definition at line 123 of file NuReweight.cxx.

123  {
124  // Don't delete fWcalc here. The GENIEReweight parent class handles it.
125  }

Member Function Documentation

double rwgt::NuReweight::CalcWeight ( const simb::MCTruth truth,
const simb::GTruth gtruth 
) const

Definition at line 127 of file NuReweight.cxx.

127  {
128  genie::EventRecord evr = this->RetrieveGHEP(truth, gtruth);
129  double wgt = this->CalculateWeight(evr);
130  //mf::LogVerbatim("GENIEReweight") << "New Event Weight is: " << wgt;
131  return wgt;
132  }
double CalculateWeight(const genie::EventRecord &evr) const
Generated Event Record. It is a GHepRecord object that can accept / be visited by EventRecordVisitorI...
Definition: EventRecord.h:38
genie::EventRecord RetrieveGHEP(const simb::MCTruth &truth, const simb::GTruth &gtruth) const
Definition: NuReweight.cxx:134
genie::EventRecord rwgt::NuReweight::RetrieveGHEP ( const simb::MCTruth truth,
const simb::GTruth gtruth 
) const
private

Definition at line 134 of file NuReweight.cxx.

134  {
135  genie::EventRecord newEvent;
136  newEvent.SetWeight(gtruth.fweight);
137  newEvent.SetProbability(gtruth.fprobability);
138  newEvent.SetXSec(gtruth.fXsec);
139 
142  newEvent.SetDiffXSec(gtruth.fDiffXsec,space);
143 
144  TLorentzVector vtx = gtruth.fVertex;
145  newEvent.SetVertex(vtx);
146 
147  for(int i = 0; i < truth.NParticles(); i++) {
148  simb::MCParticle mcpart = truth.GetParticle(i);
149 
150  int gmid = mcpart.PdgCode();
152  int gmmo = mcpart.Mother();
153  int ndaughters = mcpart.NumberDaughters();
154  //find the track ID of the first and last daughter particles
155  int fdtrkid = 0;
156  int ldtrkid = 0;
157  if(ndaughters !=0) {
158  fdtrkid = mcpart.Daughter(0);
159  if(ndaughters == 1) {
160  ldtrkid = 1;
161  }
162  else if(ndaughters >1) {
163  fdtrkid = mcpart.Daughter(ndaughters-1);
164  }
165  }
166  int gmfd = -1;
167  int gmld = -1;
168  //Genie uses the index in the particle array to reference the daughter particles.
169  //MCTruth keeps the particles in the same order so use the track ID to find the proper index.
170  for(int j = 0; j < truth.NParticles(); j++) {
171  simb::MCParticle temp = truth.GetParticle(i);
172  if(temp.TrackId() == fdtrkid) {
173  gmfd = j;
174  }
175  if(temp.TrackId() == ldtrkid) {
176  gmld = j;
177  }
178  }
179 
180  double gmpx = mcpart.Px(0);
181  double gmpy = mcpart.Py(0);
182  double gmpz = mcpart.Pz(0);
183  double gme = mcpart.E(0);
184  double gmvx = mcpart.Gvx();
185  double gmvy = mcpart.Gvy();
186  double gmvz = mcpart.Gvz();
187  double gmvt = mcpart.Gvt();
188  int gmri = mcpart.Rescatter();
189 
190  genie::GHepParticle gpart(gmid, gmst, gmmo, -1, gmfd, gmld, gmpx, gmpy, gmpz, gme, gmvx, gmvy, gmvz, gmvt);
191  gpart.SetRescatterCode(gmri);
192  TVector3 polz = mcpart.Polarization();
193  if(polz.x() !=0 || polz.y() !=0 || polz.z() !=0) {
194  gpart.SetPolarization(polz);
195  }
196  newEvent.AddParticle(gpart);
197 
198  }
199 
200  genie::ProcessInfo proc_info;
203 
204  proc_info.Set(gscty,ginty);
205 
206  genie::XclsTag gxt;
207 
208  //Set Exclusive Final State particle numbers
210  gxt.SetResonance(gres);
211  gxt.SetDecayMode(gtruth.fDecayMode);
212  gxt.SetNPions(gtruth.fNumPiPlus, gtruth.fNumPi0, gtruth.fNumPiMinus);
213  gxt.SetNNucleons(gtruth.fNumProton, gtruth.fNumNeutron);
214 
215  if (gtruth.fIsCharm) {
216  gxt.SetCharm(gtruth.fCharmHadronPdg);
217  } else {
218  gxt.UnsetCharm();
219  }
220 
221  if (gtruth.fIsStrange) {
222  gxt.SetStrange(gtruth.fStrangeHadronPdg);
223  } else {
224  gxt.UnsetStrange();
225  }
226 
227  //Set the GENIE kinematic info
228  genie::Kinematics gkin;
229  gkin.Setx(gtruth.fgX, true);
230  gkin.Sety(gtruth.fgY, true);
231  gkin.Sett(gtruth.fgT, true);
232  gkin.SetW(gtruth.fgW, true);
233  gkin.SetQ2(gtruth.fgQ2, true);
234  gkin.Setq2(gtruth.fgq2, true);
235  simb::MCNeutrino nu = truth.GetNeutrino();
236  simb::MCParticle lep = nu.Lepton();
237  gkin.SetFSLeptonP4(lep.Px(), lep.Py(), lep.Pz(), lep.E());
238  gkin.SetHadSystP4(gtruth.fFShadSystP4.Px(), gtruth.fFShadSystP4.Py(), gtruth.fFShadSystP4.Pz(), gtruth.fFShadSystP4.E());
239 
240  //Set the GENIE final state interaction info
241  genie::Interaction * p_gint = new genie::Interaction;
242  genie::InitialState * p_ginstate = p_gint->InitStatePtr();
243 
244  p_ginstate->SetPdgs(gtruth.ftgtPDG, gtruth.fProbePDG);
245 
246  int targetNucleon = nu.HitNuc();
247  int struckQuark = nu.HitQuark();
248 
249  genie::Target* target123 = p_ginstate->TgtPtr();
250  target123->SetHitNucPdg(targetNucleon);
251  target123->SetHitQrkPdg(struckQuark);
252  target123->SetHitSeaQrk(gtruth.fIsSeaQuark);
253  target123->SetHitNucPosition(gtruth.fHitNucPos);
254 
255  if(newEvent.HitNucleonPosition()> 0) {
256  genie::GHepParticle * hitnucleon = newEvent.HitNucleon();
257  std::unique_ptr<TLorentzVector> p4hitnucleon(hitnucleon->GetP4());
258  target123->SetHitNucP4(*p4hitnucleon);
259  }
260  else {
261  TLorentzVector dummy(0.,0.,0.,0.);
262  target123->SetHitNucP4(dummy);
263  }
264 
265  genie::GHepParticle * probe = newEvent.Probe();
266  std::unique_ptr<TLorentzVector> p4probe(probe->GetP4());
267  p_ginstate->SetProbeP4(*p4probe);
268  if(newEvent.TargetNucleusPosition()> 0) {
270  std::unique_ptr<TLorentzVector> p4target(target->GetP4());
271  p_ginstate->SetTgtP4(*p4target);
272  } else {
273  TLorentzVector dummy(0.,0.,0.,0.);
274  p_ginstate->SetTgtP4(dummy);
275  }
276  p_gint->SetProcInfo(proc_info);
277  p_gint->SetKine(gkin);
278  p_gint->SetExclTag(gxt);
279  newEvent.AttachSummary(p_gint);
280 
281  /*
282  //For temporary debugging purposes
283  genie::Interaction *inter = newEvent.Summary();
284  const genie::InitialState &initState = inter->InitState();
285  const genie::Target &tgt = initState.Tgt();
286 
287  std::cout << "TargetPDG as Recorded: " << gtruth.ftgtPDG << std::endl;
288  std::cout << "TargetZ as Recorded: " << gtruth.ftgtZ << std::endl;
289  std::cout << "TargetA as Recorded: " << gtruth.ftgtA << std::endl;
290  std::cout << "TargetPDG as Recreated: " << tgt.Pdg() << std::endl;
291  std::cout << "TargetZ as Recreated: " << tgt.Z() << std::endl;
292  std::cout << "TargetA as Recreated: " << tgt.A() << std::endl;
293  */
294 
295  return newEvent;
296 
297  }
double fgW
Definition: GTruth.h:63
double E(const int i=0) const
Definition: MCParticle.h:232
int fGint
interaction code
Definition: GTruth.h:56
const TVector3 & Polarization() const
Definition: MCParticle.h:213
void SetTgtP4(const TLorentzVector &P4)
virtual void SetXSec(double xsec)
Definition: GHepRecord.h:133
void SetNPions(int npi_plus, int npi_0, int npi_minus)
Definition: XclsTag.cxx:97
int PdgCode() const
Definition: MCParticle.h:211
virtual void SetWeight(double wght)
Definition: GHepRecord.h:131
double fgq2
Definition: GTruth.h:62
const simb::MCNeutrino & GetNeutrino() const
Definition: MCTruth.h:77
double fgX
Definition: GTruth.h:65
double Py(const int i=0) const
Definition: MCParticle.h:230
double Gvz() const
Definition: MCParticle.h:247
void SetProbeP4(const TLorentzVector &P4)
int HitQuark() const
Definition: MCNeutrino.h:153
void SetQ2(double Q2, bool selected=false)
Definition: Kinematics.cxx:265
void Setq2(double q2, bool selected=false)
Definition: Kinematics.cxx:277
int Mother() const
Definition: MCParticle.h:212
enum genie::EGHepStatus GHepStatus_t
virtual void SetProbability(double prob)
Definition: GHepRecord.h:132
void UnsetStrange(void)
Definition: XclsTag.cxx:91
int fGPhaseSpace
phase space system of DiffXSec
Definition: GTruth.h:32
double Px(const int i=0) const
Definition: MCParticle.h:229
int fNumNeutron
number of neutrons after reaction, before FSI
Definition: GTruth.h:76
void SetHitNucP4(const TLorentzVector &p4)
Definition: Target.cxx:206
Generated/set kinematical variables for an event.
Definition: Kinematics.h:40
int HitNuc() const
Definition: MCNeutrino.h:152
void SetHitQrkPdg(int pdgc)
Definition: Target.cxx:201
virtual int HitNucleonPosition(void) const
Definition: GHepRecord.cxx:454
double fXsec
cross section of interaction
Definition: GTruth.h:30
void SetHitNucPosition(double r)
Definition: Target.cxx:227
int fNumPiPlus
number of pi pluses after reaction, before FSI
Definition: GTruth.h:78
int StatusCode() const
Definition: MCParticle.h:210
int fNumPiMinus
number of pi minuses after reaction, before FSI
Definition: GTruth.h:79
double Gvx() const
Definition: MCParticle.h:245
int NParticles() const
Definition: MCTruth.h:75
void SetKine(const Kinematics &kine)
double Gvy() const
Definition: MCParticle.h:246
enum genie::EKinePhaseSpace KinePhaseSpace_t
enum genie::EResonance Resonance_t
void SetResonance(Resonance_t res)
Definition: XclsTag.cxx:137
void Set(ScatteringType_t sc_type, InteractionType_t int_type)
bool fIsStrange
strange production // added version 13
Definition: GTruth.h:73
int NumberDaughters() const
Definition: MCParticle.h:216
virtual void AttachSummary(Interaction *interaction)
Definition: GHepRecord.cxx:143
Contains minimal information for tagging exclusive processes.
Definition: XclsTag.h:39
int TrackId() const
Definition: MCParticle.h:209
int Daughter(const int i) const
Definition: MCParticle.cxx:112
void SetCharm(int charm_pdgc=0)
Definition: XclsTag.cxx:68
virtual GHepParticle * Probe(void) const
Definition: GHepRecord.cxx:321
int fResNum
resonance number
Definition: GTruth.h:80
Summary information for an interaction.
Definition: Interaction.h:55
int fNumProton
number of protons after reaction, before FSI
Definition: GTruth.h:75
int fDecayMode
Definition: GTruth.h:81
double fprobability
interaction probability
Definition: GTruth.h:29
const simb::MCParticle & Lepton() const
Definition: MCNeutrino.h:147
int fProbePDG
Definition: GTruth.h:39
void SetFSLeptonP4(const TLorentzVector &p4)
Definition: Kinematics.cxx:307
int fGscatter
neutrino scattering code
Definition: GTruth.h:55
void SetStrange(int strange_pdgc=0)
Definition: XclsTag.cxx:85
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:43
int fNumPi0
number of pi0 after reaction, before FSI
Definition: GTruth.h:77
int fCharmHadronPdg
Definition: GTruth.h:72
TLorentzVector * GetP4(void) const
A Neutrino Interaction Target. Is a transparent encapsulation of quite different physical systems suc...
Definition: Target.h:41
void Sett(double t, bool selected=false)
Definition: Kinematics.cxx:301
bool fIsCharm
did the interaction produce a charmed hadron?
Definition: GTruth.h:71
double fweight
event interaction weight (genie internal)
Definition: GTruth.h:28
virtual GHepParticle * TargetNucleus(void) const
Definition: GHepRecord.cxx:330
void SetPdgs(int tgt_pdgc, int probe_pdgc)
void SetDecayMode(int decay_mode)
Definition: XclsTag.cxx:142
const simb::MCParticle & GetParticle(int i) const
Definition: MCTruth.h:76
Generated Event Record. It is a GHepRecord object that can accept / be visited by EventRecordVisitorI...
Definition: EventRecord.h:38
void Setx(double x, bool selected=false)
Definition: Kinematics.cxx:241
void SetExclTag(const XclsTag &xcls)
void SetW(double W, bool selected=false)
Definition: Kinematics.cxx:289
void SetNNucleons(int np, int nn)
Definition: XclsTag.cxx:104
double fHitNucPos
Definition: GTruth.h:52
enum genie::EScatteringType ScatteringType_t
int ftgtPDG
PDG of Target Nucleus, nucleon only if free.
Definition: GTruth.h:47
double fgQ2
Definition: GTruth.h:61
void SetHitNucPdg(int pdgc)
Definition: Target.cxx:188
void Sety(double y, bool selected=false)
Definition: Kinematics.cxx:253
virtual GHepParticle * HitNucleon(void) const
Definition: GHepRecord.cxx:350
Target * TgtPtr(void) const
Definition: InitialState.h:68
virtual void SetVertex(double x, double y, double z, double t)
Definition: GHepRecord.cxx:863
double Gvt() const
Definition: MCParticle.h:248
TLorentzVector fFShadSystP4
generated final state hadronic system (LAB frame)
Definition: GTruth.h:68
double Pz(const int i=0) const
Definition: MCParticle.h:231
int Rescatter() const
Definition: MCParticle.h:251
void UnsetCharm(void)
Definition: XclsTag.cxx:74
void SetHadSystP4(const TLorentzVector &p4)
Definition: Kinematics.cxx:317
cet::LibraryManager dummy("noplugin")
InitialState * InitStatePtr(void) const
Definition: Interaction.h:73
virtual void AddParticle(const GHepParticle &p)
Definition: GHepRecord.cxx:535
double fgT
Definition: GTruth.h:64
void SetHitSeaQrk(bool tf)
Definition: Target.cxx:212
Event generator information.
Definition: MCNeutrino.h:18
bool fIsSeaQuark
Definition: GTruth.h:50
TLorentzVector fVertex
Definition: GTruth.h:26
enum genie::EInteractionType InteractionType_t
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
double fgY
Definition: GTruth.h:66
double fDiffXsec
differential cross section of interaction
Definition: GTruth.h:31
void SetProcInfo(const ProcessInfo &proc)
virtual int TargetNucleusPosition(void) const
Definition: GHepRecord.cxx:406
virtual void SetDiffXSec(double xsec, KinePhaseSpace_t ps)
Definition: GHepRecord.h:134
Initial State information.
Definition: InitialState.h:49
int fStrangeHadronPdg
Definition: GTruth.h:74

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