Namespaces | Functions
genie::utils::kinematics Namespace Reference

Kinematical utilities. More...

Namespaces

 electromagnetic
 

Functions

double PhaseSpaceVolume (const Interaction *const i, KinePhaseSpace_t ps)
 
double Jacobian (const Interaction *const i, KinePhaseSpace_t f, KinePhaseSpace_t t)
 
bool TransformMatched (KinePhaseSpace_t ia, KinePhaseSpace_t ib, KinePhaseSpace_t a, KinePhaseSpace_t b, bool &fwd)
 
Range1D_t InelWLim (double Ev, double M, double ml)
 
Range1D_t InelQ2Lim_W (double Ev, double M, double ml, double W, double Q2min_cut=controls::kMinQ2Limit)
 
Range1D_t Inelq2Lim_W (double Ev, double M, double ml, double W, double q2min_cut=-1 *controls::kMinQ2Limit)
 
Range1D_t InelQ2Lim (double Ev, double M, double ml, double Q2min_cut=controls::kMinQ2Limit)
 
Range1D_t Inelq2Lim (double Ev, double M, double ml, double q2min_cut=-1 *controls::kMinQ2Limit)
 
Range1D_t InelXLim (double Ev, double M, double ml)
 
Range1D_t InelYLim (double Ev, double M, double ml)
 
Range1D_t InelYLim_X (double Ev, double M, double ml, double x)
 
Range1D_t CohW2Lim (double Mn, double m_produced, double mlep, double Ev, double Q2)
 
Range1D_t CohNuLim (double W2min, double W2max, double Q2, double Mn, double xsi)
 
Range1D_t CohYLim (double Mn, double m_produced, double mlep, double Ev, double Q2, double xsi)
 
Range1D_t CohYLim (double EvL, double ml)
 
Range1D_t CohXLim (void)
 
Range1D_t CohQ2Lim (double Mn, double m_produced, double mlep, double Ev)
 
Range1D_t Cohq2Lim (double Mn, double m_produced, double mlep, double Ev)
 
Range1D_t CEvNSQ2Lim (double Ev)
 
Range1D_t DarkWLim (double Ev, double M, double ml)
 
Range1D_t DarkQ2Lim_W (double Ev, double M, double ml, double W, double Q2min_cut=controls::kMinQ2Limit)
 
Range1D_t Darkq2Lim_W (double Ev, double M, double ml, double W, double q2min_cut=-1 *controls::kMinQ2Limit)
 
Range1D_t DarkQ2Lim (double Ev, double M, double ml, double Q2min_cut=controls::kMinQ2Limit)
 
Range1D_t Darkq2Lim (double Ev, double M, double ml, double q2min_cut=-1 *controls::kMinQ2Limit)
 
Range1D_t DarkXLim (double Ev, double M, double ml)
 
Range1D_t DarkYLim (double Ev, double M, double ml)
 
Range1D_t DarkYLim_X (double Ev, double M, double ml, double x)
 
double CohW2Min (double Mn, double m_produced)
 
double QD2toQ2 (double QD2)
 
double Q2toQD2 (double Q2)
 
void WQ2toXY (double Ev, double M, double W, double Q2, double &x, double &y)
 
void XYtoWQ2 (double Ev, double M, double &W, double &Q2, double x, double y)
 
void XQ2toWY (double Ev, double M, double &W, double Q2, double x, double &y)
 
double XYtoW (double Ev, double M, double x, double y)
 
double XYtoQ2 (double Ev, double M, double x, double y)
 
double Q2YtoX (double Ev, double M, double Q2, double y)
 
void UpdateWQ2FromXY (const Interaction *in)
 
void UpdateXYFromWQ2 (const Interaction *in)
 
void UpdateWYFromXQ2 (const Interaction *in)
 
void UpdateXFromQ2Y (const Interaction *in)
 
void ApplyCutsToKineLimits (Range1D_t &r, double min, double max)
 
double Q2 (const Interaction *const i)
 
double W (const Interaction *const i)
 
bool IsAboveCharmThreshold (double x, double Q2, double M, double mc)
 
double SlowRescalingVar (double x, double Q2, double M, double mc)
 
double RESImportanceSamplingEnvelope (double *x, double *par)
 
double DISImportanceSamplingEnvelope (double *x, double *par)
 
double COHImportanceSamplingEnvelope (double *x, double *par)
 

Detailed Description

Kinematical utilities.

Author
Costas Andreopoulos <constantinos.andreopoulos cern.ch> University of Liverpool & STFC Rutherford Appleton Laboratory

Changes required to implement the GENIE Boosted Dark Matter module were installed by Josh Berger (Univ. of Wisconsin)

November 26, 2004

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

Function Documentation

void genie::utils::kinematics::ApplyCutsToKineLimits ( Range1D_t r,
double  min,
double  max 
)

Definition at line 1258 of file KineUtils.cxx.

1260 {
1261  // if the min,max are within the existing limits, the cut can be applied
1262  // by narrowing down the xisting limits
1263  if ( utils::math::IsWithinLimits(min_cut, range ) ) range.min = min_cut;
1264  if ( utils::math::IsWithinLimits(max_cut, range ) ) range.max = max_cut;
1265 
1266  // if the min-cut is above the existing max-limit or
1267  // if the max-cut is below the existing min-limit then
1268  // the range should be invalidated
1269 
1270  if (min_cut > range.max || max_cut < range.min) {
1271 
1272  range.min = 0;
1273  range.max = 0;
1274  }
1275 }
bool IsWithinLimits(double x, Range1D_t range)
Definition: MathUtils.cxx:258
Range1D_t genie::utils::kinematics::CEvNSQ2Lim ( double  Ev)

Definition at line 873 of file KineUtils.cxx.

874 {
876  return Q2;
877 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
A simple [min,max] interval for doubles.
Definition: Range1.h:42
static const double kAVerySmallNum
Definition: Controls.h:39
double genie::utils::kinematics::COHImportanceSamplingEnvelope ( double *  x,
double *  par 
)

Definition at line 1453 of file KineUtils.cxx.

1455 {
1456  //-- inputs
1457  double xb = x[0]; // x
1458  double yb = x[1]; // y
1459 
1460  //-- parameters
1461  double xsmax = 3*par[0]; // safety factor * max cross section in (x,y)
1462  double Ev = par[1]; // neutrino energy;
1463 
1464  if(yb<0.|| yb>1.) return 0.;
1465  if(xb<0.|| xb>1.) return 0.;
1466 
1467  if(Ev<1) return xsmax;
1468  if(xb/Ev<1E-4 && yb>0.95) return 5*xsmax;
1469 
1470  double func = 0;
1471  double xp = 0.1;
1472  double yp = (Ev>2.5) ? 2.5/Ev : 1;
1473 
1474  if(xb>xp) {
1475  double xs0=0;
1476  if(yb<yp) {
1477  xs0 = xsmax;
1478  } else {
1479  xs0 = xsmax - (yb-yp)*xsmax;
1480  }
1481  double d = TMath::Power( (xb-xp)/0.075, 2);
1482  func = xs0/(1 + d);
1483  } else {
1484  if(yb>yp) {
1485  func = xsmax - (yb-yp)*xsmax;
1486  } else {
1487  func = xsmax;
1488  }
1489  }
1490  return func;
1491 }
def func()
Definition: docstring.py:7
list x
Definition: train.py:276
Range1D_t genie::utils::kinematics::CohNuLim ( double  W2min,
double  W2max,
double  Q2,
double  Mn,
double  xsi 
)

Definition at line 810 of file KineUtils.cxx.

812 {
813  Range1D_t nul;
814  nul.min = -1;
815  nul.max = -1;
816 
817  double nu_min = (W2min + Q2 - Mn * Mn) / (2.0 * Mn);
818  double nu_max = (W2max + Q2 - Mn * Mn) / (2.0 * Mn);
819  double xsiQ = xsi * TMath::Sqrt(Q2);
820 
821  nul.min = (xsiQ > nu_min) ? xsiQ : nu_min;
822  nul.max = nu_max;
823 
824  return nul;
825 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
A simple [min,max] interval for doubles.
Definition: Range1.h:42
double max
Definition: Range1.h:53
double min
Definition: Range1.h:52
Range1D_t genie::utils::kinematics::CohQ2Lim ( double  Mn,
double  m_produced,
double  mlep,
double  Ev 
)

Definition at line 730 of file KineUtils.cxx.

731 {
732  // The expressions for Q^2 min appears in PRD 74, 054007 (2006) by
733  // Kartavtsev, Paschos, and Gounaris
734 
735  // That expression is specified for the case of the pion.
736  // GENIE has also Coherent interaction with Single gamma production and Rho
737  // so that formula will be adapted for a generic m_produced
738 
739  Range1D_t Q2;
740  Q2.min = 0.0;
741  Q2.max = std::numeric_limits<double>::max(); // Value must be overriden in user options
742 
743  double Mn2 = Mn * Mn;
744  double mlep2 = mlep * mlep;
745  double s = Mn2 + 2.0 * Mn * Ev;
746  double W2min = CohW2Min(Mn, m_produced);
747 
748  // Looks like Q2min = A * B - C, where A, B, and C are complicated
749  double a = 1.0;
750  double b = mlep2 / s;
751  double c = W2min / s;
752  double lambda = a * a + b * b + c * c - 2.0 * a * b - 2.0 * a * c - 2.0 * b * c;
753  if (lambda > 0) {
754  double A = (s - Mn * Mn) / 2.0;
755  double B = 1 - TMath::Sqrt(lambda);
756  double C = 0.5 * (W2min + mlep2 - Mn2 * (W2min - mlep2) / s );
757  if (A * B - C < 0) {
758  SLOG("KineLimits", pERROR)
759  << "Q2 kinematic limits calculation failed for CohQ2Lim. "
760  << "Assuming Q2min = 0.0";
761  }
762  Q2.min = TMath::Max(0., A * B - C);
763  } else {
764  SLOG("KineLimits", pERROR)
765  << "Q2 kinematic limits calculation failed for CohQ2Lim. "
766  << "Assuming Q2min = 0.0";
767  }
768 
769  return Q2;
770 }
#define pERROR
Definition: Messenger.h:59
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
A simple [min,max] interval for doubles.
Definition: Range1.h:42
const double a
static int max(int a, int b)
double max
Definition: Range1.h:53
double CohW2Min(double Mn, double m_produced)
Definition: KineUtils.cxx:861
#define A
Definition: memgrp.cpp:38
static bool * b
Definition: config.cpp:1043
double min
Definition: Range1.h:52
#define SLOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a short string (using the FUNCTION and...
Definition: Messenger.h:84
static QCString * s
Definition: config.cpp:1042
Range1D_t genie::utils::kinematics::Cohq2Lim ( double  Mn,
double  m_produced,
double  mlep,
double  Ev 
)

Definition at line 772 of file KineUtils.cxx.

773 {
774  Range1D_t Q2 = utils::kinematics::CohQ2Lim(Mn, m_produced, mlep, Ev);
775  Range1D_t q2;
776  q2.min = - Q2.max;
777  q2.max = - Q2.min;
778  return q2;
779 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
A simple [min,max] interval for doubles.
Definition: Range1.h:42
double max
Definition: Range1.h:53
Range1D_t CohQ2Lim(double Mn, double m_produced, double mlep, double Ev)
Definition: KineUtils.cxx:730
double min
Definition: Range1.h:52
Range1D_t genie::utils::kinematics::CohW2Lim ( double  Mn,
double  m_produced,
double  mlep,
double  Ev,
double  Q2 
)

Definition at line 781 of file KineUtils.cxx.

783 {
784  // These expressions for W^2 min and max appear in PRD 74, 054007 (2006) by
785  // Kartavtsev, Paschos, and Gounaris
786 
787  // That expression is specified for the case of the pion.
788  // GENIE has also Coherent interaction with Single gamma production and Rho
789  // so that formula will be adapted for a generic m_produced
790 
791  Range1D_t W2l;
792  W2l.min = -1;
793  W2l.max = -1;
794 
795  double s = Mn * Mn + 2.0 * Mn * Ev;
796  double Mnterm = 1 - Mn * Mn / s;
797  double Mlterm = 1 - mlep * mlep / s;
798  // Here T1, T2 are generically "term 1" and "term 2" in a long expression
799  double T1 = 0.25 * s * s * Mnterm * Mnterm * Mlterm;
800  double T2 = Q2 - (0.5 * s * Mnterm) + (0.5 * mlep * mlep * Mnterm);
801 
802  W2l.min = CohW2Min(Mn, m_produced);
803  W2l.max = (T1 - T2 * T2 ) *
804  (1.0 / Mnterm) *
805  (1.0 / (Q2 + mlep * mlep));
806 
807  return W2l;
808 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
A simple [min,max] interval for doubles.
Definition: Range1.h:42
Definition: 013_class.h:14
Definition: 013_class.h:10
double max
Definition: Range1.h:53
double CohW2Min(double Mn, double m_produced)
Definition: KineUtils.cxx:861
double min
Definition: Range1.h:52
static QCString * s
Definition: config.cpp:1042
double genie::utils::kinematics::CohW2Min ( double  Mn,
double  m_produced 
)

Definition at line 861 of file KineUtils.cxx.

862 {
863  // These expressions for W^2 min and max appear in PRD 74, 054007 (2006) by
864  // Kartavtsev, Paschos, and Gounaris
865 
866  // That expression is specified for the case of the pion.
867  // GENIE has also Coherent interaction with Single gamma production and Rho
868  // so that formula will be adapted for a generic m_produced
869 
870  return (Mn + m_produced) * (Mn + m_produced);
871 }
Range1D_t genie::utils::kinematics::CohXLim ( void  )

Definition at line 722 of file KineUtils.cxx.

723 {
724 // Computes x limits for coherent v interactions
725 
727  return x;
728 }
A simple [min,max] interval for doubles.
Definition: Range1.h:42
static const double kASmallNum
Definition: Controls.h:40
list x
Definition: train.py:276
Range1D_t genie::utils::kinematics::CohYLim ( double  Mn,
double  m_produced,
double  mlep,
double  Ev,
double  Q2,
double  xsi 
)

Definition at line 827 of file KineUtils.cxx.

829 {
830  Range1D_t ylim;
831  ylim.min = -1;
832  ylim.max = -1;
833 
834  Range1D_t W2lim = genie::utils::kinematics::CohW2Lim(Mn, m_produced, mlep, Ev, Q2);
835  if (W2lim.min > W2lim.max) {
836  LOG("KineLimits", pDEBUG)
837  << "Kinematically forbidden region in CohYLim. W2min = " << W2lim.min
838  << "; W2max =" << W2lim.max;
839  LOG("KineLimits", pDEBUG)
840  << " Mn = " << Mn << "; m_had_sys = " << m_produced << "; mlep = "
841  << mlep << "; Ev = " << Ev << "; Q2 = " << Q2;
842  return ylim;
843  }
845  Q2, Mn, xsi);
846  ylim.min = nulim.min / Ev;
847  ylim.max = nulim.max / Ev;
848 
849  return ylim;
850 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
Range1D_t CohW2Lim(double Mn, double m_produced, double mlep, double Ev, double Q2)
Definition: KineUtils.cxx:781
A simple [min,max] interval for doubles.
Definition: Range1.h:42
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
Range1D_t CohNuLim(double W2min, double W2max, double Q2, double Mn, double xsi)
Definition: KineUtils.cxx:810
double max
Definition: Range1.h:53
double min
Definition: Range1.h:52
#define pDEBUG
Definition: Messenger.h:63
Range1D_t genie::utils::kinematics::CohYLim ( double  EvL,
double  ml 
)

Definition at line 852 of file KineUtils.cxx.

853 {
854 // Computes y limits for coherent v interactions
855 
857  1.-ml/EvL - controls::kASmallNum);
858  return y;
859 }
A simple [min,max] interval for doubles.
Definition: Range1.h:42
static const double kASmallNum
Definition: Controls.h:40
static const double kPionMass
Definition: Constants.h:73
Range1D_t genie::utils::kinematics::DarkQ2Lim ( double  Ev,
double  M,
double  ml,
double  Q2min_cut = controls::kMinQ2Limit 
)

Definition at line 960 of file KineUtils.cxx.

962 {
963 // Computes Q2 (>0) limits irrespective of W for inelastic v interactions
964 
965  Range1D_t Q2;
966  Q2.min = -1;
967  Q2.max = -1;
968 
970  if(W.min<0) return Q2;
971 
972  Q2 = utils::kinematics::DarkQ2Lim_W(Ev,M,ml,W.min,Q2min_cut);
973  return Q2;
974 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
A simple [min,max] interval for doubles.
Definition: Range1.h:42
Range1D_t DarkQ2Lim_W(double Ev, double M, double ml, double W, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:902
double max
Definition: Range1.h:53
double min
Definition: Range1.h:52
Range1D_t DarkWLim(double Ev, double M, double ml)
Definition: KineUtils.cxx:879
Range1D_t genie::utils::kinematics::Darkq2Lim ( double  Ev,
double  M,
double  ml,
double  q2min_cut = -1*controls::kMinQ2Limit 
)

Definition at line 976 of file KineUtils.cxx.

978 {
979 // Computes Q2 (>0) limits irrespective of W for inelastic v interactions
980 
981  Range1D_t Q2 = utils::kinematics::DarkQ2Lim(Ev,M,ml,-1.*q2min_cut);
982  Range1D_t q2;
983  q2.min = - Q2.max;
984  q2.max = - Q2.min;
985  return q2;
986 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
A simple [min,max] interval for doubles.
Definition: Range1.h:42
Range1D_t DarkQ2Lim(double Ev, double M, double ml, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:960
double max
Definition: Range1.h:53
double min
Definition: Range1.h:52
Range1D_t genie::utils::kinematics::DarkQ2Lim_W ( double  Ev,
double  M,
double  ml,
double  W,
double  Q2min_cut = controls::kMinQ2Limit 
)

Definition at line 902 of file KineUtils.cxx.

904 {
905 // Computes Q2 limits (>0) @ the input W for inelastic v interactions
906 
907  Range1D_t Q2;
908  Q2.min = -1;
909  Q2.max = -1;
910 
911  double M2 = TMath::Power(M, 2.);
912  double ml2 = TMath::Power(ml, 2.);
913  double W2 = TMath::Power(W, 2.);
914  double s = M2 + 2*M*Ev + ml2;
915  assert(s > 0);
916  double sqs = TMath::Sqrt(s);
917  double E1CM = (s + ml2 - M2) / (2.*sqs);
918  double p1CM = TMath::Max(0., E1CM*E1CM - ml2);
919  p1CM = TMath::Sqrt(p1CM);
920  double E3CM = (s + ml2 - W2) / (2.*sqs);
921  double p3CM = TMath::Max(0., E3CM*E3CM - ml2);
922  p3CM = TMath::Sqrt(p3CM);
923 
924  SLOG("KineLimits", pDEBUG) << "s = " << s;
925  SLOG("KineLimits", pDEBUG) << "Ev = " << Ev;
926  SLOG("KineLimits", pDEBUG) << "M = " << M;
927  SLOG("KineLimits", pDEBUG) << "W = " << W;
928  SLOG("KineLimits", pDEBUG) << "E1_CM = " << E1CM;
929  SLOG("KineLimits", pDEBUG) << "p1_CM = " << p1CM;
930  SLOG("KineLimits", pDEBUG) << "E3_CM = " << E3CM;
931  SLOG("KineLimits", pDEBUG) << "p3_CM = " << p3CM;
932 
933  Q2.min = TMath::Power(p3CM - p1CM,2) - TMath::Power((W2 - M2) / (2.*sqs),2);
934  Q2.max = TMath::Power(p3CM + p1CM,2) - TMath::Power((W2 - M2) / (2.*sqs),2);
935 
936  SLOG("KineLimits", pDEBUG) << "Nominal Q^2 limits: " << Q2.min << " , " << Q2.max;
937  // guard against overflows
938  Q2.max = TMath::Max(0., Q2.max);
939  Q2.min = TMath::Max(0., Q2.min);
940 
941  // limit the minimum Q2
942  if(Q2.min < Q2min_cut) {Q2.min = Q2min_cut; }
943  if(Q2.max < Q2.min ) {Q2.min = -1; Q2.max = -1;}
944 
945  return Q2;
946 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
A simple [min,max] interval for doubles.
Definition: Range1.h:42
double max
Definition: Range1.h:53
double min
Definition: Range1.h:52
#define SLOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a short string (using the FUNCTION and...
Definition: Messenger.h:84
static QCString * s
Definition: config.cpp:1042
#define pDEBUG
Definition: Messenger.h:63
Range1D_t genie::utils::kinematics::Darkq2Lim_W ( double  Ev,
double  M,
double  ml,
double  W,
double  q2min_cut = -1*controls::kMinQ2Limit 
)

Definition at line 948 of file KineUtils.cxx.

950 {
951 // Computes q2 (<0) limits @ the input W for inelastic v interactions
952 
953  Range1D_t Q2 = utils::kinematics::DarkQ2Lim_W(Ev,M,ml,W,-1.*q2min_cut);
954  Range1D_t q2;
955  q2.min = - Q2.max;
956  q2.max = - Q2.min;
957  return q2;
958 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
A simple [min,max] interval for doubles.
Definition: Range1.h:42
Range1D_t DarkQ2Lim_W(double Ev, double M, double ml, double W, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:902
double max
Definition: Range1.h:53
double min
Definition: Range1.h:52
Range1D_t genie::utils::kinematics::DarkWLim ( double  Ev,
double  M,
double  ml 
)

Definition at line 879 of file KineUtils.cxx.

880 {
881 // Computes W limits for inelastic v interactions
882 //
883  double M2 = TMath::Power(M,2);
884  double ml2 = TMath::Power(ml,2);
885  double s = M2 + 2*M*Ev + ml2;
886  assert (s>0);
887 
888  Range1D_t W;
890 // W.min = kNeutronMass + kPionMass;
891  W.max = TMath::Sqrt(s) - ml;
892  if(W.max<=W.min) {
893  W.min = -1;
894  W.max = -1;
895  return W;
896  }
899  return W;
900 }
A simple [min,max] interval for doubles.
Definition: Range1.h:42
static const double kPhotontest
Definition: Constants.h:79
static const double kNeutronMass
Definition: Constants.h:76
static const double kASmallNum
Definition: Controls.h:40
double max
Definition: Range1.h:53
double min
Definition: Range1.h:52
static QCString * s
Definition: config.cpp:1042
Range1D_t genie::utils::kinematics::DarkXLim ( double  Ev,
double  M,
double  ml 
)

Definition at line 988 of file KineUtils.cxx.

990 {
991 // Computes Bjorken x limits for inelastic interactions
992 // For the dark matter case, it is relatively straightforward
993  Range1D_t Wl = utils::kinematics::DarkWLim(Ev, M, ml);
994  double Wmin = Wl.min;
995  double W2min = Wmin*Wmin;
996  SLOG("KineLimits", pDEBUG) << "W^2_min = " << W2min;
997  Range1D_t Q2l = utils::kinematics::DarkQ2Lim_W(Ev, M, ml, Wmin);
998  SLOG("KineLimits", pDEBUG) << "Q^2 range : " << Q2l.min << " , " << Q2l.max;
999  double M2 = M*M;
1000  Range1D_t x;
1001  x.min = Q2l.min / (Q2l.min + W2min - M2);
1002  x.max = Q2l.max / (Q2l.max + W2min - M2);
1003 
1004  SLOG("KineLimits", pDEBUG) << "x = [" << x.min << ", " << x.max << "]";
1005  return x;
1006 }
A simple [min,max] interval for doubles.
Definition: Range1.h:42
Range1D_t DarkQ2Lim_W(double Ev, double M, double ml, double W, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:902
double max
Definition: Range1.h:53
list x
Definition: train.py:276
double min
Definition: Range1.h:52
#define SLOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a short string (using the FUNCTION and...
Definition: Messenger.h:84
Range1D_t DarkWLim(double Ev, double M, double ml)
Definition: KineUtils.cxx:879
#define pDEBUG
Definition: Messenger.h:63
Range1D_t genie::utils::kinematics::DarkYLim ( double  Ev,
double  M,
double  ml 
)

Definition at line 1008 of file KineUtils.cxx.

1009 {
1010  // For dark inelastic scattering, can compute exactly and is much simpler
1011  Range1D_t Wl = utils::kinematics::DarkWLim(Ev, M, ml);
1012  double Wmin = Wl.min;
1013  double W2min = Wmin*Wmin;
1014  Range1D_t Q2l = utils::kinematics::DarkQ2Lim_W(Ev, M, ml, Wmin);
1015  double M2 = M*M;
1016  Range1D_t y;
1017  y.min = (Q2l.min + W2min - M2) / (2*Ev*M);
1018  y.max = (Q2l.max + W2min - M2) / (2*Ev*M);
1019 
1020  SLOG("KineLimits", pDEBUG) << "y = [" << y.min << ", " << y.max << "]";
1021  return y;
1022 }
A simple [min,max] interval for doubles.
Definition: Range1.h:42
Range1D_t DarkQ2Lim_W(double Ev, double M, double ml, double W, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:902
double max
Definition: Range1.h:53
double min
Definition: Range1.h:52
#define SLOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a short string (using the FUNCTION and...
Definition: Messenger.h:84
Range1D_t DarkWLim(double Ev, double M, double ml)
Definition: KineUtils.cxx:879
#define pDEBUG
Definition: Messenger.h:63
Range1D_t genie::utils::kinematics::DarkYLim_X ( double  Ev,
double  M,
double  ml,
double  x 
)

Definition at line 1024 of file KineUtils.cxx.

1027 {
1028  // Computes y limits @ the input x for inelastic interactions
1029  // We hit y_min when W = W_min and y_max when Q^2 = Q_min^2 or Q^2 = Q_max^2
1030 
1031  Range1D_t y;
1032  y.min = -1;
1033  y.max = -1;
1034 
1035  Range1D_t Wl = utils::kinematics::DarkWLim(Ev, M, ml);
1036  double Wmin = Wl.min;
1037  double W2min = Wmin*Wmin;
1038  double M2 = M*M;
1039  double ml2 = ml*ml;
1040  y.min = (W2min - M2) / (1.-x) / (2*Ev*M);
1041  y.max = 2.* M * x *(Ev*Ev - ml2) / Ev / (2. * M * Ev * x + M2 * x * x + ml2);
1042 
1043  return y;
1044 }
A simple [min,max] interval for doubles.
Definition: Range1.h:42
double max
Definition: Range1.h:53
list x
Definition: train.py:276
double min
Definition: Range1.h:52
Range1D_t DarkWLim(double Ev, double M, double ml)
Definition: KineUtils.cxx:879
double genie::utils::kinematics::DISImportanceSamplingEnvelope ( double *  x,
double *  par 
)

Definition at line 1421 of file KineUtils.cxx.

1423 {
1424  //-- inputs
1425  double xb = x[0]; // scaling variable x brorken
1426  //double y = x[1]; // inelasticity y
1427 
1428  //-- parameters
1429  double xpeak = par[0]; // x peak
1430  double xmin = par[1]; // min x
1431  double xmax = 1.; // max x
1432  double xsmax = par[2]; // safety factor * max cross section in (x,y)
1433 
1434  double func = 0;
1435 
1436  if(xb < xpeak/20.) {
1437  //low-x falling edge
1438  double plateau_edge = xpeak/20.;
1439  double slope = xsmax/(xmin - plateau_edge);
1440  func = xsmax - slope * (xb - plateau_edge);
1441  } else if (xb > 2*xpeak) {
1442  //high-x falling edge
1443  double plateau_edge = 2*xpeak;
1444  double slope = xsmax/(xmax - plateau_edge);
1445  func = xsmax - slope * (xb - plateau_edge);
1446  } else {
1447  // plateau
1448  func = xsmax;
1449  }
1450  return func;
1451 }
def func()
Definition: docstring.py:7
list x
Definition: train.py:276
Range1D_t genie::utils::kinematics::InelQ2Lim ( double  Ev,
double  M,
double  ml,
double  Q2min_cut = controls::kMinQ2Limit 
)

Definition at line 416 of file KineUtils.cxx.

418 {
419 // Computes Q2 (>0) limits irrespective of W for inelastic v interactions
420 
421  Range1D_t Q2;
422  Q2.min = -1;
423  Q2.max = -1;
424 
426  if(W.min<0) return Q2;
427 
428  Q2 = utils::kinematics::InelQ2Lim_W(Ev,M,ml,W.min,Q2min_cut);
429  return Q2;
430 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
Range1D_t InelWLim(double Ev, double M, double ml)
Definition: KineUtils.cxx:345
A simple [min,max] interval for doubles.
Definition: Range1.h:42
Range1D_t InelQ2Lim_W(double Ev, double M, double ml, double W, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:366
double max
Definition: Range1.h:53
double min
Definition: Range1.h:52
Range1D_t genie::utils::kinematics::Inelq2Lim ( double  Ev,
double  M,
double  ml,
double  q2min_cut = -1*controls::kMinQ2Limit 
)

Definition at line 432 of file KineUtils.cxx.

434 {
435 // Computes Q2 (>0) limits irrespective of W for inelastic v interactions
436 
437  Range1D_t Q2 = utils::kinematics::InelQ2Lim(Ev,M,ml,-1.*q2min_cut);
438  Range1D_t q2;
439  q2.min = - Q2.max;
440  q2.max = - Q2.min;
441  return q2;
442 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
A simple [min,max] interval for doubles.
Definition: Range1.h:42
Range1D_t InelQ2Lim(double Ev, double M, double ml, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:416
double max
Definition: Range1.h:53
double min
Definition: Range1.h:52
Range1D_t genie::utils::kinematics::InelQ2Lim_W ( double  Ev,
double  M,
double  ml,
double  W,
double  Q2min_cut = controls::kMinQ2Limit 
)

Definition at line 366 of file KineUtils.cxx.

368 {
369 // Computes Q2 limits (>0) @ the input W for inelastic v interactions
370 
371  Range1D_t Q2;
372  Q2.min = -1;
373  Q2.max = -1;
374 
375  double M2 = TMath::Power(M, 2.);
376  double ml2 = TMath::Power(ml, 2.);
377  double W2 = TMath::Power(W, 2.);
378  double s = M2 + 2*M*Ev;
379 
380  SLOG("KineLimits", pDEBUG) << "s = " << s;
381  SLOG("KineLimits", pDEBUG) << "Ev = " << Ev;
382  assert (s>0);
383 
384  double auxC = 0.5*(s-M2)/s;
385  double aux1 = s + ml2 - W2;
386  double aux2 = aux1*aux1 - 4*s*ml2;
387 
388  (aux2 < 0) ? ( aux2 = 0 ) : ( aux2 = TMath::Sqrt(aux2) );
389 
390  Q2.max = -ml2 + auxC * (aux1 + aux2); // => 0
391  Q2.min = -ml2 + auxC * (aux1 - aux2); // => 0
392 
393  // guard against overflows
394  Q2.max = TMath::Max(0., Q2.max);
395  Q2.min = TMath::Max(0., Q2.min);
396 
397  // limit the minimum Q2
398  if(Q2.min < Q2min_cut) {Q2.min = Q2min_cut; }
399  if(Q2.max < Q2.min ) {Q2.min = -1; Q2.max = -1;}
400 
401  return Q2;
402 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
A simple [min,max] interval for doubles.
Definition: Range1.h:42
double max
Definition: Range1.h:53
double min
Definition: Range1.h:52
#define SLOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a short string (using the FUNCTION and...
Definition: Messenger.h:84
static QCString * s
Definition: config.cpp:1042
#define pDEBUG
Definition: Messenger.h:63
Range1D_t genie::utils::kinematics::Inelq2Lim_W ( double  Ev,
double  M,
double  ml,
double  W,
double  q2min_cut = -1*controls::kMinQ2Limit 
)

Definition at line 404 of file KineUtils.cxx.

406 {
407 // Computes q2 (<0) limits @ the input W for inelastic v interactions
408 
409  Range1D_t Q2 = utils::kinematics::InelQ2Lim_W(Ev,M,ml,W,-1.*q2min_cut);
410  Range1D_t q2;
411  q2.min = - Q2.max;
412  q2.max = - Q2.min;
413  return q2;
414 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
A simple [min,max] interval for doubles.
Definition: Range1.h:42
Range1D_t InelQ2Lim_W(double Ev, double M, double ml, double W, double Q2min_cut=controls::kMinQ2Limit)
Definition: KineUtils.cxx:366
double max
Definition: Range1.h:53
double min
Definition: Range1.h:52
Range1D_t genie::utils::kinematics::InelWLim ( double  Ev,
double  M,
double  ml 
)

Definition at line 345 of file KineUtils.cxx.

346 {
347 // Computes W limits for inelastic v interactions
348 //
349  double M2 = TMath::Power(M,2);
350  double s = M2 + 2*M*Ev;
351  assert (s>0);
352 
353  Range1D_t W;
355  W.max = TMath::Sqrt(s) - ml;
356  if(W.max<=W.min) {
357  W.min = -1;
358  W.max = -1;
359  return W;
360  }
363  return W;
364 }
A simple [min,max] interval for doubles.
Definition: Range1.h:42
static const double kPhotontest
Definition: Constants.h:79
static const double kNeutronMass
Definition: Constants.h:76
static const double kASmallNum
Definition: Controls.h:40
double max
Definition: Range1.h:53
double min
Definition: Range1.h:52
static QCString * s
Definition: config.cpp:1042
Range1D_t genie::utils::kinematics::InelXLim ( double  Ev,
double  M,
double  ml 
)

Definition at line 444 of file KineUtils.cxx.

446 {
447 // Computes Bjorken x limits for inelastic v interactions
448 
449  double M2 = TMath::Power(M, 2.);
450  double ml2 = TMath::Power(ml,2.);
451  double s = M2 + 2*M*Ev;
452 
453  SLOG("KineLimits", pDEBUG) << "s = " << s;
454  SLOG("KineLimits", pDEBUG) << "Ev = " << Ev;
455  assert (s>M2);
456 
457  Range1D_t x;
458  x.min = ml2/(s-M2) + controls::kAVerySmallNum;
459  x.max = 1. - controls::kAVerySmallNum;
460 
461  return x;
462 }
A simple [min,max] interval for doubles.
Definition: Range1.h:42
double max
Definition: Range1.h:53
static const double kAVerySmallNum
Definition: Controls.h:39
list x
Definition: train.py:276
double min
Definition: Range1.h:52
#define SLOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a short string (using the FUNCTION and...
Definition: Messenger.h:84
static QCString * s
Definition: config.cpp:1042
#define pDEBUG
Definition: Messenger.h:63
Range1D_t genie::utils::kinematics::InelYLim ( double  Ev,
double  M,
double  ml 
)

Definition at line 464 of file KineUtils.cxx.

465 {
466 // Computes y limits for inelastic v interactions
467 
468  Range1D_t y;
469  y.min = 999;
470  y.max = -999;
471 
472  Range1D_t xl = kinematics::InelXLim(Ev,M,ml);
473  assert(xl.min>0 && xl.max>0);
474 
475  const unsigned int N=100;
476  const double logxmin = TMath::Log10(xl.min);
477  const double logxmax = TMath::Log10(xl.max);
478  const double dlogx = (logxmax-logxmin) / (double)(N-1);
479 
480  for(unsigned int i=0; i<N; i++) {
481  double x = TMath::Power(10, logxmin + i*dlogx);
482 
483  Range1D_t y_x = kinematics::InelYLim_X(Ev,M,ml,x);
484  if(y_x.min>=0 && y_x.min<=1) y.min = TMath::Min(y.min, y_x.min);
485  if(y_x.max>=0 && y_x.max<=1) y.max = TMath::Max(y.max, y_x.max);
486  }
487 
488  if(y.max >= 0 && y.max <= 1 && y.min >= 0 && y.min <= 1) {
489  y.min = TMath::Max(y.min, controls::kAVerySmallNum);
490  y.max = TMath::Min(y.max, 1 - controls::kAVerySmallNum);
491  } else {
492  y.min = -1;
493  y.max = -1;
494  }
495  SLOG("KineLimits", pDEBUG) << "y = [" << y.min << ", " << y.max << "]";
496  return y;
497 }
A simple [min,max] interval for doubles.
Definition: Range1.h:42
Range1D_t InelXLim(double Ev, double M, double ml)
Definition: KineUtils.cxx:444
double max
Definition: Range1.h:53
Range1D_t InelYLim_X(double Ev, double M, double ml, double x)
Definition: KineUtils.cxx:499
static const double kAVerySmallNum
Definition: Controls.h:39
list x
Definition: train.py:276
double min
Definition: Range1.h:52
#define SLOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a short string (using the FUNCTION and...
Definition: Messenger.h:84
#define pDEBUG
Definition: Messenger.h:63
Range1D_t genie::utils::kinematics::InelYLim_X ( double  Ev,
double  M,
double  ml,
double  x 
)

Definition at line 499 of file KineUtils.cxx.

502 {
503 // Computes y limits @ the input x for inelastic v interactions
504 
505  Range1D_t y;
506  y.min = -1;
507  y.max = -1;
508 
509  double Ev2 = TMath::Power(Ev,2);
510  double ml2 = TMath::Power(ml,2);
511 
512  SLOG("KineLimits", pDEBUG) << "x = " << x;
513  SLOG("KineLimits", pDEBUG) << "Ev = " << Ev;
514 
515  assert (Ev>0);
516  assert (x>0&&x<1);
517 
518  double a = 0.5 * ml2/(M*Ev*x);
519  double b = ml2/Ev2;
520  double c = 1 + 0.5*x*M/Ev;
521  double d = TMath::Max(0., TMath::Power(1-a,2.) - b);
522 
523  double A = 0.5 * (1-a-0.5*b)/c;
524  double B = 0.5 * TMath::Sqrt(d)/c;
525 
526  y.min = TMath::Max(0., A-B) + controls::kAVerySmallNum;
527  y.max = TMath::Min(1., A+B) - controls::kAVerySmallNum;
528 
529  return y;
530 }
A simple [min,max] interval for doubles.
Definition: Range1.h:42
const double a
double max
Definition: Range1.h:53
static const double kAVerySmallNum
Definition: Controls.h:39
#define A
Definition: memgrp.cpp:38
static bool * b
Definition: config.cpp:1043
list x
Definition: train.py:276
double min
Definition: Range1.h:52
#define SLOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a short string (using the FUNCTION and...
Definition: Messenger.h:84
#define pDEBUG
Definition: Messenger.h:63
bool genie::utils::kinematics::IsAboveCharmThreshold ( double  x,
double  Q2,
double  M,
double  mc 
)

Definition at line 1225 of file KineUtils.cxx.

1227 {
1228 // x : scaling variable (plain or modified)
1229 // Q2: momentum transfer
1230 // M : hit nucleon "mass" (nucleon can be off the mass shell)
1231 // mc: charm mass
1232 //
1233  double M2 = TMath::Power(M,2);
1234  double v = 0.5*Q2/(M*x);
1235  double W2 = TMath::Max(0., M2+2*M*v-Q2);
1236  double W = TMath::Sqrt(W2);
1237  double Wmin = M + kLightestChmHad;
1238  double xc = utils::kinematics::SlowRescalingVar(x,Q2,M,mc);
1239 
1240  if(xc>=1 || W<=Wmin) return false;
1241  else return true;
1242 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
static const double kLightestChmHad
Definition: Constants.h:78
double SlowRescalingVar(double x, double Q2, double M, double mc)
Definition: KineUtils.cxx:1244
list x
Definition: train.py:276
double genie::utils::kinematics::Jacobian ( const Interaction *const  i,
KinePhaseSpace_t  f,
KinePhaseSpace_t  t 
)

Definition at line 130 of file KineUtils.cxx.

132 {
133 // Returns the Jacobian for a kinematical transformation:
134 // from_ps (x,y,z,...) -> to_ps (u,v,w,...).
135 //
136 // Note on the convention used here:
137 // The differential cross-section d^n sigma / dxdydz..., is transformed to
138 // d^n sigma / dudvdw... using the Jacobian, computed in this function, as:
139 // as d^n sigma / dudvdw... = J * d^n sigma / dxdydz...
140 // where:
141 //
142 // dxdxdz... = J * dudvdw...
143 //
144 // | dx/du dx/dv dx/dw ... |
145 // | |
146 // d {x,y,z,...} | dy/du dy/dv dy/dw ... |
147 // J = ------------ = | |
148 // d {u,v,w,...} | dz/du dz/dv dz/dw ... |
149 // | |
150 // | ... ... ... ... |
151 //
152  SLOG("KineLimits", pDEBUG)
153  << "Computing Jacobian for transformation: "
154  << KinePhaseSpace::AsString(fromps) << " --> "
155  << KinePhaseSpace::AsString(tops);
156 
157  double J=0;
158  bool forward;
159  const Kinematics & kine = i->Kine();
160 
161  // cover the simple case
162  if ( fromps == tops )
163  {
164  forward = true;
165  J = 1.;
166  }
167  //
168  // transformation: {Q2}|E -> {lnQ2}|E
169  //
170  else
171  if ( TransformMatched(fromps,tops,kPSQ2fE,kPSlogQ2fE,forward) )
172  {
173  J = 1. / kine.Q2();
174  }
175 
176  //
177  // transformation: {QD2}|E -> {Q2}|E
178  //
179  else
180  if ( TransformMatched(fromps,tops,kPSQD2fE,kPSQ2fE,forward) )
181  {
182  J = TMath::Power(1+kine.Q2()/controls::kMQD2,-2)/controls::kMQD2;
183  }
184 
185  //
186  // transformation: {x,y}|E -> {lnx,lny}|E
187  //
188  else
189  if ( TransformMatched(fromps,tops,kPSxyfE,kPSlogxlogyfE,forward) )
190  {
191  J = 1. / (kine.x() * kine.y());
192  }
193 
194  //
195  // transformation: {log10x,log10Q2}|E -> {x,Q2}|E
196  //
197  else
198  if ( TransformMatched(fromps,tops,kPSxQ2fE,kPSlog10xlog10Q2fE,forward) )
199  {
200  J = TMath::Log(10.)*kine.x() * TMath::Log(10.)*kine.Q2();
201  }
202 
203  //
204  // transformation: {Q2,y}|E -> {lnQ2,lny}|E
205  //
206  else
207  if ( TransformMatched(fromps,tops,kPSQ2yfE,kPSlogQ2logyfE,forward) )
208  {
209  J = 1. / (kine.Q2() * kine.y());
210  }
211 
212  //
213  // transformation: {Q2,y}|E -> {x,y}|E
214  //
215  else
216  if ( TransformMatched(fromps,tops,kPSQ2yfE,kPSxyfE,forward) )
217  {
218  const InitialState & init_state = i->InitState();
219  double Ev = init_state.ProbeE(kRfHitNucRest);
220  double M = init_state.Tgt().HitNucP4Ptr()->M();
221  double y = kine.y();
222  J = 2*y*Ev*M;
223  }
224 
225  //
226  // transformation: {W,Q2}|E -> {W,lnQ2}|E
227  //
228  else if ( TransformMatched(fromps,tops,kPSWQ2fE,kPSWlogQ2fE,forward) )
229  {
230  J = 1. / kine.Q2();
231  }
232 
233  //
234  // transformation: {W,QD2}|E -> {W,Q2}|E
235  //
236  else
237  if ( TransformMatched(fromps,tops,kPSWQD2fE,kPSWQ2fE,forward) )
238  {
239  J = TMath::Power(1+kine.Q2()/controls::kMQD2,-2)/controls::kMQD2;
240  }
241 
242  //
243  // transformation: {W2,Q2}|E --> {x,y}|E
244  //
245  else
246  if ( TransformMatched(fromps,tops,kPSW2Q2fE,kPSxyfE,forward) )
247  {
248  const InitialState & init_state = i->InitState();
249  double Ev = init_state.ProbeE(kRfHitNucRest);
250  double M = init_state.Tgt().HitNucP4Ptr()->M();
251  double y = kine.y();
252  J = TMath::Power(2*M*Ev,2) * y;
253  }
254 
255  //
256  // transformation: {W,Q2}|E -> {x,y}|E
257  //
258  else
259  if ( TransformMatched(fromps,tops,kPSWQ2fE,kPSxyfE,forward) )
260  {
261  const InitialState & init_state = i->InitState();
262  double Ev = init_state.ProbeE(kRfHitNucRest);
263  double M = init_state.Tgt().HitNucP4Ptr()->M();
264  double y = kine.y();
265  double W = kine.W();
266  J = 2*TMath::Power(M*Ev,2) * y/W;
267  }
268 
269  // Transformation: {Omegalep,Omegapi}|E -> {Omegalep,Thetapi}|E
270  else if ( TransformMatched(fromps,tops,kPSElOlOpifE,kPSElOlTpifE,forward) ) {
271  // Use symmetry to turn 4d integral into 3d * 2pi
272  J = 2*constants::kPi;
273  }
274 
275  // Transformation: {Tl,ctl} -> {W,Q2}|E
276  // IMPORTANT NOTE: This transformation can't be done exactly in two
277  // dimensions due to Fermi motion. For fixed {W,Q2}, the azimuthal angle phi
278  // is uniformly distributed in the hit nucleon rest frame, while for fixed
279  // {Tl,ctl} it is uniform in the lab frame (i.e., the target nucleus rest
280  // frame). A 3D Jacobian is needed in order to account for the relationship
281  // between these two azimuthal angles. In the implementation below, I choose
282  // to neglect Fermi motion. Under this approximation, the lab and hit nucleon
283  // rest frames are the same. Doing things this way is a stopgap solution for
284  // the new XSecShape_CCMEC reweighting dial developed for MicroBooNE. A full
285  // solution should use the Jacobian that connects the 3D phase spaces which
286  // each include phi. - S. Gardiner, 29 July 2020
287  else if ( TransformMatched(fromps, tops, kPSTlctl, kPSWQ2fE, forward) )
288  {
289  // Probe properties (mass, energy, momentum)
290  const InitialState& init_state = i->InitState();
291  double mv = init_state.Probe()->Mass();
292  double Ev = init_state.ProbeE( kRfLab );
293  double pv = std::sqrt( std::max(0., Ev*Ev - mv*mv) );
294 
295  // Invariant mass of the initial hit nucleon
296  const TLorentzVector& hit_nuc_P4 = init_state.Tgt().HitNucP4();
297  double M = hit_nuc_P4.M();
298 
299  // Outgoing lepton mass
300  double ml = i->FSPrimLepton()->Mass();
301 
302  double W = i->Kine().GetKV( kKVW );
303  double Q2 = i->Kine().GetKV( kKVQ2 );
304  double El = Ev - ( (W*W + Q2 - M*M) / (2.*M) );
305  double pl = std::sqrt( std::max(0., El*El - ml*ml) );
306 
307  // Compute the Jacobian for {Tl, ctl} --> {W, Q2}
308  // (it will be inverted below for the inverse transformation)
309  J = W / ( 2. * pv * pl * M );
310  }
311 
312  else {
313  std::ostringstream msg;
314  msg << "Can not compute Jacobian for transforming: "
315  << KinePhaseSpace::AsString(fromps) << " --> "
316  << KinePhaseSpace::AsString(tops);
317  SLOG("KineLimits", pFATAL) << "*** " << msg.str();
319  //exit(1);
320  }
321 
322  // if any of the above transforms was reverse, invert the jacobian
323  if(!forward) J = 1./J;
324 
325  return J;
326 }
double W(bool selected=false) const
Definition: Kinematics.cxx:157
bool TransformMatched(KinePhaseSpace_t ia, KinePhaseSpace_t ib, KinePhaseSpace_t a, KinePhaseSpace_t b, bool &fwd)
Definition: KineUtils.cxx:328
double J(double q0, double q3, double Enu, double ml)
Definition: MECUtils.cxx:147
void msg(const char *fmt,...)
Definition: message.cpp:107
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
#define pFATAL
Definition: Messenger.h:56
Generated/set kinematical variables for an event.
Definition: Kinematics.h:39
TParticlePDG * Probe(void) const
double x(bool selected=false) const
Definition: Kinematics.cxx:99
double y(bool selected=false) const
Definition: Kinematics.cxx:112
const TLorentzVector & HitNucP4(void) const
Definition: Target.h:91
Exception used inside Interaction classes.
const Kinematics & Kine(void) const
Definition: Interaction.h:71
double GetKV(KineVar_t kv) const
Definition: Kinematics.cxx:323
static int max(int a, int b)
TParticlePDG * FSPrimLepton(void) const
final state primary lepton
static const double kMQD2
Definition: Controls.h:65
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:247
const InitialState & InitState(void) const
Definition: Interaction.h:69
const char * AsString(Resonance_t res)
resonance id -> string
double Q2(bool selected=false) const
Definition: Kinematics.cxx:125
const Target & Tgt(void) const
Definition: InitialState.h:66
#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
static const double kPi
Definition: Constants.h:37
Initial State information.
Definition: InitialState.h:48
#define pDEBUG
Definition: Messenger.h:63
double genie::utils::kinematics::PhaseSpaceVolume ( const Interaction *const  i,
KinePhaseSpace_t  ps 
)

Definition at line 36 of file KineUtils.cxx.

38 {
39  double vol = 0;
40 
41  const KPhaseSpace & phase_space = in->PhaseSpace();
42 
43  switch(ps) {
44 
45  case(kPSQ2fE):
46  {
47  Range1D_t Q2 = phase_space.Limits(kKVQ2);
48  vol = Q2.max - Q2.min;
49  return vol;
50  break;
51  }
52  case(kPSq2fE):
53  {
54  Range1D_t q2 = phase_space.Limits(kKVq2);
55  vol = q2.min - q2.max;
56  return vol;
57  break;
58  }
59  case(kPSWfE):
60  {
61  Range1D_t W = phase_space.Limits(kKVW);
62  vol = W.max - W.min;
63  return vol;
64  break;
65  }
66  case(kPSWQ2fE):
67  {
68  Range1D_t W = phase_space.Limits(kKVW);
69  if(W.max<0) return 0;
70  const int kNW = 100;
71  double dW = (W.max-W.min)/(kNW-1);
73  for(int iw=0; iw<kNW; iw++) {
74  interaction.KinePtr()->SetW(W.min + iw*dW);
75  Range1D_t Q2 = interaction.PhaseSpace().Q2Lim_W();
76  double dQ2 = (Q2.max-Q2.min);
77  vol += (dW*dQ2);
78  }
79  return vol;
80  break;
81  }
82  case(kPSxyfE):
83  {
84  Range1D_t W = phase_space.Limits(kKVW);
85  if(W.max<0) return 0;
86 
87  const InitialState & init_state = in->InitState();
88  double Ev = init_state.ProbeE(kRfHitNucRest);
89  double M = init_state.Tgt().HitNucP4Ptr()->M();
90 
91  const int kNx = 100;
92  const int kNy = 100;
93  const double kminx = controls::kMinX;
94  const double kminy = controls::kMinY;
95  const double kdx = (controls::kMaxX - kminx) / (kNx-1);
96  const double kdy = (controls::kMaxY - kminy) / (kNy-1);
97  const double kdV = kdx*kdy;
98 
99  double cW=-1, cQ2 = -1;
100 
102 
103  for(int ix=0; ix<kNx; ix++) {
104  double x = kminx+ix*kdx;
105  for(int iy=0; iy<kNy; iy++) {
106  double y = kminy+iy*kdy;
107 
108  XYtoWQ2(Ev, M, cW, cQ2, x, y);
109  if(!math::IsWithinLimits(cW, W)) continue;
110 
111  interaction.KinePtr()->SetW(cW);
112  Range1D_t Q2 = interaction.PhaseSpace().Q2Lim_W();
113  if(!math::IsWithinLimits(cQ2, Q2)) continue;
114 
115  vol += kdV;
116  }
117  }
118  return vol;
119  break;
120  }
121  default:
122  SLOG("KineLimits", pERROR)
123  << "Couldn't compute phase space volume for "
124  << KinePhaseSpace::AsString(ps) << "using \n" << *in;
125  return 0;
126  }
127  return 0;
128 }
static const double kMaxX
Definition: Controls.h:44
#define pERROR
Definition: Messenger.h:59
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
A simple [min,max] interval for doubles.
Definition: Range1.h:42
static const double kMaxY
Definition: Controls.h:46
bool IsWithinLimits(double x, Range1D_t range)
Definition: MathUtils.cxx:258
static const double kMinY
Definition: Controls.h:45
Range1D_t Limits(KineVar_t kvar) const
Return the kinematical variable limits.
Summary information for an interaction.
Definition: Interaction.h:56
void XYtoWQ2(double Ev, double M, double &W, double &Q2, double x, double y)
Definition: KineUtils.cxx:1142
Kinematical phase space.
Definition: KPhaseSpace.h:33
static const double kMinX
Definition: Controls.h:43
static constexpr double ps
Definition: Units.h:99
double max
Definition: Range1.h:53
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:247
const char * AsString(Resonance_t res)
resonance id -> string
list x
Definition: train.py:276
double min
Definition: Range1.h:52
const Target & Tgt(void) const
Definition: InitialState.h:66
#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
Initial State information.
Definition: InitialState.h:48
double genie::utils::kinematics::Q2 ( const Interaction *const  i)

Definition at line 1064 of file KineUtils.cxx.

1065 {
1066 // Get Q^2 from kinematics object
1067 // If Q^2 is not set but x,y are, then compute Q2 from x,y
1068 
1069  const Kinematics & kinematics = interaction->Kine();
1070 
1071  if (kinematics.KVSet(kKVQ2) || kinematics.KVSet(kKVq2)) {
1072  double Q2 = kinematics.Q2();
1073  return Q2;
1074  }
1075  if (kinematics.KVSet(kKVy)) {
1076  const InitialState & init_state = interaction->InitState();
1077  double Mn = init_state.Tgt().HitNucP4Ptr()->M(); // can be off m/shell
1078  double x = kinematics.x();
1079  double y = kinematics.y();
1080  double Ev = init_state.ProbeE(kRfHitNucRest);
1081  double Q2 = 2*Mn*Ev*x*y;
1082  return Q2;
1083  }
1084  SLOG("KineLimits", pERROR) << "Couldn't compute Q^2 for \n"<< *interaction;
1085  return 0;
1086 }
#define pERROR
Definition: Messenger.h:59
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
Generated/set kinematical variables for an event.
Definition: Kinematics.h:39
double x(bool selected=false) const
Definition: Kinematics.cxx:99
double y(bool selected=false) const
Definition: Kinematics.cxx:112
bool KVSet(KineVar_t kv) const
Definition: Kinematics.cxx:317
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:247
list x
Definition: train.py:276
double Q2(bool selected=false) const
Definition: Kinematics.cxx:125
const Target & Tgt(void) const
Definition: InitialState.h:66
#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
Initial State information.
Definition: InitialState.h:48
double genie::utils::kinematics::Q2toQD2 ( double  Q2)

Definition at line 1048 of file KineUtils.cxx.

1049 {
1050 // Q2 -> QD2 transformation. QD2 takes out the dipole form of the form factors
1051 // making the differential cross section to be flatter and speeding up the
1052 // kinematical selection.
1053 
1054  assert(Q2>0);
1055  return TMath::Power(1+Q2/controls::kMQD2, -1);
1056 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
static const double kMQD2
Definition: Controls.h:65
double genie::utils::kinematics::Q2YtoX ( double  Ev,
double  M,
double  Q2,
double  y 
)

Definition at line 1209 of file KineUtils.cxx.

1211 {
1212 // Converts (Q^2,y) => x
1213 // Ev is the neutrino energy at the struck nucleon rest frame
1214 // M is the nucleon mass - it does not need to be on the mass shell
1215  assert(Ev > 0. && M > 0. && Q2 > 0. && y > 0.);
1216 
1217  double x = Q2 / (2. * y * M * Ev);
1218 
1219  LOG("KineLimits", pDEBUG) << "(Ev=" << Ev << ",Q2=" << Q2
1220  << ",y=" << y << ",M=" << M << ") => x=" << x;
1221 
1222  return x;
1223 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
list x
Definition: train.py:276
#define pDEBUG
Definition: Messenger.h:63
double genie::utils::kinematics::QD2toQ2 ( double  QD2)

Definition at line 1058 of file KineUtils.cxx.

1059 {
1060  assert(QD2>0);
1061  return controls::kMQD2*(1/QD2-1);
1062 }
static const double kMQD2
Definition: Controls.h:65
double genie::utils::kinematics::RESImportanceSamplingEnvelope ( double *  x,
double *  par 
)

Definition at line 1359 of file KineUtils.cxx.

1361 {
1362  //-- inputs
1363  double QD2 = x[0]; // QD2 (Q2 transformed to take out the dipole form)
1364  double W = x[1]; // invariant mass
1365 
1366  //-- parameters
1367  double mD = par[0]; // resonance mass
1368  double gD = par[1]; // resonance width
1369  double xsmax = par[2]; // safety factor * max cross section in (W,Q2)
1370  double Wmax = par[3]; // kinematically allowed Wmax
1371  //double E = par[4]; // neutrino energy
1372 
1373 // return 3*xsmax; ////////////////NOTE
1374 
1375  xsmax*=5;
1376 
1377  double func = 0;
1378 
1379  if(Wmax > mD) {
1380  // -- if the resonance mass is within the kinematical range then
1381  // -- the envelope is defined as a plateau above the resonance peak,
1382  // -- a steeply falling leading edge (low-W side) and a more slowly
1383  // -- falling trailing edge (high-W side)
1384 
1385  double hwfe = mD+2*gD; // high W falling edge
1386  double lwfe = mD-gD/2; // low W falling edge
1387 
1388  if(W < lwfe) {
1389  //low-W falling edge
1390  func = xsmax / (1 + 5* TMath::Power((W-lwfe)/gD,2));
1391  } else if (W > hwfe) {
1392  //high-W falling edge
1393  func = xsmax / (1 + TMath::Power((W-hwfe)/gD,2));
1394  } else {
1395  // plateau
1396  func = xsmax;
1397  }
1398  } else {
1399  // -- if the resonance mass is above the kinematical range then the
1400  // -- envelope is a small plateau just bellow Wmax and a falling edge
1401  // -- at lower W
1402 
1403  double plateau_edge = Wmax-0.1;
1404 
1405  if (W > plateau_edge) {
1406  // plateau
1407  func = xsmax;
1408  } else {
1409  //low-W falling edge
1410  func = xsmax / (1 + TMath::Power((W-plateau_edge)/gD,2));
1411  }
1412  }
1413 
1414  if(QD2<0.5) {
1415  func *= (1 - (0.5-QD2));
1416  }
1417 
1418  return func;
1419 }
def func()
Definition: docstring.py:7
list x
Definition: train.py:276
double genie::utils::kinematics::SlowRescalingVar ( double  x,
double  Q2,
double  M,
double  mc 
)

Definition at line 1244 of file KineUtils.cxx.

1246 {
1247 // x : scaling variable (plain or modified)
1248 // Q2: momentum transfer
1249 // M : hit nucleon "mass" (nucleon can be off the mass shell)
1250 // mc: charm mass
1251 
1252  double mc2 = TMath::Power(mc,2);
1253  double v = 0.5*Q2/(M*x);
1254  double xc = x + 0.5*mc2/(M*v);
1255  return xc;
1256 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
list x
Definition: train.py:276
bool genie::utils::kinematics::TransformMatched ( KinePhaseSpace_t  ia,
KinePhaseSpace_t  ib,
KinePhaseSpace_t  a,
KinePhaseSpace_t  b,
bool fwd 
)

Definition at line 328 of file KineUtils.cxx.

331 {
332 
333  // match a->b or b->a
334  bool matched = ( (a==inpa&&b==inpb) || (a==inpb&&b==inpa) );
335 
336  if(matched) {
337  if(a==inpa&&b==inpb) fwd = true; // matched forward transform: a->b
338  else fwd = false; // matched reverse transform: b->a
339  }
340  return matched;
341 }
const double a
static bool * b
Definition: config.cpp:1043
void genie::utils::kinematics::UpdateWQ2FromXY ( const Interaction in)

Definition at line 1277 of file KineUtils.cxx.

1278 {
1279  Kinematics * kine = in->KinePtr();
1280 
1281  if(kine->KVSet(kKVx) && kine->KVSet(kKVy)) {
1282  const InitialState & init_state = in->InitState();
1283  double Ev = init_state.ProbeE(kRfHitNucRest);
1284  double M = init_state.Tgt().HitNucP4Ptr()->M(); // can be off mass shell
1285  double x = kine->x();
1286  double y = kine->y();
1287 
1288  double Q2=-1,W=-1;
1289  kinematics::XYtoWQ2(Ev,M,W,Q2,x,y);
1290  kine->SetQ2(Q2);
1291  kine->SetW(W);
1292  }
1293 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
void SetQ2(double Q2, bool selected=false)
Definition: Kinematics.cxx:255
Kinematics * KinePtr(void) const
Definition: Interaction.h:76
Generated/set kinematical variables for an event.
Definition: Kinematics.h:39
double x(bool selected=false) const
Definition: Kinematics.cxx:99
double y(bool selected=false) const
Definition: Kinematics.cxx:112
bool KVSet(KineVar_t kv) const
Definition: Kinematics.cxx:317
void XYtoWQ2(double Ev, double M, double &W, double &Q2, double x, double y)
Definition: KineUtils.cxx:1142
void SetW(double W, bool selected=false)
Definition: Kinematics.cxx:279
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:247
const InitialState & InitState(void) const
Definition: Interaction.h:69
list x
Definition: train.py:276
const Target & Tgt(void) const
Definition: InitialState.h:66
double ProbeE(RefFrame_t rf) const
Initial State information.
Definition: InitialState.h:48
void genie::utils::kinematics::UpdateWYFromXQ2 ( const Interaction in)

Definition at line 1313 of file KineUtils.cxx.

1314 {
1315  Kinematics * kine = in->KinePtr();
1316 
1317  if(kine->KVSet(kKVx) && kine->KVSet(kKVQ2)) {
1318  const InitialState & init_state = in->InitState();
1319  double Ev = init_state.ProbeE(kRfHitNucRest);
1320  double M = init_state.Tgt().HitNucP4Ptr()->M(); // can be off mass shell
1321  double x = kine->x();
1322  double Q2 = kine->Q2();
1323 
1324  double W=-1,y=-1;
1325  kinematics::XQ2toWY(Ev,M,W,Q2,x,y);
1326  kine->SetW(W);
1327  kine->Sety(y);
1328  }
1329 }
void XQ2toWY(double Ev, double M, double &W, double Q2, double x, double &y)
Definition: KineUtils.cxx:1160
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
Kinematics * KinePtr(void) const
Definition: Interaction.h:76
Generated/set kinematical variables for an event.
Definition: Kinematics.h:39
double x(bool selected=false) const
Definition: Kinematics.cxx:99
bool KVSet(KineVar_t kv) const
Definition: Kinematics.cxx:317
void SetW(double W, bool selected=false)
Definition: Kinematics.cxx:279
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:247
void Sety(double y, bool selected=false)
Definition: Kinematics.cxx:243
const InitialState & InitState(void) const
Definition: Interaction.h:69
list x
Definition: train.py:276
double Q2(bool selected=false) const
Definition: Kinematics.cxx:125
const Target & Tgt(void) const
Definition: InitialState.h:66
double ProbeE(RefFrame_t rf) const
Initial State information.
Definition: InitialState.h:48
void genie::utils::kinematics::UpdateXFromQ2Y ( const Interaction in)

Definition at line 1331 of file KineUtils.cxx.

1332 {
1333  Kinematics * kine = in->KinePtr();
1334 
1335  if(kine->KVSet(kKVy) && kine->KVSet(kKVQ2)) {
1336 
1337  const ProcessInfo & pi = in->ProcInfo();
1338  const InitialState & init_state = in->InitState();
1339  double M = 0.0;
1340  double Ev = 0.0;
1341 
1342  if (pi.IsCoherentProduction()) {
1343  M = in->InitState().Tgt().Mass(); // nucleus mass
1344  Ev = init_state.ProbeE(kRfLab);
1345  }
1346  else {
1347  M = in->InitState().Tgt().HitNucP4Ptr()->M(); //can be off m/shell
1348  Ev = init_state.ProbeE(kRfHitNucRest);
1349  }
1350 
1351  double y = kine->y();
1352  double Q2 = kine->Q2();
1353 
1354  double x = kinematics::Q2YtoX(Ev,M,Q2,y);
1355  kine->Setx(x);
1356  }
1357 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
Kinematics * KinePtr(void) const
Definition: Interaction.h:76
Generated/set kinematical variables for an event.
Definition: Kinematics.h:39
bool IsCoherentProduction(void) const
Definition: ProcessInfo.cxx:99
double Mass(void) const
Definition: Target.cxx:224
double y(bool selected=false) const
Definition: Kinematics.cxx:112
bool KVSet(KineVar_t kv) const
Definition: Kinematics.cxx:317
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:46
void Setx(double x, bool selected=false)
Definition: Kinematics.cxx:231
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:247
double Q2YtoX(double Ev, double M, double Q2, double y)
Definition: KineUtils.cxx:1209
const InitialState & InitState(void) const
Definition: Interaction.h:69
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
list x
Definition: train.py:276
float pi
Definition: units.py:11
double Q2(bool selected=false) const
Definition: Kinematics.cxx:125
const Target & Tgt(void) const
Definition: InitialState.h:66
double ProbeE(RefFrame_t rf) const
Initial State information.
Definition: InitialState.h:48
void genie::utils::kinematics::UpdateXYFromWQ2 ( const Interaction in)

Definition at line 1295 of file KineUtils.cxx.

1296 {
1297  Kinematics * kine = in->KinePtr();
1298 
1299  if(kine->KVSet(kKVW) && kine->KVSet(kKVQ2)) {
1300  const InitialState & init_state = in->InitState();
1301  double Ev = init_state.ProbeE(kRfHitNucRest);
1302  double M = init_state.Tgt().HitNucP4Ptr()->M(); // can be off mass shell
1303  double W = kine->W();
1304  double Q2 = kine->Q2();
1305 
1306  double x=-1,y=-1;
1307  kinematics::WQ2toXY(Ev,M,W,Q2,x,y);
1308  kine->Setx(x);
1309  kine->Sety(y);
1310  }
1311 }
double W(bool selected=false) const
Definition: Kinematics.cxx:157
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
Kinematics * KinePtr(void) const
Definition: Interaction.h:76
Generated/set kinematical variables for an event.
Definition: Kinematics.h:39
bool KVSet(KineVar_t kv) const
Definition: Kinematics.cxx:317
void WQ2toXY(double Ev, double M, double W, double Q2, double &x, double &y)
Definition: KineUtils.cxx:1119
void Setx(double x, bool selected=false)
Definition: Kinematics.cxx:231
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:247
void Sety(double y, bool selected=false)
Definition: Kinematics.cxx:243
const InitialState & InitState(void) const
Definition: Interaction.h:69
list x
Definition: train.py:276
double Q2(bool selected=false) const
Definition: Kinematics.cxx:125
const Target & Tgt(void) const
Definition: InitialState.h:66
double ProbeE(RefFrame_t rf) const
Initial State information.
Definition: InitialState.h:48
double genie::utils::kinematics::W ( const Interaction *const  i)

Definition at line 1088 of file KineUtils.cxx.

1089 {
1090  const ProcessInfo & process_info = interaction->ProcInfo();
1091 
1092  if(process_info.IsQuasiElastic()) {
1093  // hadronic inv. mass is equal to the recoil nucleon on-shell mass
1094  int rpdgc = interaction->RecoilNucleonPdg();
1095  double M = PDGLibrary::Instance()->Find(rpdgc)->Mass();
1096  return M;
1097  }
1098 
1099  const Kinematics & kinematics = interaction->Kine();
1100  if(kinematics.KVSet(kKVW)) {
1101  double W = kinematics.W();
1102  return W;
1103  }
1104  if(kinematics.KVSet(kKVx) && kinematics.KVSet(kKVy)) {
1105  const InitialState & init_state = interaction->InitState();
1106  double Ev = init_state.ProbeE(kRfHitNucRest);
1107  double M = init_state.Tgt().HitNucP4Ptr()->M();
1108  double M2 = M*M;
1109  double x = kinematics.x();
1110  double y = kinematics.y();
1111  double W2 = TMath::Max(0., M2 + 2*Ev*M*y*(1-x));
1112  double W = TMath::Sqrt(W2);
1113  return W;
1114  }
1115  SLOG("KineLimits", pERROR) << "Couldn't compute W for \n"<< *interaction;
1116  return 0;
1117 }
double W(bool selected=false) const
Definition: Kinematics.cxx:157
#define pERROR
Definition: Messenger.h:59
bool IsQuasiElastic(void) const
Definition: ProcessInfo.cxx:69
Generated/set kinematical variables for an event.
Definition: Kinematics.h:39
double x(bool selected=false) const
Definition: Kinematics.cxx:99
double y(bool selected=false) const
Definition: Kinematics.cxx:112
bool KVSet(KineVar_t kv) const
Definition: Kinematics.cxx:317
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:46
TLorentzVector * HitNucP4Ptr(void) const
Definition: Target.cxx:247
list x
Definition: train.py:276
const Target & Tgt(void) const
Definition: InitialState.h:66
#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
Initial State information.
Definition: InitialState.h:48
void genie::utils::kinematics::WQ2toXY ( double  Ev,
double  M,
double  W,
double  Q2,
double &  x,
double &  y 
)

Definition at line 1119 of file KineUtils.cxx.

1121 {
1122 // Converts (W,Q2) => (x,y)
1123 // Uses the system: a) W^2 - M^2 = 2*Ev*M*y*(1-x) and b) Q^2 = 2*x*y*M*Ev
1124 // Ev is the neutrino energy at the struck nucleon rest frame
1125 // M is the nucleon mass - it does not need to be on the mass shell
1126 
1127  double M2 = TMath::Power(M,2);
1128  double W2 = TMath::Power(W,2);
1129 
1130  x = Q2 / (W2-M2+Q2);
1131  y = (W2-M2+Q2) / (2*M*Ev);
1132 
1133  x = TMath::Min(1.,x);
1134  y = TMath::Min(1.,y);
1135  x = TMath::Max(0.,x);
1136  y = TMath::Max(0.,y);
1137 
1138  LOG("KineLimits", pDEBUG)
1139  << "(W=" << W << ",Q2=" << Q2 << ") => (x="<< x << ", y=" << y<< ")";
1140 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
list x
Definition: train.py:276
#define pDEBUG
Definition: Messenger.h:63
void genie::utils::kinematics::XQ2toWY ( double  Ev,
double  M,
double &  W,
double  Q2,
double  x,
double &  y 
)

Definition at line 1160 of file KineUtils.cxx.

1162 {
1163 // Converts (x,Q2) => (W,Y)
1164 // Uses the system: a) W^2 - M^2 = 2*Ev*M*y*(1-x) and b) Q^2 = 2*x*y*M*Ev
1165 // Ev is the neutrino energy at the struck nucleon rest frame
1166 // M is the nucleon mass - it does not need to be on the mass shell
1167 
1168  double M2 = TMath::Power(M,2);
1169  y = Q2 / (2 * x * M * Ev);
1170  y = TMath::Min(1.,y);
1171 
1172  double W2 = M2 + 2*Ev*M*y*(1-x);
1173  W = TMath::Sqrt(TMath::Max(0., W2));
1174 
1175  LOG("KineLimits", pDEBUG)
1176  << "(x=" << x << ",Q2=" << Q2 << " => (W=" << W << ",Y=" << y << ")";
1177 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
list x
Definition: train.py:276
#define pDEBUG
Definition: Messenger.h:63
double genie::utils::kinematics::XYtoQ2 ( double  Ev,
double  M,
double  x,
double  y 
)

Definition at line 1195 of file KineUtils.cxx.

1197 {
1198 // Converts (x,y) => Q2
1199 // Ev is the neutrino energy at the struck nucleon rest frame
1200 // M is the nucleon mass - it does not need to be on the mass shell
1201 
1202  double Q2 = 2*x*y*M*Ev;
1203 
1204  LOG("KineLimits", pDEBUG) << "(x=" << x << ",y=" << y << ") => Q2=" << Q2;
1205 
1206  return Q2;
1207 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
list x
Definition: train.py:276
#define pDEBUG
Definition: Messenger.h:63
double genie::utils::kinematics::XYtoW ( double  Ev,
double  M,
double  x,
double  y 
)

Definition at line 1179 of file KineUtils.cxx.

1181 {
1182 // Converts (x,y) => W
1183 // Ev is the neutrino energy at the struck nucleon rest frame
1184 // M is the nucleon mass - it does not need to be on the mass shell
1185 
1186  double M2 = TMath::Power(M,2);
1187  double W2 = M2 + 2*Ev*M*y*(1-x);
1188  double W = TMath::Sqrt(TMath::Max(0., W2));
1189 
1190  LOG("KineLimits", pDEBUG) << "(x=" << x << ",y=" << y << ") => W=" << W;
1191 
1192  return W;
1193 }
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
list x
Definition: train.py:276
#define pDEBUG
Definition: Messenger.h:63
void genie::utils::kinematics::XYtoWQ2 ( double  Ev,
double  M,
double &  W,
double &  Q2,
double  x,
double  y 
)

Definition at line 1142 of file KineUtils.cxx.

1144 {
1145 // Converts (x,y) => (W,Q2)
1146 // Uses the system: a) W^2 - M^2 = 2*Ev*M*y*(1-x) and b) Q^2 = 2*x*y*M*Ev
1147 // Ev is the neutrino energy at the struck nucleon rest frame
1148 // M is the nucleon mass - it does not need to be on the mass shell
1149 
1150  double M2 = TMath::Power(M,2);
1151  double W2 = M2 + 2*Ev*M*y*(1-x);
1152 
1153  W = TMath::Sqrt(TMath::Max(0., W2));
1154  Q2 = 2*x*y*M*Ev;
1155 
1156  LOG("KineLimits", pDEBUG)
1157  << "(x=" << x << ",y=" << y << " => (W=" << W << ",Q2=" << Q2 << ")";
1158 }
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1064
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
list x
Definition: train.py:276
#define pDEBUG
Definition: Messenger.h:63