gambit is hosted by Hepforge, IPPP Durham
GAMBIT  v1.5.0-252-gf9a3f78
a Global And Modular Bsm Inference Tool
DarkBit_standalone_ScalarSingletDM_Z2.cpp
Go to the documentation of this file.
1 // GAMBIT: Global and Modular BSM Inference Tool
2 // *********************************************
20 
27 
28 using namespace DarkBit::Functown; // Functors wrapping the module's actual module functions
29 using namespace BackendIniBit::Functown; // Functors wrapping the backend initialisation functions
30 
33 QUICK_FUNCTION(DarkBit, cascadeMC_gammaSpectra, OLD_CAPABILITY, CMC_dummy, DarkBit::stringFunkMap, ())
34 
35 
36 namespace Gambit
37 {
38 
39  namespace DarkBit
40  {
41 
42  // Dummy functor for returning empty cascadeMC result spectra
43  void CMC_dummy(DarkBit::stringFunkMap& result)
44  {
46  result = sfm;
47  }
48 
49  // Create spectrum object from SLHA file SM.slha and ScalarSingletDM_Z2 model parameters
50  void createSpectrum(Spectrum& outSpec)
51  {
52  using namespace Pipes::createSpectrum;
53  std::string inputFileName = "DarkBit/data/SM.slha";
54 
56  singletmodel.HiggsPoleMass = 125.;
57  singletmodel.HiggsVEV = 246.;
58  singletmodel.SingletPoleMass = *Param["mS"];
59  singletmodel.SingletLambda = *Param["lambda_hS"];
60 
61  SLHAstruct slhaea = read_SLHA(inputFileName);
62  outSpec = spectrum_from_SLHAea<Models::ScalarSingletDM_Z2SimpleSpec, Models::ScalarSingletDM_Z2Model>(singletmodel, slhaea, Spectrum::mc_info(), Spectrum::mr_info());
63  }
64 
65  // Create decay object from SLHA file decays.slha
66  void createDecays(DecayTable& outDecays)
67  {
68  using namespace Pipes::createDecays;
69 
70  std::string filename = "DarkBit/data/decays.slha";
71  outDecays = DecayTable(filename);
72  }
73  }
74 }
75 
76 int main()
77 {
78 
79  try
80  {
81 
82  std::cout << std::endl
83  << "Start DarkBit Scalar Singlet DM standalone example!" << std::endl;
84  std::cout << "---------------------------------------------------" << std::endl;
85  std::cout << std::endl;
86  std::cout << "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!" << std::endl;
87  std::cout << "This program needs DarkBit/data/SM.slha for SM parameters and " << std::endl;
88  std::cout << "DarkBit/data/decays.slha for the Higgs width and branching fraction. If " << std::endl;
89  std::cout << "these are not present, it dies!" << std::endl;
90  std::cout << "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!" << std::endl;
91  std::cout << std::endl;
92 
93  // ---- Initialise logging and exceptions ----
94 
95  initialise_standalone_logs("runs/DarkBit_standalone_ScalarSingletDM_Z2/logs/");
96  logger()<<"Running DarkBit standalone example"<<LogTags::info<<EOM;
97  model_warning().set_fatal(true);
98 
99 
100  // ---- Check that required backends are present ----
101 
102  if (not Backends::backendInfo().works["DarkSUSY_generic_wimp6.2.2"]) backend_error().raise(LOCAL_INFO, "DarkSUSY 6.2.2 for a generic WIMP is missing!");
103  if (not Backends::backendInfo().works["MicrOmegas_ScalarSingletDM_Z23.6.9.2"]) backend_error().raise(LOCAL_INFO, "MicrOmegas 3.6.9.2 for ScalarSingletDM_Z2 is missing!");
104  if (not Backends::backendInfo().works["gamLike1.0.1"]) backend_error().raise(LOCAL_INFO, "gamLike 1.0.1 is missing!");
105  if (not Backends::backendInfo().works["DDCalc2.2.0"]) backend_error().raise(LOCAL_INFO, "DDCalc 2.2.0 is missing!");
106  if (not Backends::backendInfo().works["nulike1.0.9"]) backend_error().raise(LOCAL_INFO, "nulike 1.0.9 is missing!");
107 
108 
109  // ---- Initialize models ----
110 
111  // Initialize ScalarSingletDM_Z2 model -- Adjust the model parameters here:
112  ModelParameters* SingletDM_primary_parameters = Models::ScalarSingletDM_Z2::Functown::primary_parameters.getcontentsPtr();
113  SingletDM_primary_parameters->setValue("mS", 1000.);
114  SingletDM_primary_parameters->setValue("lambda_hS", 1.0);
115 
116  // Initialize halo model
117  ModelParameters* Halo_primary_parameters = Models::Halo_Einasto::Functown::primary_parameters.getcontentsPtr();
118  Halo_primary_parameters->setValue("rho0", 0.4);
119  Halo_primary_parameters->setValue("rhos", 0.08);
120  Halo_primary_parameters->setValue("vrot", 235.);
121  Halo_primary_parameters->setValue("v0", 235.);
122  Halo_primary_parameters->setValue("vesc", 550.);
123  Halo_primary_parameters->setValue("rs", 20.);
124  Halo_primary_parameters->setValue("r_sun", 8.5);
125  Halo_primary_parameters->setValue("alpha", 0.17);
126 
127 
128  // --- Resolve halo dependencies ---
129  ExtractLocalMaxwellianHalo.notifyOfModel("Halo_Einasto");
130  ExtractLocalMaxwellianHalo.resolveDependency(&Models::Halo_Einasto::Functown::primary_parameters);
131  ExtractLocalMaxwellianHalo.reset_and_calculate();
132 
133  GalacticHalo_Einasto.notifyOfModel("Halo_Einasto");
134  GalacticHalo_Einasto.resolveDependency(&Models::Halo_Einasto::Functown::primary_parameters);
135  GalacticHalo_Einasto.reset_and_calculate();
136 
137  // Initialize nuclear_params_fnq model
138  ModelParameters* nuclear_params_fnq = Models::nuclear_params_fnq::Functown::primary_parameters.getcontentsPtr();
139  nuclear_params_fnq->setValue("fpd", 0.034);
140  nuclear_params_fnq->setValue("fpu", 0.023);
141  nuclear_params_fnq->setValue("fps", 0.14);
142  nuclear_params_fnq->setValue("fnd", 0.041);
143  nuclear_params_fnq->setValue("fnu", 0.019);
144  nuclear_params_fnq->setValue("fns", 0.14);
145  nuclear_params_fnq->setValue("deltad", -0.40);
146  nuclear_params_fnq->setValue("deltau", 0.74);
147  nuclear_params_fnq->setValue("deltas", -0.12);
148 
149 
150  // ---- Initialize spectrum and decays ---
151 
152  createSpectrum.notifyOfModel("ScalarSingletDM_Z2");
153  createSpectrum.resolveDependency(&Models::ScalarSingletDM_Z2::Functown::primary_parameters);
154  createSpectrum.reset_and_calculate();
155 
156  createDecays.notifyOfModel("ScalarSingletDM_Z2");
157  createDecays.reset_and_calculate();
158 
159 
160  // ---- Set up basic internal structures for direct & indirect detection ----
161 
162  // Set identifier for DM particle
163  DarkMatter_ID_ScalarSingletDM.notifyOfModel("ScalarSingletDM_Z2");
164  DarkMatter_ID_ScalarSingletDM.reset_and_calculate();
165 
166  // Set up process catalog
167  TH_ProcessCatalog_ScalarSingletDM_Z2.notifyOfModel("ScalarSingletDM_Z2");
170  TH_ProcessCatalog_ScalarSingletDM_Z2.reset_and_calculate();
171 
172  // Assume for direct and indirect detection likelihoods that dark matter
173  // density is always the measured one (despite relic density results)
174  RD_fraction_one.reset_and_calculate();
175 
176  // Set generic WIMP mass object
178  mwimp_generic.resolveDependency(&DarkMatter_ID_ScalarSingletDM);
179  mwimp_generic.reset_and_calculate();
180 
181  // Set generic annihilation rate in late universe (v->0 limit)
184  sigmav_late_universe.reset_and_calculate();
185 
186  // ---- Initialize backends ----
187 
188  // Initialize nulike backend
189  Backends::nulike_1_0_9::Functown::nulike_bounds.setStatus(2);
190  nulike_1_0_9_init.reset_and_calculate();
191 
192  // Initialize gamLike backend
193  gamLike_1_0_1_init.reset_and_calculate();
194 
195  // Initialize MicrOmegas backend (specific for ScalarSingletDM_Z2)
196  MicrOmegas_ScalarSingletDM_Z2_3_6_9_2_init.notifyOfModel("ScalarSingletDM_Z2");
197  MicrOmegas_ScalarSingletDM_Z2_3_6_9_2_init.resolveDependency(&createSpectrum);
198  MicrOmegas_ScalarSingletDM_Z2_3_6_9_2_init.resolveDependency(&createDecays);
199  MicrOmegas_ScalarSingletDM_Z2_3_6_9_2_init.reset_and_calculate();
200  // For the below VXdecay = 0 - no 3 body final states via virtual X
201  // 1 - annihilations to 3 body final states via virtual X
202  // 2 - (co)annihilations to 3 body final states via virtual X
203  MicrOmegas_ScalarSingletDM_Z2_3_6_9_2_init.setOption<int>("VZdecay", 1);
204  MicrOmegas_ScalarSingletDM_Z2_3_6_9_2_init.setOption<int>("VWdecay", 1);
205  MicrOmegas_ScalarSingletDM_Z2_3_6_9_2_init.reset_and_calculate();
206 
207  // Initialize DarkSUSY backend
208  DarkSUSY_generic_wimp_6_2_2_init.reset_and_calculate();
209 
210  // Initialize DarkSUSY Local Halo Model
213  DarkSUSY_PointInit_LocalHalo_func.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::dshmcom);
214  DarkSUSY_PointInit_LocalHalo_func.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::dshmisodf);
215  DarkSUSY_PointInit_LocalHalo_func.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::dshmframevelcom);
216  DarkSUSY_PointInit_LocalHalo_func.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::dshmnoclue);
217  DarkSUSY_PointInit_LocalHalo_func.reset_and_calculate();
218 
219  // ---- Relic density ----
220 
221  // Relic density calculation with MicrOmegas
222  RD_oh2_Xf_MicrOmegas.notifyOfModel("ScalarSingletDM_Z2");
223  RD_oh2_Xf_MicrOmegas.resolveBackendReq(&Backends::MicrOmegas_ScalarSingletDM_Z2_3_6_9_2::Functown::darkOmega);
224  RD_oh2_Xf_MicrOmegas.reset_and_calculate();
225  RD_oh2_MicrOmegas.resolveDependency(&RD_oh2_Xf_MicrOmegas);
226  RD_oh2_MicrOmegas.reset_and_calculate();
227 
228  // Retrieve and print MicrOmegas result
229  logger() << "Omega h^2 from MicrOmegas: " << RD_oh2_MicrOmegas(0) << LogTags::info << EOM;
230 
231  // Relic density calculation with GAMBIT (DarkSUSY Boltzmann solver)
234  RD_spectrum_from_ProcessCatalog.reset_and_calculate();
235 
236  RD_eff_annrate_from_ProcessCatalog.resolveDependency(&TH_ProcessCatalog_ScalarSingletDM_Z2);
237  RD_eff_annrate_from_ProcessCatalog.resolveDependency(&DarkMatter_ID_ScalarSingletDM);
238  RD_eff_annrate_from_ProcessCatalog.reset_and_calculate();
239 
241  RD_spectrum_ordered_func.reset_and_calculate();
242 
243  RD_oh2_DS_general.resolveDependency(&RD_spectrum_ordered_func);
244  RD_oh2_DS_general.resolveDependency(&RD_eff_annrate_from_ProcessCatalog);
245  RD_oh2_DS_general.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::rdpars);
246  RD_oh2_DS_general.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::rdtime);
247  RD_oh2_DS_general.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::dsrdcom);
248  RD_oh2_DS_general.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::dsrdstart);
249  RD_oh2_DS_general.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::dsrdens);
250  RD_oh2_DS_general.setOption<int>("fast", 1); // 0: normal; 1: fast; 2: dirty
251  RD_oh2_DS_general.reset_and_calculate();
252 
253  // Retrieve and print GAMBIT result
254  logger() << "Omega h^2 from GAMBIT: " << RD_oh2_DS_general(0) << LogTags::info << EOM;
255 
256  // Calculate WMAP likelihoods
257  // Uncomment one of the following two line to choose which Omega h^2 value to use
258  //lnL_oh2_Simple.resolveDependency(&RD_oh2_MicrOmegas);
259  lnL_oh2_Simple.resolveDependency(&RD_oh2_DS_general);
260  lnL_oh2_Simple.reset_and_calculate();
261 
262  logger() << "Relic density lnL: " << lnL_oh2_Simple((0)) << LogTags::info << EOM;
263 
264  // ---- Direct detection -----
265 
266  // Calculate DD couplings with Micromegas
267 
268  DD_couplings_MicrOmegas.notifyOfModel("ScalarSingletDM_Z2");
269  DD_couplings_MicrOmegas.notifyOfModel("nuclear_params_fnq");
270  DD_couplings_MicrOmegas.resolveDependency(&Models::nuclear_params_fnq::Functown::primary_parameters);
271  DD_couplings_MicrOmegas.resolveBackendReq(&Backends::MicrOmegas_ScalarSingletDM_Z2_3_6_9_2::Functown::nucleonAmplitudes);
272  DD_couplings_MicrOmegas.resolveBackendReq(&Backends::MicrOmegas_ScalarSingletDM_Z2_3_6_9_2::Functown::FeScLoop);
273  DD_couplings_MicrOmegas.resolveBackendReq(&Backends::MicrOmegas_ScalarSingletDM_Z2_3_6_9_2::Functown::mocommon_);
274  DD_couplings_MicrOmegas.reset_and_calculate();
275 
276  // Calculate DD couplings with GAMBIT
277 
278  DD_couplings_ScalarSingletDM_Z2.notifyOfModel("nuclear_params_fnq");
279  DD_couplings_ScalarSingletDM_Z2.notifyOfModel("ScalarSingletDM_Z2");
280  DD_couplings_ScalarSingletDM_Z2.resolveDependency(&Models::nuclear_params_fnq::Functown::primary_parameters);
282  DD_couplings_ScalarSingletDM_Z2.reset_and_calculate();
283 
284  // Set generic scattering cross-sections for later use
285  double sigma_SI_p_GB, sigma_SI_p_MO;
286 
287  sigma_SI_p_simple.resolveDependency(&DD_couplings_MicrOmegas);
288  sigma_SI_p_simple.resolveDependency(&mwimp_generic);
289  sigma_SI_p_simple.reset_and_calculate();
290  sigma_SI_p_MO = sigma_SI_p_simple(0);
291 
292  sigma_SD_p_simple.resolveDependency(&DD_couplings_MicrOmegas);
293  sigma_SD_p_simple.resolveDependency(&mwimp_generic);
294  sigma_SD_p_simple.reset_and_calculate();
295 
296  logger() << "sigma_SI,p with MicrOmegas: " << sigma_SI_p_simple(0) << LogTags::info << EOM;
297 
298  // Set generic scattering cross-sections for later use
300  sigma_SI_p_simple.reset_and_calculate();
301  sigma_SI_p_GB = sigma_SI_p_simple(0);
302 
304  sigma_SD_p_simple.reset_and_calculate();
305 
306  logger() << "sigma_SI,p with GAMBIT: " << sigma_SI_p_simple(0) << LogTags::info << EOM;
307 
308  // Initialize DDCalc backend
309  Backends::DDCalc_2_2_0::Functown::DDCalc_CalcRates_simple.setStatus(2);
310  Backends::DDCalc_2_2_0::Functown::DDCalc_Experiment.setStatus(2);
311  Backends::DDCalc_2_2_0::Functown::DDCalc_LogLikelihood.setStatus(2);
312 
313  DDCalc_2_2_0_init.resolveDependency(&ExtractLocalMaxwellianHalo);
314  DDCalc_2_2_0_init.resolveDependency(&RD_fraction_one);
315  DDCalc_2_2_0_init.resolveDependency(&mwimp_generic);
316  // Choose one of the two below lines to determine where the couplings used in the likelihood
317  // calculation come from
318  DDCalc_2_2_0_init.resolveDependency(&DD_couplings_ScalarSingletDM_Z2);
319  //DDCalc_2_2_0_init.resolveDependency(&DD_couplings_MicrOmegas);
320  DDCalc_2_2_0_init.reset_and_calculate();
321 
322  // Calculate direct detection rates for LUX 2016
323  LUX_2016_Calc.resolveBackendReq(&Backends::DDCalc_2_2_0::Functown::DDCalc_Experiment);
324  LUX_2016_Calc.resolveBackendReq(&Backends::DDCalc_2_2_0::Functown::DDCalc_CalcRates_simple);
325  LUX_2016_Calc.reset_and_calculate();
326 
327  // Calculate direct detection likelihood for LUX 2016
328  LUX_2016_GetLogLikelihood.resolveDependency(&LUX_2016_Calc);
329  LUX_2016_GetLogLikelihood.resolveBackendReq(&Backends::DDCalc_2_2_0::Functown::DDCalc_Experiment);
330  LUX_2016_GetLogLikelihood.resolveBackendReq(&Backends::DDCalc_2_2_0::Functown::DDCalc_LogLikelihood);
331  LUX_2016_GetLogLikelihood.reset_and_calculate();
332 
333  logger() << "LUX_2016 lnL: " << LUX_2016_GetLogLikelihood(0) << LogTags::info << EOM;
334 
335  // ---- Gamma-ray yields ----
336 
337  // Initialize tabulated gamma-ray yields
338  SimYieldTable_DarkSUSY.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::dsanyield_sim);
339  SimYieldTable_DarkSUSY.reset_and_calculate();
340 
341  // Collect missing final states for simulation in cascade MC
343  GA_missingFinalStates.resolveDependency(&SimYieldTable_DarkSUSY);
345  GA_missingFinalStates.reset_and_calculate();
346 
347 
348  // Infer for which type of final states particles MC should be performed
349  cascadeMC_FinalStates.setOption<std::vector<std::string>>("cMC_finalStates", daFunk::vec<std::string>("gamma"));
350  cascadeMC_FinalStates.reset_and_calculate();
351 
352  // Collect decay information for cascade MC
354  cascadeMC_DecayTable.resolveDependency(&SimYieldTable_DarkSUSY);
355  cascadeMC_DecayTable.reset_and_calculate();
356 
357  // Set up MC loop manager for cascade MC
358  cascadeMC_LoopManager.resolveDependency(&GA_missingFinalStates);
359  std::vector<functor*> nested_functions = initVector<functor*>(
361  cascadeMC_LoopManager.setNestedList(nested_functions);
362 
363  // Set up initial state for cascade MC step
364  cascadeMC_InitialState.resolveDependency(&GA_missingFinalStates);
365  cascadeMC_InitialState.resolveLoopManager(&cascadeMC_LoopManager);
366 
367  // Perform MC step for cascade MC
368  cascadeMC_GenerateChain.resolveDependency(&cascadeMC_InitialState);
369  cascadeMC_GenerateChain.resolveDependency(&cascadeMC_DecayTable);
370  cascadeMC_GenerateChain.resolveLoopManager(&cascadeMC_LoopManager);
371 
372  // Generate histogram for cascade MC
373  cascadeMC_Histograms.resolveDependency(&cascadeMC_InitialState);
374  cascadeMC_Histograms.resolveDependency(&cascadeMC_GenerateChain);
375  cascadeMC_Histograms.resolveDependency(&TH_ProcessCatalog_ScalarSingletDM_Z2);
376  cascadeMC_Histograms.resolveDependency(&SimYieldTable_DarkSUSY);
377  cascadeMC_Histograms.resolveDependency(&cascadeMC_FinalStates);
378  cascadeMC_Histograms.resolveLoopManager(&cascadeMC_LoopManager);
379 
380  // Check convergence of cascade MC
381  cascadeMC_EventCount.resolveDependency(&cascadeMC_InitialState);
382  cascadeMC_EventCount.resolveLoopManager(&cascadeMC_LoopManager);
383 
384  // Start cascade MC loop
385  cascadeMC_LoopManager.reset_and_calculate();
386 
387  // Infer gamma-ray spectra for recorded MC results
388  cascadeMC_gammaSpectra.resolveDependency(&GA_missingFinalStates);
389  cascadeMC_gammaSpectra.resolveDependency(&cascadeMC_FinalStates);
390  cascadeMC_gammaSpectra.resolveDependency(&cascadeMC_Histograms);
391  cascadeMC_gammaSpectra.resolveDependency(&cascadeMC_EventCount);
392  cascadeMC_gammaSpectra.reset_and_calculate();
393 
394  // Calculate total gamma-ray yield (cascade MC + tabulated results)
396  GA_AnnYield_General.resolveDependency(&SimYieldTable_DarkSUSY);
398  GA_AnnYield_General.resolveDependency(&cascadeMC_gammaSpectra);
399  GA_AnnYield_General.reset_and_calculate();
400 
401  // Calculate Fermi LAT dwarf likelihood
403  lnL_FermiLATdwarfs_gamLike.resolveDependency(&RD_fraction_one);
404  lnL_FermiLATdwarfs_gamLike.resolveBackendReq(&Backends::gamLike_1_0_1::Functown::lnL);
405  lnL_FermiLATdwarfs_gamLike.reset_and_calculate();
406 
407  logger() << "Fermi LAT dwarf spheroidal lnL: " << lnL_FermiLATdwarfs_gamLike(0) << LogTags::info << EOM;
408 
409  // ---- IceCube limits ----
410 
411  // Infer WIMP capture rate in Sun
412  capture_rate_Sun_const_xsec.resolveDependency(&mwimp_generic);
415  capture_rate_Sun_const_xsec.resolveDependency(&RD_fraction_one);
416  capture_rate_Sun_const_xsec.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::dssenu_capsuntab);
419  capture_rate_Sun_const_xsec.reset_and_calculate();
420 
421  // Infer WIMP equilibration time in Sun
424  equilibration_time_Sun.resolveDependency(&mwimp_generic);
426  equilibration_time_Sun.reset_and_calculate();
427 
428  // Infer WIMP annihilation rate in Sun
429  annihilation_rate_Sun.resolveDependency(&equilibration_time_Sun);
431  annihilation_rate_Sun.reset_and_calculate();
432 
433  // Infer neutrino yield from Sun
435  nuyield_from_DS.resolveDependency(&mwimp_generic);
436  nuyield_from_DS.resolveDependency(&sigmav_late_universe);
438  nuyield_from_DS.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::dsgenericwimp_nusetup);
439  nuyield_from_DS.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::neutrino_yield);
440  nuyield_from_DS.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::DS_neutral_h_decay_channels);
441  nuyield_from_DS.resolveBackendReq(&Backends::DarkSUSY_generic_wimp_6_2_2::Functown::DS_charged_h_decay_channels);
442  nuyield_from_DS.reset_and_calculate();
443 
444  // Calculate number of events at IceCube
445  IC79WH_full.resolveDependency(&mwimp_generic);
446  IC79WH_full.resolveDependency(&annihilation_rate_Sun);
447  IC79WH_full.resolveDependency(&nuyield_from_DS);
448  IC79WH_full.resolveBackendReq(&Backends::nulike_1_0_9::Functown::nulike_bounds);
449  IC79WH_full.reset_and_calculate();
450  IC79WL_full.resolveDependency(&mwimp_generic);
451  IC79WL_full.resolveDependency(&annihilation_rate_Sun);
452  IC79WL_full.resolveDependency(&nuyield_from_DS);
453  IC79WL_full.resolveBackendReq(&Backends::nulike_1_0_9::Functown::nulike_bounds);
454  IC79WL_full.reset_and_calculate();
455  IC79SL_full.resolveDependency(&mwimp_generic);
456  IC79SL_full.resolveDependency(&annihilation_rate_Sun);
457  IC79SL_full.resolveDependency(&nuyield_from_DS);
458  IC79SL_full.resolveBackendReq(&Backends::nulike_1_0_9::Functown::nulike_bounds);
459  IC79SL_full.reset_and_calculate();
460 
461  // Calculate IceCube likelihood
462  IC79WH_bgloglike.resolveDependency(&IC79WH_full);
463  IC79WH_bgloglike.reset_and_calculate();
464  IC79WH_loglike.resolveDependency(&IC79WH_full);
465  IC79WH_loglike.reset_and_calculate();
466  IC79WL_bgloglike.resolveDependency(&IC79WL_full);
467  IC79WL_bgloglike.reset_and_calculate();
468  IC79WL_loglike.resolveDependency(&IC79WL_full);
469  IC79WL_loglike.reset_and_calculate();
470  IC79SL_bgloglike.resolveDependency(&IC79SL_full);
471  IC79SL_bgloglike.reset_and_calculate();
472  IC79SL_loglike.resolveDependency(&IC79SL_full);
473  IC79SL_loglike.reset_and_calculate();
474  IC79_loglike.resolveDependency(&IC79WH_bgloglike);
475  IC79_loglike.resolveDependency(&IC79WH_loglike);
476  IC79_loglike.resolveDependency(&IC79WL_bgloglike);
477  IC79_loglike.resolveDependency(&IC79WL_loglike);
478  IC79_loglike.resolveDependency(&IC79SL_bgloglike);
479  IC79_loglike.resolveDependency(&IC79SL_loglike);
480  IC79_loglike.reset_and_calculate();
481 
482  logger() << "IceCube 79 lnL: " << IC79_loglike(0) << LogTags::info << EOM;
483 
484  // ---- Dump results on screen ----
485 
486  cout << endl;
487  cout << "Omega h^2 from MicrOmegas: " << RD_oh2_MicrOmegas(0) << endl;
488  cout << "Omega h^2 from DarkSUSY (via DarkBit): " << RD_oh2_DS_general(0) << endl;
489  cout << "Relic density lnL: " << lnL_oh2_Simple(0) << endl;
490  cout << "sigma_SI,p with MicrOmegas: " << sigma_SI_p_MO << endl;
491  cout << "sigma_SI,p with DarkBit: " << sigma_SI_p_GB << endl;
492  cout << "LUX_2016 lnL: " << LUX_2016_GetLogLikelihood(0) << endl;
493  cout << "Fermi LAT dwarf spheroidal lnL: " << lnL_FermiLATdwarfs_gamLike(0) << endl;
494  cout << "IceCube 79 lnL: " << IC79_loglike(0) << endl;
495  }
496 
497  catch (std::exception& e)
498  {
499  std::cout << "DarkBit_standalone_ScalarSingletDM_Z2 has exited with fatal exception: " << e.what() << std::endl;
500  }
501 
502  return 0;
503 }
void cascadeMC_FinalStates(std::vector< std::string > &list)
Function for retrieving list of final states for cascade decays.
Definition: Cascades.cpp:41
void GA_missingFinalStates(std::vector< std::string > &result)
Identification of final states that are not yet tabulated.
Definition: GamYields.cpp:61
void IC79WL_full(nudata &result)
79-string IceCube WL sample: predicted signal and background counts, observed counts and likelihoods...
void DD_couplings_ScalarSingletDM_Z2(DM_nucleon_couplings &result)
Direct detection couplings for Z2 scalar singlet DM.
void equilibration_time_Sun(double &result)
Equilibration time for capture and annihilation of dark matter in the Sun (s)
void sigmav_late_universe(double &result)
Retrieve the total thermally-averaged annihilation cross-section for indirect detection (cm^3 / s)...
Definition: DarkBit.cpp:70
void DarkSUSY_PointInit_LocalHalo_func(bool &result)
Function to set Local Halo Parameters in DarkSUSY (DS 6)
void IC79SL_bgloglike(double &result)
#define NEW_CAPABILITY
void sigma_SI_p_simple(double &result)
Simple calculator of the spin-independent WIMP-proton cross-section.
Includes everything needed to use a GAMBIT module as a standalone analysis code.
void IC79WL_bgloglike(double &result)
void RD_fraction_one(double &result)
#define LOCAL_INFO
Definition: local_info.hpp:34
void SimYieldTable_DarkSUSY(SimYieldTable &result)
SimYieldTable based on DarkSUSY6 tabulated results.
Definition: GamYields.cpp:558
void cascadeMC_InitialState(std::string &pID)
Function selecting initial state for decay chain.
Definition: Cascades.cpp:140
SLHAstruct read_SLHA(str slha)
Read an SLHA file in to an SLHAea object with some error-checking.
void IC79WL_loglike(double &result)
void RD_oh2_Xf_MicrOmegas(ddpair &result)
Relic density directly from a call of initialized MicrOmegas.
void createDecays(DecayTable &outDecays)
Helpers for making simple spectra from SLHAea objects.
void TH_ProcessCatalog_ScalarSingletDM_Z2(DarkBit::TH_ProcessCatalog &result)
Set up process catalog for Z2 scalar singlet DM.
void nuyield_from_DS(nuyield_info &result)
Neutrino yield function pointer and setup.
void annihilation_rate_Sun(double &result)
Annihilation rate of dark matter in the Sun (s^-1)
void IC79WH_loglike(double &result)
void IC79SL_loglike(double &result)
warning & model_warning()
Model warnings.
Routines to help translate between SLHA2 sfermions and SLHA1 (or similar) sfermions.
void RD_spectrum_from_ProcessCatalog(RD_spectrum_type &result)
Collects information about resonances and threshold energies directly from the ProcessCatalog [NB: th...
std::map< str, daFunk::Funk > stringFunkMap
Definition: SimpleHist.hpp:101
General small utility functions.
void lnL_oh2_Simple(double &result)
Likelihood for cosmological relic density constraints.
void DD_couplings_MicrOmegas(DM_nucleon_couplings &result)
Get direct detection couplings from initialized MicrOmegas.
error & backend_error()
Backend errors.
SLHAea::Coll SLHAstruct
Less confusing name for SLHAea container class.
void cascadeMC_EventCount(std::map< std::string, int > &counts)
Event counter for cascade decays.
Definition: Cascades.cpp:175
void setValue(std::string const &inkey, double const &value)
Set single parameter value.
void IC79SL_full(nudata &result)
79-string IceCube SL sample: predicted signal and background counts, observed counts and likelihoods...
void cascadeMC_GenerateChain(DarkBit::DecayChain::ChainContainer &chain)
Function for generating decay chains.
Definition: Cascades.cpp:198
void DarkMatter_ID_ScalarSingletDM(std::string &result)
void capture_rate_Sun_const_xsec(double &result)
Capture rate of regular dark matter in the Sun (no v-dependent or q-dependent cross-sections) (s^-1)...
void GA_AnnYield_General(daFunk::Funk &result)
General routine to derive annihilation yield.
Definition: GamYields.cpp:191
void initialise_standalone_logs(str)
Logger setup standalone utility function.
void lnL_FermiLATdwarfs_gamLike(double &result)
Fermi LAT dwarf likelihoods, using gamLike backend.
const Logging::endofmessage EOM
Explicit const instance of the end of message struct in Gambit namespace.
Definition: logger.hpp:99
Logging::LogMaster & logger()
Function to retrieve a reference to the Gambit global log object.
Definition: logger.cpp:95
void IC79WH_full(nudata &result)
79-string IceCube WH sample: predicted signal and background counts, observed counts and likelihoods...
void IC79_loglike(double &result)
Composite IceCube 79-string likelihood function.
void ExtractLocalMaxwellianHalo(LocalMaxwellianHalo &result)
Module function providing local density and velocity dispersion parameters.
Definition: DarkBit.cpp:152
QUICK_FUNCTION(DarkBit, decay_rates, NEW_CAPABILITY, createDecays, DecayTable,()) QUICK_FUNCTION(DarkBit
Rollcall header for module DarkBit.
Spectrum ScalarSingletDM_Z2
void IC79WH_bgloglike(double &result)
GAMBIT native decay table class.
Definition: decay_table.hpp:35
std::vector< YAML::sdd > mc_info
Typedefs for making it easier to manipulate mass cut and mass ratio cut info.
Definition: spectrum.hpp:119
void mwimp_generic(double &result)
Retrieve the DM mass in GeV for generic models (GeV)
Definition: DarkBit.cpp:60
void sigma_SD_p_simple(double &result)
Simple calculator of the spin-dependent WIMP-proton cross-section.
void cascadeMC_LoopManager()
Loop manager for cascade decays.
Definition: Cascades.cpp:80
Simple extension of the SMHiggsSimpleSpec "model object" to include scalar singlet DM parameters We c...
void cascadeMC_Histograms(std::map< std::string, std::map< std::string, SimpleHist > > &result)
Function responsible for histogramming, and evaluating end conditions for event loop.
Definition: Cascades.cpp:363
void RD_oh2_MicrOmegas(double &result)
void cascadeMC_DecayTable(DarkBit::DecayChain::DecayTable &table)
Function setting up the decay table used in decay chains.
Definition: Cascades.cpp:58
void RD_oh2_DS_general(double &result)
General routine for calculation of relic density, using DarkSUSY 6+ Boltzmann solver.
std::vector< YAML::ssdd > mr_info
Definition: spectrum.hpp:120
DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry DecayTable::Entry decay_rates
void GalacticHalo_Einasto(GalacticHaloProperties &result)
Module function to generate GalacticHaloProperties for Einasto profile.
Definition: DarkBit.cpp:140
void RD_spectrum_ordered_func(RD_spectrum_type &result)
Order RD_spectrum object and derive coannihilation thresholds.
TODO: see if we can use this one:
Definition: Analysis.hpp:33
"Standard Model" (low-energy) plus high-energy model container class
Definition: spectrum.hpp:110