12auto checkH1 = [](
TList *out) {
16 std::cout <<
"checkH1 >>> Test failure: output list not found\n";
21 auto hdmd =
dynamic_cast<TH1F *
>(out->FindObject(
"hdmd"));
23 std::cout <<
"checkH1 >>> Test failure: 'hdmd' histo not found\n";
26 if ((
Int_t)(hdmd->GetEntries()) != 7525) {
27 std::cout <<
"checkH1 >>> Test failure: 'hdmd' histo: wrong number"
29 << (
Int_t)(hdmd->GetEntries()) <<
": expected 7525) \n";
32 if (
TMath::Abs((hdmd->GetMean() - 0.15512023) / 0.15512023) > 0.001) {
33 std::cout <<
"checkH1 >>> Test failure: 'hdmd' histo: wrong mean (" << hdmd->
GetMean()
34 <<
": expected 0.15512023) \n";
38 auto h2 =
dynamic_cast<TH2F *
>(out->FindObject(
"h2"));
40 std::cout <<
"checkH1 >>> Test failure: 'h2' histo not found\n";
43 if ((
Int_t)(h2->GetEntries()) != 7525) {
44 std::cout <<
"checkH1 >>> Test failure: 'h2' histo: wrong number"
46 << (
Int_t)(h2->GetEntries()) <<
": expected 7525) \n";
49 if (
TMath::Abs((h2->GetMean() - 0.15245688) / 0.15245688) > 0.001) {
50 std::cout <<
"checkH1 >>> Test failure: 'h2' histo: wrong mean (" << h2->
GetMean() <<
": expected 0.15245688) \n";
59auto doFit = [](
TList *out,
const char *lfn = 0) ->
Int_t {
69 if (hdmd == 0 || h2 == 0) {
70 std::cout <<
"doFit: hdmd = " << hdmd <<
" , h2 = " << h2 <<
"\n";
79 c1->SetBottomMargin(0.15);
80 hdmd->GetXaxis()->SetTitle(
"m_{K#pi#pi} - m_{K#pi}[GeV/c^{2}]");
81 hdmd->GetXaxis()->SetTitleOffset(1.4);
84 if (
gROOT->GetListOfFunctions()->FindObject(
"f5"))
85 delete gROOT->GetFunction(
"f5");
94 par[2] / 2.5066 / par[4] *
TMath::Exp(-xp3 / 2 / par[4] / par[4]));
98 auto f5 =
new TF1(
"f5",
fdm5, 0.139, 0.17, 5);
99 f5->SetParameters(1000000, .25, 2000, .1454, .001);
100 hdmd->Fit(
"f5",
"lr");
103 Double_t ref_f5[4] = {959915.0, 0.351114, 1185.03, 0.145569};
104 for (
int i : {0, 1, 2, 3}) {
105 if ((
TMath::Abs((f5->GetParameters())[i] - ref_f5[i]) / ref_f5[i]) > 0.001) {
106 std::cout <<
"\n >>> Test failure: fit to 'f5': parameter '" << f5->GetParName(i) <<
"' has wrong value ("
107 << (f5->GetParameters())[i] <<
": expected" << ref_f5[i] <<
") \n";
117 auto c2 =
new TCanvas(
"c2",
"tauD0", 100, 100, 800, 600);
119 c2->SetBottomMargin(0.15);
125 if (
gROOT->GetListOfFunctions()->FindObject(
"f2"))
126 delete gROOT->GetFunction(
"f2");
129 const auto dxbin = (0.17 - 0.13) / 40;
130 const auto sigma = 0.0012;
134 auto xp3 = (
x - 0.1454) * (
x - 0.1454);
140 auto f2 =
new TF1(
"f2",
fdm2, 0.139, 0.17, 2);
141 f2->SetParameters(10000, 10);
144 std::cout <<
"doFit: restricting fit to two bins only in this example...\n";
146 h2->FitSlicesX(f2, 10, 20, 10,
"g5 l");
149 Double_t ref_f2[2] = {52432.2, 105.481};
150 for (
int i : {0, 1}) {
151 if ((
TMath::Abs((f2->GetParameters())[i] - ref_f2[i]) / ref_f2[i]) > 0.001) {
152 std::cout <<
"\n >>> Test failure: fit to 'f2': parameter '" << f2->GetParName(i) <<
"' has wrong value ("
153 << (f2->GetParameters())[i] <<
": expected" << ref_f2[i] <<
") \n";
161 h2_1->GetXaxis()->SetTitle(
"#tau[ps]");
162 h2_1->SetMarkerStyle(21);
170 auto psdmd = (
TPaveStats *)hdmd->GetListOfFunctions()->FindObject(
"stats");
171 psdmd->SetOptStat(1110);
184 auto hdmd =
new TH1F(
"hdmd",
"Dm_d", 40, 0.13, 0.17);
185 auto h2 =
new TH2F(
"h2",
"ptD0 vs Dm_d", 30, 0.135, 0.165, 30, -3, 6);
203 while (reader.Next()) {
215 if (fNhitrp.At(*fIk) * fNhitrp.At(*fIpi) <= 1)
218 if (fRend.At(*fIk) - fRstart.At(*fIk) <= 22)
220 if (fRend.At(*fIpi) - fRstart.At(*fIpi) <= 22)
222 if (fNlhk.At(*fIk) <= 0.1)
224 if (fNlhpi.At(*fIpi) <= 0.1)
227 if (fNlhpi.At(*fIpis) <= 0.1)
234 h2->Fill(*fDm_d, *fRpd0_t / 0.029979 * 1.8646 / *fPtd0_d);
250 auto elist =
new TEntryList(
"elist",
"H1 selection from Cut");
265 while (reader.Next()) {
277 if (fNhitrp.At(*fIk) * fNhitrp.At(*fIpi) <= 1)
280 if (fRend.At(*fIk) - fRstart.At(*fIk) <= 22)
282 if (fRend.At(*fIpi) - fRstart.At(*fIpi) <= 22)
284 if (fNlhk.At(*fIk) <= 0.1)
286 if (fNlhpi.At(*fIpi) <= 0.1)
289 if (fNlhpi.At(*fIpis) <= 0.1)
295 elist->Enter(reader.GetCurrentEntry(), reader.GetTree());
305 auto hdmd =
new TH1F(
"hdmd",
"Dm_d", 40, 0.13, 0.17);
306 auto h2 =
new TH2F(
"h2",
"ptD0 vs Dm_d", 30, 0.135, 0.165, 30, -3, 6);
312 while (reader.Next()) {
315 h2->Fill(*fDm_d, *fRpd0_t / 0.029979 * 1.8646 / *fPtd0_d);
R__EXTERN TStyle * gStyle
R__EXTERN TSystem * gSystem
A List of entry numbers in a TTree or TChain.
1-D histogram with a double per channel (see TH1 documentation)}
1-D histogram with a float per channel (see TH1 documentation)}
virtual Double_t GetMean(Int_t axis=1) const
For axis = 1,2 or 3 returns the mean value of the histogram along X,Y or Z axis.
2-D histogram with a float per channel (see TH1 documentation)}
Use the TLine constructor to create a simple line.
TObject * FindObject(const char *name) const override
Find an object in this list using its name.
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
The histogram statistics painter class.
void SetOptStat(Int_t stat=1)
The type of information printed in the histogram statistics box can be selected via the parameter mod...
void SetOptFit(Int_t fit=1)
The type of information about fit parameters printed in the histogram statistics box can be selected ...
virtual Int_t RedirectOutput(const char *name, const char *mode="a", RedirectHandle_t *h=nullptr)
Redirect standard output (stdout, stderr) to the specified file.
An interface for reading collections stored in ROOT columnar datasets.
An interface for reading values stored in ROOT columnar datasets.
A simple, robust and fast interface to read values from ROOT columnar datasets such as TTree,...
Double_t fdm5(Double_t *xx, Double_t *par)
Double_t fdm2(Double_t *xx, Double_t *par)
Double_t Exp(Double_t x)
Returns the base-e exponential function of x, which is e raised to the power x.
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Returns x raised to the power y.
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.