hpstr
The Heavy Photon Search Toolkit for Reconstruction (hpstr) provides an interface to physics data from the HPS experiment saved in the LCIO format and converts it into an ROOT based format.
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
BumpHunter.h
Go to the documentation of this file.
1
10#ifndef __BUMP_HUNTER_H__
11#define __BUMP_HUNTER_H__
12
13//----------------//
14// C++ StdLib //
15//----------------//
16#include <cmath>
17#include <cstdio>
18#include <exception>
19#include <fstream>
20#include <map>
21#include <vector>
22
23//----------//
24// ROOT //
25//----------//
26#include <TF1.h>
27#include <TFitResult.h>
28#include <TFitResultPtr.h>
29#include <TH1.h>
30#include <TMath.h>
31#include <TCanvas.h>
32#include <TFile.h>
33#include <TRandom.h>
34#include <TString.h>
35
36#include <Math/ProbFunc.h>
37#include <Math/MinimizerOptions.h>
38
39//---//
40#include "HpsFitResult.h"
42#include "LegendreFitFunction.h"
43
50 public:
61 BumpHunter(FitFunction::BkgModel model, int poly_order, int toy_poly_order,
62 int res_factor, double res_scale = 1.00, bool asymptotic_limit = true);
63
66
74 HpsFitResult* performSearch(TH1* histogram, double mass_hypothesis, bool skip_bkg_fit,
75 bool skip_ul);
76
89 void initialize(TH1* histogram, double &mass_hypothesis);
90
96 void calculatePValue(HpsFitResult* result);
97
103
110 void setBounds(double low_bound, double high_bound);
111
117 void enableDebug(bool debug = true) { this->debug = debug; };
118
124 void runBatchMode(bool batch = true) { _batch = batch; };
125
131 void writeResults(bool write_results = true) { _write_results = write_results; };
132
139 void getUpperLimit(TH1* histogram, HpsFitResult* result);
140
147 void getUpperLimitAsymptotic(TH1* histogram, HpsFitResult* result);
148
155 void getUpperLimitAsymCLs(TH1* histogram, HpsFitResult* result);
156
163 void getUpperLimitPower(TH1* histogram, HpsFitResult* result);
164
170 void setResolutionScale(double res_scale) { res_scale_ = res_scale; }
171
177 void setWindowSizeUsesResScale(bool window_use_res_scale) {
178 window_use_res_scale_ = window_use_res_scale;
179 }
180
192 std::vector<TH1*> generateToys(TH1* histogram, double n_toys, int seed, int toy_sig_samples,
193 int bkg_mult = 1, TH1* signal_hist = nullptr);
194
202 double getMassResolution(double mass) {
203 // Omar's 2015 mass resolution.
204 //return res_scale_ * ((-6.2 * mass * mass * mass) + (0.91 * mass * mass) - (0.00297 * mass) + 0.000579);
205 //return 0.0389938364847*mass - 0.0000713783511061; // ideal
206 //return 0.0501460737193*mass - 0.0000917925595224; // scaled to moller mass from data
207 //return 0.0532190838657*mass - 0.0000922283032152; // scaled to moller mass + sys
208 //return res_scale_*(0.000955 - 0.004198 * mass + 0.2367 * mass * mass - 0.7009 * mass * mass * mass);
209
210 // Rafo's 2016 mass resolution.
211 //return res_scale_ * (0.000379509 + (0.0416842 * mass) - (0.271364 * mass * mass) + (3.49537 * mass * mass * mass) - (11.1153 * mass * mass * mass * mass));
212 return getMassResolution(mass, res_scale_);
213 };
214
215 private:
223 inline double getMassResolution(double mass, double res_scale) {
224 return res_scale * (0.000379509 + (0.0416842 * mass) - (0.271364 * mass * mass) + (3.49537 * mass * mass * mass) - (11.1153 * mass * mass * mass * mass));
225 }
226
233 double correctMass(double mass);
234
240 void printDebug(std::string message);
241
250 void getChi2Prob(double min_nll_null, double min_nll, double &q0, double &p_value);
251
254
256 std::ofstream* ofs;
257
258 //
259 // Variable definitions
260 //
262 double corr_mass_{0};
263
265 double lower_bound_{0.016};
266
268 double upper_bound_{0.115};
269
271 double integral_{0};
272
275
281
286 double res_factor_{13};
287
289 double window_size_{0};
290
292 int bins_{0};
293
295 double bin_width_{0.0};
296
299
302
304 //double res_scale_{1.56};
305 double res_scale_{1.00};
306
311 bool _batch{false};
312
314 bool debug{false};
315
317 bool _write_results{false};
318
320 double window_start_{0};
321
323 double window_end_{0};
324
327
330
333};
334
335#endif // __BUMP_HUNTER_H__
description
Definition BumpHunter.h:49
void getChi2Prob(double min_nll_null, double min_nll, double &q0, double &p_value)
description
void setWindowSizeUsesResScale(bool window_use_res_scale)
Sets whether the window size is scaled according to the resolution scaling factor.
Definition BumpHunter.h:177
double bin_width_
Definition BumpHunter.h:295
bool _write_results
Definition BumpHunter.h:317
int toy_poly_order_
Definition BumpHunter.h:301
void printDebug(std::string message)
Print debug statement.
double correctMass(double mass)
description
double window_end_
Definition BumpHunter.h:323
double mass_hypothesis_
Definition BumpHunter.h:326
void calculatePValue(HpsFitResult *result)
description
double getMassResolution(double mass, double res_scale)
Get Mass Resolution.
Definition BumpHunter.h:223
double window_size_
Definition BumpHunter.h:289
double window_start_
Definition BumpHunter.h:320
void getUpperLimitPower(TH1 *histogram, HpsFitResult *result)
description
void getUpperLimitAsymptotic(TH1 *histogram, HpsFitResult *result)
description
void getUpperLimitAsymCLs(TH1 *histogram, HpsFitResult *result)
description
HpsFitResult * bkg_only_result_
Definition BumpHunter.h:253
void initialize(TH1 *histogram, double &mass_hypothesis)
Given the mass of interest, setup the window parameters and initialize the fit parameters.
double getMassResolution(double mass)
Definition BumpHunter.h:202
void runBatchMode(bool batch=true)
Enable batch running.
Definition BumpHunter.h:124
void setResolutionScale(double res_scale)
Set the resolution after instantiation.
Definition BumpHunter.h:170
int poly_order_
Definition BumpHunter.h:298
double res_scale_
Definition BumpHunter.h:305
double integral_
Definition BumpHunter.h:271
double res_factor_
Definition BumpHunter.h:286
void enableDebug(bool debug=true)
Enable/disable debug.
Definition BumpHunter.h:117
double lower_bound_
Definition BumpHunter.h:265
bool window_use_res_scale_
Definition BumpHunter.h:332
void writeResults(bool write_results=true)
Write the fit results to a text file.
Definition BumpHunter.h:131
double upper_bound_
Definition BumpHunter.h:268
FitFunction::BkgModel bkg_model_
Definition BumpHunter.h:274
void fitBkgOnly()
Fit using a background only model.
void setBounds(double low_bound, double high_bound)
Set the histogram bounds.
std::vector< TH1 * > generateToys(TH1 *histogram, double n_toys, int seed, int toy_sig_samples, int bkg_mult=1, TH1 *signal_hist=nullptr)
description
void getUpperLimit(TH1 *histogram, HpsFitResult *result)
Get the signal upper limit.
bool asymptotic_limit_
Definition BumpHunter.h:280
HpsFitResult * performSearch(TH1 *histogram, double mass_hypothesis, bool skip_bkg_fit, bool skip_ul)
Perform a search for a resonance at the given mass hypothesis.
double mass_resolution_
Definition BumpHunter.h:329
std::ofstream * ofs
Definition BumpHunter.h:256
double corr_mass_
Definition BumpHunter.h:262
BkgModel
description
Definition FitFunction.h:39
description