9 std::cout <<
"[SimpEquations]::Parameters " << std::endl;
12 std::cout <<
"ratio of dark pion mass to dark pion decay constant: " <<
ratio_mPi_to_fPi_ << std::endl;
13 std::cout <<
"Dark Vector decay lepton mass: " <<
m_l_ << std::endl;
14 std::cout <<
"Alpha Dark: " <<
alpha_dark_ << std::endl;
22 std::cout <<
"[SimpEquations]::Parameters Loaded From Json Config " << std::endl;
25 std::cout <<
"ratio of dark pion mass to dark pion decay constant: " <<
ratio_mPi_to_fPi_ << std::endl;
26 std::cout <<
"Dark Vector decay lepton mass: " <<
m_l_ << std::endl;
27 std::cout <<
"Alpha Dark: " <<
alpha_dark_ << std::endl;
34 std::ifstream i_file(paramsConfigFile);
37 std::cout << el.key() <<
" : " << el.value() <<
"\n";
47 m_l_ = param.value().at(
"lepton_mass");
53 std::cout <<
"[SimpEquations]::ERROR::REQUIRED PARAMETER CONFIGURATIONS NOT FOUND IN "
54 << paramsConfigFile << std::endl;
59 double coeff = (2.0*alpha_dark/3.0) * m_Ap;
60 double pow1 = std::pow((1-(4*m_pi*m_pi/(m_Ap*m_Ap))),3/2.);
61 double pow2 = std::pow(((m_V*m_V)/((m_Ap*m_Ap)-(m_V*m_V))),2);
62 return coeff*pow1*pow2;
66 double alpha_dark,
double f_pi){
67 double x = m_pi / m_Ap;
68 double y = m_V / m_Ap;
70 double coeff = alpha_dark*
Tv/(192.*std::pow(M_PI,4));
71 return coeff * std::pow((m_Ap/m_pi),2) * std::pow(m_V/m_pi,2) * std::pow((m_pi/f_pi),4) * m_Ap*std::pow(
Beta(x,y),3./2.);
75 double alpha_dark,
double f_pi){
76 double x = m_pi / m_Ap;
77 double y = m_V / m_Ap;
79 double coeff = alpha_dark*
Tv/(192.*std::pow(M_PI,4));
80 return coeff * std::pow((m_Ap/m_pi),2) * std::pow(m_V/m_pi,2) * std::pow((m_pi/f_pi),4) * m_Ap*std::pow(
Beta(x,y),3./2.);
83 double alpha_dark,
double f_pi){
84 double x = m_pi / m_Ap;
85 double y = m_V / m_Ap;
86 double Tv = 18.0-((3.0/2.0) +(3.0/4.0)) ;
87 double coeff = alpha_dark*
Tv/(192.*std::pow(M_PI,4));
88 return coeff * std::pow((m_Ap/m_pi),2) * std::pow(m_V/m_pi,2) * std::pow((m_pi/f_pi),4) * m_Ap*std::pow(
Beta(x,y),3./2.);
93 double total_rate =
rate_Vrho_pi(m_Ap, m_pi, m_V, alpha_dark, f_pi) +
96 rate_2pi(m_Ap, m_pi, m_V, alpha_dark);
98 total_rate +
rate_2V(m_Ap, m_V, alpha_dark);
99 return rate_2pi(m_Ap, m_pi, m_V, alpha_dark)/total_rate;
104 double total_rate =
rate_Vrho_pi(m_Ap, m_pi, m_V, alpha_dark, f_pi) +
107 rate_2pi(m_Ap, m_pi, m_V, alpha_dark);
109 total_rate +
rate_2V(m_Ap, m_V, alpha_dark);
110 return rate_Vrho_pi(m_Ap, m_pi, m_V, alpha_dark, f_pi)/total_rate;
115 double total_rate =
rate_Vrho_pi(m_Ap, m_pi, m_V, alpha_dark, f_pi) +
118 rate_2pi(m_Ap, m_pi, m_V, alpha_dark);
120 total_rate +
rate_2V(m_Ap, m_V, alpha_dark);
121 return rate_Vphi_pi(m_Ap, m_pi, m_V, alpha_dark, f_pi)/total_rate;
125 double alpha_dark,
double f_pi){
126 double total_rate =
rate_Vrho_pi(m_Ap, m_pi, m_V, alpha_dark, f_pi) +
129 rate_2pi(m_Ap, m_pi, m_V, alpha_dark);
131 total_rate +
rate_2V(m_Ap, m_V, alpha_dark);
135double SimpEquations::br_2V(
double m_Ap,
double m_pi,
double m_V,
double alpha_dark,
double f_pi,
double rho,
double phi){
136 double total_rate =
rate_Vrho_pi(m_Ap, m_pi, m_V, alpha_dark, f_pi) +
139 rate_2pi(m_Ap, m_pi, m_V, alpha_dark);
141 total_rate +
rate_2V(m_Ap, m_V, alpha_dark);
144 return rate_2V(m_Ap,m_V,alpha_dark)/total_rate;
135double SimpEquations::br_2V(
double m_Ap,
double m_pi,
double m_V,
double alpha_dark,
double f_pi,
double rho,
double phi) {
…}
148 if(rho)
return 3./4.;
149 else if(phi)
return 3./2.;
154 return (1+std::pow(y,2)-std::pow(x,2)-2*y)*(1+std::pow(y,2)-std::pow(x,2)+2*y);
159 return alpha_dark/6. * m_Ap *
f(r);
163 double num = 1 + 16*std::pow(r,2) - 68*std::pow(r,4) - 48*std::pow(r,6);
164 double den = std::pow(1-std::pow(r,2),2);
165 return num/den * std::pow((1-4*std::pow(r,2)),0.5);
168double SimpEquations::rate_2l(
double m_Ap,
double m_pi,
double m_V,
double eps,
double alpha_dark,
double f_pi,
double m_l,
bool rho){
169 double alpha = 1./137.0;
170 double coeff = 16*M_PI*alpha_dark*alpha*std::pow(eps,2)*std::pow(f_pi,2)/(3*std::pow(m_V,2));
171 double term1 = std::pow((std::pow(m_V,2)/(std::pow(m_Ap,2) - std::pow(m_V,2))),2);
172 double term2 = std::pow((1-(4*std::pow(m_l,2)/std::pow(m_V,2))),0.5);
173 double term3 = 1+(2*std::pow(m_l,2)/std::pow(m_V,2));
175 if(rho) constant = 2;
176 return coeff * term1 * term2 * term3 * m_V * constant;
168double SimpEquations::rate_2l(
double m_Ap,
double m_pi,
double m_V,
double eps,
double alpha_dark,
double f_pi,
double m_l,
bool rho) {
…}
179double SimpEquations::getCtau(
double m_Ap,
double m_pi,
double m_V,
double eps,
double alpha_dark,
double f_pi,
double m_l,
bool rho){
181 double hbar = 6.58e-22;
182 double rate =
rate_2l(m_Ap,m_pi,m_V,eps,alpha_dark,f_pi,m_l,rho);
183 double tau = hbar/rate;
179double SimpEquations::getCtau(
double m_Ap,
double m_pi,
double m_V,
double eps,
double alpha_dark,
double f_pi,
double m_l,
bool rho) {
…}
189 double gamma = E_V/m_V;
195 double E_V, TEfficiency* effCalc_h,
double dNdm,
double radFrac,
double radAcc,
double target_pos,
double zcut){
204 double gcTau = ctau *
gamma(m_V/1000.0, E_V);
205 std::cout <<
"gcTau: " << gcTau << std::endl;
206 std::cout <<
"radFrac: " << radFrac << std::endl;
207 std::cout <<
"radAcc: " << radAcc << std::endl;
208 std::cout <<
"dNdm:" << dNdm << std::endl;
212 for(
int zbin = 0; zbin < (effCalc_h->GetTotalHistogram()->GetNbinsX())+1; zbin++){
213 double zz = effCalc_h->GetTotalHistogram()->GetBinLowEdge(zbin);
214 if(zz < zcut)
continue;
215 effVtx += (TMath::Exp((target_pos-zz)/gcTau)/gcTau)*
216 (effCalc_h->GetEfficiency(zbin) - effCalc_h->GetEfficiencyErrorLow(zbin))*
217 effCalc_h->GetTotalHistogram()->GetBinWidth(zbin);
221 double apProduction = (3.*137/2.)*3.14159*(m_Ap*eps*eps*radFrac*dNdm)/radAcc;
230 std::cout <<
"Branching ratio is " << br_VPi << std::endl;
232 double br_V_ee = 1.0;
235 double expSignal = apProduction * effVtx * br_VPi * br_V_ee;
241 double E_V, TEfficiency* effCalc_h,
double target_pos,
double zcut){
250 double gcTau = ctau *
gamma(m_V/1000.0, E_V);
254 for(
int zbin = 0; zbin < (effCalc_h->GetTotalHistogram()->GetNbinsX())+1; zbin++){
255 double zz = effCalc_h->GetTotalHistogram()->GetBinLowEdge(zbin);
256 if(zz < zcut)
continue;
257 effVtx += (TMath::Exp((target_pos-zz)/gcTau)/gcTau)*
258 (effCalc_h->GetEfficiency(zbin) - effCalc_h->GetEfficiencyErrorLow(zbin))*
259 effCalc_h->GetTotalHistogram()->GetBinWidth(zbin);
275 double br_V_ee = 1.0;
278 double expSignal = apProduction * effVtx * br_VPi * br_V_ee;
284 double f_pi,
double m_l,
bool rho,
double E_V, TEfficiency* effCalc_h,
double target_pos,
double zcut){
287 double ctau =
getCtau(m_Ap,m_pi, m_V,eps,alpha_dark,f_pi,m_l,rho);
288 double gcTau = ctau *
gamma(m_V/1000.0, E_V);
292 for(
int zbin = 0; zbin < (effCalc_h->GetTotalHistogram()->GetNbinsX())+1; zbin++){
293 double zz = effCalc_h->GetTotalHistogram()->GetBinLowEdge(zbin);
294 if(zz < zcut)
continue;
295 effVtx += (TMath::Exp((target_pos-zz)/gcTau)/gcTau)*
296 (effCalc_h->GetEfficiency(zbin) - effCalc_h->GetEfficiencyErrorLow(zbin))*
297 effCalc_h->GetTotalHistogram()->GetBinWidth(zbin);
313 double br_V_ee = 1.0;
316 double expSignal = apProduction * effVtx * br_VPi * br_V_ee;
322 double radFrac = 0.0;
324 radFrac = ( -1.92497e-01 + 1.47144e-02*m_Ap + -2.91966e-04*std::pow(m_Ap,2) + 2.65603e-06*std::pow(m_Ap,3) + -1.12471e-8*std::pow(m_Ap,4) + 1.74765e-11*std::pow(m_Ap,5) + 2.235718e-15*std::pow(m_Ap,6));
333 radAcc = ( -7.93151e-01 + 1.04324e-01*m_Ap + -5.55225e-03*std::pow(m_Ap,2) + 1.55480e-04*std::pow(m_Ap,3) + -2.53281e-06*std::pow(m_Ap,4) + 2.54558e-08*std::pow(m_Ap,5) + -1.60877e-10*std::pow(m_Ap,6) + 6.24627e-13*std::pow(m_Ap,7) + -1.36375e-15*std::pow(m_Ap,8) + 1.28312e-18*std::pow(m_Ap,9) );
340 double massRes = 0.0;
342 massRes = 9.73217e-01 + 3.63659e-02*m_V + -7.32046e-05*m_V*m_V;
358 std::cout <<
"[SimpEquations]::WARNING! USING TEMPORARY BACKGROUND MODEL" << std::endl;
362 dNdm = 204656 + -2.08037e7*std::pow(m_Ap,1) + 8.17706e8*std::pow(m_Ap,2) + -1.62741e10*std::pow(m_Ap,3)
363 + 1.85484e11*std::pow(m_Ap,4) + -1.26987e12*std::pow(m_Ap,5) + 5.17974e12*std::pow(m_Ap,6)
364 + -1.16277e13*std::pow(m_Ap,7) + 1.10758e13*std::pow(m_Ap, 8);
367 std::cout <<
"LOOK: Background Rate taken from Control Region InvMass is " << dNdm << std::endl;
double Tv(bool rho, bool phi)
double rate_Vrho_pi(double m_Ap, double m_pi, double m_V, double alpha_dark, double f_pi)
double expectedSignalCalculation(double m_V, double eps, bool rho, double E_V, TEfficiency *effCalc_h, double target_pos, double zcut)
double getCtau(double m_Ap, double m_pi, double m_V, double eps, double alpha_dark, double f_pi, double m_l, bool rho)
double rate_Vphi_pi(double m_Ap, double m_pi, double m_V, double alpha_dark, double f_pi)
double rate_2V(double m_Ap, double m_V, double alpha_dark)
double controlRegionBackgroundRate(double m_Ap)
double rate_Vcharged_pi(double m_Ap, double m_pi, double m_V, double alpha_dark, double f_pi)
double radiativeFraction(double m_Ap)
void loadParametersConfig(const std::string paramsConfigFile)
double Beta(double x, double y)
double br_Vrho_pi(double m_Ap, double m_pi, double m_V, double alpha_dark, double f_pi)
int year_
year (used to specify polynomial choices)
double rate_2pi(double m_Ap, double m_pi, double m_V, double alpha_dark)
double radiativeAcceptance(double m_Ap)
double br_Vphi_pi(double m_Ap, double m_pi, double m_V, double alpha_dark, double f_pi)
double m_l_
default lepton mass (ele/pos only)
double gamma(double m_V, double E_V)
double br_Vcharged_pi(double m_Ap, double m_pi, double m_V, double alpha_dark, double f_pi)
double br_2V(double m_Ap, double m_pi, double m_V, double alpha_dark, double f_pi, double rho, double phi)
double mass_ratio_Ap_to_Pid_
default Ap to dark pion mass ratio
double alpha_dark_
default A' to dark meson coupling
double rate_2l(double m_Ap, double m_pi, double m_V, double eps, double alpha_dark, double f_pi, double m_l, bool rho)
double br_2pi(double m_Ap, double m_pi, double m_V, double alpha_dark, double f_pi)
double ratio_mPi_to_fPi_
defualt dark pion mass to decay constant ratio
json params_config_
read in simp params
double mass_ratio_Ap_to_Vd_
default Ap to dark vector mass ratio
double massResolution(double m_V)