gambit is hosted by Hepforge, IPPP Durham
GAMBIT  v1.5.0-252-gf9a3f78
a Global And Modular Bsm Inference Tool
ColliderBit_Higgs.cpp
Go to the documentation of this file.
1 // GAMBIT: Global and Modular BSM Inference Tool
2 // *********************************************
27 
28 #include <cmath>
29 #include <string>
30 #include <iostream>
31 #include <fstream>
32 #include <memory>
33 #include <numeric>
34 #include <sstream>
35 #include <vector>
36 
39 
40 //#define COLLIDERBIT_DEBUG
41 
42 
43 namespace Gambit
44 {
45 
46  namespace ColliderBit
47  {
48 
50  void set_CS(hb_ModelParameters &result, const HiggsCouplingsTable& couplings, int n_neutral_higgses)
51  {
52  for(int i = 0; i < n_neutral_higgses; i++)
53  {
54  result.CS_bg_hjb_ratio[i] = couplings.C_bb2[i];
55  result.CS_bb_hj_ratio[i] = couplings.C_bb2[i];
56  result.CS_lep_bbhj_ratio[i] = couplings.C_bb2[i];
57 
58  result.CS_lep_tautauhj_ratio[i] = couplings.C_tautau2[i];
59 
60  result.CS_lep_hjZ_ratio[i] = couplings.C_ZZ2[i];
61  result.CS_gg_hjZ_ratio[i] = 0.;
62  result.CS_dd_hjZ_ratio[i] = couplings.C_ZZ2[i];
63  result.CS_uu_hjZ_ratio[i] = couplings.C_ZZ2[i];
64  result.CS_ss_hjZ_ratio[i] = couplings.C_ZZ2[i];
65  result.CS_cc_hjZ_ratio[i] = couplings.C_ZZ2[i];
66  result.CS_bb_hjZ_ratio[i] = couplings.C_ZZ2[i];
67 
68  result.CS_ud_hjWp_ratio[i] = couplings.C_WW2[i];
69  result.CS_cs_hjWp_ratio[i] = couplings.C_WW2[i];
70  result.CS_ud_hjWm_ratio[i] = couplings.C_WW2[i];
71  result.CS_cs_hjWm_ratio[i] = couplings.C_WW2[i];
72 
73  result.CS_tev_vbf_ratio[i] = couplings.C_WW2[i];
74  result.CS_lhc7_vbf_ratio[i] = couplings.C_WW2[i];
75  result.CS_lhc8_vbf_ratio[i] = couplings.C_WW2[i];
76 
77  result.CS_gg_hj_ratio[i] = couplings.C_gg2[i];
78 
79  result.CS_tev_tthj_ratio[i] = couplings.C_tt2[i];
80  result.CS_lhc7_tthj_ratio[i] = couplings.C_tt2[i];
81  result.CS_lhc8_tthj_ratio[i] = couplings.C_tt2[i];
82 
83  for(int j = 0; j < n_neutral_higgses; j++)
84  {
85  result.CS_lep_hjhi_ratio[i][j] = couplings.C_hiZ2[i][j];
86  }
87  }
88  // LEP H+ H- x-section ratio
89  result.CS_lep_HpjHmi_ratio[0] = 1.;
90  }
91 
93  void set_SMLikeHiggs_ModelParameters(const SubSpectrum& spec, const HiggsCouplingsTable& couplings, hb_ModelParameters &result)
94  {
95  // Retrieve the decays
96  const DecayTable::Entry& decays = couplings.get_neutral_decays(0);
97 
98  // Set the CP of the lone higgs
99  result.CP[0] = couplings.CP[0];
100 
101  // Set h mass and uncertainty
102  result.Mh[0] = spec.get(Par::Pole_Mass,25,0);
103  bool has_high_err = spec.has(Par::Pole_Mass_1srd_high, 25, 0);
104  bool has_low_err = spec.has(Par::Pole_Mass_1srd_low, 25, 0);
105  if (has_high_err and has_low_err)
106  {
107  double upper = spec.get(Par::Pole_Mass_1srd_high, 25, 0);
108  double lower = spec.get(Par::Pole_Mass_1srd_low, 25, 0);
109  result.deltaMh[0] = result.Mh[0] * std::max(upper,lower);
110  }
111  else
112  {
113  result.deltaMh[0] = 0.;
114  }
115 
116  // Set the total h width
117  result.hGammaTot[0] = decays.width_in_GeV;
118 
119  // Set the branching fractions
120  result.BR_hjss[0] = decays.BF("s", "sbar");
121  result.BR_hjcc[0] = decays.BF("c", "cbar");
122  result.BR_hjbb[0] = decays.BF("b", "bbar");
123  result.BR_hjmumu[0] = decays.BF("mu+", "mu-");
124  result.BR_hjtautau[0] = decays.BF("tau+", "tau-");
125  result.BR_hjWW[0] = decays.BF("W+", "W-");
126  result.BR_hjZZ[0] = decays.BF("Z0", "Z0");
127  result.BR_hjZga[0] = decays.BF("gamma", "Z0");
128  result.BR_hjgaga[0] = decays.BF("gamma", "gamma");
129  result.BR_hjgg[0] = decays.BF("g", "g");
130 
131  // Add the invisibles
132  result.BR_hjinvisible[0] = 0.;
133  for (auto it = couplings.invisibles.begin(); it != couplings.invisibles.end(); ++it)
134  {
135  result.BR_hjinvisible[0] += decays.BF(*it, *it);
136  }
137 
138  // Retrieve cross-section ratios from the HiggsCouplingsTable
139  set_CS(result, couplings, 1);
140 
141  // Zero all heavy neutral higgs masses, widths and effective couplings
142  for(int i = 1; i < 3; i++)
143  {
144  result.Mh[i] = 0.;
145  result.deltaMh[i] = 0.;
146  result.hGammaTot[i] = 0.;
147  result.CP[i] = 0.;
148  result.BR_hjss[i] = 0.;
149  result.BR_hjcc[i] = 0.;
150  result.BR_hjbb[i] = 0.;
151  result.BR_hjmumu[i] = 0.;
152  result.BR_hjtautau[i] = 0.;
153  result.BR_hjWW[i] = 0.;
154  result.BR_hjZZ[i] = 0.;
155  result.BR_hjZga[i] = 0.;
156  result.BR_hjgaga[i] = 0.;
157  result.BR_hjgg[i] = 0.;
158  result.BR_hjinvisible[i] = 0.;
159  result.CS_lep_hjZ_ratio[i] = 0.;
160  result.CS_lep_bbhj_ratio[i] = 0.;
161  result.CS_lep_tautauhj_ratio[i] = 0.;
162  result.CS_gg_hj_ratio[i] = 0.;
163  result.CS_bb_hj_ratio[i] = 0.;
164  result.CS_bg_hjb_ratio[i] = 0.;
165  result.CS_ud_hjWp_ratio[i] = 0.;
166  result.CS_cs_hjWp_ratio[i] = 0.;
167  result.CS_ud_hjWm_ratio[i] = 0.;
168  result.CS_cs_hjWm_ratio[i] = 0.;
169  result.CS_gg_hjZ_ratio[i] = 0.;
170  result.CS_dd_hjZ_ratio[i] = 0.;
171  result.CS_uu_hjZ_ratio[i] = 0.;
172  result.CS_ss_hjZ_ratio[i] = 0.;
173  result.CS_cc_hjZ_ratio[i] = 0.;
174  result.CS_bb_hjZ_ratio[i] = 0.;
175  result.CS_tev_vbf_ratio[i] = 0.;
176  result.CS_tev_tthj_ratio[i] = 0.;
177  result.CS_lhc7_vbf_ratio[i] = 0.;
178  result.CS_lhc7_tthj_ratio[i] = 0.;
179  result.CS_lhc8_vbf_ratio[i] = 0.;
180  result.CS_lhc8_tthj_ratio[i] = 0.;
181  for(int j = 0; j < 3; j++) result.BR_hjhihi[i][j] = 0.;
182  for(int j = 0; j < 3; j++) result.CS_lep_hjhi_ratio[i][j] = 0.;
183  }
184 
185  // Zero all H+ masses, widths and effective couplings
186  result.MHplus[0] = 0.;
187  result.deltaMHplus[0] = 0.;
188  result.HpGammaTot[0] = 0.;
189  result.BR_tWpb = 0.;
190  result.BR_tHpjb[0] = 0.;
191  result.BR_Hpjcs[0] = 0.;
192  result.BR_Hpjcb[0] = 0.;
193  result.BR_Hptaunu[0] = 0.;
194  result.CS_lep_HpjHmi_ratio[0] = 0.;
195  }
196 
198  void SMHiggs_ModelParameters(hb_ModelParameters &result)
199  {
200  using namespace Pipes::SMHiggs_ModelParameters;
201  set_SMLikeHiggs_ModelParameters(Dep::SM_spectrum->get_HE(), *Dep::Higgs_Couplings, result);
202  }
203 
205  void SMLikeHiggs_ModelParameters(hb_ModelParameters &result)
206  {
207  using namespace Pipes::SMLikeHiggs_ModelParameters;
208  dep_bucket<Spectrum>* spectrum_dependency = nullptr;
209  if (ModelInUse("ScalarSingletDM_Z2") or ModelInUse("ScalarSingletDM_Z2_running")) spectrum_dependency = &Dep::ScalarSingletDM_Z2_spectrum;
210  else if (ModelInUse("ScalarSingletDM_Z3") or ModelInUse("ScalarSingletDM_Z3_running")) spectrum_dependency = &Dep::ScalarSingletDM_Z3_spectrum;
211  else ColliderBit_error().raise(LOCAL_INFO, "No valid model for SMLikeHiggs_ModelParameters.");
212  const SubSpectrum& spec = (*spectrum_dependency)->get_HE();
213  set_SMLikeHiggs_ModelParameters(spec, *Dep::Higgs_Couplings, result);
214  }
215 
217  void MSSMHiggs_ModelParameters(hb_ModelParameters &result)
218  {
219  using namespace Pipes::MSSMHiggs_ModelParameters;
220 
221  // Set up neutral Higgses
222  static const std::vector<str> sHneut = initVector<str>("h0_1", "h0_2", "A0");
223 
224  // Set the CP of the Higgs states.
225  for (int i = 0; i < 3; i++) result.CP[i] = Dep::Higgs_Couplings->CP[i];
226 
227  // Retrieve higgs partial widths
228  const HiggsCouplingsTable::h0_decay_array_type& h0_widths = Dep::Higgs_Couplings->get_neutral_decays_array(3);
229  const DecayTable::Entry& H_plus_widths = Dep::Higgs_Couplings->get_charged_decays(0);
230  const DecayTable::Entry& t_widths = Dep::Higgs_Couplings->get_t_decays();
231 
232  // Retrieve masses
233  const Spectrum& fullspectrum = *Dep::MSSM_spectrum;
234  const SubSpectrum& spec = fullspectrum.get_HE();
235 
236  // Neutral higgs masses and errors
237  for(int i = 0; i < 3; i++)
238  {
239  result.Mh[i] = spec.get(Par::Pole_Mass,sHneut[i]);
240  double upper = spec.get(Par::Pole_Mass_1srd_high,sHneut[i]);
241  double lower = spec.get(Par::Pole_Mass_1srd_low,sHneut[i]);
242  result.deltaMh[i] = result.Mh[i] * std::max(upper,lower);
243  }
244 
245  // Loop over all neutral Higgses, setting their branching fractions and total widths.
246  for(int i = 0; i < 3; i++)
247  {
248  result.hGammaTot[i] = h0_widths[i]->width_in_GeV;
249  result.BR_hjss[i] = h0_widths[i]->BF("s", "sbar");
250  result.BR_hjcc[i] = h0_widths[i]->BF("c", "cbar");
251  result.BR_hjbb[i] = h0_widths[i]->BF("b", "bbar");
252  result.BR_hjmumu[i] = h0_widths[i]->BF("mu+", "mu-");
253  result.BR_hjtautau[i] = h0_widths[i]->BF("tau+", "tau-");
254  result.BR_hjWW[i] = h0_widths[i]->BF("W+", "W-");
255  result.BR_hjZZ[i] = h0_widths[i]->BF("Z0", "Z0");
256  result.BR_hjZga[i] = h0_widths[i]->BF("gamma", "Z0");
257  result.BR_hjgaga[i] = h0_widths[i]->BF("gamma", "gamma");
258  result.BR_hjgg[i] = h0_widths[i]->BF("g", "g");
259  // Do decays to invisibles
260  result.BR_hjinvisible[i] = 0.;
261  for (auto it = Dep::Higgs_Couplings->invisibles.begin(); it != Dep::Higgs_Couplings->invisibles.end(); ++it)
262  {
263  result.BR_hjinvisible[i] += h0_widths[i]->BF(*it, *it);
264  }
265  // Do decays to other neutral higgses
266  for (int j = 0; j < 3; j++)
267  {
268  if (2.*result.Mh[j] < result.Mh[i] and h0_widths[i]->has_channel(sHneut[j],sHneut[j]))
269  {
270  result.BR_hjhihi[i][j] = h0_widths[i]->BF(sHneut[j],sHneut[j]);
271  }
272  else
273  {
274  result.BR_hjhihi[i][j] = 0.;
275  }
276  }
277  }
278 
279  // Charged higgs masses and errors
280  result.MHplus[0] = spec.get(Par::Pole_Mass,"H+");
281  double upper = spec.get(Par::Pole_Mass_1srd_high,"H+");
282  double lower = spec.get(Par::Pole_Mass_1srd_low,"H+");
283  result.deltaMHplus[0] = result.MHplus[0] * std::max(upper,lower);
284 
285  // Set charged Higgs branching fractions and total width.
286  result.HpGammaTot[0] = H_plus_widths.width_in_GeV;
287  result.BR_Hpjcs[0] = H_plus_widths.BF("c", "sbar");
288  result.BR_Hpjcb[0] = H_plus_widths.BF("c", "bbar");
289  result.BR_Hptaunu[0] = H_plus_widths.BF("tau+", "nu_tau");
290 
291  // Set top branching fractions
292  result.BR_tWpb = t_widths.BF("W+", "b");
293  result.BR_tHpjb[0] = t_widths.has_channel("H+", "b") ? t_widths.BF("H+", "b") : 0.0;
294 
295  // Retrieve cross-section ratios from the HiggsCouplingsTable
296  set_CS(result, *Dep::Higgs_Couplings, 3);
297  }
298 
299 
301  void calc_HB_LEP_LogLike(double &result)
302  {
303  using namespace Pipes::calc_HB_LEP_LogLike;
304 
305  hb_ModelParameters ModelParam = *Dep::HB_ModelParameters;
306 
307  Farray<double, 1,3, 1,3> CS_lep_hjhi_ratio;
308  Farray<double, 1,3, 1,3> BR_hjhihi;
309  for(int i = 0; i < 3; i++) for(int j = 0; j < 3; j++)
310  {
311  CS_lep_hjhi_ratio(i+1,j+1) = ModelParam.CS_lep_hjhi_ratio[i][j];
312  BR_hjhihi(i+1,j+1) = ModelParam.BR_hjhihi[i][j];
313  }
314 
315  BEreq::HiggsBounds_neutral_input_part(&ModelParam.Mh[0], &ModelParam.hGammaTot[0], &ModelParam.CP[0],
316  &ModelParam.CS_lep_hjZ_ratio[0], &ModelParam.CS_lep_bbhj_ratio[0],
317  &ModelParam.CS_lep_tautauhj_ratio[0], CS_lep_hjhi_ratio,
318  &ModelParam.CS_gg_hj_ratio[0], &ModelParam.CS_bb_hj_ratio[0],
319  &ModelParam.CS_bg_hjb_ratio[0], &ModelParam.CS_ud_hjWp_ratio[0],
320  &ModelParam.CS_cs_hjWp_ratio[0], &ModelParam.CS_ud_hjWm_ratio[0],
321  &ModelParam.CS_cs_hjWm_ratio[0], &ModelParam.CS_gg_hjZ_ratio[0],
322  &ModelParam.CS_dd_hjZ_ratio[0], &ModelParam.CS_uu_hjZ_ratio[0],
323  &ModelParam.CS_ss_hjZ_ratio[0], &ModelParam.CS_cc_hjZ_ratio[0],
324  &ModelParam.CS_bb_hjZ_ratio[0], &ModelParam.CS_tev_vbf_ratio[0],
325  &ModelParam.CS_tev_tthj_ratio[0], &ModelParam.CS_lhc7_vbf_ratio[0],
326  &ModelParam.CS_lhc7_tthj_ratio[0], &ModelParam.CS_lhc8_vbf_ratio[0],
327  &ModelParam.CS_lhc8_tthj_ratio[0], &ModelParam.BR_hjss[0],
328  &ModelParam.BR_hjcc[0], &ModelParam.BR_hjbb[0],
329  &ModelParam.BR_hjmumu[0], &ModelParam.BR_hjtautau[0],
330  &ModelParam.BR_hjWW[0], &ModelParam.BR_hjZZ[0],
331  &ModelParam.BR_hjZga[0], &ModelParam.BR_hjgaga[0],
332  &ModelParam.BR_hjgg[0], &ModelParam.BR_hjinvisible[0], BR_hjhihi);
333 
334  BEreq::HiggsBounds_charged_input(&ModelParam.MHplus[0], &ModelParam.HpGammaTot[0], &ModelParam.CS_lep_HpjHmi_ratio[0],
335  &ModelParam.BR_tWpb, &ModelParam.BR_tHpjb[0], &ModelParam.BR_Hpjcs[0],
336  &ModelParam.BR_Hpjcb[0], &ModelParam.BR_Hptaunu[0]);
337 
338  BEreq::HiggsBounds_set_mass_uncertainties(&ModelParam.deltaMh[0],&ModelParam.deltaMHplus[0]);
339 
340  // run Higgs bounds 'classic'
341  double obsratio;
342  int HBresult, chan, ncombined;
343  BEreq::run_HiggsBounds_classic(HBresult,chan,obsratio,ncombined);
344 
345  // extract the LEP chisq
346  double chisq_withouttheory,chisq_withtheory;
347  int chan2;
348  double theor_unc = 1.5; // theory uncertainty
349  BEreq::HB_calc_stats(theor_unc,chisq_withouttheory,chisq_withtheory,chan2);
350 
351  // Catch HiggsBound's error value, chisq = -999
352  if( fabs(chisq_withouttheory - (-999.)) < 1e-6)
353  {
354  ColliderBit_warning().raise(LOCAL_INFO, "Got chisq=-999 from HB_calc_stats in HiggsBounds, indicating a cross-section outside tabulated range. Will use chisq=0.");
355  chisq_withouttheory = 0.0;
356  }
357 
358  result = -0.5*chisq_withouttheory;
359  }
360 
362  void calc_HS_LHC_LogLike(double &result)
363  {
364  using namespace Pipes::calc_HS_LHC_LogLike;
365 
366  hb_ModelParameters ModelParam = *Dep::HB_ModelParameters;
367 
368  Farray<double, 1,3, 1,3> CS_lep_hjhi_ratio;
369  Farray<double, 1,3, 1,3> BR_hjhihi;
370  for(int i = 0; i < 3; i++) for(int j = 0; j < 3; j++)
371  {
372  CS_lep_hjhi_ratio(i+1,j+1) = ModelParam.CS_lep_hjhi_ratio[i][j];
373  BR_hjhihi(i+1,j+1) = ModelParam.BR_hjhihi[i][j];
374  }
375 
376  BEreq::HiggsBounds_neutral_input_part_HS(&ModelParam.Mh[0], &ModelParam.hGammaTot[0], &ModelParam.CP[0],
377  &ModelParam.CS_lep_hjZ_ratio[0], &ModelParam.CS_lep_bbhj_ratio[0],
378  &ModelParam.CS_lep_tautauhj_ratio[0], CS_lep_hjhi_ratio,
379  &ModelParam.CS_gg_hj_ratio[0], &ModelParam.CS_bb_hj_ratio[0],
380  &ModelParam.CS_bg_hjb_ratio[0], &ModelParam.CS_ud_hjWp_ratio[0],
381  &ModelParam.CS_cs_hjWp_ratio[0], &ModelParam.CS_ud_hjWm_ratio[0],
382  &ModelParam.CS_cs_hjWm_ratio[0], &ModelParam.CS_gg_hjZ_ratio[0],
383  &ModelParam.CS_dd_hjZ_ratio[0], &ModelParam.CS_uu_hjZ_ratio[0],
384  &ModelParam.CS_ss_hjZ_ratio[0], &ModelParam.CS_cc_hjZ_ratio[0],
385  &ModelParam.CS_bb_hjZ_ratio[0], &ModelParam.CS_tev_vbf_ratio[0],
386  &ModelParam.CS_tev_tthj_ratio[0], &ModelParam.CS_lhc7_vbf_ratio[0],
387  &ModelParam.CS_lhc7_tthj_ratio[0], &ModelParam.CS_lhc8_vbf_ratio[0],
388  &ModelParam.CS_lhc8_tthj_ratio[0], &ModelParam.BR_hjss[0],
389  &ModelParam.BR_hjcc[0], &ModelParam.BR_hjbb[0],
390  &ModelParam.BR_hjmumu[0], &ModelParam.BR_hjtautau[0],
391  &ModelParam.BR_hjWW[0], &ModelParam.BR_hjZZ[0],
392  &ModelParam.BR_hjZga[0], &ModelParam.BR_hjgaga[0],
393  &ModelParam.BR_hjgg[0], &ModelParam.BR_hjinvisible[0], BR_hjhihi);
394 
395  BEreq::HiggsBounds_charged_input_HS(&ModelParam.MHplus[0], &ModelParam.HpGammaTot[0], &ModelParam.CS_lep_HpjHmi_ratio[0],
396  &ModelParam.BR_tWpb, &ModelParam.BR_tHpjb[0], &ModelParam.BR_Hpjcs[0],
397  &ModelParam.BR_Hpjcb[0], &ModelParam.BR_Hptaunu[0]);
398 
399  BEreq::HiggsSignals_neutral_input_MassUncertainty(&ModelParam.deltaMh[0]);
400 
401  // add uncertainties to cross-sections and branching ratios
402  // double dCS[5] = {0.,0.,0.,0.,0.};
403  // double dBR[5] = {0.,0.,0.,0.,0.};
404  // BEreq::setup_rate_uncertainties(dCS,dBR);
405 
406  // run HiggsSignals
407  int mode = 1; // 1- peak-centered chi2 method (recommended)
408  double csqmu, csqmh, csqtot, Pvalue;
409  int nobs;
410  BEreq::run_HiggsSignals(mode, csqmu, csqmh, csqtot, nobs, Pvalue);
411 
412  result = -0.5*csqtot;
413 
414  #ifdef COLLIDERBIT_DEBUG
415  std::ofstream f;
416  f.open ("HB_ModelParameters_contents.dat");
417  f<<"LHC log-likleihood";
418  for (int i = 0; i < 3; i++) f<<
419  " higgs index" <<
420  " "<<i<<":CP"<<
421  " "<<i<<":Mh"<<
422  " "<<i<<":hGammaTot"<<
423  " "<<i<<":CS_lep_hjZ_ratio"<<
424  " "<<i<<":CS_tev_vbf_ratio"<<
425  " "<<i<<":CS_lep_bbhj_ratio"<<
426  " "<<i<<":CS_lep_tautauhj_ratio"<<
427  " "<<i<<":CS_gg_hj_ratio"<<
428  " "<<i<<":CS_tev_tthj_ratio"<<
429  " "<<i<<":CS_lhc7_tthj_ratio"<<
430  " "<<i<<":CS_lhc8_tthj_ratio"<<
431  " "<<i<<":CS_lep_hjhi_ratio[0]"<<
432  " "<<i<<":CS_lep_hjhi_ratio[1]"<<
433  " "<<i<<":CS_lep_hjhi_ratio[2]"<<
434  " "<<i<<":BR_ss"<<
435  " "<<i<<":BR_cc"<<
436  " "<<i<<":BR_bb"<<
437  " "<<i<<":BR_mumu"<<
438  " "<<i<<":BR_tautau"<<
439  " "<<i<<":BR_WW"<<
440  " "<<i<<":BR_ZZ"<<
441  " "<<i<<":BR_Zga"<<
442  " "<<i<<":BR_gamgam"<<
443  " "<<i<<":BR_gg"<<
444  " "<<i<<":BR_invisible"<<
445  " "<<i<<":BR_hihi[0]"<<
446  " "<<i<<":BR_hihi[1]"<<
447  " "<<i<<":BR_hihi[2]";
448  f<<
449  " higgs index" <<
450  " "<<4<<"MHplus"<<
451  " "<<4<<":HpGammaTot"<<
452  " "<<4<<":CS_lep_HpjHmi_ratio"<<
453  " "<<4<<":BR_H+->cs"<<
454  " "<<4<<":BR_H+->cb"<<
455  " "<<4<<":BR_H+->taunu"<<
456  " "<<4<<":BR_t->W+b"<<
457  " "<<4<<":BR_t->H+b";
458  f << endl << std::setw(18) << result;
459  const int w = 24;
460  for (int i = 0; i < 3; i++)
461  {
462  f << std::setw(w) << i << std::setw(w) <<
463  ModelParam.CP[i] << std::setw(w) <<
464  ModelParam.Mh[i] << std::setw(w) <<
465  ModelParam.hGammaTot[i] << std::setw(w) <<
466  ModelParam.CS_lep_hjZ_ratio[i] << std::setw(w) <<
467  ModelParam.CS_tev_vbf_ratio[i] << std::setw(w) <<
468  ModelParam.CS_lep_bbhj_ratio[i] << std::setw(w) <<
469  ModelParam.CS_lep_tautauhj_ratio[i] << std::setw(w) <<
470  ModelParam.CS_gg_hj_ratio[i] << std::setw(w) <<
471  ModelParam.CS_tev_tthj_ratio[i] << std::setw(w) <<
472  ModelParam.CS_lhc7_tthj_ratio[i] << std::setw(w) <<
473  ModelParam.CS_lhc8_tthj_ratio[i];
474  for (int j = 0; j < 3; j++) f << std::setw(w) << ModelParam.CS_lep_hjhi_ratio[i][j];
475  f << std::setw(w) <<
476  ModelParam.BR_hjss[i] << std::setw(w) <<
477  ModelParam.BR_hjcc[i] << std::setw(w) <<
478  ModelParam.BR_hjbb[i] << std::setw(w) <<
479  ModelParam.BR_hjmumu[i] << std::setw(w) <<
480  ModelParam.BR_hjtautau[i] << std::setw(w) <<
481  ModelParam.BR_hjWW[i] << std::setw(w) <<
482  ModelParam.BR_hjZZ[i] << std::setw(w) <<
483  ModelParam.BR_hjZga[i] << std::setw(w) <<
484  ModelParam.BR_hjgaga[i] << std::setw(w) <<
485  ModelParam.BR_hjgg[i] << std::setw(w) <<
486  ModelParam.BR_hjinvisible[i];
487  for (int j = 0; j < 3; j++) f << std::setw(w) << ModelParam.BR_hjhihi[i][j];
488  }
489  f << std::setw(w) << 4 << std::setw(w) <<
490  ModelParam.MHplus[0] << std::setw(w) <<
491  ModelParam.HpGammaTot[0] << std::setw(w) <<
492  ModelParam.CS_lep_HpjHmi_ratio[0] << std::setw(w) <<
493  ModelParam.BR_Hpjcs[0] << std::setw(w) <<
494  ModelParam.BR_Hpjcb[0] << std::setw(w) <<
495  ModelParam.BR_Hptaunu[0] << std::setw(w) <<
496  ModelParam.BR_tWpb << std::setw(w) <<
497  ModelParam.BR_tHpjb[0];
498  f.close();
499  #endif
500 
501  }
502 
504  void FH_HiggsProd(fh_HiggsProd &result)
505  {
506  using namespace Pipes::FH_HiggsProd;
507 
508  Farray<fh_real, 1,52> prodxs;
509 
510  fh_HiggsProd HiggsProd;
511  int error;
512  fh_real sqrts;
513 
514  // Tevatron
515  sqrts = 2.;
516  error = 1;
517  BEreq::FHHiggsProd(error, sqrts, prodxs);
518  if (error != 0)
519  {
520  std::ostringstream err;
521  err << "BEreq::FHHiggsProd raised error flag for Tevatron: " << error << ".";
522  invalid_point().raise(err.str());
523  }
524  for(int i = 0; i < 52; i++) HiggsProd.prodxs_Tev[i] = prodxs(i+1);
525  // LHC7
526  sqrts = 7.;
527  error = 1;
528  BEreq::FHHiggsProd(error, sqrts, prodxs);
529  if (error != 0)
530  {
531  std::ostringstream err;
532  err << "BEreq::FHHiggsProd raised error flag for LHC7: " << error << ".";
533  invalid_point().raise(err.str());
534  }
535  for(int i = 0; i < 52; i++) HiggsProd.prodxs_LHC7[i] = prodxs(i+1);
536  // LHC8
537  sqrts = 8.;
538  error = 1;
539  BEreq::FHHiggsProd(error, sqrts, prodxs);
540  if (error != 0)
541  {
542  std::ostringstream err;
543  err << "BEreq::FHHiggsProd raised error flag for LHC8: " << error << ".";
544  invalid_point().raise(err.str());
545  }
546  for(int i = 0; i < 52; i++) HiggsProd.prodxs_LHC8[i] = prodxs(i+1);
547 
548  // The ttbar production cross-sections for the (BSM,SM) model can be found at (prodxs_X[h+27], prodxs_X[h+30]),
549  // where h is the higgs index (0 = h0_1, 1 = h0_2, 2 = A0) and X is one of Tev, LHC7 or LHC8.
550  result = HiggsProd;
551 
552  }
553 
554 
555  }
556 }
Rollcall header for ColliderBit module.
Array class that matches the memory structure and functionality of arrays in Fortran codes Syntax: Fa...
Definition: util_types.hpp:306
DecayTable entry class. Holds the info on all decays of a given particle.
Definition: decay_table.hpp:96
double C_ZZ2[max_neutral_higgses]
double BF(const std::vector< std::pair< int, int > > &) const
Retrieve branching fraction for decay to a given final state.
std::vector< str > invisibles
Particles that higgses can decay invisibly to.
double C_tt2[max_neutral_higgses]
void MSSMHiggs_ModelParameters(hb_ModelParameters &result)
MSSM Higgs model parameters.
double CP[max_neutral_higgses]
CP of neutral higgses.
double C_WW2[max_neutral_higgses]
Effective couplings for neutral higgses.
const DecayTable::Entry & get_neutral_decays(int) const
Retrieve decays of a specific neutral Higgs, by index.
void set_CS(hb_ModelParameters &result, const HiggsCouplingsTable &couplings, int n_neutral_higgses)
Helper function to set HiggsBounds/Signals parameters cross-section ratios from a GAMBIT HiggsCouplin...
void FH_HiggsProd(fh_HiggsProd &result)
Higgs production cross-sections from FeynHiggs.
#define LOCAL_INFO
Definition: local_info.hpp:34
GAMBIT native higgs coupling table class.
void SMLikeHiggs_ModelParameters(hb_ModelParameters &result)
SM-like (SM + possible invisibles) Higgs model parameters for HiggsBounds/Signals.
Spectrum Spectrum Spectrum Spectrum Spectrum Spectrum Spectrum ScalarSingletDM_Z3_spectrum
GAMBIT error class.
Definition: exceptions.hpp:136
void calc_HB_LEP_LogLike(double &result)
Get a LEP chisq from HiggsBounds.
virtual void raise(const std::string &)
Raise the exception, i.e. throw it.
Definition: exceptions.cpp:422
virtual double get(const Par::Tags, const str &, const SpecOverrideOptions=use_overrides, const SafeBool check_antiparticle=SafeBool(true)) const =0
double C_bb2[max_neutral_higgses]
double C_tautau2[max_neutral_higgses]
double width_in_GeV
Total particle width (in GeV)
bool has_channel(const std::vector< std::pair< int, int > > &) const
Check if a given final state exists in this DecayTable::Entry.
An interface class for module dependencies.
Header file that includes all GAMBIT headers required for a module source file.
const DecayTable::Entry * h0_decay_array_type[max_neutral_higgses]
Types to make returning decay arrays easier.
Virtual base class for interacting with spectrum generator output.
Definition: subspectrum.hpp:87
void SMHiggs_ModelParameters(hb_ModelParameters &result)
SM Higgs model parameters for HiggsBounds/Signals.
Spectrum Spectrum Spectrum Spectrum Spectrum Spectrum SM_spectrum
invalid_point_exception & invalid_point()
Invalid point exceptions.
double C_gg2[max_neutral_higgses]
void set_SMLikeHiggs_ModelParameters(const SubSpectrum &spec, const HiggsCouplingsTable &couplings, hb_ModelParameters &result)
Helper function for populating a HiggsBounds/Signals ModelParameters object for SM-like Higgs...
TODO: see if we can use this one:
Definition: Analysis.hpp:33
virtual bool has(const Par::Tags, const str &, const SpecOverrideOptions=use_overrides, const SafeBool check_antiparticle=SafeBool(true)) const =0
Getters/Setters etc.
SubSpectrum & get_HE()
Definition: spectrum.cpp:225
"Standard Model" (low-energy) plus high-energy model container class
Definition: spectrum.hpp:110
double C_hiZ2[max_neutral_higgses][max_neutral_higgses]
void calc_HS_LHC_LogLike(double &result)
Get an LHC chisq from HiggsSignals.