All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Classes | Enumerations | Functions
evgen::ldm Namespace Reference

Classes

class  MeVPrtlGen
 
class  MeVPrtlTestRayTrace
 
class  Constants
 
class  EmptyKaonGen
 EmptyKaonGen class definiton. More...
 
class  HiggsMakeDecay
 HiggsMakeDecay class definiton. More...
 
class  Kaon2HiggsFlux
 Kaon2HiggsFlux class definiton Implementation of model taken from: https://arxiv.org/abs/1909.11670. More...
 
class  MonoEnergyHiggsFlux
 MonoEnergyHiggsFlux class definiton. More...
 
class  HNLMakeDecay
 HNLMakeDecay class definiton. More...
 
class  Kaon2HNLFlux
 Kaon2HNLFlux class definiton. More...
 
class  IMesonGen
 IMesonGen interface class definiton. More...
 
class  IMeVPrtlDecay
 IMeVPrtlDecay interface class definiton. More...
 
class  IMeVPrtlFlux
 IMeVPrtlFlux interface class definiton. More...
 
class  IMeVPrtlStage
 IMeVPrtlStage interface class definiton. General interface behind each stage. Provides random number generation. More...
 
class  IRayTrace
 IRayTrace interface class definiton. More...
 
struct  RayWeightInfo
 
class  MixedWeightRayTraceBox
 MixedWeightRayTraceBox class definiton. More...
 
class  NuMiKaonGen
 NuMiKaonGen class definiton. More...
 
class  RayTraceBox
 RayTraceBox class definiton. More...
 
class  ReThrowRayTraceBox
 ReThrowRayTraceBox class definiton. More...
 
class  WeightedRayTraceBox
 WeightedRayTraceBox class definiton. More...
 
class  KaonParent
 
class  MeVPrtlDecay
 
class  MeVPrtlFlux
 
class  MeVPrtlTruth
 

Enumerations

enum  Generator { kUnknown =-1, kDissonantHiggs =0, kHNL =1 }
 

Functions

double twobody_momentum (double parent_mass, double childA_mass, double childB_mass)
 
int calcPrtlRayWgt (double rest_frame_p, double M, TVector3 dir, TVector3 boost, double rand, double &lab_frame_p_out, double &costh_rest_out, double &wgt)
 
double minKinematicCosTheta (double parentM, double secM, double prtlM, double parentE)
 
double forwardPrtlEnergy (double parentM, double secM, double prtlM, double parentE)
 
double secPDG2Mass (int pdg)
 
double flat_to_exp_rand (double x, double mean, double a, double b)
 
double forcedecay_weight (double mean, double a, double b)
 
double higgs_momentum (double kaon_mass, double pion_mass, double higs_mass)
 
double LeptonPartialWidth (double lep_mass, double higs_mass, double mixing)
 
double ElectronPartialWidth (double higs_mass, double mixing)
 
double MuonPartialWidth (double higs_mass, double mixing)
 
double PionPartialWidth (double pion_mass, double higs_mass, double mixing)
 
double PiPlusPartialWidth (double higs_mass, double mixing)
 
double PiZeroPartialWidth (double higs_mass, double mixing)
 
double KaonPlusBranchingRatio (double higs_mass, double mixing)
 
double KaonLongBranchingRatio (double higs_mass, double mixing)
 
double SMKaonBR (int kaon_pdg)
 
int PionPdg (int kaon_pdg)
 
double HNLNuDiLepLNVDalitz (TLorentzVector K, TLorentzVector LA, TLorentzVector LB, TLorentzVector LBD, TLorentzVector NU)
 
double HNLLepPiLNVDalitz (TLorentzVector K, TLorentzVector LA, TLorentzVector N, TLorentzVector PI, TLorentzVector LB)
 
double HNLLepPiLNCDalitz (TLorentzVector K, TLorentzVector LA, TLorentzVector N, TLorentzVector PI, TLorentzVector LB)
 
double HNLLepPiDalitzMax (double mK, double mA, double mN, double mP, double mB)
 
double lambda (double a, double b, double c)
 
double I1_integrand (double s, void *param)
 
double I2_integrand (double s, void *param)
 
double hnl_momentum (double kaon_mass, double lep_mass, double hnl_mass)
 
double BranchingRatio (double hnl_mass, double u4, bool is_muon)
 
std::pair< double, bool > Branch (double hnl_mass, double ue4, double um4, double rand)
 
double QEstimator (unsigned nsuccess, unsigned nfail, unsigned r)
 

Enumeration Type Documentation

Enumerator
kUnknown 
kDissonantHiggs 
kHNL 

Definition at line 12 of file MeVPrtlTruth.h.

Function Documentation

std::pair<double, bool> evgen::ldm::Branch ( double  hnl_mass,
double  ue4,
double  um4,
double  rand 
)

Definition at line 137 of file Kaon2HNLFlux_tool.cc.

137  {
138  double kplus_mass = Constants::Instance().kplus_mass;
139 
140  double br_muon = (um4 > 0. && hnl_mass < kplus_mass - Constants::Instance().muon_mass) ? BranchingRatio(hnl_mass, um4, true) : 0.;
141  double br_elec = (ue4 > 0. && hnl_mass < kplus_mass - Constants::Instance().elec_mass) ? BranchingRatio(hnl_mass, ue4, false): 0.;
142  if (br_muon == 0. && br_elec == 0.) return {0., false};
143 
144  double br_weight = br_muon + br_elec;
145 
146  bool is_muon = rand < (br_muon / br_weight);
147 
148  return {br_weight, is_muon};
149 }
double BranchingRatio(double hnl_mass, double u4, bool is_muon)
double evgen::ldm::BranchingRatio ( double  hnl_mass,
double  u4,
bool  is_muon 
)

Definition at line 120 of file Kaon2HNLFlux_tool.cc.

120  {
121  double lep_mass = is_muon ? Constants::Instance().muon_mass : Constants::Instance().elec_mass;
122  double kplus_mass = Constants::Instance().kplus_mass;
123 
124  if (hnl_mass > kplus_mass - lep_mass) return 0.;
125 
126  double smbr = is_muon ? Constants::Instance().kaonp_mup_numu : Constants::Instance().kaonp_ep_nue;
127  double lep_ratio = (lep_mass * lep_mass) / (kplus_mass * kplus_mass);
128  double hnl_ratio = (hnl_mass * hnl_mass) / (kplus_mass * kplus_mass);
129  double kinematic_factor = (lep_ratio + hnl_ratio - (lep_ratio - hnl_ratio) * (lep_ratio - hnl_ratio)) \
130  * sqrt(1 + hnl_ratio * hnl_ratio + lep_ratio * lep_ratio - 2*(hnl_ratio + lep_ratio + hnl_ratio*lep_ratio)) \
131  / (lep_ratio * (1. - lep_ratio) * (1. - lep_ratio));
132 
133  // scale the branching ratio
134  return smbr * (u4 / (1. - u4)) * kinematic_factor;
135 }
int evgen::ldm::calcPrtlRayWgt ( double  rest_frame_p,
double  M,
TVector3  dir,
TVector3  boost,
double  rand,
double &  lab_frame_p_out,
double &  costh_rest_out,
double &  wgt 
)

Definition at line 102 of file Constants.cpp.

104 {
105  // preset values
106  lab_frame_p_out = 0.;
107  costh_rest_out = 0.;
108  wgt = 0.;
109 
110  double beta = boost.Mag();
111  double gamma = 1. / sqrt(1 - beta*beta);
112  double E_rest = sqrt(rest_frame_p*rest_frame_p + M*M);
113 
114  // cos angle between kaon direction and daughter direction (in lab frame)
115  double costh_lab = boost.Unit().Dot(dir.Unit());
116  double sinth_lab_sq = 1 - costh_lab*costh_lab;
117 
118  // The problem: we are given the rest frame momentum and lab frame angle. Given this,
119  // determine the lab frame momentum and rest frame angle.
120  //
121  // Unlike in the massless case, there are 1 OR 0 OR 2 possible solutions to this problem.
122  auto costh_rest = [costh_lab, sinth_lab_sq, beta, gamma, M, rest_frame_p, E_rest](int SIGN) { return -(-SIGN * costh_lab*\
123  sqrt(-M*M + (E_rest*E_rest)/(gamma*gamma) + (M*M)*(beta*beta)*(costh_lab*costh_lab)) \
124  +E_rest*beta*gamma*sinth_lab_sq) / \
125  (rest_frame_p*gamma*(1 - beta*beta*costh_lab*costh_lab)); };
126 
127  double costh_rest_plus = costh_rest(1);
128  double costh_rest_minus = costh_rest(-1);
129 
130  auto lab_frame_p = [costh_lab, beta, gamma, M, E_rest](int SIGN) { return (1. / (1 - beta*beta * costh_lab*costh_lab)) *\
131  (E_rest * beta * costh_lab / gamma \
132  + SIGN * sqrt(-M*M + (E_rest*E_rest)/(gamma*gamma) + M*M * beta*beta * costh_lab*costh_lab)); };
133 
134  double lab_frame_p_plus = lab_frame_p(1);
135  double lab_frame_p_minus = lab_frame_p(-1);
136 
137  bool plus_valid = !isnan(lab_frame_p_plus) && lab_frame_p_plus > 0.;
138  bool minus_valid = !isnan(lab_frame_p_minus) && lab_frame_p_minus > 0.;
139 
140  // Now compute the weight
141  auto wgtf = [costh_lab, rest_frame_p, E_rest, M, beta, gamma](double p_lab, int SIGN) {
142  double E_lab = sqrt(p_lab*p_lab + M*M);
143 
144  double dp_labdp_rest = (rest_frame_p / (E_rest * gamma)) * (beta*costh_lab +\
145  SIGN*E_rest / (gamma * sqrt(E_rest*E_rest / (gamma*gamma) - M*M + M*M * beta*beta * costh_lab*costh_lab))) /\
146  (1 - beta*beta * costh_lab*costh_lab);
147  return (E_rest / E_lab) * (p_lab*p_lab / (rest_frame_p*rest_frame_p)) * abs(dp_labdp_rest);
148  };
149 
150  double plus_wgt = plus_valid ? wgtf(lab_frame_p_plus, 1) : 0.;
151  double minus_wgt = minus_valid ? wgtf(lab_frame_p_minus, -1) : 0.;
152 
153  //double threshold = (plus_valid || minus_valid) ? minus_valid / (plus_valid + minus_valid) : -1;
154  double threshold = 0.5;
155 
156  // Select the solution to use
157  bool select_plus = plus_valid && (!minus_valid || (rand >= threshold));
158  bool select_minus = minus_valid && (!plus_valid || (rand < threshold));
159 
160  // Prints out stuff
161  std::cout << "COSTH LAB: " << costh_lab << std::endl;
162  std::cout << "PREST: " << rest_frame_p << std::endl;
163  std::cout << "MASS: " << M << std::endl;
164  std::cout << "BETA: " << beta << std::endl;
165  std::cout << "GAMMA: " << gamma << std::endl;
166 
167  std::cout << "COSTH REST PLUS: " << costh_rest_plus << std::endl;
168  std::cout << "COSTH REST MINUS: " << costh_rest_minus << std::endl;
169  std::cout << "PLAB PLUS: " << lab_frame_p_plus << std::endl;
170  std::cout << "PLAB MINUS: " << lab_frame_p_minus << std::endl;
171  std::cout << "WGT PLUS: " << plus_wgt << std::endl;
172  std::cout << "WGT MINUS: " << minus_wgt << std::endl;
173 
174  if (select_plus) std::cout << "SELECTED PLUS" << std::endl;
175  else if (select_minus) std::cout << "SELECTED MINUS" << std::endl;
176  else std::cout << "SELECTED NONE" << std::endl;
177 
178  if (select_plus) {
179  costh_rest_out = costh_rest_plus;
180  lab_frame_p_out = lab_frame_p_plus;
181  wgt = plus_wgt * (plus_valid + minus_valid);
182  }
183  else if (select_minus) {
184  costh_rest_out = costh_rest_minus;
185  lab_frame_p_out = lab_frame_p_minus;
186  wgt = minus_wgt * (plus_valid + minus_valid);
187  }
188  else { // No solution!
189  return -1;
190  }
191 
192  return 0;
193 }
T abs(T value)
tuple dir
Definition: dropbox.py:28
BEGIN_PROLOG could also be cout
double evgen::ldm::ElectronPartialWidth ( double  higs_mass,
double  mixing 
)

Definition at line 121 of file HiggsMakeDecay_tool.cc.

121  {
122  return LeptonPartialWidth(Constants::Instance().elec_mass, higs_mass, mixing);
123 }
double LeptonPartialWidth(double lep_mass, double higs_mass, double mixing)
double evgen::ldm::flat_to_exp_rand ( double  x,
double  mean,
double  a,
double  b 
)

Definition at line 91 of file HiggsMakeDecay_tool.cc.

91  {
92  double A = (1. - exp(-(b-a)/mean));
93  return - mean * log(1 - x * A) + a;
94 }
process_name opflash particleana ie x
process_name gaushit a
double mean(const std::vector< short > &wf, size_t start, size_t nsample)
Definition: UtilFunc.cxx:13
float A
Definition: dedx.py:137
double evgen::ldm::forcedecay_weight ( double  mean,
double  a,
double  b 
)

Definition at line 97 of file HiggsMakeDecay_tool.cc.

97  {
98  return exp(-a/mean) - exp(-b/mean);
99 }
process_name gaushit a
double mean(const std::vector< short > &wf, size_t start, size_t nsample)
Definition: UtilFunc.cxx:13
double evgen::ldm::forwardPrtlEnergy ( double  parentM,
double  secM,
double  prtlM,
double  parentE 
)

Definition at line 212 of file Constants.cpp.

212  {
213  // look up the rest frame prtl momentum
214  double rest_prtlP = twobody_momentum(parentM, secM, prtlM);
215 
216  // set the dir and boost in the same direction
217  double gamma = parentE / parentM;
218  double beta = sqrt(1. - 1. / (gamma*gamma));
219  TVector3 dir(0., 0., 1.);
220  TVector3 boost(0., 0., beta);
221 
222  // make sure we get the positive (more energetic) solution
223  double rand = 1.;
224 
225  double lab_frame_p_out, costh_rest_out, wgt;
226  int err = calcPrtlRayWgt(rest_prtlP, prtlM, dir, boost, rand,
227  lab_frame_p_out, costh_rest_out, wgt);
228  assert(!err);
229  (void) err;
230 
231  double lab_prtlE = sqrt(lab_frame_p_out*lab_frame_p_out + prtlM*prtlM);
232  return lab_prtlE;
233 }
EResult err(const char *call)
double twobody_momentum(double parent_mass, double childA_mass, double childB_mass)
Definition: Constants.cpp:73
j template void())
Definition: json.hpp:3108
tuple dir
Definition: dropbox.py:28
int calcPrtlRayWgt(double rest_frame_p, double M, TVector3 dir, TVector3 boost, double rand, double &lab_frame_p_out, double &costh_rest_out, double &wgt)
Definition: Constants.cpp:102
double evgen::ldm::higgs_momentum ( double  kaon_mass,
double  pion_mass,
double  higs_mass 
)

Definition at line 101 of file HiggsMakeDecay_tool.cc.

101  {
102  return sqrt(kaon_mass * kaon_mass * kaon_mass * kaon_mass
103  -2 * kaon_mass * kaon_mass * pion_mass * pion_mass
104  -2 * kaon_mass * kaon_mass * higs_mass * higs_mass
105  + pion_mass * pion_mass * pion_mass * pion_mass
106  + higs_mass * higs_mass * higs_mass * higs_mass
107  -2 * pion_mass * pion_mass * higs_mass * higs_mass) / ( 2 * kaon_mass );
108 }
double evgen::ldm::hnl_momentum ( double  kaon_mass,
double  lep_mass,
double  hnl_mass 
)

Definition at line 88 of file Kaon2HNLFlux_tool.cc.

88  {
89  if (kaon_mass - lep_mass < hnl_mass) return -1.;
90 
91  return sqrt(kaon_mass * kaon_mass * kaon_mass * kaon_mass
92  -2 * kaon_mass * kaon_mass * lep_mass * lep_mass
93  -2 * kaon_mass * kaon_mass * hnl_mass * hnl_mass
94  + lep_mass * lep_mass * lep_mass * lep_mass
95  + hnl_mass * hnl_mass * hnl_mass * hnl_mass
96  -2 * lep_mass * lep_mass * hnl_mass * hnl_mass) / ( 2 * kaon_mass );
97 }
double evgen::ldm::HNLLepPiDalitzMax ( double  mK,
double  mA,
double  mN,
double  mP,
double  mB 
)

Definition at line 55 of file HNLDecayDalitz.cpp.

55  {
56  double pNA = evgen::ldm::twobody_momentum(mK, mA, mN);
57  double pPB = evgen::ldm::twobody_momentum(mN, mP, mB);
58 
59  TLorentzVector K(TVector3(0, 0, pNA * mK / mN), mK * sqrt(1 + pNA * pNA / (mN*mN)));
60  TLorentzVector LA(TVector3(0, 0, pNA * mK / mN), sqrt(mA*mA + pNA * pNA * mK * mK / (mN * mN)));
61  TLorentzVector N(TVector3(0, 0, 0), mN);
62  // In the LNC case, if mN > mA, then the direction of lB wants to be aligned with lA (+z)
63  // (because of angular momentum conservation)
64  // Flipping the mass ordering or LNC/LNV flips the sign
65  double sign = mN > mA ? 1 : -1;
66  TLorentzVector LB(TVector3(0, 0, sign * pPB), sqrt(mB * mB + pPB * pPB));
67  TLorentzVector PI(TVector3(0, 0, -sign * pPB), sqrt(mP * mP + pPB * pPB));
68 
69  return HNLLepPiLNCDalitz(K, LA, N, PI, LB);
70 }
double twobody_momentum(double parent_mass, double childA_mass, double childB_mass)
Definition: Constants.cpp:73
constexpr double PI
int sign(double val)
Definition: UtilFunc.cxx:104
process_name largeant stream1 can override from command line with o or output physics producers generator N
double HNLLepPiLNCDalitz(TLorentzVector K, TLorentzVector LA, TLorentzVector N, TLorentzVector PI, TLorentzVector LB)
double evgen::ldm::HNLLepPiLNCDalitz ( TLorentzVector  K,
TLorentzVector  LA,
TLorentzVector  N,
TLorentzVector  PI,
TLorentzVector  LB 
)

Definition at line 23 of file HNLDecayDalitz.cpp.

23  {
24  return 8 * K.Dot(LA) * PI.Dot(LB) * K.Dot(N) * PI.Dot(N) \
25  -4 * K.Dot(LA) * PI.Dot(LB) * K.Dot(PI) * N.Dot(N) \
26  -4 * K.Dot(K) * PI.Dot(LB) * LA.Dot(N) * PI.Dot(N) \
27  +2 * K.Dot(K) * PI.Dot(LB) * LA.Dot(PI) * N.Dot(N) \
28  -4 * K.Dot(LA) * PI.Dot(PI) * K.Dot(N) * LB.Dot(N) \
29  +2 * K.Dot(LA) * PI.Dot(PI) * K.Dot(LB) * N.Dot(N) \
30  +2 * K.Dot(K) * PI.Dot(PI) * LA.Dot(N) * LB.Dot(N) \
31  -1 * K.Dot(K) * PI.Dot(PI) * LA.Dot(LB) * N.Dot(N);
32 }
constexpr double PI
process_name largeant stream1 can override from command line with o or output physics producers generator N
double evgen::ldm::HNLLepPiLNVDalitz ( TLorentzVector  K,
TLorentzVector  LA,
TLorentzVector  N,
TLorentzVector  PI,
TLorentzVector  LB 
)

Definition at line 34 of file HNLDecayDalitz.cpp.

34  {
35  return N.Dot(N) * (
36  4 * K.Dot(LA) * K.Dot(PI) * LB.Dot(PI) \
37  -2 * K.Dot(K) * LA.Dot(PI) * LB.Dot(PI) \
38  -2 * K.Dot(LA) * K.Dot(LB) * PI.Dot(PI) \
39  +1 * K.Dot(K) * LA.Dot(LB) * PI.Dot(PI)
40  );
41 }
constexpr double PI
process_name largeant stream1 can override from command line with o or output physics producers generator N
double evgen::ldm::HNLNuDiLepLNVDalitz ( TLorentzVector  K,
TLorentzVector  LA,
TLorentzVector  LB,
TLorentzVector  LBD,
TLorentzVector  NU 
)

Definition at line 6 of file HNLDecayDalitz.cpp.

6  {
7 
8  auto pcross = [K, LA](TLorentzVector F) {
9  return 2. * LA.Dot(K) * F.Dot(K) - LA.Dot(F) * K.Dot(K);
10  };
11 
12  double A1 = 0.;
13  double A2 = 0.;
14  double A3 = 0.;
15 
16  return A1*2*NU.Dot(LB)*pcross(LBD) + A2*2*NU.Dot(LBD)*pcross(LB) + A3*LB.Dot(LB)*pcross(NU);
17 }
double evgen::ldm::I1_integrand ( double  s,
void *  param 
)

Definition at line 203 of file HNLMakeDecay_tool.cc.

203  {
204  double *xyz = (double *)param;
205  double x = xyz[0];
206  double y = xyz[1];
207  double z = xyz[2];
208 
209  return 12.*(s - x*x - y*y)*(1 + z*z - s)*sqrt(lambda(s,x*x,y*y))*sqrt(lambda(1.,s,z*z))/s;
210 }
process_name opflash particleana ie ie ie z
process_name opflash particleana ie x
double lambda(double a, double b, double c)
process_name opflash particleana ie ie y
then echo File list $list not found else cat $list while read file do echo $file sed s
Definition: file_to_url.sh:60
double evgen::ldm::I2_integrand ( double  s,
void *  param 
)

Definition at line 212 of file HNLMakeDecay_tool.cc.

212  {
213  double *xyz = (double *)param;
214  double x = xyz[0];
215  double y = xyz[1];
216  double z = xyz[2];
217 
218  return 24*y*z*(1. + x*x - s)*sqrt(lambda(s,y*y,z*z))*sqrt(lambda(1.,s,z*z))/s;
219 }
process_name opflash particleana ie ie ie z
process_name opflash particleana ie x
double lambda(double a, double b, double c)
process_name opflash particleana ie ie y
then echo File list $list not found else cat $list while read file do echo $file sed s
Definition: file_to_url.sh:60
double evgen::ldm::KaonLongBranchingRatio ( double  higs_mass,
double  mixing 
)

Definition at line 128 of file Kaon2HiggsFlux_tool.cc.

128  {
129  double klong_mass = Constants::Instance().klong_mass;
130  double pizero_mass = Constants::Instance().pizero_mass;
131  double tquark_mass = Constants::Instance().tquark_mass;
132  double higgs_vev = Constants::Instance().higgs_vev;
133  double rel_VtsVtd_squared = Constants::Instance().rel_VtsVtd_squared;
134 
135  double M_KL = (1. / 2.) * (3. / (16. * M_PI * M_PI * higgs_vev * higgs_vev * higgs_vev)) * (klong_mass * klong_mass) * (tquark_mass * tquark_mass) * mixing;
136  double M_KL2 = M_KL * M_KL * rel_VtsVtd_squared;
137 
138  double klong_width = (2 * higgs_momentum(klong_mass, pizero_mass, higs_mass) / klong_mass) * M_KL2 / (16 * M_PI * klong_mass);
139 
140  double klong_2higgs_lifetime = Constants::Instance().hbar / klong_width;
141 
142  return Constants::Instance().klong_lifetime / klong_2higgs_lifetime;
143 
144 }
double higgs_momentum(double kaon_mass, double pion_mass, double higs_mass)
double evgen::ldm::KaonPlusBranchingRatio ( double  higs_mass,
double  mixing 
)

Definition at line 104 of file Kaon2HiggsFlux_tool.cc.

104  {
105  double kplus_mass = Constants::Instance().kplus_mass;
106  double piplus_mass = Constants::Instance().piplus_mass;
107  double tquark_mass = Constants::Instance().tquark_mass;
108  double higgs_vev = Constants::Instance().higgs_vev;
109  double abs_VtsVtd_squared = Constants::Instance().abs_VtsVtd_squared;
110 
111  // Kplus width
112  //
113  // matrix element for kplus
114  double M_KP = (1. / 2.) * ( 3. / (16. * M_PI * M_PI * higgs_vev * higgs_vev * higgs_vev)) * (kplus_mass * kplus_mass) * (tquark_mass * tquark_mass) * mixing;
115  double M_KP2 = M_KP * M_KP * abs_VtsVtd_squared;
116 
117  double kplus_width = (2 * higgs_momentum(kplus_mass, piplus_mass, higs_mass)/kplus_mass) * M_KP2 / (16 * M_PI * kplus_mass);
118 
119  // convert to partial lifetime
120  double kplus_2higgs_lifetime = Constants::Instance().hbar / kplus_width;
121 
122  // and branching ratio
123  //
124  // (this higgs decay would make a negligible contribution to the overall lifetime)
125  return Constants::Instance().kplus_lifetime / kplus_2higgs_lifetime;
126 }
double higgs_momentum(double kaon_mass, double pion_mass, double higs_mass)
double evgen::ldm::lambda ( double  a,
double  b,
double  c 
)

Definition at line 134 of file HNLMakeDecay_tool.cc.

134  {
135  return a*a + b*b + c*c - 2*a*b - 2*b*c - 2*c*a;
136 }
process_name gaushit a
double evgen::ldm::LeptonPartialWidth ( double  lep_mass,
double  higs_mass,
double  mixing 
)

Definition at line 112 of file HiggsMakeDecay_tool.cc.

112  {
113  if (lep_mass * 2 >= higs_mass) return 0.;
114 
115  double higgs_vev = Constants::Instance().higgs_vev;
116 
117  double width = (mixing * mixing * lep_mass * lep_mass * higs_mass / (8 * M_PI * higgs_vev * higgs_vev)) * pow(1 - 4 * lep_mass * lep_mass / (higs_mass * higs_mass), 3. / 2.);
118  return width;
119 }
double evgen::ldm::minKinematicCosTheta ( double  parentM,
double  secM,
double  prtlM,
double  parentE 
)

Definition at line 195 of file Constants.cpp.

195  {
196  // guard against parent-at-rest
197  if (parentE - parentM < 1e-4) return 0.;
198 
199  // look up the rest frame prtl momentum
200  double rest_prtlP = twobody_momentum(parentM, secM, prtlM);
201 
202  // set the dir and boost in the same direction
203  double gamma = parentE / parentM;
204  double beta = sqrt(1. - 1. / (gamma*gamma));
205 
206  double rest_prtlE = sqrt(prtlM*prtlM + rest_prtlP*rest_prtlP);
207 
208  return sqrt(prtlM*prtlM - rest_prtlE * rest_prtlE / (gamma*gamma)) / (prtlM*beta);
209 }
double twobody_momentum(double parent_mass, double childA_mass, double childB_mass)
Definition: Constants.cpp:73
do i e
double evgen::ldm::MuonPartialWidth ( double  higs_mass,
double  mixing 
)

Definition at line 125 of file HiggsMakeDecay_tool.cc.

125  {
126  return LeptonPartialWidth(Constants::Instance().muon_mass, higs_mass, mixing);
127 }
double LeptonPartialWidth(double lep_mass, double higs_mass, double mixing)
double evgen::ldm::PionPartialWidth ( double  pion_mass,
double  higs_mass,
double  mixing 
)

Definition at line 129 of file HiggsMakeDecay_tool.cc.

129  {
130  if (pion_mass * 2 >= higs_mass) return 0.;
131 
132  double higgs_vev = Constants::Instance().higgs_vev;
133 
134  double form_factor = (2. / 9.) * higs_mass * higs_mass + (11. / 9.) * pion_mass * pion_mass;
135 
136  double width = (mixing * mixing * 3 * form_factor * form_factor / (32 * M_PI * higgs_vev * higgs_vev * higs_mass)) * pow(1- 4. * pion_mass * pion_mass / (higs_mass * higs_mass), 1. / 2.);
137 
138  return width;
139 }
int evgen::ldm::PionPdg ( int  kaon_pdg)

Definition at line 167 of file Kaon2HiggsFlux_tool.cc.

167  {
168  switch (kaon_pdg) {
169  case 321:
170  return 211;
171  case -321:
172  return -211;
173  case 130:
174  return 111;
175  default:
176  return -1;
177  }
178 }
double evgen::ldm::PiPlusPartialWidth ( double  higs_mass,
double  mixing 
)

Definition at line 141 of file HiggsMakeDecay_tool.cc.

141  {
142  return PionPartialWidth(Constants::Instance().piplus_mass, higs_mass, mixing);
143 }
double PionPartialWidth(double pion_mass, double higs_mass, double mixing)
double evgen::ldm::PiZeroPartialWidth ( double  higs_mass,
double  mixing 
)

Definition at line 145 of file HiggsMakeDecay_tool.cc.

145  {
146  return PionPartialWidth(Constants::Instance().pizero_mass, higs_mass, mixing);
147 }
double PionPartialWidth(double pion_mass, double higs_mass, double mixing)
double evgen::ldm::QEstimator ( unsigned  nsuccess,
unsigned  nfail,
unsigned  r 
)

Definition at line 96 of file MixedWeightRayTraceBox_tool.cc.

96  {
97  if (nsuccess < r) return 0.;
98 
99  return ((double)(r - 1) / (r+nfail-1));
100 }
esac echo uname r
double evgen::ldm::secPDG2Mass ( int  pdg)

Definition at line 236 of file Constants.cpp.

236  {
237  switch (pdg) {
238  case 11:
239  case -11:
240  return Constants::Instance().elec_mass;
241  break;
242  case 13:
243  case -13:
244  return Constants::Instance().muon_mass;
245  break;
246  case 211:
247  case -211:
248  return Constants::Instance().piplus_mass;
249  break;
250  default:
251  std::cerr << "BAD secondary pdg: " << pdg << std::endl;
252  std::cout << "BAD secondary pdg: " << pdg << std::endl;
253  return 0.;
254  }
255 }
var pdg
Definition: selectors.fcl:14
BEGIN_PROLOG could also be cerr
BEGIN_PROLOG could also be cout
double evgen::ldm::SMKaonBR ( int  kaon_pdg)

Definition at line 146 of file Kaon2HiggsFlux_tool.cc.

146  {
147  // The Kaons in Dk2nu file only include those that decay to neutrinos.
148  //
149  // We want all kaons -- in order to account for this, we divide by the
150  // branching-ratio of kaons to neutrinos
151  //
152  // Taken from:
153  // /cvmfs/minerva.opensciencegrid.org/minerva/beamsim/x86_64/geant4/source/particles/hadrons/mesons/src/G4KaonPlus.cc
154  // /cvmfs/minerva.opensciencegrid.org/minerva/beamsim/x86_64/geant4/source/particles/hadrons/mesons/src/G4KaonZerLong.cc
155  switch (kaon_pdg) {
156  case 321:
157  return 0.6339 /* 5 */ + 0.0559 /* 6 */ + 0.0330 /* 7 */;
158  case -321:
159  return 0.6339 /* 8 */ + 0.0559 /* 9 */ + 0.0330 /* 10 */;
160  case 130:
161  return 0.2020 /* 1 */ + 0.2020 /* 2 */ + 0.1348 /* 3 */ + 0.1348 /* 4 */;
162  default:
163  return -1;
164  }
165 }
double evgen::ldm::twobody_momentum ( double  parent_mass,
double  childA_mass,
double  childB_mass 
)

Definition at line 73 of file Constants.cpp.

73  {
74  if (parent_mass < childA_mass + childB_mass) return -1.;
75 
76  return sqrt(parent_mass * parent_mass * parent_mass * parent_mass
77  -2 * parent_mass * parent_mass * childA_mass * childA_mass
78  -2 * parent_mass * parent_mass * childB_mass * childB_mass
79  + childA_mass * childA_mass * childA_mass * childA_mass
80  + childB_mass * childB_mass * childB_mass * childB_mass
81  -2 * childA_mass * childA_mass * childB_mass * childB_mass) / ( 2 * parent_mass );
82 
83 }