RESKinematicsGenerator.cxx
Go to the documentation of this file.
1 //____________________________________________________________________________
2 /*
3  Copyright (c) 2003-2019, The GENIE Collaboration
4  For the full text of the license visit http://copyright.genie-mc.org
5  or see $GENIE/LICENSE
6 
7  Author: Costas Andreopoulos <costas.andreopoulos \at stfc.ac.uk>
8  University of Liverpool & STFC Rutherford Appleton Lab
9 
10  For the class documentation see the corresponding header file.
11 
12  Important revisions after version 2.0.0 :
13  @ Mar 03, 2009 - CA
14  Moved into the new RES package from its previous location (EVGModules)
15  @ Sep 21, 2009 - CA
16  In ComputeMaxXSec() protect against NW=1 when calculating the step size dW
17  as pointed out by Robert.
18  @ Oct 20, 2009 - CA
19  For charged lepton scattering do a more brute force kinematical selection
20  and do not use the importance sampling envelope used for v scattering
21  (was very inefficient)
22  @ Feb 06, 2013 - CA
23  When the value of the differential cross-section for the selected kinematics
24  is set to the event, set the corresponding KinePhaseSpace_t value too.
25  @ Jul 26, 2018 - IL (Afroditi Papadopoulou, Adi Ashkenazi - Massachusetts Institute of Technology)
26  Included importance sampling envelop both for neutrino and electron scattering
27 */
28 //____________________________________________________________________________
29 
30 #include <TMath.h>
31 #include <TF2.h>
32 #include <TROOT.h>
33 
51 
52 using namespace genie;
53 using namespace genie::controls;
54 using namespace genie::utils;
55 
56 //___________________________________________________________________________
58 KineGeneratorWithCache("genie::RESKinematicsGenerator")
59 {
60  fEnvelope = 0;
61 }
62 //___________________________________________________________________________
64 KineGeneratorWithCache("genie::RESKinematicsGenerator", config)
65 {
66  fEnvelope = 0;
67 }
68 //___________________________________________________________________________
70 {
71  if(fEnvelope) delete fEnvelope;
72 }
73 //___________________________________________________________________________
75 {
76  if(fGenerateUniformly) {
77  LOG("RESKinematics", pNOTICE)
78  << "Generating kinematics uniformly over the allowed phase space";
79  }
80 
81  //-- Get the random number generators
83 
84  //-- Access cross section algorithm for running thread
86  const EventGeneratorI * evg = rtinfo->RunningThread();
87  fXSecModel = evg->CrossSectionAlg();
88 
89  //-- Get the interaction from the GHEP record
90  Interaction * interaction = evrec->Summary();
91  interaction->SetBit(kISkipProcessChk);
92 
93  //-- EM or CC/NC process (different importance sampling methods)
94  bool is_em = interaction->ProcInfo().IsEM();
95 
96  //-- Compute the W limits
97  // (the physically allowed W's, unless an external cut is imposed)
98  const KPhaseSpace & kps = interaction->PhaseSpace();
99  Range1D_t W = kps.Limits(kKVW);
100  //costas says this is bad if(W.max>1.7) W.max=1.7;
101 //assert(W.min>=0. && W.min<W.max);
102 
103  if(W.max <=0 || W.min>=W.max) {
104  LOG("RESKinematics", pWARN) << "No available phase space";
105  evrec->EventFlags()->SetBitNumber(kKineGenErr, true);
107  exception.SetReason("No available phase space");
108  exception.SwitchOnFastForward();
109  throw exception;
110  }
111 
112  const InitialState & init_state = interaction -> InitState();
113  double E = init_state.ProbeE(kRfHitNucRest);
114  // double M = init_state.Tgt().HitNucP4().M();
115  // double ml = interaction->FSPrimLepton()->Mass();
116 
117  //-- For the subsequent kinematic selection with the rejection method:
118  // Calculate the max differential cross section or retrieve it from the
119  // cache. Throw an exception and quit the evg thread if a non-positive
120  // value is found.
121  // If the kinematics are generated uniformly over the allowed phase
122  // space the max xsec is irrelevant
123  double xsec_max = (fGenerateUniformly) ? -1 : this->MaxXSec(evrec);
124 
125  //-- Try to select a valid W, Q2 pair using the rejection method
126  double dW = W.max - W.min;
127  double xsec = -1;
128 
129  unsigned int iter = 0;
130  bool accept = false;
131  while(1) {
132  iter++;
133  if(iter > kRjMaxIterations) {
134  LOG("RESKinematics", pWARN)
135  << "*** Could not select a valid (W,Q^2) pair after "
136  << iter << " iterations";
137  evrec->EventFlags()->SetBitNumber(kKineGenErr, true);
139  exception.SetReason("Couldn't select kinematics");
140  exception.SwitchOnFastForward();
141  throw exception;
142  }
143 
144  double gW = 0; // current hadronic invariant mass
145  double gQ2 = 0; // current momentum transfer
146  double gQD2 = 0; // tranformed Q2 to take out dipole form
147 
148  if(fGenerateUniformly) {
149  //-- Generate a W uniformly in the kinematically allowed range.
150  // For the generated W, compute the Q2 range and generate a value
151  // uniformly over that range
152  gW = W.min + dW * rnd->RndKine().Rndm();
153  Range1D_t Q2 = kps.Q2Lim_W();
154  if(Q2.max<=0. || Q2.min>=Q2.max) continue;
155  gQ2 = Q2.min + (Q2.max-Q2.min) * rnd->RndKine().Rndm();
156 
157  interaction->SetBit(kISkipKinematicChk);
158 
159  } else {
160 
161 
162  // > neutrino scattering
163  // Selecting unweighted event kinematics using an importance sampling
164  // method. Q2 with be transformed to QD2 to take out the dipole form.
165  // An importance sampling envelope will be constructed for W.
166  // first pass, configure the sampling envelope
167  if(iter==1) {
168  LOG("RESKinematics", pINFO) << "Initializing the sampling envelope";
169  if(!fEnvelope) {
170  LOG("RESKinematics", pFATAL) << "Null sampling envelope!";
171  exit(1);
172  }
173  interaction->KinePtr()->SetW(W.min);
174  Range1D_t Q2 = kps.Q2Lim_W();
175  double Q2min = -99.;
176  if (is_em) { Q2min = Q2.min + kASmallNum; }
177  else { Q2min = 0 + kASmallNum; }
178  double Q2max = Q2.max - kASmallNum;
179 
180  // In unweighted mode - use transform that takes out the dipole form
181  double QD2min = utils::kinematics::Q2toQD2(Q2max);
182  double QD2max = utils::kinematics::Q2toQD2(Q2min);
183 
184 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
185  LOG("RESKinematics", pDEBUG)
186  << "Q^2: [" << Q2min << ", " << Q2max << "] => "
187  << "QD^2: [" << QD2min << ", " << QD2max << "]";
188 #endif
189  double mR, gR;
190  if(!interaction->ExclTag().KnownResonance()) {
191  mR = 1.2;
192  gR = 0.6;
193  } else {
194  Resonance_t res = interaction->ExclTag().Resonance();
195  mR = res::Mass(res);
196  gR = (E>mR) ? 0.220 : 0.400;
197  }
198 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
199  LOG("RESKinematics", pDEBUG)
200  << "(m,g) = (" << mR << ", " << gR
201  << "), max(xsec,W) = (" << xsec_max << ", " << W.max << ")";
202 #endif
203  fEnvelope->SetRange(QD2min,W.min,QD2max,W.max); // range
204  fEnvelope->SetParameter(0, mR); // resonance mass
205  fEnvelope->SetParameter(1, gR); // resonance width
206  fEnvelope->SetParameter(2, xsec_max); // max differential xsec
207  fEnvelope->SetParameter(3, W.max); // kinematically allowed Wmax
208  }// first pass
209 
210  // Generate W,QD2 using the 2-D envelope as PDF
211  fEnvelope->GetRandom2(gQD2,gW);
212 
213  // QD2 -> Q2
214  gQ2 = utils::kinematics::QD2toQ2(gQD2);
215  } // uniformly over phase space?
216 
217  LOG("RESKinematics", pINFO) << "Trying: W = " << gW << ", Q2 = " << gQ2;
218 
219  //-- Set kinematics for current trial
220  interaction->KinePtr()->SetW(gW);
221  interaction->KinePtr()->SetQ2(gQ2);
222 
223  //-- Computing cross section for the current kinematics
224  xsec = fXSecModel->XSec(interaction, kPSWQ2fE);
225 
226  //-- Decide whether to accept the current kinematics
227  if(!fGenerateUniformly) {
228 
229  // unified neutrino / electron scattering
230  double max = fEnvelope->Eval(gQD2, gW);
231  double t = max * rnd->RndKine().Rndm();
232  double J = kinematics::Jacobian(interaction,kPSWQ2fE,kPSWQD2fE);
233 
234  this->AssertXSecLimits(interaction, xsec, max);
235 
236 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
237  LOG("RESKinematics", pDEBUG)
238  << "xsec= " << xsec << ", J= " << J << ", Rnd= " << t;
239 #endif
240  accept = (t < J*xsec);
241  } // charged lepton or neutrino scattering?
242  else {
243  accept = (xsec>0);
244  } // uniformly over phase space
245 
246  //-- If the generated kinematics are accepted, finish-up module's job
247  if(accept) {
248  LOG("RESKinematics", pINFO)
249  << "Selected: W = " << gW << ", Q2 = " << gQ2;
250  // reset 'trust' bits
251  interaction->ResetBit(kISkipProcessChk);
252  interaction->ResetBit(kISkipKinematicChk);
253 
254  // compute x,y for selected W,Q2
255  // note: hit nucleon can be off the mass-shell
256  double gx=-1, gy=-1;
257  double M = init_state.Tgt().HitNucP4().M();
258  //double M = init_state.Tgt().HitNucMass();
259  kinematics::WQ2toXY(E,M,gW,gQ2,gx,gy);
260 
261  // set the cross section for the selected kinematics
262  evrec->SetDiffXSec(xsec,kPSWQ2fE);
263 
264  // for uniform kinematics, compute an event weight as
265  // wght = (phase space volume)*(differential xsec)/(event total xsec)
266  if(fGenerateUniformly) {
267  double vol = kinematics::PhaseSpaceVolume(interaction,kPSWQ2fE);
268  double totxsec = evrec->XSec();
269  double wght = (vol/totxsec)*xsec;
270  LOG("RESKinematics", pNOTICE) << "Kinematics wght = "<< wght;
271 
272  // apply computed weight to the current event weight
273  wght *= evrec->Weight();
274  LOG("RESKinematics", pNOTICE) << "Current event wght = " << wght;
275  evrec->SetWeight(wght);
276  }
277 
278  // lock selected kinematics & clear running values
279  interaction->KinePtr()->SetQ2(gQ2, true);
280  interaction->KinePtr()->SetW (gW, true);
281  interaction->KinePtr()->Setx (gx, true);
282  interaction->KinePtr()->Sety (gy, true);
283  interaction->KinePtr()->ClearRunningValues();
284 
285  return;
286  } // accept
287  } // iterations
288 }
289 //___________________________________________________________________________
291 {
292  Algorithm::Configure(config);
293  this->LoadConfig();
294 }
295 //____________________________________________________________________________
297 {
298  Algorithm::Configure(config);
299  this->LoadConfig();
300 }
301 //____________________________________________________________________________
303 {
304  // Safety factor for the maximum differential cross section
305  this->GetParamDef("MaxXSec-SafetyFactor", fSafetyFactor, 1.25);
306 
307  // Minimum energy for which max xsec would be cached, forcing explicit
308  // calculation for lower eneries
309  this->GetParamDef("Cache-MinEnergy", fEMin, 0.5);
310 
311  // Load Wcut used in DIS/RES join scheme
312  this->GetParam("Wcut", fWcut);
313 
314  // Maximum allowed fractional cross section deviation from maxim cross
315  // section used in rejection method
316  this->GetParamDef("MaxXSec-DiffTolerance", fMaxXSecDiffTolerance, 999999.);
318 
319  // Generate kinematics uniformly over allowed phase space and compute
320  // an event weight?
321  this->GetParamDef("UniformOverPhaseSpace", fGenerateUniformly, false);
322 
323  // Envelope employed when importance sampling is used
324  // (initialize with dummy range)
325  if(fEnvelope) delete fEnvelope;
326  fEnvelope = new TF2("res-envelope",
328  // stop ROOT from deleting this object of its own volition
329  gROOT->GetListOfFunctions()->Remove(fEnvelope);
330 }
331 //____________________________________________________________________________
333  const Interaction * interaction) const
334 {
335 // Computes the maximum differential cross section in the requested phase
336 // space. This method overloads KineGeneratorWithCache::ComputeMaxXSec
337 // method and the value is cached at a circular cache branch for retrieval
338 // during subsequent event generation.
339 // The computed max differential cross section does not need to be the exact
340 // maximum. The number used in the rejection method will be scaled up by a
341 // safety factor. But this needs to be fast - do not use a very fine grid.
342 
343  double max_xsec = 0.;
344 
345  const InitialState & init_state = interaction -> InitState();
346  double E = init_state.ProbeE(kRfHitNucRest);
347  bool is_em = interaction->ProcInfo().IsEM();
349 
350 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
351  LOG("RESKinematics", pDEBUG) << "Scanning phase space for E= " << E;
352 #endif
353  //bool scan1d = (E>1.0);
354  bool scan1d = false;
355 
356  double md;
357  if(!interaction->ExclTag().KnownResonance()) md=1.23;
358  else {
359  Resonance_t res = interaction->ExclTag().Resonance();
360  md=res::Mass(res);
361  }
362 
363  if(scan1d) {
364  // ** 1-D Scan
365  //
366 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
367  LOG("RESKinematics", pDEBUG)
368  << "Will search for max{xsec} along W(=MRes) = " << md;
369 #endif
370  // Set W around the value where d^2xsec/dWdQ^2 peaks
371  interaction->KinePtr()->SetW(md);
372 
373  const KPhaseSpace & kps = interaction->PhaseSpace();
374  Range1D_t rQ2 = kps.Q2Lim_W();
375  if( rQ2.max < Q2Thres || rQ2.min <=0 ) return 0.;
376 
377  int NQ2 = 25;
378  int NQ2b = 5;
379  double logQ2min = TMath::Log(rQ2.min+kASmallNum);
380  double logQ2max = TMath::Log(rQ2.max-kASmallNum);
381  double dlogQ2 = (logQ2max - logQ2min) /(NQ2-1);
382 
383  double xseclast = -1;
384  bool increasing = true;
385 
386  for(int iq2=0; iq2<NQ2; iq2++) {
387  double Q2 = TMath::Exp(logQ2min + iq2 * dlogQ2);
388  interaction->KinePtr()->SetQ2(Q2);
389  double xsec = fXSecModel->XSec(interaction, kPSWQ2fE);
390 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
391  LOG("RESKinematics", pDEBUG)
392  << "xsec(W= " << md << ", Q2= " << Q2 << ") = " << xsec;
393 #endif
394  max_xsec = TMath::Max(xsec, max_xsec);
395  increasing = xsec-xseclast>=0;
396  xseclast=xsec;
397 
398  // once the cross section stops increasing, I reduce the step size and
399  // step backwards a little bit to handle cases that the max cross section
400  // is grossly underestimated (very peaky distribution & large step)
401  if(!increasing) {
402  dlogQ2/=NQ2b;
403  for(int iq2b=0; iq2b<NQ2b; iq2b++) {
404  Q2 = TMath::Exp(TMath::Log(Q2) - dlogQ2);
405  if(Q2 < rQ2.min) continue;
406  interaction->KinePtr()->SetQ2(Q2);
407  xsec = fXSecModel->XSec(interaction, kPSWQ2fE);
408 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
409  LOG("RESKinematics", pDEBUG)
410  << "xsec(W= " << md << ", Q2= " << Q2 << ") = " << xsec;
411 #endif
412  max_xsec = TMath::Max(xsec, max_xsec);
413  }
414  break;
415  }
416  } // Q2
417 
418  } else {
419 
420  // ** 2-D Scan
421  //
422  const KPhaseSpace & kps = interaction->PhaseSpace();
423  Range1D_t rW = kps.WLim();
424 
425  int NW = 20;
426  double Wmin = rW.min + kASmallNum;
427  double Wmax = rW.max - kASmallNum;
428 
429  Wmax = TMath::Min(Wmax,fWcut);
430 
431  Wmin = TMath::Max(Wmin, md-.3);
432  Wmax = TMath::Min(Wmax, md+.3);
433 
434  if(Wmax-Wmin<0.05) { NW=1; Wmin=Wmax; }
435 
436  double dW = (NW>1) ? (Wmax-Wmin)/(NW-1) : 0.;
437 
438  for(int iw=0; iw<NW; iw++) {
439  double W = Wmin + iw*dW;
440  interaction->KinePtr()->SetW(W);
441 
442  int NQ2 = 25;
443  int NQ2b = 4;
444 
445  Range1D_t rQ2 = kps.Q2Lim_W();
446  if( rQ2.max < Q2Thres || rQ2.min <=0 ) continue;
447  if( rQ2.max-rQ2.min<0.02 ) {NQ2=5; NQ2b=3;}
448 
449  double logQ2min = TMath::Log(rQ2.min+kASmallNum);
450  double logQ2max = TMath::Log(rQ2.max-kASmallNum);
451  double dlogQ2 = (logQ2max - logQ2min) /(NQ2-1);
452  double xseclast = -1;
453  bool increasing = true;
454 
455  for(int iq2=0; iq2<NQ2; iq2++) {
456  double Q2 = TMath::Exp(logQ2min + iq2 * dlogQ2);
457  interaction->KinePtr()->SetQ2(Q2);
458  double xsec = fXSecModel->XSec(interaction, kPSWQ2fE);
459  LOG("RESKinematics", pDEBUG)
460  << "xsec(W= " << W << ", Q2= " << Q2 << ") = " << xsec;
461  max_xsec = TMath::Max(xsec, max_xsec);
462  increasing = xsec-xseclast>=0;
463  xseclast=xsec;
464 
465  // once the cross section stops increasing, I reduce the step size and
466  // step backwards a little bit to handle cases that the max cross section
467  // is grossly underestimated (very peaky distribution & large step)
468  if(!increasing) {
469  dlogQ2/=NQ2b;
470  for(int iq2b=0; iq2b<NQ2b; iq2b++) {
471  Q2 = TMath::Exp(TMath::Log(Q2) - dlogQ2);
472  if(Q2 < rQ2.min) continue;
473  interaction->KinePtr()->SetQ2(Q2);
474  xsec = fXSecModel->XSec(interaction, kPSWQ2fE);
475  LOG("RESKinematics", pDEBUG)
476  << "xsec(W= " << W << ", Q2= " << Q2 << ") = " << xsec;
477  max_xsec = TMath::Max(xsec, max_xsec);
478  }
479  break;
480  }
481  } // Q2
482  }//W
483  }//2d scan
484 
485  // Apply safety factor, since value retrieved from the cache might
486  // correspond to a slightly different energy
487  // Apply larger safety factor for smaller energies.
488  max_xsec *= ( (E<md) ? 2. : fSafetyFactor);
489 
490 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
491  LOG("RESKinematics", pDEBUG) << interaction->AsString();
492  LOG("RESKinematics", pDEBUG) << "Max xsec in phase space = " << max_xsec;
493  LOG("RESKinematics", pDEBUG) << "Computed using " << fXSecModel->Id();
494 #endif
495 
496  return max_xsec;
497 }
498 //___________________________________________________________________________
double fWcut
Wcut parameter in DIS/RES join scheme.
virtual double MaxXSec(GHepRecord *evrec) const
const KPhaseSpace & PhaseSpace(void) const
Definition: Interaction.h:73
virtual void SetWeight(double wght)
Definition: GHepRecord.h:131
bool fGenerateUniformly
uniform over allowed phase space + event weight?
double J(double q0, double q3, double Enu, double ml)
Definition: MECUtils.cxx:141
THE MAIN GENIE PROJECT NAMESPACE
Definition: GeneratorBase.h:8
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:991
virtual Interaction * Summary(void) const
Definition: GHepRecord.cxx:135
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:79
void SetQ2(double Q2, bool selected=false)
Definition: Kinematics.cxx:265
Kinematics * KinePtr(void) const
Definition: Interaction.h:76
double fSafetyFactor
maxxsec -> maxxsec * safety_factor
virtual double XSec(const Interaction *i, KinePhaseSpace_t k=kPSfE) const =0
Compute the cross section for the input interaction.
A simple [min,max] interval for doubles.
Definition: Range1.h:43
virtual void AssertXSecLimits(const Interaction *in, double xsec, double xsec_max) const
bool KnownResonance(void) const
Definition: XclsTag.h:61
#define pFATAL
Definition: Messenger.h:57
double fMaxXSecDiffTolerance
max{100*(xsec-maxxsec)/.5*(xsec+maxxsec)} if xsec>maxxsec
Defines the EventGeneratorI interface.
double Wmax[kNWBins]
virtual double Weight(void) const
Definition: GHepRecord.h:125
::xsd::cxx::tree::exception< char > exception
Definition: Database.h:225
double Mass(Resonance_t res)
resonance mass (GeV)
Definition: config.py:1
Range1D_t Q2Lim_W(void) const
Q2 limits @ fixed W.
Abstract class. Provides a data caching mechanism for for concrete implementations of the EventRecord...
enum genie::EResonance Resonance_t
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:30
Range1D_t Limits(KineVar_t kvar) const
Return the kinematical variable limits.
string AsString(void) const
static const double kMinQ2Limit
Definition: Controls.h:41
double PhaseSpaceVolume(const Interaction *const i, KinePhaseSpace_t ps)
Definition: KineUtils.cxx:34
Summary information for an interaction.
Definition: Interaction.h:56
const TLorentzVector & HitNucP4(void) const
Definition: Target.h:92
An exception thrown by EventRecordVisitorI when the normal processing sequence has to be disrupted (f...
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
double QD2toQ2(double QD2)
Definition: KineUtils.cxx:985
Float_t E
Definition: plot.C:20
void WQ2toXY(double Ev, double M, double W, double Q2, double &x, double &y)
Definition: KineUtils.cxx:1046
virtual void Configure(const Registry &config)
Definition: Algorithm.cxx:70
Kinematical phase space.
Definition: KPhaseSpace.h:34
static const double kASmallNum
Definition: Controls.h:40
#define pINFO
Definition: Messenger.h:63
TRandom3 & RndKine(void) const
rnd number generator used by kinematics generators
Definition: RandomGen.h:51
TF2 * fEnvelope
2-D envelope used for importance sampling
Resonance_t Resonance(void) const
Definition: XclsTag.h:62
double Q2toQD2(double Q2)
Definition: KineUtils.cxx:975
Misc GENIE control constants.
Double_t xsec[nknots]
Definition: testXsec.C:47
#define pWARN
Definition: Messenger.h:61
const UInt_t kISkipKinematicChk
if set, skip kinematic validity checks
Definition: Interaction.h:48
virtual const XSecAlgorithmI * CrossSectionAlg(void) const =0
bool IsEM(void) const
void Setx(double x, bool selected=false)
Definition: Kinematics.cxx:241
static RunningThreadInfo * Instance(void)
virtual const AlgId & Id(void) const
Get algorithm ID.
Definition: Algorithm.h:98
double max
Definition: Range1.h:54
void SetW(double W, bool selected=false)
Definition: Kinematics.cxx:289
virtual TBits * EventFlags(void) const
Definition: GHepRecord.h:118
A registry. Provides the container for algorithm configuration parameters.
Definition: Registry.h:66
void Sety(double y, bool selected=false)
Definition: Kinematics.cxx:253
const XclsTag & ExclTag(void) const
Definition: Interaction.h:72
double fEMin
min E for which maxxsec is cached - forcing explicit calc.
void ProcessEventRecord(GHepRecord *event_rec) const
exit(0)
assert(nhit_max >=nhit_nbins)
virtual double XSec(void) const
Definition: GHepRecord.h:127
double Jacobian(const Interaction *const i, KinePhaseSpace_t f, KinePhaseSpace_t t)
Definition: KineUtils.cxx:128
double gQ2
Definition: gtestDISSF.cxx:55
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
double Wmin[kNWBins]
double min
Definition: Range1.h:53
#define pNOTICE
Definition: Messenger.h:62
void ClearRunningValues(void)
Definition: Kinematics.cxx:357
bool GetParamDef(const RgKey &name, T &p, const T &def) const
bool GetParam(const RgKey &name, T &p, bool is_top_call=true) const
const Target & Tgt(void) const
Definition: InitialState.h:67
static const unsigned int kRjMaxIterations
Definition: Controls.h:27
void Configure(const Registry &config)
T max(sqlite3 *const db, std::string const &table_name, std::string const &column_name)
Definition: statistics.h:68
const EventGeneratorI * RunningThread(void)
#define W(x)
double ProbeE(RefFrame_t rf) const
GENIE&#39;s GHEP MC event record.
Definition: GHepRecord.h:46
Keep info on the event generation thread currently on charge. This is used so that event generation m...
double ComputeMaxXSec(const Interaction *interaction) const
Root of GENIE utility namespaces.
Range1D_t WLim(void) const
W limits.
const UInt_t kISkipProcessChk
if set, skip process validity checks
Definition: Interaction.h:47
virtual void SetDiffXSec(double xsec, KinePhaseSpace_t ps)
Definition: GHepRecord.h:134
Initial State information.
Definition: InitialState.h:49
#define pDEBUG
Definition: Messenger.h:64
double RESImportanceSamplingEnvelope(double *x, double *par)
Definition: KineUtils.cxx:1249