Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
hist041_TProfile2Poly_realistic.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_hist
3/// \notebook
4/// Different charges depending on region
5///
6/// \macro_image
7/// \macro_code
8///
9/// \date November 2022
10/// \author Filip Ilic
11
12#include <iostream>
13#include <fstream>
14#include <TProfile2Poly.h>
15#include <TProfile2D.h>
16#include <TCanvas.h>
17#include <TRandom.h>
18#include <TROOT.h>
19
21{
22 int NUM_LS = 8;
23 TCanvas *c1 = new TCanvas("c1", "moving charge", 900, 400);
24 TCanvas *c2 = new TCanvas("c2", "Merge Individual moving charge plots", 800, 400);
25
26 c1->Divide(NUM_LS, 3);
27 c2->Divide(3, 1);
28
29 // -------------------- Construct Reference plot bins ------------------------
30 auto new_avg = new TProfile2Poly();
31
35 auto tot_merge = new TProfile2Poly();
36 auto det_avg_merge = new TProfile2Poly();
37 auto det_err_merge = new TProfile2Poly();
38
39 float minx = -15;
40 float maxx = 15;
41 float miny = -15;
42 float maxy = 15;
43 float binsz = 0.5;
44
45 for (float i = minx; i < maxx; i += binsz) {
46 for (float j = miny; j < maxy; j += binsz) {
47 tot_merge->AddBin(i, j, i + binsz, j + binsz);
48 for (int l = 0; l < NUM_LS; ++l) {
49 tot_avg_ls[l].AddBin(i, j, i + binsz, j + binsz);
50 }
51 }
52 }
53
54 // -------------------- Construct detector bins ------------------------
55 auto h2p = new TH2Poly();
56 auto tp2p = new TProfile2Poly();
57 std::ifstream infile;
58 TString dir = gROOT->GetTutorialDir();
59 dir.Append("/hist/data/tprofile2poly_tutorial.data");
60 infile.open(dir.Data());
61
62 if (!infile) // Verify that the file was open successfully
63 {
64 std::cerr << dir.Data() << std::endl; // Report error
65 std::cerr << "Error code: " << std::strerror(errno) << std::endl; // Get some info as to why
66 return;
67 }
68 std::cout << " WE ARE AFTER LOADING DATA " << std::endl;
69
70 std::vector<std::pair<Double_t, Double_t>> allCoords;
71 Double_t a, b;
72 while (infile >> a >> b)
73 allCoords.emplace_back(a, b);
74
75 if (allCoords.size() % 3 != 0) {
76 cout << "[ERROR] Bad file" << endl;
77 return;
78 }
79
80 Double_t x[3], y[3];
81 for (int i = 0; i < allCoords.size(); i += 3) {
82 x[0] = allCoords[i + 0].first;
83 y[0] = allCoords[i + 0].second;
84 x[1] = allCoords[i + 1].first;
85 y[1] = allCoords[i + 1].second;
86 x[2] = allCoords[i + 2].first;
87 y[2] = allCoords[i + 2].second;
88
89 det_avg_merge->AddBin(3, x, y);
90 det_err_merge->AddBin(3, x, y);
91
92 for (int l = 0; l < NUM_LS; ++l) {
93 det_avg_ls[l].AddBin(3, x, y);
94 det_err_ls[l].AddBin(3, x, y);
95 }
96 }
97
98 std::cout << " WE ARE AFTER ADDING BINS " << std::endl;
99
100 // -------------------- Simulate particles ------------------------
101 TRandom ran;
102
103 // moving error
104 Double_t xoffset1 = 0;
105 Double_t yoffset1 = 0;
106 Double_t xoffset2 = 0;
107 Double_t yoffset2 = 0;
108
109 for (int i = 0; i <= NUM_LS - 1; ++i) { // LumiSection
110 std::cout << "[In Progress] LumiSection " << i << std::endl;
111 for (int j = 0; j < numEvents; ++j) { // Events
112 Double_t r1 = ran.Gaus(0, 10);
113 Double_t r2 = ran.Gaus(0, 8);
114 Double_t rok = ran.Gaus(10, 1);
115 Double_t rbad1 = ran.Gaus(8, 5);
116 Double_t rbad2 = ran.Gaus(-8, 5);
117
118 Double_t val = rok;
119
120 xoffset1 += 0.00002;
121 yoffset1 += 0.00002;
122
123 xoffset2 += 0.00003;
124 yoffset2 += 0.00004;
125
126 if (r2 > 3. - yoffset1 && r2 < 8. - yoffset1 && r1 > 1. + xoffset1 && r1 < 5. + xoffset1) {
127 val -= rbad1;
128 }
129
130 if (r2 > -10 + yoffset2 && r2 < -8 + yoffset2 && r1 > -6 + xoffset2 && r1 < 8 + xoffset2) {
131 val -= rbad2;
132 }
133
134 tot_avg_ls[i].Fill(r1, r2, val);
135 det_avg_ls[i].Fill(r1, r2, val);
136 det_err_ls[i].Fill(r1, r2, val);
137 }
138
139 std::string title;
140
141 c1->cd(i + 1);
142 title = "Global View: Avg in LS " + std::to_string(i);
143 tot_avg_ls[i].SetTitle(title.c_str());
144 tot_avg_ls[i].SetStats(false);
145 tot_avg_ls[i].Draw("COLZ");
146 c1->Update();
147
148 c1->cd((i + 1) + NUM_LS);
149 title = "Detector View: Avg in LS " + std::to_string(i);
150 det_avg_ls[i].SetTitle(title.c_str());
151 det_avg_ls[i].SetStats(false);
152 det_avg_ls[i].Draw("COLZ");
153 c1->Update();
154
155 c1->cd((i + 1) + (NUM_LS * 2));
156 title = "Detector View: Error in LS " + std::to_string(i);
157 det_err_ls[i].SetTitle(title.c_str());
158 det_err_ls[i].SetStats(false);
159 det_err_ls[i].SetContentToError();
160 det_err_ls[i].Draw("COLZ");
161 c1->Update();
162 }
163
164 std::vector<TProfile2Poly *> tot_avg_v;
165 std::vector<TProfile2Poly *> det_avg_v;
166 for (int t = 0; t < NUM_LS; t++) {
167 tot_avg_v.push_back(&tot_avg_ls[t]);
168 det_avg_v.push_back(&det_avg_ls[t]);
169 }
170
171 std::cout << "[In Progress] Merging" << std::endl;
172
173 tot_merge->Merge(tot_avg_v);
174 c2->cd(1);
175 tot_merge->SetTitle("Total average merge");
176 tot_merge->Draw("COLZ");
177
178 det_avg_merge->Merge(det_avg_v);
179 c2->cd(2);
180 det_avg_merge->SetTitle("Detector average merge");
181 det_avg_merge->SetContentToAverage(); // implicit
182 det_avg_merge->Draw("COLZ");
183
184 det_err_merge->Merge(det_avg_v);
185 c2->cd(3);
186 det_err_merge->SetTitle("Detector error merge");
187 det_err_merge->SetContentToError();
188 det_err_merge->Draw("COLZ");
189}
#define b(i)
Definition RSha256.hxx:100
#define a(i)
Definition RSha256.hxx:99
int Int_t
Definition RtypesCore.h:45
double Double_t
Definition RtypesCore.h:59
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
#define gROOT
Definition TROOT.h:406
The Canvas class.
Definition TCanvas.h:23
2D Histogram with Polygonal Bins
Definition TH2Poly.h:66
2D Profile Histogram with Polygonal Bins.
This is the base class for the ROOT Random number generators.
Definition TRandom.h:27
Basic string class.
Definition TString.h:139
const char * Data() const
Definition TString.h:376
TString & Append(const char *cs)
Definition TString.h:572
Double_t y[n]
Definition legend1.C:17
return c1
Definition legend1.C:41
Double_t x[n]
Definition legend1.C:17
return c2
Definition legend2.C:14
TLine l
Definition textangle.C:4