GENIEGenerator
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
DMELKinematicsGenerator.cxx
Go to the documentation of this file.
1 //____________________________________________________________________________
2 /*
3  Copyright (c) 2003-2024, The GENIE Collaboration
4  For the full text of the license visit http://copyright.genie-mc.org
5 
6 
7  Author: Joshua Berger <jberger \at physics.wisc.edu>
8  University of Wisconsin-Madison
9 
10  Costas Andreopoulos <c.andreopoulos \at cern.ch>
11  University of Liverpool
12 */
13 //____________________________________________________________________________
14 
15 #include <TMath.h>
16 
17 #include "Framework/Conventions/GBuild.h"
34 
35 using namespace genie;
36 using namespace genie::controls;
37 using namespace genie::constants;
38 using namespace genie::utils;
39 
40 //___________________________________________________________________________
42 KineGeneratorWithCache("genie::DMELKinematicsGenerator")
43 {
44 
45 }
46 //___________________________________________________________________________
48 KineGeneratorWithCache("genie::DMELKinematicsGenerator", config)
49 {
50 
51 }
52 //___________________________________________________________________________
54 {
55 
56 }
57 //___________________________________________________________________________
59 {
60  if(fGenerateUniformly) {
61  LOG("DMELKinematics", pNOTICE)
62  << "Generating kinematics uniformly over the allowed phase space";
63  }
64 
65  //-- Get the random number generators
67 
68  //-- Access cross section algorithm for running thread
70  const EventGeneratorI * evg = rtinfo->RunningThread();
71  fXSecModel = evg->CrossSectionAlg();
72 
73  //-- Get the interaction and set the 'trust' bits
74  Interaction * interaction = evrec->Summary();
75  interaction->SetBit(kISkipProcessChk);
76  interaction->SetBit(kISkipKinematicChk);
77 
78  // store the struck nucleon position for use by the xsec method
79  double hitNucPos = evrec->HitNucleon()->X4()->Vect().Mag();
80  interaction->InitStatePtr()->TgtPtr()->SetHitNucPosition(hitNucPos);
81 
82  //-- Note: The kinematic generator would be using the free nucleon cross
83  // section (even for nuclear targets) so as not to double-count nuclear
84  // suppression. This assumes that a) the nuclear suppression was turned
85  // on when computing the cross sections for selecting the current event
86  // and that b) if the event turns out to be unphysical (Pauli-blocked)
87  // the next attempted event will be forced to DM again.
88  // (discussion with Hugh - GENIE/NeuGEN integration workshop - 07APR2006
89  interaction->SetBit(kIAssumeFreeNucleon);
90 
91  //-- Get the limits for the generated Q2
92  const KPhaseSpace & kps = interaction->PhaseSpace();
93  Range1D_t Q2 = kps.Limits(kKVQ2);
94 
95  if(Q2.max <=0 || Q2.min>=Q2.max) {
96  LOG("DMELKinematics", pWARN) << "No available phase space";
97  evrec->EventFlags()->SetBitNumber(kKineGenErr, true);
99  exception.SetReason("No available phase space");
100  exception.SwitchOnFastForward();
101  throw exception;
102  }
103 
104  //-- For the subsequent kinematic selection with the rejection method:
105  // Calculate the max differential cross section or retrieve it from the
106  // cache. Throw an exception and quit the evg thread if a non-positive
107  // value is found.
108  // If the kinematics are generated uniformly over the allowed phase
109  // space the max xsec is irrelevant
110  LOG("DMELKinematics",pNOTICE) << "Setting max XSec";
111  double xsec_max = (fGenerateUniformly) ? -1 : this->MaxXSec(evrec);
112  LOG("DMELKinematics",pNOTICE) << "Set max XSec to " << xsec_max;
113 
114  //-- Try to select a valid Q2 using the rejection method
115 
116  // kinematical limits
117  double Q2min = Q2.min+kASmallNum;
118  double Q2max = Q2.max-kASmallNum;
119 //double QD2min = utils::kinematics::Q2toQD2(Q2min);
120 //double QD2max = utils::kinematics::Q2toQD2(Q2max);
121  double xsec = -1.;
122  double gQ2 = 0.;
123 
124  unsigned int iter = 0;
125  bool accept = false;
126  while(1) {
127  iter++;
128  if(iter > kRjMaxIterations) {
129  LOG("DMELKinematics", pWARN)
130  << "Couldn't select a valid Q^2 after " << iter << " iterations";
131  evrec->EventFlags()->SetBitNumber(kKineGenErr, true);
133  exception.SetReason("Couldn't select kinematics");
134  exception.SwitchOnFastForward();
135  throw exception;
136  }
137 
138  //-- Generate a Q2 value within the allowed phase space
139 /*
140  if(fGenerateUniformly) {
141  gQ2 = Q2min + (Q2max-Q2min) * rnd->RndKine().Rndm();
142  } else {
143  // In unweighted mode - use transform that takes out the dipole form
144  double gQD2 = QD2min + (QD2max-QD2min) * rnd->RndKine().Rndm();
145  gQ2 = utils::kinematics::QD2toQ2(gQD2);
146  }
147 */
148  LOG("DMELKinematics",pNOTICE) << "Attempting to sample";
149  gQ2 = Q2min + (Q2max-Q2min) * rnd->RndKine().Rndm();
150  interaction->KinePtr()->SetQ2(gQ2);
151  LOG("DMELKinematics", pINFO) << "Trying: Q^2 = " << gQ2;
152 
153  //-- Computing cross section for the current kinematics
154  xsec = fXSecModel->XSec(interaction, kPSQ2fE);
155 
156  //-- Decide whether to accept the current kinematics
157  if(!fGenerateUniformly) {
158  this->AssertXSecLimits(interaction, xsec, xsec_max);
159 
160  double t = xsec_max * rnd->RndKine().Rndm();
161  //double J = kinematics::Jacobian(interaction,kPSQ2fE,kPSQD2fE);
162  double J = 1.;
163 
164 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
165  LOG("DMELKinematics", pDEBUG)
166  << "xsec= " << xsec << ", J= " << J << ", Rnd= " << t;
167 #endif
168  accept = (t < J*xsec);
169  } else {
170  accept = (xsec>0);
171  }
172 
173  //-- If the generated kinematics are accepted, finish-up module's job
174  if(accept) {
175  LOG("DMELKinematics", pINFO) << "Selected: Q^2 = " << gQ2;
176 
177  // reset bits
178  interaction->ResetBit(kISkipProcessChk);
179  interaction->ResetBit(kISkipKinematicChk);
180  interaction->ResetBit(kIAssumeFreeNucleon);
181 
182  // compute the rest of the kinematical variables
183 
184  // get neutrino energy at struck nucleon rest frame and the
185  // struck nucleon mass (can be off the mass shell)
186  const InitialState & init_state = interaction->InitState();
187  double E = init_state.ProbeE(kRfHitNucRest);
188  double M = init_state.Tgt().HitNucP4().M();
189 
190  LOG("DMELKinematics", pNOTICE) << "E = " << E << ", M = "<< M;
191 
192  // The hadronic inv. mass is equal to the recoil nucleon on-shell mass.
193  // For DM/Charm events it is set to be equal to the on-shell mass of
194  // the generated charm baryon (Lamda_c+, Sigma_c+ or Sigma_c++)
195  // Similarly for strange baryons
196  //
197  const XclsTag & xcls = interaction->ExclTag();
198  int rpdgc = 0;
199  if(xcls.IsCharmEvent()) { rpdgc = xcls.CharmHadronPdg(); }
200  else if(xcls.IsStrangeEvent()) { rpdgc = xcls.StrangeHadronPdg(); }
201  else { rpdgc = interaction->RecoilNucleonPdg(); }
202  assert(rpdgc);
203  double gW = PDGLibrary::Instance()->Find(rpdgc)->Mass();
204 
205  LOG("DMELKinematics", pNOTICE) << "Selected: W = "<< gW;
206 
207  // (W,Q2) -> (x,y)
208  double gx=0, gy=0;
209  kinematics::WQ2toXY(E,M,gW,gQ2,gx,gy);
210 
211  // set the cross section for the selected kinematics
212  evrec->SetDiffXSec(xsec,kPSQ2fE);
213 
214  // for uniform kinematics, compute an event weight as
215  // wght = (phase space volume)*(differential xsec)/(event total xsec)
216  if(fGenerateUniformly) {
217  double vol = kinematics::PhaseSpaceVolume(interaction,kPSQ2fE);
218  double totxsec = evrec->XSec();
219  double wght = (vol/totxsec)*xsec;
220  LOG("DMELKinematics", pNOTICE) << "Kinematics wght = "<< wght;
221 
222  // apply computed weight to the current event weight
223  wght *= evrec->Weight();
224  LOG("DMELKinematics", pNOTICE) << "Current event wght = " << wght;
225  evrec->SetWeight(wght);
226  }
227 
228  // lock selected kinematics & clear running values
229  interaction->KinePtr()->SetQ2(gQ2, true);
230  interaction->KinePtr()->SetW (gW, true);
231  interaction->KinePtr()->Setx (gx, true);
232  interaction->KinePtr()->Sety (gy, true);
233  interaction->KinePtr()->ClearRunningValues();
234 
235  return;
236  }
237  }// iterations
238 }
239 //___________________________________________________________________________
241  GHepRecord * evrec) const
242 {
243  //-- Get the random number generators
244  RandomGen * rnd = RandomGen::Instance();
245 
246  //-- Access cross section algorithm for running thread
248  const EventGeneratorI * evg = rtinfo->RunningThread();
249  fXSecModel = evg->CrossSectionAlg();
250 
251  //-- Get the interaction and set the 'trust' bits
252  Interaction * interaction = new Interaction(*evrec->Summary());
253  interaction->SetBit(kISkipProcessChk);
254  interaction->SetBit(kISkipKinematicChk);
255 
256  // store the struck nucleon position for use by the xsec method
257  double hitNucPos = evrec->HitNucleon()->X4()->Vect().Mag();
258  interaction->InitStatePtr()->TgtPtr()->SetHitNucPosition(hitNucPos);
259 
260  //-- Note: The kinematic generator would be using the free nucleon cross
261  // section (even for nuclear targets) so as not to double-count nuclear
262  // suppression. This assumes that a) the nuclear suppression was turned
263  // on when computing the cross sections for selecting the current event
264  // and that b) if the event turns out to be unphysical (Pauli-blocked)
265  // the next attempted event will be forced to DM again.
266  // (discussion with Hugh - GENIE/NeuGEN integration workshop - 07APR2006
267  interaction->SetBit(kIAssumeFreeNucleon);
268 
269  //-- Assume scattering off a nucleon on the mass shell (PWIA prescription)
270  double Mn = interaction->InitState().Tgt().HitNucMass(); // PDG mass, take it to be on-shell
271  double pxn = interaction->InitState().Tgt().HitNucP4().Px();
272  double pyn = interaction->InitState().Tgt().HitNucP4().Py();
273  double pzn = interaction->InitState().Tgt().HitNucP4().Pz();
274  double En = interaction->InitState().Tgt().HitNucP4().Energy();
275  double En0 = TMath::Sqrt(pxn*pxn + pyn*pyn + pzn*pzn + Mn*Mn);
276  double Eb = En0 - En;
277  interaction->InitStatePtr()->TgtPtr()->HitNucP4Ptr()->SetE(En0);
278  double ml = interaction->FSPrimLepton()->Mass();
279 
280  //-- Get the limits for the generated Q2
281  const KPhaseSpace & kps = interaction->PhaseSpace();
282  Range1D_t Q2 = kps.Limits(kKVQ2);
283 
284  if(Q2.max <=0 || Q2.min>=Q2.max) {
285  LOG("DMELKinematics", pWARN) << "No available phase space";
286  evrec->EventFlags()->SetBitNumber(kKineGenErr, true);
288  exception.SetReason("No available phase space");
289  exception.SwitchOnFastForward();
290  throw exception;
291  }
292 
293  //-- For the subsequent kinematic selection with the rejection method:
294  // Calculate the max differential cross section or retrieve it from the
295  // cache. Throw an exception and quit the evg thread if a non-positive
296  // value is found.
297  // If the kinematics are generated uniformly over the allowed phase
298  // space the max xsec is irrelevant
299 // double xsec_max = (fGenerateUniformly) ? -1 : this->MaxXSec(evrec);
300  double xsec_max = this->MaxXSec(evrec);
301 
302  // get neutrino energy at struck nucleon rest frame and the
303  // struck nucleon mass (can be off the mass shell)
304  const InitialState & init_state = interaction->InitState();
305  double E = init_state.ProbeE(kRfHitNucRest);
306 
307  LOG("DMELKinematics", pNOTICE) << "E = " << E << ", M = "<< Mn;
308 
309  //-- Try to select a valid Q2 using the rejection method
310 
311  // kinematical limits
312  double Q2min = Q2.min+kASmallNum;
313  double Q2max = Q2.max-kASmallNum;
314  double xsec = -1.;
315  double gQ2 = 0.;
316  double gW = 0.;
317  double gx = 0.;
318  double gy = 0.;
319 
320  unsigned int iter = 0;
321  bool accept = false;
322  while(1) {
323  iter++;
324  if(iter > kRjMaxIterations) {
325  LOG("DMELKinematics", pWARN)
326  << "Couldn't select a valid Q^2 after " << iter << " iterations";
327  evrec->EventFlags()->SetBitNumber(kKineGenErr, true);
329  exception.SetReason("Couldn't select kinematics");
330  exception.SwitchOnFastForward();
331  throw exception;
332  }
333 
334  //-- Generate a Q2 value within the allowed phase space
335  gQ2 = Q2min + (Q2max-Q2min) * rnd->RndKine().Rndm();
336  LOG("DMELKinematics", pNOTICE) << "Trying: Q^2 = " << gQ2;
337 
338  // The hadronic inv. mass is equal to the recoil nucleon on-shell mass.
339  // For DM/Charm events it is set to be equal to the on-shell mass of
340  // the generated charm baryon (Lamda_c+, Sigma_c+ or Sigma_c++)
341  //
342  const XclsTag & xcls = interaction->ExclTag();
343  int rpdgc = 0;
344  if(xcls.IsCharmEvent()) { rpdgc = xcls.CharmHadronPdg(); }
345  else { rpdgc = interaction->RecoilNucleonPdg(); }
346  assert(rpdgc);
347  gW = PDGLibrary::Instance()->Find(rpdgc)->Mass();
348 
349  // (W,Q2) -> (x,y)
350  kinematics::WQ2toXY(E,Mn,gW,gQ2,gx,gy);
351 
352  LOG("DMELKinematics", pNOTICE) << "W = "<< gW;
353  LOG("DMELKinematics", pNOTICE) << "x = "<< gx;
354  LOG("DMELKinematics", pNOTICE) << "y = "<< gy;
355 
356  // v
357  double gv = gy * E;
358  double gv2 = gv*gv;
359 
360  LOG("DMELKinematics", pNOTICE) << "v = "<< gv;
361 
362  // v -> v~
363  double gvtilde = gv + Mn - Eb - TMath::Sqrt(Mn*Mn+pxn*pxn+pyn*pyn+pzn*pzn);
364  double gvtilde2 = gvtilde*gvtilde;
365 
366  LOG("DMELKinematics", pNOTICE) << "v~ = "<< gvtilde;
367 
368  // Q~^2
369  double gQ2tilde = gQ2 - gv2 + gvtilde2;
370 
371  LOG("DMELKinematics", pNOTICE) << "Q~^2 = "<< gQ2tilde;
372 
373  // Set updated Q2
374  interaction->KinePtr()->SetQ2(gQ2tilde);
375 
376  //-- Computing cross section for the current kinematics
377  xsec = fXSecModel->XSec(interaction, kPSQ2fE);
378 
379  //-- Decide whether to accept the current kinematics
380 // if(!fGenerateUniformly) {
381  this->AssertXSecLimits(interaction, xsec, xsec_max);
382 
383  double t = xsec_max * rnd->RndKine().Rndm();
384 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
385  LOG("DMELKinematics", pDEBUG)
386  << "xsec= " << xsec << ", Rnd= " << t;
387 #endif
388  accept = (t < xsec);
389 // } else {
390 // accept = (xsec>0);
391 // }
392 
393  //-- If the generated kinematics are accepted, finish-up module's job
394  if(accept) {
395  LOG("DMELKinematics", pNOTICE) << "Selected: Q^2 = " << gQ2;
396 
397  // reset bits
398 // interaction->ResetBit(kISkipProcessChk);
399 // interaction->ResetBit(kISkipKinematicChk);
400 // interaction->ResetBit(kIAssumeFreeNucleon);
401 
402  // set the cross section for the selected kinematics
403  evrec->SetDiffXSec(xsec,kPSQ2fE);
404 
405  // for uniform kinematics, compute an event weight as
406  // wght = (phase space volume)*(differential xsec)/(event total xsec)
407 // if(fGenerateUniformly) {
408 // double vol = kinematics::PhaseSpaceVolume(interaction,kPSQ2fE);
409 // double totxsec = evrec->XSec();
410 // double wght = (vol/totxsec)*xsec;
411 // LOG("DMELKinematics", pNOTICE) << "Kinematics wght = "<< wght;
412 
413  // apply computed weight to the current event weight
414 // wght *= evrec->Weight();
415 // LOG("DMELKinematics", pNOTICE) << "Current event wght = " << wght;
416 // evrec->SetWeight(wght);
417 // }
418 
419  // lock selected kinematics & clear running values
420 // interaction->KinePtr()->SetQ2(gQ2, true);
421 // interaction->KinePtr()->SetW (gW, true);
422 // interaction->KinePtr()->Setx (gx, true);
423 // interaction->KinePtr()->Sety (gy, true);
424 // interaction->KinePtr()->ClearRunningValues();
425 
426  evrec->Summary()->KinePtr()->SetQ2(gQ2, true);
427  evrec->Summary()->KinePtr()->SetW (gW, true);
428  evrec->Summary()->KinePtr()->Setx (gx, true);
429  evrec->Summary()->KinePtr()->Sety (gy, true);
430  evrec->Summary()->KinePtr()->ClearRunningValues();
431  delete interaction;
432 
433  return;
434  }
435  }// iterations
436 }
437 //___________________________________________________________________________
439 {
440  Algorithm::Configure(config);
441  this->LoadConfig();
442 }
443 //____________________________________________________________________________
445 {
446  Algorithm::Configure(config);
447  this->LoadConfig();
448 }
449 //____________________________________________________________________________
451 {
452 // Load sub-algorithms and config data to reduce the number of registry
453 // lookups
454 
455  //-- Safety factor for the maximum differential cross section
456  GetParamDef( "MaxXSec-SafetyFactor", fSafetyFactor , 1.25 ) ;
457 
458  //-- Minimum energy for which max xsec would be cached, forcing explicit
459  // calculation for lower eneries
460  GetParamDef( "Cache-MinEnergy", fEMin, 1.00 ) ;
461 
462  //-- Maximum allowed fractional cross section deviation from maxim cross
463  // section used in rejection method
464  GetParamDef( "MaxXSec-DiffTolerance", fMaxXSecDiffTolerance, 999999. ) ;
465  assert(fMaxXSecDiffTolerance>=0);
466 
467  //-- Generate kinematics uniformly over allowed phase space and compute
468  // an event weight?
469  GetParamDef( "UniformOverPhaseSpace", fGenerateUniformly, false ) ;
470 
471 }
472 //____________________________________________________________________________
474  const Interaction * interaction) const
475 {
476 // Computes the maximum differential cross section in the requested phase
477 // space. This method overloads KineGeneratorWithCache::ComputeMaxXSec
478 // method and the value is cached at a circular cache branch for retrieval
479 // during subsequent event generation.
480 // The computed max differential cross section does not need to be the exact
481 // maximum. The number used in the rejection method will be scaled up by a
482 // safety factor. But it needs to be fast - do not use a very small dQ2 step.
483 
484  double max_xsec = 0.0;
485 
486  const KPhaseSpace & kps = interaction->PhaseSpace();
487  Range1D_t rQ2 = kps.Limits(kKVQ2);
488  LOG("DMELKinematics", pDEBUG) << "Range of Q^2: " << rQ2.min << " to " << rQ2.max;
489  if(rQ2.min <=0 || rQ2.max <= rQ2.min) return 0.;
490 
491  const double logQ2min = TMath::Log(rQ2.min + kASmallNum);
492  const double logQ2max = TMath::Log(rQ2.max - kASmallNum);
493 
494  const int N = 15;
495  const int Nb = 10;
496 
497  double dlogQ2 = (logQ2max - logQ2min) /(N-1);
498  double xseclast = -1;
499  bool increasing;
500 
501  for(int i=0; i<N; i++) {
502  double Q2 = TMath::Exp(logQ2min + i * dlogQ2);
503  interaction->KinePtr()->SetQ2(Q2);
504  double xsec = fXSecModel->XSec(interaction, kPSQ2fE);
505 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
506  LOG("DMELKinematics", pDEBUG) << "xsec(Q2= " << Q2 << ") = " << xsec;
507 #endif
508  max_xsec = TMath::Max(xsec, max_xsec);
509  increasing = xsec-xseclast>=0;
510  xseclast = xsec;
511 
512  // once the cross section stops increasing, I reduce the step size and
513  // step backwards a little bit to handle cases that the max cross section
514  // is grossly underestimated (very peaky distribution & large step)
515  if(!increasing) {
516  dlogQ2/=(Nb+1);
517  for(int ib=0; ib<Nb; ib++) {
518  Q2 = TMath::Exp(TMath::Log(Q2) - dlogQ2);
519  if(Q2 < rQ2.min) continue;
520  interaction->KinePtr()->SetQ2(Q2);
521  xsec = fXSecModel->XSec(interaction, kPSQ2fE);
522 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
523  LOG("DMELKinematics", pDEBUG) << "xsec(Q2= " << Q2 << ") = " << xsec;
524 #endif
525  max_xsec = TMath::Max(xsec, max_xsec);
526  }
527  break;
528  }
529  }//Q^2
530 
531  // Apply safety factor, since value retrieved from the cache might
532  // correspond to a slightly different energy
533  max_xsec *= fSafetyFactor;
534 
535 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
536  SLOG("DMELKinematics", pDEBUG) << interaction->AsString();
537  SLOG("DMELKinematics", pDEBUG) << "Max xsec in phase space = " << max_xsec;
538  SLOG("DMELKinematics", pDEBUG) << "Computed using alg = " << *fXSecModel;
539 #endif
540 
541  return max_xsec;
542 }
543 //___________________________________________________________________________
const KPhaseSpace & PhaseSpace(void) const
Definition: Interaction.h:73
virtual void SetWeight(double wght)
Definition: GHepRecord.h:130
bool fGenerateUniformly
uniform over allowed phase space + event weight?
double J(double q0, double q3, double Enu, double ml)
Definition: MECUtils.cxx:147
void ProcessEventRecord(GHepRecord *event_rec) const
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1077
virtual Interaction * Summary(void) const
Definition: GHepRecord.cxx:91
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:71
void SetQ2(double Q2, bool selected=false)
Definition: Kinematics.cxx:255
Kinematics * KinePtr(void) const
Definition: Interaction.h:76
double fSafetyFactor
ComputeMaxXSec -&gt; ComputeMaxXSec * fSafetyFactor.
int RecoilNucleonPdg(void) const
recoil nucleon pdg
A simple [min,max] interval for doubles.
Definition: Range1.h:42
virtual void AssertXSecLimits(const Interaction *in, double xsec, double xsec_max) const
double HitNucMass(void) const
Definition: Target.cxx:233
int CharmHadronPdg(void) const
Definition: XclsTag.h:52
bool IsStrangeEvent(void) const
Definition: XclsTag.h:53
double fMaxXSecDiffTolerance
max{100*(xsec-maxxsec)/.5*(xsec+maxxsec)} if xsec&gt;maxxsec
Defines the EventGeneratorI interface.
virtual double MaxXSec(GHepRecord *evrec, const int nkey=0) const
void Configure(const Registry &config)
virtual double Weight(void) const
Definition: GHepRecord.h:124
void SetHitNucPosition(double r)
Definition: Target.cxx:210
Abstract class. Provides a data caching mechanism for for concrete implementations of the EventRecord...
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:29
Range1D_t Limits(KineVar_t kvar) const
Return the kinematical variable limits.
string AsString(void) const
Contains minimal information for tagging exclusive processes.
Definition: XclsTag.h:39
double PhaseSpaceVolume(const Interaction *const i, KinePhaseSpace_t ps)
Definition: KineUtils.cxx:36
bool IsCharmEvent(void) const
Definition: XclsTag.h:50
Summary information for an interaction.
Definition: Interaction.h:56
const TLorentzVector & HitNucP4(void) const
Definition: Target.h:91
An exception thrown by EventRecordVisitorI when the normal processing sequence has to be disrupted (f...
int StrangeHadronPdg(void) const
Definition: XclsTag.h:55
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
void WQ2toXY(double Ev, double M, double W, double Q2, double &x, double &y)
Definition: KineUtils.cxx:1132
virtual void Configure(const Registry &config)
Definition: Algorithm.cxx:62
Kinematical phase space.
Definition: KPhaseSpace.h:33
static const double kASmallNum
Definition: Controls.h:40
#define pINFO
Definition: Messenger.h:62
TRandom3 & RndKine(void) const
rnd number generator used by kinematics generators
Definition: RandomGen.h:50
#define pWARN
Definition: Messenger.h:60
const UInt_t kISkipKinematicChk
if set, skip kinematic validity checks
Definition: Interaction.h:48
virtual double XSec(const Interaction *i, KinePhaseSpace_t k=kPSfE) const =0
Compute the cross section for the input interaction.
TParticlePDG * FSPrimLepton(void) const
final state primary lepton
double ComputeMaxXSec(const Interaction *in) const
void Setx(double x, bool selected=false)
Definition: Kinematics.cxx:231
static RunningThreadInfo * Instance(void)
double max
Definition: Range1.h:53
void SetW(double W, bool selected=false)
Definition: Kinematics.cxx:279
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:247
static PDGLibrary * Instance(void)
Definition: PDGLibrary.cxx:68
virtual TBits * EventFlags(void) const
Definition: GHepRecord.h:117
const TLorentzVector * X4(void) const
Definition: GHepParticle.h:79
A registry. Provides the container for algorithm configuration parameters.
Definition: Registry.h:65
void Sety(double y, bool selected=false)
Definition: Kinematics.cxx:243
const UInt_t kIAssumeFreeNucleon
Definition: Interaction.h:49
const XclsTag & ExclTag(void) const
Definition: Interaction.h:72
virtual GHepParticle * HitNucleon(void) const
Definition: GHepRecord.cxx:313
Target * TgtPtr(void) const
Definition: InitialState.h:67
double fEMin
min E for which maxxsec is cached - forcing explicit calc.
virtual const XSecAlgorithmI * CrossSectionAlg(void) const =0
virtual double XSec(void) const
Definition: GHepRecord.h:126
double gQ2
Definition: gtestDISSF.cxx:55
InitialState * InitStatePtr(void) const
Definition: Interaction.h:74
const InitialState & InitState(void) const
Definition: Interaction.h:69
double min
Definition: Range1.h:52
void SpectralFuncExperimentalCode(GHepRecord *event_rec) const
TParticlePDG * Find(int pdgc, bool must_exist=true)
Definition: PDGLibrary.cxx:86
#define pNOTICE
Definition: Messenger.h:61
void ClearRunningValues(void)
Definition: Kinematics.cxx:347
bool GetParamDef(const RgKey &name, T &p, const T &def) const
const Target & Tgt(void) const
Definition: InitialState.h:66
static const unsigned int kRjMaxIterations
Definition: Controls.h:26
const EventGeneratorI * RunningThread(void)
#define SLOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a short string (using the FUNCTION and...
Definition: Messenger.h:84
double ProbeE(RefFrame_t rf) const
GENIE&#39;s GHEP MC event record.
Definition: GHepRecord.h:45
Keep info on the event generation thread currently on charge. This is used so that event generation m...
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:133
Initial State information.
Definition: InitialState.h:48
#define pDEBUG
Definition: Messenger.h:63