Public Member Functions | Private Member Functions | Private Attributes | List of all members
cafrwgt::CAFReweight Class Reference

#include "/cvmfs/nova-development.opensciencegrid.org/novasoft/releases/N20-11-25/CAFReweight/CAFReweight.h"

Public Member Functions

 CAFReweight ()
 
 ~CAFReweight ()
 
double CalcWeight (const std::vector< caf::SRGenerator > &gen, unsigned int i=0)
 
std::vector< double > CalcWeightVector (const std::vector< caf::SRGenerator > &gen)
 
void WgtMultiNeutrino ()
 
void WgtBestNeutrino ()
 
void AddReweightValue (rwgt::ReweightLabel_t rLabel, double value, unsigned int i=0)
 
rwgt::GENIEReweight * GenieReweight (int i)
 
void Configure (int i=-1)
 
int NReweights () const
 

Private Member Functions

genie::EventRecord RetrieveGHEP (const caf::SRGenerator &gen)
 

Private Attributes

bool fWgtMultiNus
 For accounting for multiple neutrino interactions in the slice. More...
 
std::vector< rwgt::GENIEReweight * > fGRwgt
 GENIE reweighting. More...
 

Detailed Description

Definition at line 24 of file CAFReweight.h.

Constructor & Destructor Documentation

cafrwgt::CAFReweight::CAFReweight ( )

Definition at line 54 of file CAFReweight.cxx.

54  :
55  fWgtMultiNus(false) {
56  }
bool fWgtMultiNus
For accounting for multiple neutrino interactions in the slice.
Definition: CAFReweight.h:46
cafrwgt::CAFReweight::~CAFReweight ( )

Definition at line 58 of file CAFReweight.cxx.

References fGRwgt, and MECModelEnuComparisons::i.

58  {
59  for(unsigned int i = 0; i < fGRwgt.size(); ++i) delete fGRwgt[i];
60  }
std::vector< rwgt::GENIEReweight * > fGRwgt
GENIE reweighting.
Definition: CAFReweight.h:47

Member Function Documentation

void cafrwgt::CAFReweight::AddReweightValue ( rwgt::ReweightLabel_t  rLabel,
double  value,
unsigned int  i = 0 
)

Definition at line 66 of file CAFReweight.cxx.

References fGRwgt, and MECModelEnuComparisons::i.

Referenced by CAFReweightExample().

66  {
67  if(fGRwgt.size() <= i) {
68  fGRwgt.push_back(new rwgt::GENIEReweight);
69  }
70 
71  fGRwgt[i]->AddReweightValue(rLabel, value);
72  }
const XML_Char int const XML_Char * value
Definition: expat.h:331
std::vector< rwgt::GENIEReweight * > fGRwgt
GENIE reweighting.
Definition: CAFReweight.h:47
double cafrwgt::CAFReweight::CalcWeight ( const std::vector< caf::SRGenerator > &  gen,
unsigned int  i = 0 
)

Definition at line 133 of file CAFReweight.cxx.

References fGRwgt, fWgtMultiNus, MECModelEnuComparisons::i, RetrieveGHEP(), and wgt.

Referenced by CAFReweightExample().

133  {
134  double wgt = 1.0;
135  if(genvec.size()==0) {
136  return wgt;
137  }
138  else {
139  if(fWgtMultiNus==false) {
140  genie::EventRecord evr = this->RetrieveGHEP(genvec[0]);
141  return fGRwgt[i]->CalculateWeight(evr);
142  }
143  else {
144  wgt = 0;
145  for(unsigned int i = 0; i < genvec.size(); i++) {
146  genie::EventRecord evr = this->RetrieveGHEP(genvec[i]);
147  double nu_wgt = genvec[i].fVisEfraction;
148  wgt += nu_wgt*(fGRwgt[i]->CalculateWeight(evr));
149  }
150  return wgt;
151  }
152  }
153  }
bool fWgtMultiNus
For accounting for multiple neutrino interactions in the slice.
Definition: CAFReweight.h:46
const ana::Var wgt
Generated Event Record. It is a GHepRecord object that can accept / be visited by EventRecordVisitorI...
Definition: EventRecord.h:38
genie::EventRecord RetrieveGHEP(const caf::SRGenerator &gen)
std::vector< rwgt::GENIEReweight * > fGRwgt
GENIE reweighting.
Definition: CAFReweight.h:47
std::vector< double > cafrwgt::CAFReweight::CalcWeightVector ( const std::vector< caf::SRGenerator > &  gen)

Definition at line 85 of file CAFReweight.cxx.

References fGRwgt, fWgtMultiNus, MECModelEnuComparisons::i, calib::j, RetrieveGHEP(), and wgt.

85  {
86 
87  std::vector<double> wgts;
88 
89  //Initialize wgts vector.
90  for(unsigned int i = 0; i < fGRwgt.size(); i++) {
91  wgts.push_back(1.0);
92  }
93 
94  //Catch the case where the generator vector in the caf is empty
95  if(genvec.size()==0) {
96  return wgts;
97  }
98  else { //If the genvec is not empty calcuate actual weights
99  //This is the case where you are not doing multi nu weighting
100  if(fWgtMultiNus==false) {
101  //Store the genie event record for the circumstances that
102  //the multiple weights with different configurations are being
103  //requested.
104  genie::EventRecord evr = this->RetrieveGHEP(genvec[0]);
105  for(unsigned int i = 0; i < fGRwgt.size(); i++) { //loop over the vector of configurations. But use the same event record
106  double wgt = fGRwgt[i]->CalculateWeight(evr);
107  wgts[i] = wgt;
108  }
109  return wgts;
110  }
111  else { //Now for the multi nu case
112  std::vector<genie::EventRecord> evr_vec; //make vector of event records.
113  for(unsigned int i = 0; i < genvec.size(); i++) {
114  genie::EventRecord evr = this->RetrieveGHEP(genvec[i]);
115  evr_vec.push_back(evr);
116  }
117  //Loop over the configurations
118  for(unsigned int i = 0; i < fGRwgt.size(); i++) {
119  double wgt = 0.0;
120  //Loop over the event records
121  for(unsigned int j = 0; j < evr_vec.size(); j++) {
122  double nu_wgt = genvec[j].fVisEfraction; //This should have the same size as evr_vec
123  wgt += nu_wgt*(fGRwgt[i]->CalculateWeight(evr_vec[j]));
124  }
125  wgts[i]=wgt;
126  }
127  return wgts; //Return the vector of weights.
128  }
129  }
130  }
bool fWgtMultiNus
For accounting for multiple neutrino interactions in the slice.
Definition: CAFReweight.h:46
const double j
Definition: BetheBloch.cxx:29
const ana::Var wgt
Generated Event Record. It is a GHepRecord object that can accept / be visited by EventRecordVisitorI...
Definition: EventRecord.h:38
genie::EventRecord RetrieveGHEP(const caf::SRGenerator &gen)
std::vector< rwgt::GENIEReweight * > fGRwgt
GENIE reweighting.
Definition: CAFReweight.h:47
void cafrwgt::CAFReweight::Configure ( int  i = -1)

Definition at line 74 of file CAFReweight.cxx.

References fGRwgt, MECModelEnuComparisons::i, and calib::j.

Referenced by CAFReweightExample().

74  {
75  if(i >= 0){
76  fGRwgt[i]->Configure();
77  }
78  else{
79  for(unsigned int j = 0; j < fGRwgt.size(); j++) {
80  fGRwgt[j]->Configure();
81  }
82  }
83  }
const double j
Definition: BetheBloch.cxx:29
std::vector< rwgt::GENIEReweight * > fGRwgt
GENIE reweighting.
Definition: CAFReweight.h:47
rwgt::GENIEReweight * cafrwgt::CAFReweight::GenieReweight ( int  i)

Definition at line 62 of file CAFReweight.cxx.

References fGRwgt, and MECModelEnuComparisons::i.

62  {
63  return fGRwgt[i];
64  }
std::vector< rwgt::GENIEReweight * > fGRwgt
GENIE reweighting.
Definition: CAFReweight.h:47
int cafrwgt::CAFReweight::NReweights ( ) const
inline

Definition at line 41 of file CAFReweight.h.

41 {return fGRwgt.size();}
std::vector< rwgt::GENIEReweight * > fGRwgt
GENIE reweighting.
Definition: CAFReweight.h:47
genie::EventRecord cafrwgt::CAFReweight::RetrieveGHEP ( const caf::SRGenerator &  gen)
private

Definition at line 155 of file CAFReweight.cxx.

References genie::GHepRecord::AttachSummary(), genie::GHepRecord::HitNucleon(), genie::GHepRecord::HitNucleonPosition(), MECModelEnuComparisons::i, genie::kPSNull, LOG_WARN, make_static_page::new, genie::GHepParticle::P4(), genie::GHepRecord::Probe(), genie::ProcessInfo::Set(), genie::XclsTag::SetCharm(), genie::GHepRecord::SetDiffXSec(), genie::Interaction::SetExclTag(), genie::Kinematics::SetFSLeptonP4(), genie::Kinematics::SetHadSystP4(), genie::Target::SetHitNucP4(), genie::Target::SetHitNucPdg(), genie::Target::SetHitQrkPdg(), genie::Target::SetHitSeaQrk(), genie::Interaction::SetKine(), genie::XclsTag::SetNNucleons(), genie::XclsTag::SetNPions(), genie::GHepRecord::SetProbability(), genie::InitialState::SetProbeP4(), genie::Kinematics::SetQ2(), genie::Kinematics::Setq2(), genie::XclsTag::SetResonance(), genie::Kinematics::Sett(), genie::InitialState::SetTgtP4(), genie::GHepRecord::SetVertex(), genie::Kinematics::SetW(), genie::GHepRecord::SetWeight(), genie::Kinematics::Setx(), genie::GHepRecord::SetXSec(), genie::Kinematics::Sety(), make_syst_table_plots::space, genie::GHepRecord::TargetNucleus(), genie::GHepRecord::TargetNucleusPosition(), genie::InitialState::TgtPtr(), and genie::XclsTag::UnsetCharm().

Referenced by CalcWeight(), and CalcWeightVector().

155  {
156 
157  genie::EventRecord newEvent;
158  newEvent.SetWeight(gen.fweight);
159  newEvent.SetProbability(gen.fprobability);
160  newEvent.SetXSec(gen.fXsec);
161 
162 #ifndef SETDIFFXSEC_1ARG
163  genie::KinePhaseSpace_t space = genie::kPSNull; // kPSQ2fE; // ????
164  // dsig/dQ2, dsig/dQ2dW, dsig/dxdy ...
165 
166  newEvent.SetDiffXSec(gen.fDiffXsec,space);
167 
168  // TODO: we don't know currently know what to use ...
169  // for now just to get things to compile ... I (Nate) needs to look at this
170  static int nmsg = 10;
171  if ( nmsg > 0 ) {
172  LOG_WARN("CAFReweight") << "RetrieveGHEP(simb::MCTruth,simb::GTruth) is not correctly setting KinePhaseSpace_t in SetDiffXSec()\n"
173  << "At the time of the conversion to R-2_8_0 (2013-05-01) this is not critical\n"
174  << "But it should be fixed";
175  --nmsg;
176  if ( nmsg == 0 ) LOG_WARN("CAFReweight") << "... last of such messages";
177  }
178  //assert(0);
179 #else
180  newEvent.SetDiffXSec(gen.fDiffXsec);
181 #endif
182 
183  TLorentzVector vtx = gen.fVertex;
184  newEvent.SetVertex(vtx);
185 
186  int nParticles = gen.fParticleList.size();
187  // Reserve space for the particles we're about to add
188  newEvent.Expand(nParticles);
189  for(int i = 0; i < nParticles; i++) {
190  int gmid = gen.fParticleList[i].fPdgCode;
191  int gmmo = gen.fParticleList[i].fMother;
192  int gmfd = gen.fParticleList[i].fFirstDaughter;
193  int gmld = gen.fParticleList[i].fLastDaughter;
194 
195  int gmri = gen.fParticleList[i].fRescatterCode;
196 
197  double gmvx = gen.fParticleList[i].fNuclVtx.X();
198  double gmvy = gen.fParticleList[i].fNuclVtx.Y();
199  double gmvz = gen.fParticleList[i].fNuclVtx.Z();
200  double gmvt = gen.fParticleList[i].fNuclVtx.T();
201 
202  double gmpx = gen.fParticleList[i].fP4.Px();
203  double gmpy = gen.fParticleList[i].fP4.Py();
204  double gmpz = gen.fParticleList[i].fP4.Pz();
205  double gme = gen.fParticleList[i].fP4.E();
206 
207  genie::GHepStatus_t gmst = (genie::GHepStatus_t)gen.fParticleList[i].fStatus;
208 
209  TVector3 polz = gen.fParticleList[i].fPolz;
210 
211  // Normally one would create a GHepParticle and add it to the event with
212  // AddParticle().
213  //
214  // But AddParticle() calls UpdateDaughterLists() every time, which is
215  // slow. And unnecessary to be called even once, because we set the
216  // daughters by ourselves above.
217  //
218  // Instead, create the particle in-place in the EventRecord by hand. This
219  // is significantly faster in tests.
220 
221  genie::GHepParticle* gpart = new (newEvent[i]) genie::GHepParticle(gmid, gmst, gmmo, -1, gmfd, gmld, gmpx, gmpy, gmpz, gme, gmvx, gmvy, gmvz, gmvt);
222  gpart->SetRescatterCode(gmri);
223 
224  if(polz.x() !=0 || polz.y() !=0 || polz.z() !=0) {
225  gpart->SetPolarization(polz);
226  }
227  }
228 
229  genie::ProcessInfo proc_info;
230  genie::ScatteringType_t gscty = (genie::ScatteringType_t)gen.fGscatter;
232 
233  proc_info.Set(gscty, ginty);
234 
235  genie::XclsTag gxt;
236 
237  //Set Exclusive Final State particle numbers
238  genie::Resonance_t gres = (genie::Resonance_t)gen.fResNum;
239  gxt.SetResonance(gres);
240  gxt.SetNPions(gen.fNumPiPlus, gen.fNumPi0, gen.fNumPiMinus);
241  gxt.SetNNucleons(gen.fNumProton, gen.fNumNeutron);
242 
243  if(gen.fIsCharm) {
244  gxt.SetCharm(0);
245  }
246  else {
247  gxt.UnsetCharm();
248  }
249 
250  //Set the GENIE kinematic info
251  genie::Kinematics gkin;
252  gkin.Setx(gen.fgX, true);
253  gkin.Sety(gen.fgY, true);
254  gkin.Sett(gen.fgT, true);
255  gkin.SetW(gen.fgW, true);
256  gkin.SetQ2(gen.fgQ2, true);
257  gkin.Setq2(gen.fgq2, true);
258  TLorentzVector lep_p4 = gen.fParticleList[gen.fFSLeptonIndx].fP4;
259 
260  gkin.SetFSLeptonP4(lep_p4.Px(), lep_p4.Py(), lep_p4.Pz(), lep_p4.E());
261  gkin.SetHadSystP4(gen.fFShadSystP4.Px(), gen.fFShadSystP4.Py(), gen.fFShadSystP4.Pz(), gen.fFShadSystP4.E());
262 
263  //Set the GENIE final state interaction info
264  genie::InitialState ginstate(gen.ftgtPDG, gen.fProbePDG);
265 
266  // Set up the target, except for PDG set above
267  genie::Target* target123 = ginstate.TgtPtr();
268 
269  //int Z = gtruth.ftgtZ;
270  //int A = gtruth.ftgtA;
271 
272  //target123->SetId(Z,A);
273 
274  //int pdg_code = pdg::IonPdgCode(A, Z);
275  //TParticlePDG * p = PDGLibrary::Instance()->Find(pdg_code);
276 
277  //mf::LogWarning("GENIEReweight") << "Setting Target Z and A";
278  //mf::LogWarning("GENIEReweight") << "Saved PDG: " << gtruth.ftgtPDG;
279  //mf::LogWarning("GENIEReweight") << "Target PDG: " << target123->Pdg();
280  target123->SetHitNucPdg(gen.fHitNucleon);
281  target123->SetHitQrkPdg(gen.fHitQuark);
282  target123->SetHitSeaQrk(gen.fIsSeaQuark);
283 
284  if(newEvent.HitNucleonPosition()> 0) {
285  target123->SetHitNucP4(*newEvent.HitNucleon()->P4());
286  }
287  else {
288  target123->SetHitNucP4(TLorentzVector(0.,0.,0.,0.));
289  }
290 
291  // Set probe information
292  ginstate.SetProbeP4(*newEvent.Probe()->P4());
293  if(newEvent.TargetNucleusPosition()> 0) {
294  ginstate.SetTgtP4(*newEvent.TargetNucleus()->P4());
295  }
296  else {
297  ginstate.SetTgtP4(TLorentzVector(0.,0.,0.,0.));
298  }
299 
300  // Build the final Interaction object
301  genie::Interaction* p_gint = new genie::Interaction(ginstate, proc_info);
302  p_gint->SetKine(gkin);
303  p_gint->SetExclTag(gxt);
304  newEvent.AttachSummary(p_gint);
305 
306  return newEvent;
307 
308  }
virtual void SetXSec(double xsec)
Definition: GHepRecord.h:133
void SetNPions(int npi_plus, int npi_0, int npi_minus)
Definition: XclsTag.cxx:97
virtual void SetWeight(double wght)
Definition: GHepRecord.h:131
void SetQ2(double Q2, bool selected=false)
Definition: Kinematics.cxx:265
const TLorentzVector * P4(void) const
Definition: GHepParticle.h:79
void Setq2(double q2, bool selected=false)
Definition: Kinematics.cxx:277
enum genie::EGHepStatus GHepStatus_t
virtual void SetProbability(double prob)
Definition: GHepRecord.h:132
void SetHitNucP4(const TLorentzVector &p4)
Definition: Target.cxx:206
Generated/set kinematical variables for an event.
Definition: Kinematics.h:40
void SetHitQrkPdg(int pdgc)
Definition: Target.cxx:201
virtual int HitNucleonPosition(void) const
Definition: GHepRecord.cxx:454
void SetKine(const Kinematics &kine)
enum genie::EKinePhaseSpace KinePhaseSpace_t
#define LOG_WARN(stream)
Definition: Messenger.h:134
enum genie::EResonance Resonance_t
void SetResonance(Resonance_t res)
Definition: XclsTag.cxx:123
void Set(ScatteringType_t sc_type, InteractionType_t int_type)
virtual void AttachSummary(Interaction *interaction)
Definition: GHepRecord.cxx:143
Contains minimal information for tagging exclusive processes.
Definition: XclsTag.h:37
void SetCharm(int charm_pdgc=0)
Definition: XclsTag.cxx:68
virtual GHepParticle * Probe(void) const
Definition: GHepRecord.cxx:321
Summary information for an interaction.
Definition: Interaction.h:56
void SetFSLeptonP4(const TLorentzVector &p4)
Definition: Kinematics.cxx:307
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
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
virtual GHepParticle * TargetNucleus(void) const
Definition: GHepRecord.cxx:330
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
enum genie::EScatteringType ScatteringType_t
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
virtual void SetVertex(double x, double y, double z, double t)
Definition: GHepRecord.cxx:863
void UnsetCharm(void)
Definition: XclsTag.cxx:74
void SetHadSystP4(const TLorentzVector &p4)
Definition: Kinematics.cxx:317
void SetHitSeaQrk(bool tf)
Definition: Target.cxx:212
enum genie::EInteractionType InteractionType_t
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
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
void cafrwgt::CAFReweight::WgtBestNeutrino ( )
inline

Definition at line 36 of file CAFReweight.h.

References Configure(), and MECModelEnuComparisons::i.

36 {fWgtMultiNus=false;}
bool fWgtMultiNus
For accounting for multiple neutrino interactions in the slice.
Definition: CAFReweight.h:46
void cafrwgt::CAFReweight::WgtMultiNeutrino ( )
inline

Definition at line 35 of file CAFReweight.h.

35 {fWgtMultiNus=true;}
bool fWgtMultiNus
For accounting for multiple neutrino interactions in the slice.
Definition: CAFReweight.h:46

Member Data Documentation

std::vector<rwgt::GENIEReweight*> cafrwgt::CAFReweight::fGRwgt
private

GENIE reweighting.

Definition at line 47 of file CAFReweight.h.

Referenced by AddReweightValue(), CalcWeight(), CalcWeightVector(), Configure(), GenieReweight(), and ~CAFReweight().

bool cafrwgt::CAFReweight::fWgtMultiNus
private

For accounting for multiple neutrino interactions in the slice.

Definition at line 46 of file CAFReweight.h.

Referenced by CalcWeight(), and CalcWeightVector().


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