Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes | List of all members
genie::KPhaseSpace Class Reference

Kinematical phase space. More...

#include "/cvmfs/nova.opensciencegrid.org/externals/genie/v3_00_06_p01/Linux64bit+3.10-2.17-e17-debug/GENIE-Generator/src/Framework/Interaction/KPhaseSpace.h"

Inheritance diagram for genie::KPhaseSpace:

Public Member Functions

 KPhaseSpace (void)
 
 KPhaseSpace (const Interaction *in)
 
 ~KPhaseSpace (void)
 
void UseInteraction (const Interaction *in)
 
double Threshold (void) const
 Energy threshold. More...
 
bool IsAboveThreshold (void) const
 Checks whether the interaction is above the energy threshold. More...
 
bool IsAllowed (void) const
 Check whether the current kinematics is in the allowed phase space. More...
 
Range1D_t Limits (KineVar_t kvar) const
 Return the kinematical variable limits. More...
 
double Minimum (KineVar_t kvar) const
 
double Maximum (KineVar_t kvar) const
 
Range1D_t WLim (void) const
 W limits. More...
 
Range1D_t Q2Lim_W (void) const
 Q2 limits @ fixed W. More...
 
Range1D_t q2Lim_W (void) const
 q2 limits @ fixed W More...
 
Range1D_t Q2Lim (void) const
 Q2 limits. More...
 
Range1D_t q2Lim (void) const
 q2 limits More...
 
Range1D_t XLim (void) const
 x limits More...
 
Range1D_t YLim (void) const
 y limits More...
 
Range1D_t YLim_X (void) const
 y limits @ fixed x More...
 
Range1D_t YLim (double xsi) const
 y limits (COH) More...
 
Range1D_t YLim_X (double xsi) const
 y limits @ fixed x (COH) More...
 
Range1D_t TLim (void) const
 t limits More...
 

Static Public Member Functions

static double GetTMaxDFR ()
 

Private Member Functions

void Init (void)
 

Private Attributes

const InteractionfInteraction
 

Detailed Description

Kinematical phase space.

Author
Costas Andreopoulos <costas.andreopoulos stfc.ac.uk> University of Liverpool & STFC Rutherford Appleton Lab

May 06, 2004

Copyright (c) 2003-2019, The GENIE Collaboration For the full text of the license visit http://copyright.genie-mc.org or see $GENIE/LICENSE

Definition at line 34 of file KPhaseSpace.h.

Constructor & Destructor Documentation

genie::KPhaseSpace::KPhaseSpace ( void  )
KPhaseSpace::KPhaseSpace ( const Interaction in)

Definition at line 47 of file KPhaseSpace.cxx.

References UseInteraction().

47  :
48 TObject(), fInteraction(NULL)
49 {
50  this->UseInteraction(in);
51 }
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
void UseInteraction(const Interaction *in)
Definition: KPhaseSpace.cxx:76
KPhaseSpace::~KPhaseSpace ( void  )

Definition at line 53 of file KPhaseSpace.cxx.

54 {
55 
56 }

Member Function Documentation

double KPhaseSpace::GetTMaxDFR ( )
static

Definition at line 58 of file KPhaseSpace.cxx.

References genie::AlgConfigPool::CommonList(), genie::Registry::GetDouble(), genie::AlgConfigPool::Instance(), and r().

Referenced by genie::DFRKinematicsGenerator::ProcessEventRecord(), and TLim().

59 {
60  static bool tMaxLoaded = false;
61  static double DFR_tMax = -1;
62 
63  if (!tMaxLoaded)
64  {
66  const Registry * r = confp->CommonList( "Param", "Diffractive" ) ;
67  double tmax = r->GetDouble("DFR-t-max");
68  DFR_tMax = tmax;
69  tMaxLoaded = true;
70  }
71 
72  return DFR_tMax;
73 
74 }
A singleton class holding all configuration registries built while parsing all loaded XML configurati...
Definition: AlgConfigPool.h:41
RgDbl GetDouble(RgKey key) const
Definition: Registry.cxx:489
Registry * CommonList(const string &file_id, const string &set_name) const
A registry. Provides the container for algorithm configuration parameters.
Definition: Registry.h:66
TRandom3 r(0)
static AlgConfigPool * Instance()
void genie::KPhaseSpace::Init ( void  )
private
bool KPhaseSpace::IsAboveThreshold ( void  ) const

Checks whether the interaction is above the energy threshold.

Definition at line 225 of file KPhaseSpace.cxx.

References E, fInteraction, genie::Interaction::InitState(), genie::ProcessInfo::IsAMNuGamma(), genie::ProcessInfo::IsCoherent(), genie::ProcessInfo::IsDarkMatterDeepInelastic(), genie::ProcessInfo::IsDarkMatterElastic(), genie::ProcessInfo::IsDeepInelastic(), genie::ProcessInfo::IsDiffractive(), genie::ProcessInfo::IsIMDAnnihilation(), genie::ProcessInfo::IsInverseBetaDecay(), genie::ProcessInfo::IsInverseMuDecay(), genie::ProcessInfo::IsMEC(), genie::ProcessInfo::IsNuElectronElastic(), genie::ProcessInfo::IsQuasiElastic(), genie::ProcessInfo::IsResonant(), genie::ProcessInfo::IsSingleKaon(), genie::kRfHitNucRest, genie::kRfLab, LOG, pDEBUG, python.hepunit::pi, genie::InitialState::ProbeE(), genie::Interaction::ProcInfo(), and Threshold().

Referenced by genie::COHXSecAR::Integrate(), genie::IMDXSec::Integrate(), genie::COHXSec::Integrate(), genie::RESXSec::Integrate(), genie::DISXSec::Integrate(), genie::QELXSec::Integrate(), genie::DMDISXSec::Integrate(), genie::DMELXSec::Integrate(), genie::MECXSec::Integrate(), genie::AlamSimoAtharVacasSKXSec::Integrate(), genie::DFRXSec::Integrate(), genie::NuElectronXSec::Integrate(), genie::SmithMonizQELCCXSec::Integrate(), genie::ReinSehgalRESXSec::Integrate(), genie::ReinSehgalRESXSecFast::Integrate(), genie::ReinSehgalSPPXSec::Integrate(), genie::FermiMover::KickHitNucleon(), genie::XSecAlgorithmI::ValidKinematics(), genie::StrumiaVissaniIBDPXSec::ValidKinematics(), and genie::KLVOxygenIBDPXSec::ValidKinematics().

226 {
227  double E = 0.;
228  double Ethr = this->Threshold();
229 
230  const ProcessInfo & pi = fInteraction->ProcInfo();
231  const InitialState & init_state = fInteraction->InitState();
232 
233  if (pi.IsCoherent() ||
234  pi.IsInverseMuDecay() ||
235  pi.IsIMDAnnihilation() ||
236  pi.IsNuElectronElastic() ||
237  pi.IsMEC())
238  {
239  E = init_state.ProbeE(kRfLab);
240  }
241 
242  if(pi.IsQuasiElastic() ||
243  pi.IsDarkMatterElastic() ||
244  pi.IsInverseBetaDecay() ||
245  pi.IsResonant() ||
246  pi.IsDeepInelastic() ||
248  pi.IsDiffractive() ||
249  pi.IsSingleKaon() ||
250  pi.IsAMNuGamma())
251  {
252  E = init_state.ProbeE(kRfHitNucRest);
253  }
254 
255  LOG("KPhaseSpace", pDEBUG) << "E = " << E << ", Ethr = " << Ethr;
256  return (E>Ethr);
257 }
bool IsResonant(void) const
Definition: ProcessInfo.cxx:92
double Threshold(void) const
Energy threshold.
Definition: KPhaseSpace.cxx:81
bool IsInverseMuDecay(void) const
bool IsQuasiElastic(void) const
Definition: ProcessInfo.cxx:67
bool IsInverseBetaDecay(void) const
bool IsDiffractive(void) const
bool IsIMDAnnihilation(void) const
bool IsSingleKaon(void) const
Definition: ProcessInfo.cxx:77
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
bool IsNuElectronElastic(void) const
Float_t E
Definition: plot.C:20
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
bool IsAMNuGamma(void) const
bool IsDarkMatterElastic(void) const
Definition: ProcessInfo.cxx:72
bool IsMEC(void) const
bool IsDeepInelastic(void) const
Definition: ProcessInfo.cxx:82
const InitialState & InitState(void) const
Definition: Interaction.h:69
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
bool IsDarkMatterDeepInelastic(void) const
Definition: ProcessInfo.cxx:87
double ProbeE(RefFrame_t rf) const
Initial State information.
Definition: InitialState.h:49
#define pDEBUG
Definition: Messenger.h:64
bool IsCoherent(void) const
Definition: ProcessInfo.cxx:97
bool KPhaseSpace::IsAllowed ( void  ) const

Check whether the current kinematics is in the allowed phase space.

Definition at line 259 of file KPhaseSpace.cxx.

References fInteraction, genie::ProcessInfo::IsCoherent(), genie::ProcessInfo::IsDarkMatterDeepInelastic(), genie::ProcessInfo::IsDarkMatterElastic(), genie::ProcessInfo::IsDeepInelastic(), genie::ProcessInfo::IsDiffractive(), genie::ProcessInfo::IsIMDAnnihilation(), genie::ProcessInfo::IsInverseBetaDecay(), genie::ProcessInfo::IsInverseMuDecay(), genie::ProcessInfo::IsMEC(), genie::ProcessInfo::IsNuElectronElastic(), genie::ProcessInfo::IsQuasiElastic(), genie::ProcessInfo::IsResonant(), genie::ProcessInfo::IsSingleKaon(), genie::utils::math::IsWithinLimits(), genie::Interaction::Kine(), LOG, genie::Range1D_t::max, genie::Range1D_t::min, pDEBUG, python.hepunit::pi, genie::Interaction::ProcInfo(), genie::Kinematics::Q2(), genie::utils::kinematics::Q2(), Q2Lim(), Q2Lim_W(), genie::Kinematics::t(), confusionMatrixTree::t, moon_position_table_new3::tl(), TLim(), genie::utils::kinematics::UpdateWQ2FromXY(), W, genie::Kinematics::W(), WLim(), genie::Kinematics::x(), submit_syst::x, XLim(), genie::Kinematics::y(), submit_syst::y, and YLim().

Referenced by genie::XSecAlgorithmI::ValidKinematics().

260 {
261  const ProcessInfo & pi = fInteraction->ProcInfo();
262  const Kinematics & kine = fInteraction->Kine();
263 
264  // ASK single kaon:
265  // XSec code returns zero when kinematics are not allowed
266  // Here just let kinematics always be allowed
267  if(pi.IsSingleKaon()) {
268  return true;
269  }
270 
271  // QEL:
272  // Check the running Q2 vs the Q2 limits
273  if(pi.IsQuasiElastic() || pi.IsInverseBetaDecay() || pi.IsDarkMatterElastic()) {
274  Range1D_t Q2l = this->Q2Lim();
275  double Q2 = kine.Q2();
276  bool in_phys = math::IsWithinLimits(Q2, Q2l);
277  bool allowed = in_phys;
278  return allowed;
279  }
280 
281  // RES
282  // Check the running W vs the W limits
283  // & the running Q2 vs Q2 limits for the given W
284  if(pi.IsResonant()) {
285  Range1D_t Wl = this->WLim();
286  Range1D_t Q2l = this->Q2Lim_W();
287  double W = kine.W();
288  double Q2 = kine.Q2();
289  bool in_phys = (math::IsWithinLimits(Q2, Q2l) && math::IsWithinLimits(W, Wl));
290  bool allowed = in_phys;
291  return allowed;
292  }
293 
294  // DIS
295  if(pi.IsDeepInelastic() || pi.IsDarkMatterDeepInelastic()) {
296  Range1D_t Wl = this->WLim();
297  Range1D_t Q2l = this->Q2Lim_W();
298  double W = kine.W();
299  double Q2 = kine.Q2();
300  bool in_phys = (math::IsWithinLimits(Q2, Q2l) && math::IsWithinLimits(W, Wl));
301  bool allowed = in_phys;
302  return allowed;
303  }
304 
305  //IMD
306  if(pi.IsInverseMuDecay() || pi.IsIMDAnnihilation() || pi.IsNuElectronElastic()) {
307  Range1D_t yl = this->YLim();
308  double y = kine.y();
309  bool in_phys = math::IsWithinLimits(y, yl);
310  bool allowed = in_phys;
311  return allowed;
312  }
313 
314  //COH
315  if (pi.IsCoherent()) {
316  Range1D_t xl = this->XLim();
317  Range1D_t yl = this->YLim();
318  double x = kine.x();
319  double y = kine.y();
320  bool in_phys = (math::IsWithinLimits(x, xl) && math::IsWithinLimits(y, yl));
321  bool allowed = in_phys;
322  return allowed;
323  }
324 
325  // DFR
326  if (pi.IsDiffractive()) {
327  // first two checks are the same as RES & DIS
328  Range1D_t Wl = this->WLim();
329  Range1D_t Q2l = this->Q2Lim_W();
330 
332  double W = kine.W();
333  double Q2 = kine.Q2();
334 
335  LOG("KPhaseSpace", pDEBUG) << " W = " << W << ", limits = [" << Wl.min << "," << Wl.max << "];";
336  LOG("KPhaseSpace", pDEBUG) << " Q2 = " << Q2 << ", limits = [" << Q2l.min << "," << Q2l.max << "];";
337  bool in_phys = math::IsWithinLimits(W, Wl);
338  in_phys = in_phys && math::IsWithinLimits(Q2, Q2l);
339 
340  // extra check: there's a t minimum.
341  // but only check if W, Q2 is reasonable
342  // (otherwise get NaNs in tmin)
343  if (in_phys)
344  {
345  double t = kine.t();
346  Range1D_t tl = this->TLim();
347  LOG("KPhaseSpace", pDEBUG) << " t = " << t << ", limits = [" << tl.min << "," << tl.max << "];";
348  in_phys = in_phys && math::IsWithinLimits(t, tl);
349  }
350  LOG("KPhaseSpace", pDEBUG) << " phase space point is " << ( in_phys ? "ALLOWED" : "NOT ALLOWED");
351 
352 
353  bool allowed = in_phys;
354  return allowed;
355  }
356 
357  // was MECTensor
358  if (pi.IsMEC()){
359  Range1D_t Q2l = this->Q2Lim();
360  double Q2 = kine.Q2();
361  bool in_phys = math::IsWithinLimits(Q2, Q2l);
362  bool allowed = in_phys;
363  return allowed;
364  }
365 
366  return false;
367 }
bool IsResonant(void) const
Definition: ProcessInfo.cxx:92
double W(bool selected=false) const
Definition: Kinematics.cxx:167
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:991
A simple [min,max] interval for doubles.
Definition: Range1.h:43
bool IsInverseMuDecay(void) const
bool IsQuasiElastic(void) const
Definition: ProcessInfo.cxx:67
bool IsWithinLimits(double x, Range1D_t range)
Definition: MathUtils.cxx:265
Generated/set kinematical variables for an event.
Definition: Kinematics.h:40
bool IsInverseBetaDecay(void) const
double x(bool selected=false) const
Definition: Kinematics.cxx:109
Range1D_t YLim(void) const
y limits
bool IsDiffractive(void) const
bool IsIMDAnnihilation(void) const
Range1D_t Q2Lim_W(void) const
Q2 limits @ fixed W.
Range1D_t Q2Lim(void) const
Q2 limits.
double y(bool selected=false) const
Definition: Kinematics.cxx:122
bool IsSingleKaon(void) const
Definition: ProcessInfo.cxx:77
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
bool IsNuElectronElastic(void) const
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
const Kinematics & Kine(void) const
Definition: Interaction.h:71
bool IsDarkMatterElastic(void) const
Definition: ProcessInfo.cxx:72
bool IsMEC(void) const
bool IsDeepInelastic(void) const
Definition: ProcessInfo.cxx:82
void UpdateWQ2FromXY(const Interaction *in)
Definition: KineUtils.cxx:1185
double max
Definition: Range1.h:54
Range1D_t XLim(void) const
x limits
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
double min
Definition: Range1.h:53
double t(bool selected=false) const
Definition: Kinematics.cxx:180
double Q2(bool selected=false) const
Definition: Kinematics.cxx:135
Range1D_t TLim(void) const
t limits
bool IsDarkMatterDeepInelastic(void) const
Definition: ProcessInfo.cxx:87
#define W(x)
Range1D_t WLim(void) const
W limits.
#define pDEBUG
Definition: Messenger.h:64
bool IsCoherent(void) const
Definition: ProcessInfo.cxx:97
Range1D_t KPhaseSpace::Limits ( KineVar_t  kvar) const

Return the kinematical variable limits.

Definition at line 191 of file KPhaseSpace.cxx.

References ana::assert(), genie::KineVar::AsString(), fInteraction, genie::kKVQ2, genie::kKVq2, genie::kKVt, genie::kKVW, genie::kKVx, genie::kKVy, LOG, pERROR, Q2Lim(), q2Lim(), R, TLim(), WLim(), XLim(), and YLim().

Referenced by genie::ReinSehgalRESXSecWithCache::CacheResExcitationXSec(), genie::IBDKinematicsGenerator::ComputeMaxXSec(), genie::NuEKinematicsGenerator::ComputeMaxXSec(), genie::QELKinematicsGenerator::ComputeMaxXSec(), genie::DISKinematicsGenerator::ComputeMaxXSec(), genie::DMELKinematicsGenerator::ComputeMaxXSec(), genie::DMDISKinematicsGenerator::ComputeMaxXSec(), genie::IMDXSec::Integrate(), genie::COHXSecAR::Integrate(), genie::COHXSec::Integrate(), genie::RESXSec::Integrate(), genie::QELXSec::Integrate(), genie::DMELXSec::Integrate(), genie::DFRXSec::Integrate(), genie::NuElectronXSec::Integrate(), genie::SmithMonizQELCCXSec::Integrate(), genie::ReinSehgalRESXSec::Integrate(), Maximum(), Minimum(), genie::utils::kinematics::PhaseSpaceVolume(), PrintLimits(), genie::DFRKinematicsGenerator::ProcessEventRecord(), genie::NuEKinematicsGenerator::ProcessEventRecord(), genie::QELKinematicsGenerator::ProcessEventRecord(), genie::RESKinematicsGenerator::ProcessEventRecord(), genie::IBDKinematicsGenerator::ProcessEventRecord(), genie::DMELKinematicsGenerator::ProcessEventRecord(), genie::DISKinematicsGenerator::ProcessEventRecord(), genie::DMDISKinematicsGenerator::ProcessEventRecord(), genie::QELKinematicsGenerator::SpectralFuncExperimentalCode(), and genie::DMELKinematicsGenerator::SpectralFuncExperimentalCode().

192 {
193  // Compute limits for the input kinematic variable irrespective of any other
194  // relevant kinematical variable
195  //
197 
198  switch(kvar) {
199  case(kKVW) : return this->WLim(); break;
200  case(kKVQ2) : return this->Q2Lim(); break;
201  case(kKVq2) : return this->q2Lim(); break;
202  case(kKVx) : return this->XLim(); break;
203  case(kKVy) : return this->YLim(); break;
204  case(kKVt) : return this->TLim(); break;
205  default:
206  LOG("KPhaseSpace", pERROR)
207  << "Couldn't compute limits for " << KineVar::AsString(kvar);
208  Range1D_t R(-1.,-1);
209  return R;
210  }
211 }
#define pERROR
Definition: Messenger.h:60
A simple [min,max] interval for doubles.
Definition: Range1.h:43
Range1D_t q2Lim(void) const
q2 limits
Range1D_t YLim(void) const
y limits
Range1D_t Q2Lim(void) const
Q2 limits.
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
#define R(x)
static string AsString(KineVar_t kv)
Definition: KineVar.h:68
assert(nhit_max >=nhit_nbins)
Range1D_t XLim(void) const
x limits
Range1D_t TLim(void) const
t limits
Range1D_t WLim(void) const
W limits.
double KPhaseSpace::Maximum ( KineVar_t  kvar) const

Definition at line 219 of file KPhaseSpace.cxx.

References Limits(), and genie::Range1D_t::max.

220 {
221  Range1D_t lim = this->Limits(kvar);
222  return lim.max;
223 }
A simple [min,max] interval for doubles.
Definition: Range1.h:43
Range1D_t Limits(KineVar_t kvar) const
Return the kinematical variable limits.
double max
Definition: Range1.h:54
double KPhaseSpace::Minimum ( KineVar_t  kvar) const

Definition at line 213 of file KPhaseSpace.cxx.

References Limits(), and genie::Range1D_t::min.

214 {
215  Range1D_t lim = this->Limits(kvar);
216  return lim.min;
217 }
A simple [min,max] interval for doubles.
Definition: Range1.h:43
Range1D_t Limits(KineVar_t kvar) const
Return the kinematical variable limits.
double min
Definition: Range1.h:53
Range1D_t KPhaseSpace::Q2Lim ( void  ) const

Q2 limits.

Definition at line 500 of file KPhaseSpace.cxx.

References genie::XclsTag::CharmHadronPdg(), genie::utils::kinematics::CohQ2Lim(), genie::utils::kinematics::DarkQ2Lim(), genie::utils::kinematics::DarkQ2Lim_W(), genie::Interaction::ExclTag(), genie::PDGLibrary::Find(), fInteraction, genie::Interaction::FSPrimLepton(), genie::Target::HitNucP4Ptr(), genie::utils::kinematics::InelQ2Lim(), genie::utils::kinematics::InelQ2Lim_W(), genie::Interaction::InitState(), genie::PDGLibrary::Instance(), genie::XclsTag::IsCharmEvent(), genie::ProcessInfo::IsCoherent(), genie::ProcessInfo::IsDarkMatterDeepInelastic(), genie::ProcessInfo::IsDarkMatterElastic(), genie::ProcessInfo::IsDeepInelastic(), genie::ProcessInfo::IsEM(), genie::ProcessInfo::IsInverseBetaDecay(), genie::ProcessInfo::IsMEC(), genie::ProcessInfo::IsQuasiElastic(), genie::ProcessInfo::IsResonant(), genie::XclsTag::IsStrangeEvent(), genie::ProcessInfo::IsWeakCC(), genie::controls::kMinQ2Limit_VLE, kPi0Mass, genie::constants::kPionMass, genie::kRfHitNucRest, genie::Range1D_t::max, genie::Range1D_t::min, python.hepunit::pi, genie::InitialState::ProbeE(), genie::Interaction::ProcInfo(), genie::Interaction::RecoilNucleon(), genie::XclsTag::StrangeHadronPdg(), genie::InitialState::Tgt(), and W.

Referenced by genie::DISXSec::CacheFreeNucleonXSec(), genie::DMDISXSec::CacheFreeNucleonXSec(), genie::utils::ComputeFullQELPXSec(), genie::DISXSec::Integrate(), genie::DMDISXSec::Integrate(), IsAllowed(), Limits(), genie::COHKinematicsGenerator::MaxXSec_BergerSehgal(), genie::COHKinematicsGenerator::MaxXSec_BergerSehgalFM(), q2Lim(), and Q2Lim_W().

501 {
502  // Computes momentum transfer (Q2>0) limits irrespective of the invariant mass
503  // For QEL this is identical to Q2Lim_W (since W is fixed)
504  // For RES & DIS, the calculation proceeds as in kinematics::InelQ2Lim().
505  //
506  Range1D_t Q2l;
507  Q2l.min = -1;
508  Q2l.max = -1;
509 
510  const ProcessInfo & pi = fInteraction->ProcInfo();
511 
512  bool is_em = pi.IsEM();
513  bool is_qel = pi.IsQuasiElastic() || pi.IsInverseBetaDecay();
514  bool is_inel = pi.IsDeepInelastic() || pi.IsResonant();
515  bool is_coh = pi.IsCoherent();
516  bool is_dme = pi.IsDarkMatterElastic();
517  bool is_dmdis = pi.IsDarkMatterDeepInelastic();
518 
519  if(!is_qel && !is_inel && !is_coh && !is_dme && !is_dmdis) return Q2l;
520 
521  const InitialState & init_state = fInteraction->InitState();
522  double Ev = init_state.ProbeE(kRfHitNucRest);
523  double M = init_state.Tgt().HitNucP4Ptr()->M(); // can be off m/shell
524  double ml = fInteraction->FSPrimLepton()->Mass();
525 
526  if(is_coh) {
527  bool pionIsCharged = pi.IsWeakCC();
528  double mpi = pionIsCharged ? kPionMass : kPi0Mass;
529  Q2l = kinematics::CohQ2Lim(M, mpi, ml, Ev);
530  return Q2l;
531  }
532 
533  const XclsTag & xcls = fInteraction->ExclTag();
534 
535  // quasi-elastic
536  if(is_qel) {
537  double W = fInteraction->RecoilNucleon()->Mass();
538  if(xcls.IsCharmEvent()) {
539  int charm_pdgc = xcls.CharmHadronPdg();
540  W = PDGLibrary::Instance()->Find(charm_pdgc)->Mass();
541  } else if(xcls.IsStrangeEvent()) {
542  int strange_pdgc = xcls.StrangeHadronPdg();
543  W = PDGLibrary::Instance()->Find(strange_pdgc)->Mass();
544  }
545  if (pi.IsInverseBetaDecay()) {
547  } else {
548  Q2l = is_em ? kinematics::electromagnetic::InelQ2Lim_W(Ev,ml,M,W) : kinematics::InelQ2Lim_W(Ev,M,ml,W);
549  }
550 
551  return Q2l;
552  }
553 
554  // dark mattter elastic
555  if(is_dme) {
556  double W = fInteraction->RecoilNucleon()->Mass();
557  if(xcls.IsCharmEvent()) {
558  int charm_pdgc = xcls.CharmHadronPdg();
559  W = PDGLibrary::Instance()->Find(charm_pdgc)->Mass();
560  } else if(xcls.IsStrangeEvent()) {
561  int strange_pdgc = xcls.StrangeHadronPdg();
562  W = PDGLibrary::Instance()->Find(strange_pdgc)->Mass();
563  }
564  if (pi.IsInverseBetaDecay()) {
566  } else {
567  Q2l = kinematics::DarkQ2Lim_W(Ev,M,ml,W);
568  }
569 
570 
571  return Q2l;
572  }
573 
574  // was MECTensor
575  // TODO: Q2maxConfig
576  if (pi.IsMEC()){
577  double W = fInteraction->RecoilNucleon()->Mass();
578  Q2l = is_em ? kinematics::electromagnetic::InelQ2Lim_W(Ev,ml,M,W) : kinematics::InelQ2Lim_W(Ev,M,ml,W);
579  double Q2maxConfig = 1.44; // need to pull from config file somehow?
580  if (Q2l.max > Q2maxConfig) Q2l.max = Q2maxConfig;
581  return Q2l;
582  }
583 
584  if (is_dmdis) {
585  Q2l = kinematics::DarkQ2Lim(Ev,M,ml);
586  return Q2l;
587  }
588 
589  // inelastic
590  Q2l = is_em ? kinematics::electromagnetic::InelQ2Lim(Ev,ml,M) : kinematics::InelQ2Lim(Ev,M,ml);
591  return Q2l;
592 }
bool IsResonant(void) const
Definition: ProcessInfo.cxx:92
const Var kPi0Mass
bool IsWeakCC(void) const
static const double kMinQ2Limit_VLE
Definition: Controls.h:42
TParticlePDG * RecoilNucleon(void) const
recoil nucleon
A simple [min,max] interval for doubles.
Definition: Range1.h:43
bool IsQuasiElastic(void) const
Definition: ProcessInfo.cxx:67
int CharmHadronPdg(void) const
Definition: XclsTag.h:50
bool IsStrangeEvent(void) const
Definition: XclsTag.h:51
Range1D_t DarkQ2Lim_W(double Ev, double M, double ml, double W, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:829
bool IsInverseBetaDecay(void) const
Range1D_t InelQ2Lim(double Ev, double M, double ml, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:359
Contains minimal information for tagging exclusive processes.
Definition: XclsTag.h:37
Range1D_t InelQ2Lim_W(double Ev, double M, double ml, double W, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:309
bool IsCharmEvent(void) const
Definition: XclsTag.h:48
Range1D_t DarkQ2Lim(double Ev, double M, double ml, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:887
int StrangeHadronPdg(void) const
Definition: XclsTag.h:53
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
Range1D_t CohQ2Lim(double Mn, double mpi, double mlep, double Ev)
Definition: KineUtils.cxx:673
bool IsDarkMatterElastic(void) const
Definition: ProcessInfo.cxx:72
TParticlePDG * FSPrimLepton(void) const
final state primary lepton
bool IsMEC(void) const
bool IsEM(void) const
bool IsDeepInelastic(void) const
Definition: ProcessInfo.cxx:82
double max
Definition: Range1.h:54
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:264
static PDGLibrary * Instance(void)
Definition: PDGLibrary.cxx:43
static const double kPionMass
Definition: Constants.h:74
const XclsTag & ExclTag(void) const
Definition: Interaction.h:72
TParticlePDG * Find(int pdgc)
Definition: PDGLibrary.cxx:61
const InitialState & InitState(void) const
Definition: Interaction.h:69
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
double min
Definition: Range1.h:53
const Target & Tgt(void) const
Definition: InitialState.h:67
bool IsDarkMatterDeepInelastic(void) const
Definition: ProcessInfo.cxx:87
void kinematics()
Definition: kinematics.C:10
#define W(x)
double ProbeE(RefFrame_t rf) const
Initial State information.
Definition: InitialState.h:49
bool IsCoherent(void) const
Definition: ProcessInfo.cxx:97
Range1D_t KPhaseSpace::q2Lim ( void  ) const

q2 limits

Definition at line 594 of file KPhaseSpace.cxx.

References genie::Range1D_t::max, genie::Range1D_t::min, genie::utils::kinematics::Q2(), q2, and Q2Lim().

Referenced by Limits().

595 {
596 // As Q2Lim(void) but with reversed sign (Q2 -> q2)
597 //
598  Range1D_t Q2 = this->Q2Lim();
599  Range1D_t q2;
600  q2.min = - Q2.max;
601  q2.max = - Q2.min;
602  return q2;
603 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:991
A simple [min,max] interval for doubles.
Definition: Range1.h:43
Range1D_t Q2Lim(void) const
Q2 limits.
Double_t q2[12][num]
Definition: f2_nu.C:137
double max
Definition: Range1.h:54
double min
Definition: Range1.h:53
Range1D_t KPhaseSpace::Q2Lim_W ( void  ) const

Q2 limits @ fixed W.

Definition at line 439 of file KPhaseSpace.cxx.

References genie::utils::kinematics::DarkQ2Lim_W(), fInteraction, genie::Interaction::FSPrimLepton(), genie::Target::HitNucP4Ptr(), genie::utils::kinematics::InelQ2Lim_W(), genie::Interaction::InitState(), genie::ProcessInfo::IsCoherent(), genie::ProcessInfo::IsDarkMatterDeepInelastic(), genie::ProcessInfo::IsDarkMatterElastic(), genie::ProcessInfo::IsDeepInelastic(), genie::ProcessInfo::IsDiffractive(), genie::ProcessInfo::IsEM(), genie::ProcessInfo::IsInverseBetaDecay(), genie::ProcessInfo::IsQuasiElastic(), genie::ProcessInfo::IsResonant(), genie::controls::kMinQ2Limit_VLE, genie::kRfHitNucRest, genie::Range1D_t::max, genie::Range1D_t::min, python.hepunit::pi, genie::InitialState::ProbeE(), genie::Interaction::ProcInfo(), Q2Lim(), genie::Interaction::RecoilNucleon(), genie::InitialState::Tgt(), and W.

Referenced by genie::DFRKinematicsGenerator::ComputeMaxXSec(), genie::RESKinematicsGenerator::ComputeMaxXSec(), genie::utils::gsl::d2XSecRESFast_dWQ2_E::DoEval(), IsAllowed(), genie::utils::kinematics::PhaseSpaceVolume(), genie::RESKinematicsGenerator::ProcessEventRecord(), and q2Lim_W().

440 {
441  // Computes momentum transfer (Q2>0) limits @ the input invariant mass
442  // The calculation proceeds as in kinematics::InelQ2Lim_W().
443  // For QEL, W is set to the recoil nucleon mass
444  //
445  // TODO: For now, choosing to handle Q2 at fixed W for coherent in the
446  // same way as for the general Q2 limits... but shouldn't we just use
447  // W = m_pi? - which we do in Q2Lim() anyway... seems like there are
448  // cleanup opportunities here.
449 
450  Range1D_t Q2l;
451  Q2l.min = -1;
452  Q2l.max = -1;
453 
454  const ProcessInfo & pi = fInteraction->ProcInfo();
455 
456  bool is_em = pi.IsEM();
457  bool is_qel = pi.IsQuasiElastic() || pi.IsInverseBetaDecay();
458  bool is_inel = pi.IsDeepInelastic() || pi.IsResonant() || pi.IsDiffractive();
459  bool is_coh = pi.IsCoherent();
460  bool is_dme = pi.IsDarkMatterElastic();
461  bool is_dmdis = pi.IsDarkMatterDeepInelastic();
462 
463  if(!is_qel && !is_inel && !is_coh && !is_dme && !is_dmdis) return Q2l;
464 
465  if(is_coh) {
466  return Q2Lim();
467  }
468 
469  const InitialState & init_state = fInteraction->InitState();
470  double Ev = init_state.ProbeE(kRfHitNucRest);
471  double M = init_state.Tgt().HitNucP4Ptr()->M(); // can be off m/shell
472  double ml = fInteraction->FSPrimLepton()->Mass();
473 
474  double W = 0;
475  if(is_qel || is_dme) W = fInteraction->RecoilNucleon()->Mass();
476  else W = kinematics::W(fInteraction);
477 
478  if (pi.IsInverseBetaDecay()) {
480  } else if (is_dme || is_dmdis) {
481  Q2l = kinematics::DarkQ2Lim_W(Ev,M,ml,W);
482  } else {
483  Q2l = is_em ? kinematics::electromagnetic::InelQ2Lim_W(Ev,ml,M,W) : Q2l = kinematics::InelQ2Lim_W(Ev,M,ml,W);
484  }
485 
486  return Q2l;
487 }
bool IsResonant(void) const
Definition: ProcessInfo.cxx:92
static const double kMinQ2Limit_VLE
Definition: Controls.h:42
TParticlePDG * RecoilNucleon(void) const
recoil nucleon
A simple [min,max] interval for doubles.
Definition: Range1.h:43
bool IsQuasiElastic(void) const
Definition: ProcessInfo.cxx:67
Range1D_t DarkQ2Lim_W(double Ev, double M, double ml, double W, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:829
bool IsInverseBetaDecay(void) const
bool IsDiffractive(void) const
Range1D_t Q2Lim(void) const
Q2 limits.
Range1D_t InelQ2Lim_W(double Ev, double M, double ml, double W, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:309
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
bool IsDarkMatterElastic(void) const
Definition: ProcessInfo.cxx:72
TParticlePDG * FSPrimLepton(void) const
final state primary lepton
bool IsEM(void) const
bool IsDeepInelastic(void) const
Definition: ProcessInfo.cxx:82
double max
Definition: Range1.h:54
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:264
const InitialState & InitState(void) const
Definition: Interaction.h:69
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
double min
Definition: Range1.h:53
const Target & Tgt(void) const
Definition: InitialState.h:67
bool IsDarkMatterDeepInelastic(void) const
Definition: ProcessInfo.cxx:87
void kinematics()
Definition: kinematics.C:10
#define W(x)
double ProbeE(RefFrame_t rf) const
Initial State information.
Definition: InitialState.h:49
bool IsCoherent(void) const
Definition: ProcessInfo.cxx:97
Range1D_t KPhaseSpace::q2Lim_W ( void  ) const

q2 limits @ fixed W

Definition at line 489 of file KPhaseSpace.cxx.

References genie::Range1D_t::max, genie::Range1D_t::min, genie::utils::kinematics::Q2(), q2, and Q2Lim_W().

490 {
491 // As Q2Lim_W(void) but with reversed sign (Q2 -> q2)
492 //
493  Range1D_t Q2 = this->Q2Lim_W();
494  Range1D_t q2;
495  q2.min = - Q2.max;
496  q2.max = - Q2.min;
497  return q2;
498 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:991
A simple [min,max] interval for doubles.
Definition: Range1.h:43
Range1D_t Q2Lim_W(void) const
Q2 limits @ fixed W.
Double_t q2[12][num]
Definition: f2_nu.C:137
double max
Definition: Range1.h:54
double min
Definition: Range1.h:53
double KPhaseSpace::Threshold ( void  ) const

Energy threshold.

Definition at line 81 of file KPhaseSpace.cxx.

References ana::assert(), genie::XclsTag::CharmHadronPdg(), genie::Interaction::ExclTag(), exit(), genie::PDGLibrary::Find(), fInteraction, genie::Interaction::FSPrimLepton(), genie::Target::HitNucIsSet(), genie::Target::HitNucP4Ptr(), genie::Target::HitNucPdg(), genie::Interaction::InitState(), genie::PDGLibrary::Instance(), genie::ProcessInfo::IsAMNuGamma(), genie::XclsTag::IsCharmEvent(), genie::ProcessInfo::IsCoherent(), genie::ProcessInfo::IsDarkMatterDeepInelastic(), genie::ProcessInfo::IsDarkMatterElastic(), genie::ProcessInfo::IsDeepInelastic(), genie::ProcessInfo::IsDiffractive(), genie::ProcessInfo::IsGlashowResonance(), genie::ProcessInfo::IsIMDAnnihilation(), genie::XclsTag::IsInclusiveCharm(), genie::ProcessInfo::IsInverseBetaDecay(), genie::ProcessInfo::IsInverseMuDecay(), genie::ProcessInfo::IsMEC(), genie::ProcessInfo::IsNuElectronElastic(), genie::ProcessInfo::IsQuasiElastic(), genie::ProcessInfo::IsResonant(), genie::ProcessInfo::IsSingleKaon(), genie::ProcessInfo::IsWeakCC(), genie::controls::kASmallNum, genie::constants::kElectronMass, genie::constants::kElectronMass2, genie::constants::kLightestChmHad, genie::constants::kMuonMass2, genie::constants::kNeutronMass, genie::constants::kNucleonMass, genie::constants::kPhotontest, kPi0Mass, genie::constants::kPionMass, genie::constants::kProtonMass, m, Munits::m2, genie::XclsTag::NProtons(), genie::Target::Pdg(), pERROR, python.hepunit::pi, genie::Interaction::ProcInfo(), genie::Interaction::RecoilNucleon(), SLOG, genie::XclsTag::StrangeHadronPdg(), genie::pdg::SwitchProtonNeutron(), genie::InitialState::Tgt(), and Wmin.

Referenced by genie::DISXSec::CacheFreeNucleonXSec(), genie::DMDISXSec::CacheFreeNucleonXSec(), genie::ReinSehgalRESXSecWithCache::CacheResExcitationXSec(), genie::ReinSehgalRESXSecWithCacheFast::CacheResExcitationXSec(), genie::XSecSplineList::CreateSpline(), IsAboveThreshold(), and genie::FermiMover::KickHitNucleon().

82 {
83  const ProcessInfo & pi = fInteraction->ProcInfo();
84  const InitialState & init_state = fInteraction->InitState();
85  const XclsTag & xcls = fInteraction->ExclTag();
86  const Target & tgt = init_state.Tgt();
87 
88  double ml = fInteraction->FSPrimLepton()->Mass();
89 
90  if (pi.IsSingleKaon()) {
91  int kaon_pdgc = xcls.StrangeHadronPdg();
92  double Mi = tgt.HitNucP4Ptr()->M(); // initial nucleon mass
93  // Final nucleon can be different for K0 interaction
94  double Mf = (xcls.NProtons()==1) ? kProtonMass : kNeutronMass;
95  double mk = PDGLibrary::Instance()->Find(kaon_pdgc)->Mass();
96  //double ml = PDGLibrary::Instance()->Find(fInteraction->FSPrimLeptonPdg())->Mass();
97  double mtot = ml + mk + Mf; // total mass of FS particles
98  double Ethresh = (mtot*mtot - Mi*Mi)/(2. * Mf);
99  return Ethresh;
100  }
101 
102  if (pi.IsCoherent()) {
103  int tgtpdgc = tgt.Pdg(); // nuclear target PDG code (10LZZZAAAI)
104  double mpi = pi.IsWeakCC() ? kPionMass : kPi0Mass;
105  double MA = PDGLibrary::Instance()->Find(tgtpdgc)->Mass();
106  double m = ml + mpi;
107  double m2 = TMath::Power(m,2);
108  double Ethr = m + 0.5*m2/MA;
109  return TMath::Max(0.,Ethr);
110  }
111 
112  if(pi.IsQuasiElastic() ||
113  pi.IsDarkMatterElastic() ||
114  pi.IsInverseBetaDecay() ||
115  pi.IsResonant() ||
116  pi.IsDeepInelastic() ||
118  pi.IsDiffractive())
119  {
120  assert(tgt.HitNucIsSet());
121  double Mn = tgt.HitNucP4Ptr()->M();
122  double Mn2 = TMath::Power(Mn,2);
123  double Wmin = kNucleonMass + kPionMass;
124  if ( pi.IsQuasiElastic() || pi.IsDarkMatterElastic() || pi.IsInverseBetaDecay() ) {
125  int finalNucPDG = tgt.HitNucPdg();
126  if ( pi.IsWeakCC() ) finalNucPDG = pdg::SwitchProtonNeutron( finalNucPDG );
127  Wmin = PDGLibrary::Instance()->Find( finalNucPDG )->Mass();
128  }
129  if (pi.IsResonant()) {
130  Wmin = kNucleonMass + kPhotontest;
131  }
132 
133  if(xcls.IsCharmEvent()) {
134  if(xcls.IsInclusiveCharm()) {
136  } else {
137  int cpdg = xcls.CharmHadronPdg();
138  double mchm = PDGLibrary::Instance()->Find(cpdg)->Mass();
139  if(pi.IsQuasiElastic() || pi.IsInverseBetaDecay()) {
140  Wmin = mchm + controls::kASmallNum;
141  }
142  else {
143  Wmin = kNeutronMass + mchm + controls::kASmallNum;
144  }
145  }//incl.?
146  }//charm?
147 
148  double smin = TMath::Power(Wmin+ml,2.);
149  double Ethr = 0.5*(smin-Mn2)/Mn;
150  // threshold is different for dark matter case
152  Ethr = TMath::Max(0.5*(smin-Mn2-ml*ml)/Mn,ml);
153  }
154 
155  return TMath::Max(0.,Ethr);
156  }
157 
158  if(pi.IsInverseMuDecay() || pi.IsIMDAnnihilation()) {
159  double Ethr = 0.5 * (kMuonMass2-kElectronMass2)/kElectronMass;
160  return TMath::Max(0.,Ethr);
161  }
162 
163  if(pi.IsNuElectronElastic() || pi.IsGlashowResonance() ) {
164  return 0;
165  }
166  if(pi.IsAMNuGamma()) {
167  return 0;
168  }
169  if (pi.IsMEC()) {
170  if (tgt.HitNucIsSet()) {
171  double Mn = tgt.HitNucP4Ptr()->M();
172  double Mn2 = TMath::Power(Mn,2);
173  double Wmin = fInteraction->RecoilNucleon()->Mass(); // mass of the recoil nucleon cluster
174  double smin = TMath::Power(Wmin+ml,2.);
175  double Ethr = 0.5*(smin-Mn2)/Mn;
176  return TMath::Max(0.,Ethr);
177  }
178  else {
179  // this was ... if (pi.IsMECTensor())
180  return ml;
181  }
182  }
183 
184  SLOG("KPhaseSpace", pERROR)
185  << "Can't compute threshold for \n" << *fInteraction;
186  exit(1);
187 
188  return 99999999;
189 }
bool IsResonant(void) const
Definition: ProcessInfo.cxx:92
const Var kPi0Mass
bool IsWeakCC(void) const
#define pERROR
Definition: Messenger.h:60
static const double kNucleonMass
Definition: Constants.h:78
int HitNucPdg(void) const
Definition: Target.cxx:321
TParticlePDG * RecoilNucleon(void) const
recoil nucleon
bool IsInverseMuDecay(void) const
bool IsQuasiElastic(void) const
Definition: ProcessInfo.cxx:67
int CharmHadronPdg(void) const
Definition: XclsTag.h:50
int Pdg(void) const
Definition: Target.h:72
static const double kPhotontest
Definition: Constants.h:80
int SwitchProtonNeutron(int pdgc)
Definition: PDGUtils.cxx:319
bool IsInverseBetaDecay(void) const
static const double kLightestChmHad
Definition: Constants.h:79
bool IsDiffractive(void) const
bool IsIMDAnnihilation(void) const
static const double kElectronMass
Definition: Constants.h:71
Contains minimal information for tagging exclusive processes.
Definition: XclsTag.h:37
bool IsCharmEvent(void) const
Definition: XclsTag.h:48
bool IsSingleKaon(void) const
Definition: ProcessInfo.cxx:77
int StrangeHadronPdg(void) const
Definition: XclsTag.h:53
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
static const double kElectronMass2
Definition: Constants.h:84
bool IsNuElectronElastic(void) const
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
bool IsAMNuGamma(void) const
static const double kMuonMass2
Definition: Constants.h:85
A Neutrino Interaction Target. Is a transparent encapsulation of quite different physical systems suc...
Definition: Target.h:41
static const double kNeutronMass
Definition: Constants.h:77
static const double kASmallNum
Definition: Controls.h:40
bool IsDarkMatterElastic(void) const
Definition: ProcessInfo.cxx:72
TParticlePDG * FSPrimLepton(void) const
final state primary lepton
bool IsMEC(void) const
bool IsDeepInelastic(void) const
Definition: ProcessInfo.cxx:82
static constexpr Double_t m2
Definition: Munits.h:145
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:264
static PDGLibrary * Instance(void)
Definition: PDGLibrary.cxx:43
static const double kPionMass
Definition: Constants.h:74
bool HitNucIsSet(void) const
Definition: Target.cxx:300
const XclsTag & ExclTag(void) const
Definition: Interaction.h:72
bool IsInclusiveCharm(void) const
Definition: XclsTag.cxx:63
TParticlePDG * Find(int pdgc)
Definition: PDGLibrary.cxx:61
exit(0)
assert(nhit_max >=nhit_nbins)
const InitialState & InitState(void) const
Definition: Interaction.h:69
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
double Wmin[kNWBins]
const Target & Tgt(void) const
Definition: InitialState.h:67
bool IsGlashowResonance(void) const
#define SLOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a short string (using the FUNCTION and...
Definition: Messenger.h:85
bool IsDarkMatterDeepInelastic(void) const
Definition: ProcessInfo.cxx:87
static const double kProtonMass
Definition: Constants.h:76
Initial State information.
Definition: InitialState.h:49
bool IsCoherent(void) const
Definition: ProcessInfo.cxx:97
int NProtons(void) const
Definition: XclsTag.h:54
Range1D_t KPhaseSpace::TLim ( void  ) const

t limits

Definition at line 809 of file KPhaseSpace.cxx.

References std::abs(), fInteraction, GetTMaxDFR(), genie::Target::HitNucMass(), genie::Interaction::InitState(), genie::ProcessInfo::IsCoherent(), genie::ProcessInfo::IsDiffractive(), std::isnan(), genie::ProcessInfo::IsWeakCC(), genie::Interaction::Kine(), kPi0Mass, genie::constants::kPionMass, genie::kRfHitNucRest, LOG, genie::Range1D_t::max, genie::Range1D_t::min, pERROR, python.hepunit::pi, genie::InitialState::ProbeE(), genie::Interaction::ProcInfo(), pWARN, genie::Kinematics::Q2(), genie::utils::kinematics::Q2(), std::sqrt(), genie::InitialState::Tgt(), moon_position_table_new3::tl(), genie::utils::kinematics::UpdateWQ2FromXY(), and genie::Kinematics::y().

Referenced by genie::DFRKinematicsGenerator::ComputeMaxXSec(), IsAllowed(), and Limits().

810 {
811  // t limits for Coherent pion production from
812  // Kartavtsev, Paschos, and Gounaris, PRD 74 054007, and
813  // Paschos and Schalla, PRD 80, 03305
814  // TODO: Attempt to assign t bounds for other reactions?
815  Range1D_t tl;
816  tl.min = -1;
817  tl.max = -1;
818 
819  const InitialState & init_state = fInteraction->InitState();
820  const ProcessInfo & pi = fInteraction->ProcInfo();
821  const Kinematics & kine = fInteraction->Kine();
823  double Ev = init_state.ProbeE(kRfHitNucRest);
824  double Q2 = kine.Q2();
825  double nu = Ev * kine.y();
826  bool pionIsCharged = pi.IsWeakCC();
827  double mpi = pionIsCharged ? kPionMass : kPi0Mass;
828  double mpi2 = mpi*mpi;
829 
830  //COH
831  if(pi.IsCoherent()) {
832  tl.min = 1.0 * (Q2 + mpi2)/(2.0 * nu) * (Q2 + mpi2)/(2.0 * nu);
833  tl.max = 0.05;
834  return tl;
835  }
836  // DFR
837  else if (pi.IsDiffractive()) {
838  // diffractive tmin from Nucl.Phys.B278,61 (1986), eq. 12
839  double M = init_state.Tgt().HitNucMass();
840  double M2 = M*M;
841  double nuSqPlusQ2 = nu*nu + Q2;
842  double nuOverM = nu / M;
843  double mpiQ2term = mpi2 - Q2 - 2*nu*nu;
844  double A1 = 1 + 2*nuOverM + nuOverM*nuOverM - nuSqPlusQ2/M2;
845  double A2 = (1+nuOverM) * mpiQ2term + 2*nuOverM*nuSqPlusQ2;
846  double A3 = mpiQ2term*mpiQ2term - 4*nuSqPlusQ2*(nu*nu - mpi2);
847 
848  tl.min = std::abs( (A2 + sqrt(A2*A2 - A1*A3)) / A1 ); // GENIE's convention is that t is positive
849  bool tminIsNaN;
850  // use std::isnan when C++11 is around
851 #if __cplusplus >= 201103L
852  tminIsNaN = std::isnan(tl.min);
853 #else
854  // this the old-fashioned way to check for NaN:
855  // NaN's aren't equal to anything, including themselves
856  tminIsNaN = tl.min != tl.min;
857 #endif
858  if (tminIsNaN)
859  {
860  LOG("KPhaseSpace", pERROR)
861  << "tmin for diffractive scattering is NaN "
862  << "( Enu = " << Ev << ", Q2 = " << Q2 << ", nu = " << nu << ")";
863  throw genie::exceptions::InteractionException("NaN tmin for diffractive scattering");
864  }
865  tl.max = this->GetTMaxDFR();
866 
867  return tl;
868  }
869 
870  // RES+DIS
871  // IMD
872  LOG("KPhaseSpace", pWARN) << "It is not sensible to ask for t limits for events that are not coherent or diffractive.";
873  return tl;
874 }
const Var kPi0Mass
bool IsWeakCC(void) const
#define pERROR
Definition: Messenger.h:60
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:991
A simple [min,max] interval for doubles.
Definition: Range1.h:43
double HitNucMass(void) const
Definition: Target.cxx:250
T sqrt(T number)
Definition: d0nt_math.hpp:156
Generated/set kinematical variables for an event.
Definition: Kinematics.h:40
bool IsDiffractive(void) const
float abs(float number)
Definition: d0nt_math.hpp:39
static double GetTMaxDFR()
Definition: KPhaseSpace.cxx:58
int isnan(const stan::math::var &a)
Definition: std_isnan.hpp:18
double y(bool selected=false) const
Definition: Kinematics.cxx:122
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
Exception used inside Interaction classes.
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
const Kinematics & Kine(void) const
Definition: Interaction.h:71
#define pWARN
Definition: Messenger.h:61
void UpdateWQ2FromXY(const Interaction *in)
Definition: KineUtils.cxx:1185
double max
Definition: Range1.h:54
static const double kPionMass
Definition: Constants.h:74
const InitialState & InitState(void) const
Definition: Interaction.h:69
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
double min
Definition: Range1.h:53
double Q2(bool selected=false) const
Definition: Kinematics.cxx:135
const Target & Tgt(void) const
Definition: InitialState.h:67
double ProbeE(RefFrame_t rf) const
Initial State information.
Definition: InitialState.h:49
bool IsCoherent(void) const
Definition: ProcessInfo.cxx:97
void KPhaseSpace::UseInteraction ( const Interaction in)

Definition at line 76 of file KPhaseSpace.cxx.

References fInteraction, and in.

Referenced by KPhaseSpace().

77 {
78  fInteraction = in;
79 }
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
ifstream in
Definition: comparison.C:7
Range1D_t KPhaseSpace::WLim ( void  ) const

W limits.

Definition at line 369 of file KPhaseSpace.cxx.

References genie::utils::kinematics::DarkWLim(), genie::Interaction::ExclTag(), fInteraction, genie::Interaction::FSPrimLepton(), genie::Target::HitNucP4Ptr(), genie::utils::kinematics::InelWLim(), genie::Interaction::InitState(), genie::XclsTag::IsCharmEvent(), genie::ProcessInfo::IsDarkMatterDeepInelastic(), genie::ProcessInfo::IsDarkMatterElastic(), genie::ProcessInfo::IsDeepInelastic(), genie::ProcessInfo::IsDiffractive(), genie::ProcessInfo::IsEM(), genie::ProcessInfo::IsInverseBetaDecay(), genie::ProcessInfo::IsQuasiElastic(), genie::ProcessInfo::IsResonant(), genie::constants::kLightestChmHad, genie::constants::kNeutronMass, genie::constants::kPionMass, genie::kRfHitNucRest, LOG, genie::Range1D_t::max, genie::Range1D_t::min, pDEBUG, python.hepunit::pi, genie::InitialState::ProbeE(), genie::Interaction::ProcInfo(), genie::Interaction::RecoilNucleon(), and genie::InitialState::Tgt().

Referenced by genie::DISXSec::CacheFreeNucleonXSec(), genie::DMDISXSec::CacheFreeNucleonXSec(), genie::DFRKinematicsGenerator::ComputeMaxXSec(), genie::RESKinematicsGenerator::ComputeMaxXSec(), genie::utils::gsl::d2XSecRESFast_dWQ2_E::d2XSecRESFast_dWQ2_E(), genie::DISXSec::Integrate(), genie::DMDISXSec::Integrate(), IsAllowed(), and Limits().

370 {
371 // Computes hadronic invariant mass limits.
372 // For QEL the range reduces to the recoil nucleon mass.
373 // For DIS & RES the calculation proceeds as in kinematics::InelWLim().
374 // It is not computed for other interactions
375 //
376  Range1D_t Wl;
377  Wl.min = -1;
378  Wl.max = -1;
379 
380  const ProcessInfo & pi = fInteraction->ProcInfo();
381 
382  bool is_em = pi.IsEM();
383  bool is_qel = pi.IsQuasiElastic() || pi.IsInverseBetaDecay() || pi.IsDarkMatterElastic();
384  bool is_inel = pi.IsDeepInelastic() || pi.IsResonant() || pi.IsDiffractive();
385  bool is_dmdis = pi.IsDarkMatterDeepInelastic();
386 
387  if(is_qel) {
388  double MR = fInteraction->RecoilNucleon()->Mass();
389  Wl.min = MR;
390  Wl.max = MR;
391  return Wl;
392  }
393  if(is_inel) {
394  const InitialState & init_state = fInteraction->InitState();
395  double Ev = init_state.ProbeE(kRfHitNucRest);
396  double M = init_state.Tgt().HitNucP4Ptr()->M(); //can be off m/shell
397  double ml = fInteraction->FSPrimLepton()->Mass();
398 
399  Wl = is_em ? kinematics::electromagnetic::InelWLim(Ev,ml,M) : kinematics::InelWLim(Ev,M,ml);
400 
402  //Wl.min = TMath::Max(Wl.min, kNeutronMass+kPionMass+kLightestChmHad);
403  Wl.min = TMath::Max(Wl.min, kNeutronMass+kLightestChmHad);
404  }
405  else if (fInteraction->ProcInfo().IsDiffractive())
406  Wl.min = TMath::Max(Wl.min, kNeutronMass+kPionMass);
407  else if ( fInteraction->ProcInfo().IsDeepInelastic() ) {
408  Wl.min = TMath::Max(Wl.min, kNeutronMass + kPionMass );
409  }
410 
411  // sanity check
412  if(Wl.min>Wl.max) {Wl.min=-1; Wl.max=-1;}
413 
414  return Wl;
415  }
416  if(is_dmdis) {
417  const InitialState & init_state = fInteraction->InitState();
418  double Ev = init_state.ProbeE(kRfHitNucRest);
419  double M = init_state.Tgt().HitNucP4Ptr()->M(); //can be off m/shell
420  double ml = fInteraction->FSPrimLepton()->Mass();
421  Wl = kinematics::DarkWLim(Ev,M,ml);
423  //Wl.min = TMath::Max(Wl.min, kNeutronMass+kPionMass+kLightestChmHad);
424  Wl.min = TMath::Max(Wl.min, kNeutronMass+kLightestChmHad);
425  }
426  else if (fInteraction->ProcInfo().IsDiffractive())
427  Wl.min = TMath::Max(Wl.min, kNeutronMass+kPionMass);
428 
429  LOG("KPhaseSpace", pDEBUG) << "Found nominal limits: " << Wl.min << ", " << Wl.max;
430 
431  // sanity check
432  if(Wl.min>Wl.max) {Wl.min=-1; Wl.max=-1;}
433 
434  return Wl;
435  }
436  return Wl;
437 }
bool IsResonant(void) const
Definition: ProcessInfo.cxx:92
TParticlePDG * RecoilNucleon(void) const
recoil nucleon
Range1D_t InelWLim(double Ev, double M, double ml)
Definition: KineUtils.cxx:288
A simple [min,max] interval for doubles.
Definition: Range1.h:43
bool IsQuasiElastic(void) const
Definition: ProcessInfo.cxx:67
bool IsInverseBetaDecay(void) const
static const double kLightestChmHad
Definition: Constants.h:79
bool IsDiffractive(void) const
bool IsCharmEvent(void) const
Definition: XclsTag.h:48
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
static const double kNeutronMass
Definition: Constants.h:77
bool IsDarkMatterElastic(void) const
Definition: ProcessInfo.cxx:72
TParticlePDG * FSPrimLepton(void) const
final state primary lepton
bool IsEM(void) const
bool IsDeepInelastic(void) const
Definition: ProcessInfo.cxx:82
double max
Definition: Range1.h:54
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:264
static const double kPionMass
Definition: Constants.h:74
const XclsTag & ExclTag(void) const
Definition: Interaction.h:72
const InitialState & InitState(void) const
Definition: Interaction.h:69
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
double min
Definition: Range1.h:53
const Target & Tgt(void) const
Definition: InitialState.h:67
bool IsDarkMatterDeepInelastic(void) const
Definition: ProcessInfo.cxx:87
void kinematics()
Definition: kinematics.C:10
double ProbeE(RefFrame_t rf) const
Range1D_t DarkWLim(double Ev, double M, double ml)
Definition: KineUtils.cxx:806
Initial State information.
Definition: InitialState.h:49
#define pDEBUG
Definition: Messenger.h:64
Range1D_t KPhaseSpace::XLim ( void  ) const

x limits

Definition at line 605 of file KPhaseSpace.cxx.

References genie::utils::kinematics::CohXLim(), genie::utils::kinematics::DarkXLim(), fInteraction, genie::Interaction::FSPrimLepton(), genie::Target::HitNucP4Ptr(), genie::utils::kinematics::InelXLim(), genie::Interaction::InitState(), genie::ProcessInfo::IsCoherent(), genie::ProcessInfo::IsDarkMatterDeepInelastic(), genie::ProcessInfo::IsDarkMatterElastic(), genie::ProcessInfo::IsDeepInelastic(), genie::ProcessInfo::IsDiffractive(), genie::ProcessInfo::IsEM(), genie::ProcessInfo::IsInverseBetaDecay(), genie::ProcessInfo::IsQuasiElastic(), genie::ProcessInfo::IsResonant(), genie::controls::kASmallNum, genie::kRfHitNucRest, genie::Range1D_t::max, genie::Range1D_t::min, python.hepunit::pi, genie::InitialState::ProbeE(), genie::Interaction::ProcInfo(), and genie::InitialState::Tgt().

Referenced by genie::DFRKinematicsGenerator::ComputeMaxXSec(), genie::utils::gsl::d5XSecAR::DoEval(), genie::utils::gsl::d5Xsec_dEldOmegaldOmegapi::DoEval(), genie::utils::gsl::d4Xsec_dEldThetaldOmegapi::DoEval(), genie::utils::gsl::d3Xsec_dOmegaldThetapi::DoEval(), IsAllowed(), and Limits().

606 {
607  // Computes x-limits;
608 
609  Range1D_t xl;
610  xl.min = -1;
611  xl.max = -1;
612 
613  const ProcessInfo & pi = fInteraction->ProcInfo();
614  bool is_em = pi.IsEM();
615 
616  //RES+DIS
617  bool is_inel = pi.IsDeepInelastic() || pi.IsResonant();
618  if(is_inel) {
619  const InitialState & init_state = fInteraction->InitState();
620  double Ev = init_state.ProbeE(kRfHitNucRest);
621  double M = init_state.Tgt().HitNucP4Ptr()->M(); // can be off m/shell
622  double ml = fInteraction->FSPrimLepton()->Mass();
623  xl = is_em ? kinematics::electromagnetic::InelXLim(Ev,ml,M) : kinematics::InelXLim(Ev,M,ml);
624  return xl;
625  }
626  //DMDIS
627  bool is_dmdis = pi.IsDarkMatterDeepInelastic();
628  if(is_dmdis) {
629  const InitialState & init_state = fInteraction->InitState();
630  double Ev = init_state.ProbeE(kRfHitNucRest);
631  double M = init_state.Tgt().HitNucP4Ptr()->M(); // can be off m/shell
632  double ml = fInteraction->FSPrimLepton()->Mass();
633  xl = kinematics::DarkXLim(Ev,M,ml);
634  return xl;
635  }
636  //COH
637  bool is_coh = pi.IsCoherent();
638  if(is_coh) {
639  xl = kinematics::CohXLim();
640  return xl;
641  }
642  //QEL
643  bool is_qel = pi.IsQuasiElastic() || pi.IsInverseBetaDecay() || pi.IsDarkMatterElastic();
644  if(is_qel) {
645  xl.min = 1;
646  xl.max = 1;
647  return xl;
648  }
649  bool is_dfr = pi.IsDiffractive();
650  if(is_dfr) {
652  xl.max = 1. - controls::kASmallNum;
653  return xl;
654  }
655 
656  return xl;
657 }
bool IsResonant(void) const
Definition: ProcessInfo.cxx:92
A simple [min,max] interval for doubles.
Definition: Range1.h:43
bool IsQuasiElastic(void) const
Definition: ProcessInfo.cxx:67
bool IsInverseBetaDecay(void) const
bool IsDiffractive(void) const
Range1D_t InelXLim(double Ev, double M, double ml)
Definition: KineUtils.cxx:387
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
Range1D_t CohXLim(void)
Definition: KineUtils.cxx:665
static const double kASmallNum
Definition: Controls.h:40
bool IsDarkMatterElastic(void) const
Definition: ProcessInfo.cxx:72
TParticlePDG * FSPrimLepton(void) const
final state primary lepton
bool IsEM(void) const
bool IsDeepInelastic(void) const
Definition: ProcessInfo.cxx:82
double max
Definition: Range1.h:54
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:264
const InitialState & InitState(void) const
Definition: Interaction.h:69
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
double min
Definition: Range1.h:53
const Target & Tgt(void) const
Definition: InitialState.h:67
bool IsDarkMatterDeepInelastic(void) const
Definition: ProcessInfo.cxx:87
void kinematics()
Definition: kinematics.C:10
double ProbeE(RefFrame_t rf) const
Range1D_t DarkXLim(double Ev, double M, double ml)
Definition: KineUtils.cxx:915
Initial State information.
Definition: InitialState.h:49
bool IsCoherent(void) const
Definition: ProcessInfo.cxx:97
Range1D_t KPhaseSpace::YLim ( void  ) const

y limits

Definition at line 659 of file KPhaseSpace.cxx.

References genie::utils::kinematics::CohYLim(), genie::utils::kinematics::DarkYLim(), fInteraction, genie::Interaction::FSPrimLepton(), genie::Target::HitNucP4Ptr(), genie::utils::kinematics::InelYLim(), genie::Interaction::InitState(), genie::ProcessInfo::IsCoherent(), genie::ProcessInfo::IsDarkMatterDeepInelastic(), genie::ProcessInfo::IsDeepInelastic(), genie::ProcessInfo::IsDiffractive(), genie::ProcessInfo::IsEM(), genie::ProcessInfo::IsIMDAnnihilation(), genie::ProcessInfo::IsInverseMuDecay(), genie::ProcessInfo::IsNuElectronElastic(), genie::ProcessInfo::IsResonant(), genie::controls::kASmallNum, genie::constants::kElectronMass, genie::constants::kPionMass, genie::kRfHitNucRest, genie::kRfLab, genie::Range1D_t::max, genie::Range1D_t::min, python.hepunit::pi, genie::InitialState::ProbeE(), genie::Interaction::ProcInfo(), and genie::InitialState::Tgt().

Referenced by genie::COHKinematicsGenerator::CalculateKin_BergerSehgal(), genie::COHKinematicsGenerator::CalculateKin_BergerSehgalFM(), genie::COHKinematicsGenerator::CalculateKin_ReinSehgal(), genie::DFRKinematicsGenerator::ComputeMaxXSec(), genie::COHElKinematicsGenerator::ComputeMaxXSec(), IsAllowed(), Limits(), genie::COHKinematicsGenerator::MaxXSec_AlvarezRuso(), genie::COHKinematicsGenerator::MaxXSec_BergerSehgal(), genie::COHKinematicsGenerator::MaxXSec_BergerSehgalFM(), genie::COHKinematicsGenerator::MaxXSec_ReinSehgal(), genie::COHElKinematicsGenerator::ProcessEventRecord(), YLim(), and YLim_X().

660 {
661  Range1D_t yl;
662  yl.min = -1;
663  yl.max = -1;
664 
665  const ProcessInfo & pi = fInteraction->ProcInfo();
666  bool is_em = pi.IsEM();
667 
668  //RES+DIS
669  bool is_inel = pi.IsDeepInelastic() || pi.IsResonant();
670  if(is_inel) {
671  const InitialState & init_state = fInteraction->InitState();
672  double Ev = init_state.ProbeE(kRfHitNucRest);
673  double M = init_state.Tgt().HitNucP4Ptr()->M(); // can be off m/shell
674  double ml = fInteraction->FSPrimLepton()->Mass();
675  yl = is_em ? kinematics::electromagnetic::InelYLim(Ev,ml,M) : kinematics::InelYLim(Ev,M,ml);
676  return yl;
677  }
678  //DMDIS
679  bool is_dmdis = pi.IsDarkMatterDeepInelastic();
680  if(is_dmdis) {
681  const InitialState & init_state = fInteraction->InitState();
682  double Ev = init_state.ProbeE(kRfHitNucRest);
683  double M = init_state.Tgt().HitNucP4Ptr()->M(); // can be off m/shell
684  double ml = fInteraction->FSPrimLepton()->Mass();
685  yl = kinematics::DarkYLim(Ev,M,ml);
686  return yl;
687  }
688  //COH
689  bool is_coh = pi.IsCoherent();
690  if(is_coh) {
691  const InitialState & init_state = fInteraction->InitState();
692  double EvL = init_state.ProbeE(kRfLab);
693  double ml = fInteraction->FSPrimLepton()->Mass();
694  yl = kinematics::CohYLim(EvL,ml);
695  return yl;
696  }
697  // IMD
698  if(pi.IsInverseMuDecay() || pi.IsIMDAnnihilation() || pi.IsNuElectronElastic()) {
699  const InitialState & init_state = fInteraction->InitState();
700  double Ev = init_state.ProbeE(kRfLab);
701  double ml = fInteraction->FSPrimLepton()->Mass();
702  double me = kElectronMass;
704  yl.max = 1 - (ml*ml + me*me)/(2*me*Ev) - controls::kASmallNum;
705  return yl;
706  }
707  bool is_dfr = pi.IsDiffractive();
708  if(is_dfr) {
709  const InitialState & init_state = fInteraction -> InitState();
710  double Ev = init_state.ProbeE(kRfHitNucRest);
711  double ml = fInteraction->FSPrimLepton()->Mass();
713  yl.max = 1. -ml/Ev - controls::kASmallNum;
714  return yl;
715  }
716  return yl;
717 }
bool IsResonant(void) const
Definition: ProcessInfo.cxx:92
A simple [min,max] interval for doubles.
Definition: Range1.h:43
bool IsInverseMuDecay(void) const
Range1D_t CohYLim(double Mn, double mpi, double mlep, double Ev, double Q2, double xsi)
Definition: KineUtils.cxx:762
bool IsDiffractive(void) const
bool IsIMDAnnihilation(void) const
static const double kElectronMass
Definition: Constants.h:71
Definition: NueSkimmer.h:24
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
bool IsNuElectronElastic(void) const
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
static const double kASmallNum
Definition: Controls.h:40
TParticlePDG * FSPrimLepton(void) const
final state primary lepton
bool IsEM(void) const
bool IsDeepInelastic(void) const
Definition: ProcessInfo.cxx:82
double max
Definition: Range1.h:54
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:264
static const double kPionMass
Definition: Constants.h:74
Range1D_t DarkYLim(double Ev, double M, double ml)
Definition: KineUtils.cxx:935
const InitialState & InitState(void) const
Definition: Interaction.h:69
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
double min
Definition: Range1.h:53
const Target & Tgt(void) const
Definition: InitialState.h:67
bool IsDarkMatterDeepInelastic(void) const
Definition: ProcessInfo.cxx:87
void kinematics()
Definition: kinematics.C:10
Range1D_t InelYLim(double Ev, double M, double ml)
Definition: KineUtils.cxx:407
double ProbeE(RefFrame_t rf) const
Initial State information.
Definition: InitialState.h:49
bool IsCoherent(void) const
Definition: ProcessInfo.cxx:97
Range1D_t KPhaseSpace::YLim ( double  xsi) const

y limits (COH)

Definition at line 764 of file KPhaseSpace.cxx.

References genie::utils::kinematics::CohYLim(), fInteraction, genie::Interaction::FSPrimLepton(), genie::Interaction::InitState(), genie::ProcessInfo::IsCoherent(), genie::ProcessInfo::IsWeakCC(), genie::Interaction::Kine(), kPi0Mass, genie::constants::kPionMass, genie::kRfHitNucRest, genie::Target::Mass(), genie::Range1D_t::max, genie::Range1D_t::min, python.hepunit::pi, genie::InitialState::ProbeE(), genie::Interaction::ProcInfo(), genie::Kinematics::Q2(), genie::utils::kinematics::Q2(), genie::InitialState::Tgt(), and YLim().

765 {
766  // Paschos-Schalla xsi parameter for y-limits in COH
767  // From PRD 80, 033005 (2009)
768 
769  Range1D_t yl;
770  yl.min = -1;
771  yl.max = -1;
772 
773  const ProcessInfo & pi = fInteraction->ProcInfo();
774 
775  //COH
776  bool is_coh = pi.IsCoherent();
777  if(is_coh) {
778  const InitialState & init_state = fInteraction->InitState();
779  const Kinematics & kine = fInteraction->Kine();
780  double Ev = init_state.ProbeE(kRfHitNucRest);
781  double Q2 = kine.Q2();
782  bool pionIsCharged = pi.IsWeakCC();
783  double Mn = init_state.Tgt().Mass();
784  double mpi = pionIsCharged ? kPionMass : kPi0Mass;
785  double mlep = fInteraction->FSPrimLepton()->Mass();
786  yl = kinematics::CohYLim(Mn, mpi, mlep, Ev, Q2, xsi);
787  return yl;
788  } else {
789  return this->YLim();
790  }
791 }
const Var kPi0Mass
bool IsWeakCC(void) const
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:991
A simple [min,max] interval for doubles.
Definition: Range1.h:43
Generated/set kinematical variables for an event.
Definition: Kinematics.h:40
Range1D_t CohYLim(double Mn, double mpi, double mlep, double Ev, double Q2, double xsi)
Definition: KineUtils.cxx:762
Range1D_t YLim(void) const
y limits
double Mass(void) const
Definition: Target.cxx:241
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
const Kinematics & Kine(void) const
Definition: Interaction.h:71
TParticlePDG * FSPrimLepton(void) const
final state primary lepton
double max
Definition: Range1.h:54
static const double kPionMass
Definition: Constants.h:74
const InitialState & InitState(void) const
Definition: Interaction.h:69
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
double min
Definition: Range1.h:53
double Q2(bool selected=false) const
Definition: Kinematics.cxx:135
const Target & Tgt(void) const
Definition: InitialState.h:67
double ProbeE(RefFrame_t rf) const
Initial State information.
Definition: InitialState.h:49
bool IsCoherent(void) const
Definition: ProcessInfo.cxx:97
Range1D_t KPhaseSpace::YLim_X ( void  ) const

y limits @ fixed x

Definition at line 719 of file KPhaseSpace.cxx.

References genie::utils::kinematics::CohYLim(), genie::utils::kinematics::DarkYLim_X(), fInteraction, genie::Interaction::FSPrimLepton(), genie::Target::HitNucP4Ptr(), genie::utils::kinematics::InelYLim_X(), genie::Interaction::InitState(), genie::ProcessInfo::IsCoherent(), genie::ProcessInfo::IsDarkMatterDeepInelastic(), genie::ProcessInfo::IsDeepInelastic(), genie::ProcessInfo::IsEM(), genie::ProcessInfo::IsResonant(), genie::Interaction::Kine(), genie::kRfHitNucRest, genie::kRfLab, genie::Range1D_t::max, genie::Range1D_t::min, python.hepunit::pi, genie::InitialState::ProbeE(), genie::Interaction::ProcInfo(), genie::InitialState::Tgt(), genie::Kinematics::x(), and submit_syst::x.

Referenced by YLim_X().

720 {
721 // Computes kinematical limits for y @ the input x
722 
723  Range1D_t yl;
724  yl.min = -1;
725  yl.max = -1;
726 
727  const ProcessInfo & pi = fInteraction->ProcInfo();
728  bool is_em = pi.IsEM();
729 
730  //RES+DIS
731  bool is_inel = pi.IsDeepInelastic() || pi.IsResonant();
732  if(is_inel) {
733  const InitialState & init_state = fInteraction->InitState();
734  double Ev = init_state.ProbeE(kRfHitNucRest);
735  double M = init_state.Tgt().HitNucP4Ptr()->M(); // can be off m/shell
736  double ml = fInteraction->FSPrimLepton()->Mass();
737  double x = fInteraction->Kine().x();
738  yl = is_em ? kinematics::electromagnetic::InelYLim_X(Ev,ml,M,x) : kinematics::InelYLim_X(Ev,M,ml,x);
739  return yl;
740  }
741  //DMDIS
742  bool is_dmdis = pi.IsDarkMatterDeepInelastic();
743  if(is_dmdis) {
744  const InitialState & init_state = fInteraction->InitState();
745  double Ev = init_state.ProbeE(kRfHitNucRest);
746  double M = init_state.Tgt().HitNucP4Ptr()->M(); // can be off m/shell
747  double ml = fInteraction->FSPrimLepton()->Mass();
748  double x = fInteraction->Kine().x();
749  yl = kinematics::DarkYLim_X(Ev,M,ml,x);
750  return yl;
751  }
752  //COH
753  bool is_coh = pi.IsCoherent();
754  if(is_coh) {
755  const InitialState & init_state = fInteraction->InitState();
756  double EvL = init_state.ProbeE(kRfLab);
757  double ml = fInteraction->FSPrimLepton()->Mass();
758  yl = kinematics::CohYLim(EvL,ml);
759  return yl;
760  }
761  return yl;
762 }
bool IsResonant(void) const
Definition: ProcessInfo.cxx:92
A simple [min,max] interval for doubles.
Definition: Range1.h:43
Range1D_t CohYLim(double Mn, double mpi, double mlep, double Ev, double Q2, double xsi)
Definition: KineUtils.cxx:762
double x(bool selected=false) const
Definition: Kinematics.cxx:109
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
const Kinematics & Kine(void) const
Definition: Interaction.h:71
TParticlePDG * FSPrimLepton(void) const
final state primary lepton
bool IsEM(void) const
bool IsDeepInelastic(void) const
Definition: ProcessInfo.cxx:82
double max
Definition: Range1.h:54
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:264
Range1D_t InelYLim_X(double Ev, double M, double ml, double x)
Definition: KineUtils.cxx:442
const InitialState & InitState(void) const
Definition: Interaction.h:69
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
double min
Definition: Range1.h:53
Range1D_t DarkYLim_X(double Ev, double M, double ml, double x)
Definition: KineUtils.cxx:951
const Target & Tgt(void) const
Definition: InitialState.h:67
bool IsDarkMatterDeepInelastic(void) const
Definition: ProcessInfo.cxx:87
void kinematics()
Definition: kinematics.C:10
double ProbeE(RefFrame_t rf) const
Initial State information.
Definition: InitialState.h:49
bool IsCoherent(void) const
Definition: ProcessInfo.cxx:97
Range1D_t KPhaseSpace::YLim_X ( double  xsi) const

y limits @ fixed x (COH)

Definition at line 793 of file KPhaseSpace.cxx.

References fInteraction, genie::ProcessInfo::IsCoherent(), python.hepunit::pi, genie::Interaction::ProcInfo(), YLim(), and YLim_X().

794 {
795  // Paschos-Schalla xsi parameter for y-limits in COH
796  // From PRD 80, 033005 (2009)
797 
798  const ProcessInfo & pi = fInteraction->ProcInfo();
799 
800  //COH
801  bool is_coh = pi.IsCoherent();
802  if(is_coh) {
803  return this->YLim(xsi);
804  } else {
805  return this->YLim_X();
806  }
807 }
Range1D_t YLim_X(void) const
y limits @ fixed x
Range1D_t YLim(void) const
y limits
const Interaction * fInteraction
Definition: KPhaseSpace.h:74
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
bool IsCoherent(void) const
Definition: ProcessInfo.cxx:97

Member Data Documentation

const Interaction* genie::KPhaseSpace::fInteraction
private

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