CGEM BOSS 6.6.5.h
BESIII Offline Software System
Loading...
Searching...
No Matches
CgemCosmicRayQA.cxx
Go to the documentation of this file.
1// **************************************************************************/
2// authors: L. Lavezzi (univ. of Torino & INFN, Italy)
3//
4
5//include system lib
6#include <iostream>
7#include <fstream>
8#include <iomanip>
9#include <string>
10#include <cmath>
11
12// Include files
13#include "GaudiKernel/AlgFactory.h"
14#include "GaudiKernel/DataObject.h"
15#include "GaudiKernel/IEventProcessor.h"
16
17#include "GaudiKernel/Incident.h"
18#include "GaudiKernel/IIncidentSvc.h"
19#include "GaudiKernel/Memory.h"
20
21#include <csignal>
22
23//for save digit & cluster
24#include "GaudiKernel/ISvcLocator.h"
25#include "GaudiKernel/IDataProviderSvc.h"
26#include "GaudiKernel/Bootstrap.h"
27#include "GaudiKernel/RegistryEntry.h"
28#include "GaudiKernel/MsgStream.h"
29
30#include "CgemRawEvent/CgemDigi.h"
31#include "CgemRecEvent/RecCgemCluster.h"
32#include "MdcRecEvent/RecMdcTrack.h"
33#include "Identifier/CgemID.h"
34
35#include "RawEvent/RawDataUtil.h"
36#include "RawEvent/DigiEvent.h"
37#include "ReconEvent/ReconEvent.h"
38#include "EventModel/EventHeader.h"
39#include "GaudiKernel/SmartDataPtr.h"
40
41#include "CgemLineFit/CgemCosmicRayQA.h"
42
43#include "TMath.h"
44//using namespace std;
45
46CgemCosmicRayQA::CgemCosmicRayQA(const std::string& name, ISvcLocator* pSvcLocator):
47 Algorithm(name,pSvcLocator){
48
49 declareProperty("filename", filename);
50 declareProperty("align_flag", align_flag);
51 // time window
52 declareProperty("minDigiTime", minDigiTime=-8875);
53 declareProperty("maxDigiTime", maxDigiTime=-8562);
54 // selection
55 declareProperty("select_size", select_size=1);
56 // cuts
57 declareProperty("cut_chi2", cut_chi2);
58 declareProperty("cut_ene_L1_x", cut_ene_L1_x=20); // fC
59 declareProperty("cut_ene_L1_v", cut_ene_L1_v=10); // fC
60 declareProperty("cut_ene_L2_x", cut_ene_L2_x=15); // fC
61 declareProperty("cut_ene_L2_v", cut_ene_L2_v=10); // fC
62 // declareProperty("LUTfile", lutfile = "/bes3fs/cgemCosmic/data/CGEM_cosmic_look_up_table_from_17_to_17.root");
63
64}
65
66CgemCosmicRayQA::~CgemCosmicRayQA(){
67 // delete lutreader;
68}
69
70
71bool CgemCosmicRayQA::read_file(TString name) {
72 cout << "reading " << name << endl;
73 input = new TFile(name, "READ");
74 tree = (TTree*) input->Get("tree");
75
76 cout << "tree " << tree << endl;
77 // hit --------------------------------------------
78 tree->SetBranchAddress("event", &event);
79 tree->SetBranchAddress("nhit", &nhit);
80 // from geometry
81 tree->SetBranchAddress("hit_strip", &hit_strip);
82 tree->SetBranchAddress("hit_view", &hit_view);
83 tree->SetBranchAddress("hit_layer", &hit_layer);
84 tree->SetBranchAddress("hit_sheet", &hit_sheet);
85 tree->SetBranchAddress("hit_length", &hit_length);
86 // from LUT
87 tree->SetBranchAddress("hit_channel", &hit_channel);
88 tree->SetBranchAddress("hit_roc", &hit_roc);
89 tree->SetBranchAddress("hit_feb", &hit_feb);
90 tree->SetBranchAddress("hit_tiger", &hit_tiger);
91 tree->SetBranchAddress("hit_chip", &hit_chip);
92 // from signal
93 tree->SetBranchAddress("hit_t", &hit_t);
94 tree->SetBranchAddress("hit_q", &hit_q);
95 tree->SetBranchAddress("hit_quality", &hit_quality);
96 tree->SetBranchAddress("hit_selgooddigi", &hit_selgooddigi);
97
98 // cluster ----------------------------------------
99 tree->SetBranchAddress("ncluster", &ncluster);
100 // 1d
101 tree->SetBranchAddress("ncluster_1d", &ncluster_1d);
102 tree->SetBranchAddress("ncluster_1d_L1_S1_x", &ncluster_1d_L1_S1_x);
103 tree->SetBranchAddress("ncluster_1d_L2_S1_x", &ncluster_1d_L2_S1_x);
104 tree->SetBranchAddress("ncluster_1d_L2_S2_x", &ncluster_1d_L2_S2_x);
105 tree->SetBranchAddress("ncluster_1d_L1_S1_v", &ncluster_1d_L1_S1_v);
106 tree->SetBranchAddress("ncluster_1d_L2_S1_v", &ncluster_1d_L2_S1_v);
107 tree->SetBranchAddress("ncluster_1d_L2_S2_v", &ncluster_1d_L2_S2_v);
108 tree->SetBranchAddress("cluster_1d_ID", &cluster_1d_ID);
109 tree->SetBranchAddress("cluster_1d_t", &cluster_1d_t);
110 tree->SetBranchAddress("cluster_1d_q", &cluster_1d_q);
111 tree->SetBranchAddress("cluster_1d_r", &cluster_1d_r);
112 tree->SetBranchAddress("cluster_1d_v", &cluster_1d_v);
113 tree->SetBranchAddress("cluster_1d_v_cc", &cluster_1d_v_cc);
114 tree->SetBranchAddress("cluster_1d_v_tpc", &cluster_1d_v_tpc);
115 tree->SetBranchAddress("cluster_1d_phi", &cluster_1d_phi);
116 tree->SetBranchAddress("cluster_1d_phi_cc", &cluster_1d_phi_cc);
117 tree->SetBranchAddress("cluster_1d_phi_tpc", &cluster_1d_phi_tpc);
118 tree->SetBranchAddress("cluster_1d_layerid", &cluster_1d_layerid);
119 tree->SetBranchAddress("cluster_1d_sheetid", &cluster_1d_sheetid);
120 tree->SetBranchAddress("cluster_1d_view", &cluster_1d_view);
121 tree->SetBranchAddress("cluster_1d_strip1", &cluster_1d_strip1);
122 tree->SetBranchAddress("cluster_1d_strip2", &cluster_1d_strip2);
123 tree->SetBranchAddress("cluster_1d_size", &cluster_1d_size);
124 tree->SetBranchAddress("cluster_1d_a_tpc", &cluster_1d_a_tpc);
125 tree->SetBranchAddress("cluster_1d_b_tpc", &cluster_1d_b_tpc);
126 tree->SetBranchAddress("cluster_1d_hitindex", cluster_1d_hitindex);
127
128 // 2d
129 tree->SetBranchAddress("ncluster_2d", &ncluster_2d);
130 tree->SetBranchAddress("cluster_2d_t", &cluster_2d_t);
131 tree->SetBranchAddress("cluster_2d_q", &cluster_2d_q);
132 tree->SetBranchAddress("cluster_2d_r", &cluster_2d_r);
133 tree->SetBranchAddress("cluster_2d_z", &cluster_2d_z);
134 tree->SetBranchAddress("cluster_2d_z_cc", &cluster_2d_z_cc);
135 tree->SetBranchAddress("cluster_2d_z_tpc", &cluster_2d_z_tpc);
136 tree->SetBranchAddress("cluster_2d_phi", &cluster_2d_phi);
137 tree->SetBranchAddress("cluster_2d_phi_cc", &cluster_2d_phi_cc);
138 tree->SetBranchAddress("cluster_2d_phi_tpc", &cluster_2d_phi_tpc);
139 tree->SetBranchAddress("cluster_2d_layerid", &cluster_2d_layerid);
140 tree->SetBranchAddress("cluster_2d_sheetid", &cluster_2d_sheetid);
141 tree->SetBranchAddress("cluster_2d_view", &cluster_2d_view);
142 tree->SetBranchAddress("cluster_2d_idx", &cluster_2d_idx);
143 tree->SetBranchAddress("cluster_2d_idv", &cluster_2d_idv);
144 tree->SetBranchAddress("cluster_2d_highest", &cluster_2d_highest);
145
146 // track ------------------------------------------
147 // parameters
148 tree->SetBranchAddress("track_dr", &track_dr);
149 tree->SetBranchAddress("track_phi0", &track_phi0);
150 tree->SetBranchAddress("track_dz", &track_dz);
151 tree->SetBranchAddress("track_tanL", &track_tanL);
152 tree->SetBranchAddress("track_chi2", &track_chi2);
153 // clusters
154 tree->SetBranchAddress("track_nfitpoint", &track_nfitpoint);
155 tree->SetBranchAddress("track_clusterid", &track_clusterid);
156 tree->SetBranchAddress("track_layerid", &track_layerid);
157 tree->SetBranchAddress("track_sheetid", &track_sheetid);
158 // layer/sheet under test (if none these are -1)
159 tree->SetBranchAddress("track_test_layerid", &track_test_layerid);
160 tree->SetBranchAddress("track_test_sheetid", &track_test_sheetid);
161 // poca
162 tree->SetBranchAddress("track_xpoca_glo", &track_xpoca_glo);
163 tree->SetBranchAddress("track_ypoca_glo", &track_ypoca_glo);
164 tree->SetBranchAddress("track_zpoca_glo", &track_zpoca_glo);
165 // intersections
166 // L1
167 tree->SetBranchAddress("track_x1top_glo", &track_x1top_glo);
168 tree->SetBranchAddress("track_y1top_glo", &track_y1top_glo);
169 tree->SetBranchAddress("track_z1top_glo", &track_z1top_glo);
170 tree->SetBranchAddress("track_phi1top_loc", &track_phi1top_loc);
171 tree->SetBranchAddress("track_v1top_loc", &track_v1top_loc);
172 tree->SetBranchAddress("track_x1bot_glo", &track_x1bot_glo);
173 tree->SetBranchAddress("track_y1bot_glo", &track_y1bot_glo);
174 tree->SetBranchAddress("track_z1bot_glo", &track_z1bot_glo);
175 tree->SetBranchAddress("track_phi1bot_loc", &track_phi1bot_loc);
176 tree->SetBranchAddress("track_v1bot_loc", &track_v1bot_loc);
177 // L2
178 tree->SetBranchAddress("track_x2top_glo", &track_x2top_glo);
179 tree->SetBranchAddress("track_y2top_glo", &track_y2top_glo);
180 tree->SetBranchAddress("track_z2top_glo", &track_z2top_glo);
181 tree->SetBranchAddress("track_phi2top_loc", &track_phi2top_loc);
182 tree->SetBranchAddress("track_v2top_loc", &track_v2top_loc);
183 tree->SetBranchAddress("track_x2bot_glo", &track_x2bot_glo);
184 tree->SetBranchAddress("track_y2bot_glo", &track_y2bot_glo);
185 tree->SetBranchAddress("track_z2bot_glo", &track_z2bot_glo);
186 tree->SetBranchAddress("track_phi2bot_loc", &track_phi2bot_loc);
187 tree->SetBranchAddress("track_v2bot_loc", &track_v2bot_loc);
188
189 // Intersection angles
190 // L1
191 tree->SetBranchAddress("ang_xy_L1", &ang_xy_L1);
192 tree->SetBranchAddress("ang_yz_L1", &ang_yz_L1);
193 // L2
194 tree->SetBranchAddress("ang_xy_L2", &ang_xy_L2);
195 tree->SetBranchAddress("ang_yz_L2", &ang_yz_L2);
196
197 // ------------------------------------------------
198 return true;
199}
200
201void CgemCosmicRayQA::define_output() {
202
203 output = new TFile("histo.root", "RECREATE");
204}
205
206
207void CgemCosmicRayQA::define_hit_histo() {
208 cout << "define_hit_histo" << endl;
209
210 double hit_mintime = -10000;
211 double hit_maxtime = -8000;
212 double hit_maxcharge = 70;
213
214 TString dname = "hit_histo";
215 output->mkdir(dname, "histograms for all the hits");
216 output->cd(dname);
217
218 // nof hits
219 // x
220 h_nofhit_L1_S1_x = new TH1F("h_nofhit_L1_S1_x", "number of hits in L1, S1, x view", MAXNOFHITS, 0, MAXNOFHITS);
221 h_nofhit_L1_S2_x = new TH1F("h_nofhit_L1_S2_x", "number of hits in L1, S2, x view", MAXNOFHITS, 0, MAXNOFHITS);
222 h_nofhit_L2_S1_x = new TH1F("h_nofhit_L2_S1_x", "number of hits in L2, S1, x view", MAXNOFHITS, 0, MAXNOFHITS);
223 h_nofhit_L2_S2_x = new TH1F("h_nofhit_L2_S2_x", "number of hits in L2, S2, x view", MAXNOFHITS, 0, MAXNOFHITS);
224 // v
225 h_nofhit_L1_S1_v = new TH1F("h_nofhit_L1_S1_v", "number of hits in L1, S1, v view", MAXNOFHITS, 0, MAXNOFHITS);
226 h_nofhit_L2_S1_v = new TH1F("h_nofhit_L2_S1_v", "number of hits in L2, S1, v view", MAXNOFHITS, 0, MAXNOFHITS);
227 h_nofhit_L2_S2_v = new TH1F("h_nofhit_L2_S2_v", "number of hits in L2, S2, v view", MAXNOFHITS, 0, MAXNOFHITS);
228
229 // charge
230 // x
231 h_hit_charge_L1_S1_x = new TH1F("h_hit_charge_L1_S1_x", "hit charge (fC) in L1, S1, x view", hit_maxcharge, 0, hit_maxcharge);
232 h_hit_charge_L1_S2_x = new TH1F("h_hit_charge_L1_S2_x", "hit charge (fC) in L1, S2, x view", hit_maxcharge, 0, hit_maxcharge);
233 h_hit_charge_L2_S1_x = new TH1F("h_hit_charge_L2_S1_x", "hit charge (fC) in L2, S1, x view", hit_maxcharge, 0, hit_maxcharge);
234 h_hit_charge_L2_S2_x = new TH1F("h_hit_charge_L2_S2_x", "hit charge (fC) in L2, S2, x view", hit_maxcharge, 0, hit_maxcharge);
235 // v
236 h_hit_charge_L1_S1_v = new TH1F("h_hit_charge_L1_S1_v", "hit charge (fC) in L1, S1, v view", hit_maxcharge, 0, hit_maxcharge);
237 h_hit_charge_L2_S1_v = new TH1F("h_hit_charge_L2_S1_v", "hit charge (fC) in L2, S1, v view", hit_maxcharge, 0, hit_maxcharge);
238 h_hit_charge_L2_S2_v = new TH1F("h_hit_charge_L2_S2_v", "hit charge (fC) in L2, S2, v view", hit_maxcharge, 0, hit_maxcharge);
239
240 // time
241 // x
242 h_hit_time_L1_S1_x = new TH1F("h_hit_time_L1_S1_x", "hit time (ns) in L1, S1, x view", 100, hit_mintime, hit_maxtime);
243 h_hit_time_L1_S2_x = new TH1F("h_hit_time_L1_S2_x", "hit time (ns) in L1, S2, x view", 100, hit_mintime, hit_maxtime);
244 h_hit_time_L2_S1_x = new TH1F("h_hit_time_L2_S1_x", "hit time (ns) in L2, S1, x view", 100, hit_mintime, hit_maxtime);
245 h_hit_time_L2_S2_x = new TH1F("h_hit_time_L2_S2_x", "hit time (ns) in L2, S2, x view", 100, hit_mintime, hit_maxtime);
246 // v
247 h_hit_time_L1_S1_v = new TH1F("h_hit_time_L1_S1_v", "hit time (ns) in L1, S1, v view", 100, hit_mintime, hit_maxtime);
248 h_hit_time_L2_S1_v = new TH1F("h_hit_time_L2_S1_v", "hit time (ns) in L2, S1, v view", 100, hit_mintime, hit_maxtime);
249 h_hit_time_L2_S2_v = new TH1F("h_hit_time_L2_S2_v", "hit time (ns) in L2, S2, v view", 100, hit_mintime, hit_maxtime);
250
251 // charge vs strip ID
252 // x
253 h_hit_charge_vs_strip_L1_S1_x = new TH2F("h_hit_charge_vs_strip_L1_S1_x", "hit charge (fC) vs stripID in L1, S1, x view", MAXNOFSTRIP_L1_x, 0, MAXNOFSTRIP_L1_x, hit_maxcharge, 0, hit_maxcharge);
254 h_hit_charge_vs_strip_L1_S2_x = new TH2F("h_hit_charge_vs_strip_L1_S2_x", "hit charge (fC) vs stripID in L1, S2, x view", MAXNOFSTRIP_L1_x, 0, MAXNOFSTRIP_L1_x, hit_maxcharge, 0, hit_maxcharge);
255 h_hit_charge_vs_strip_L2_S1_x = new TH2F("h_hit_charge_vs_strip_L2_S1_x", "hit charge (fC) vs stripID in L2, S1, x view", MAXNOFSTRIP_L2_x, 0, MAXNOFSTRIP_L2_x, hit_maxcharge, 0, hit_maxcharge);
256 h_hit_charge_vs_strip_L2_S2_x = new TH2F("h_hit_charge_vs_strip_L2_S2_x", "hit charge (fC) vs stripID in L2, S2, x view", MAXNOFSTRIP_L2_x, 0, MAXNOFSTRIP_L2_x, hit_maxcharge, 0, hit_maxcharge);
257 // v
258 h_hit_charge_vs_strip_L1_S1_v = new TH2F("h_hit_charge_vs_strip_L1_S1_v", "hit charge (fC) vs stripID in L1, S1, v view", MAXNOFSTRIP_L1_v, 0, MAXNOFSTRIP_L1_v, hit_maxcharge, 0, hit_maxcharge);
259 h_hit_charge_vs_strip_L2_S1_v = new TH2F("h_hit_charge_vs_strip_L2_S1_v", "hit charge (fC) vs stripID in L2, S1, v view", MAXNOFSTRIP_L2_v, 0, MAXNOFSTRIP_L2_v, hit_maxcharge, 0, hit_maxcharge);
260 h_hit_charge_vs_strip_L2_S2_v = new TH2F("h_hit_charge_vs_strip_L2_S2_v", "hit charge (fC) vs stripID in L2, S2, v view", MAXNOFSTRIP_L2_v, 0, MAXNOFSTRIP_L2_v, hit_maxcharge, 0, hit_maxcharge);
261
262 // time vs strip ID
263 // x
264 h_hit_time_vs_strip_L1_S1_x = new TH2F("h_hit_time_vs_strip_L1_S1_x", "hit time (ns) vs stripID in L1, S1, x view", MAXNOFSTRIP_L1_x, 0, MAXNOFSTRIP_L1_x, 100, hit_mintime, hit_maxtime);
265 h_hit_time_vs_strip_L1_S2_x = new TH2F("h_hit_time_vs_strip_L1_S2_x", "hit time (ns) vs stripID in L1, S2, x view", MAXNOFSTRIP_L1_x, 0, MAXNOFSTRIP_L1_x, 100, hit_mintime, hit_maxtime);
266 h_hit_time_vs_strip_L2_S1_x = new TH2F("h_hit_time_vs_strip_L2_S1_x", "hit time (ns) vs stripID in L2, S1, x view", MAXNOFSTRIP_L2_x, 0, MAXNOFSTRIP_L2_x, 100, hit_mintime, hit_maxtime);
267 h_hit_time_vs_strip_L2_S2_x = new TH2F("h_hit_time_vs_strip_L2_S2_x", "hit time (ns) vs stripID in L2, S2, x view", MAXNOFSTRIP_L2_x, 0, MAXNOFSTRIP_L2_x, 100, hit_mintime, hit_maxtime);
268 // v
269 h_hit_time_vs_strip_L1_S1_v = new TH2F("h_hit_time_vs_strip_L1_S1_v", "hit time (ns) vs stripID in L1, S1, v view", MAXNOFSTRIP_L1_v, 0, MAXNOFSTRIP_L1_v, 100, hit_mintime, hit_maxtime);
270 h_hit_time_vs_strip_L2_S1_v = new TH2F("h_hit_time_vs_strip_L2_S1_v", "hit time (ns) vs stripID in L2, S1, v view", MAXNOFSTRIP_L2_v, 0, MAXNOFSTRIP_L2_v, 100, hit_mintime, hit_maxtime);
271 h_hit_time_vs_strip_L2_S2_v = new TH2F("h_hit_time_vs_strip_L2_S2_v", "hit time (ns) vs stripID in L2, S2, v view", MAXNOFSTRIP_L2_v, 0, MAXNOFSTRIP_L2_v, 100, hit_mintime, hit_maxtime);
272
273 // charge vs time
274 // x
275 h_hit_charge_vs_time_L1_S1_x = new TH2F("h_hit_charge_vs_time_L1_S1_x", "hit charge (fC) vs time (ns) in L1, S1, x view", 100, hit_mintime, hit_maxtime, hit_maxcharge, 0, hit_maxcharge);
276 h_hit_charge_vs_time_L1_S2_x = new TH2F("h_hit_charge_vs_time_L1_S2_x", "hit charge (fC) vs time (ns) in L1, S2, x view", 100, hit_mintime, hit_maxtime, hit_maxcharge, 0, hit_maxcharge);
277 h_hit_charge_vs_time_L2_S1_x = new TH2F("h_hit_charge_vs_time_L2_S1_x", "hit charge (fC) vs time (ns) in L2, S1, x view", 100, hit_mintime, hit_maxtime, hit_maxcharge, 0, hit_maxcharge);
278 h_hit_charge_vs_time_L2_S2_x = new TH2F("h_hit_charge_vs_time_L2_S2_x", "hit charge (fC) vs time (ns) in L2, S2, x view", 100, hit_mintime, hit_maxtime, hit_maxcharge, 0, hit_maxcharge);
279 // v
280 h_hit_charge_vs_time_L1_S1_v = new TH2F("h_hit_charge_vs_time_L1_S1_v", "hit charge (fC) vs time (ns) in L1, S1, v view", 100, hit_mintime, hit_maxtime, hit_maxcharge, 0, hit_maxcharge);
281 h_hit_charge_vs_time_L2_S1_v = new TH2F("h_hit_charge_vs_time_L2_S1_v", "hit charge (fC) vs time (ns) in L2, S1, v view", 100, hit_mintime, hit_maxtime, hit_maxcharge, 0, hit_maxcharge);
282 h_hit_charge_vs_time_L2_S2_v = new TH2F("h_hit_charge_vs_time_L2_S2_v", "hit charge (fC) vs time (ns) in L2, S2, v view", 100, hit_mintime, hit_maxtime, hit_maxcharge, 0, hit_maxcharge);
283
284 // charge vs length
285 // v
286 h_hit_charge_vs_length_L1_S1_v = new TH2F("h_hit_charge_vs_length_L1_S1_v", "hit charge (fC) vs length (mm) in L1, S1, v view", 500, 0, MAXLENGTH_L1_v, hit_maxcharge, 0, hit_maxcharge);
287 h_hit_charge_vs_length_L2_S1_v = new TH2F("h_hit_charge_vs_length_L2_S1_v", "hit charge (fC) vs length (mm) in L2, S1, v view", 500, 0, MAXLENGTH_L2_v, hit_maxcharge, 0, hit_maxcharge);
288 h_hit_charge_vs_length_L2_S2_v = new TH2F("h_hit_charge_vs_length_L2_S2_v", "hit charge (fC) vs length (mm) in L2, S2, v view", 500, 0, MAXLENGTH_L2_v, hit_maxcharge, 0, hit_maxcharge);
289
290 output->cd();
291}
292
293void CgemCosmicRayQA::define_hit_in_time_histo() {
294 cout << "define_hit_in_time_histo" << endl;
295
296 double hit_mintime = -10000;
297 double hit_maxtime = -8000;
298 double hit_maxcharge = 70;
299
300 TString dname = "hit_in_time_histo";
301 output-> mkdir(dname, "histograms for the hits in time window");
302 output->cd(dname);
303
304 // charge in time
305 // x
306 h_hit_in_time_charge_L1_S1_x = new TH1F("h_hit_in_time_charge_L1_S1_x", "(in time) hit charge (fC) in L1, S1, x view", hit_maxcharge, 0, hit_maxcharge);
307 h_hit_in_time_charge_L1_S2_x = new TH1F("h_hit_in_time_charge_L1_S2_x", "(in time) hit charge (fC) in L1, S2, x view", hit_maxcharge, 0, hit_maxcharge);
308 h_hit_in_time_charge_L2_S1_x = new TH1F("h_hit_in_time_charge_L2_S1_x", "(in time) hit charge (fC) in L2, S1, x view", hit_maxcharge, 0, hit_maxcharge);
309 h_hit_in_time_charge_L2_S2_x = new TH1F("h_hit_in_time_charge_L2_S2_x", "(in time) hit charge (fC) in L2, S2, x view", hit_maxcharge, 0, hit_maxcharge);
310 // v
311 h_hit_in_time_charge_L1_S1_v = new TH1F("h_hit_in_time_charge_L1_S1_v", "(in time) hit charge (fC) in L1, S1, v view", hit_maxcharge, 0, hit_maxcharge);
312 h_hit_in_time_charge_L2_S1_v = new TH1F("h_hit_in_time_charge_L2_S1_v", "(in time) hit charge (fC) in L2, S1, v view", hit_maxcharge, 0, hit_maxcharge);
313 h_hit_in_time_charge_L2_S2_v = new TH1F("h_hit_in_time_charge_L2_S2_v", "(in time) hit charge (fC) in L2, S2, v view", hit_maxcharge, 0, hit_maxcharge);
314
315 // count on each strip
316 // x
317 h_hit_in_time_strip_L1_S1_x = new TH1F("h_hit_in_time_strip_L1_S1_x", "(in time) hit counts vs stripID in L1, S1, x view", MAXNOFSTRIP_L1_x, 0, MAXNOFSTRIP_L1_x);
318 h_hit_in_time_strip_L1_S2_x = new TH1F("h_hit_in_time_strip_L1_S2_x", "(in time) hit counts vs stripID in L1, S2, x view", MAXNOFSTRIP_L1_x, 0, MAXNOFSTRIP_L1_x);
319 h_hit_in_time_strip_L2_S1_x = new TH1F("h_hit_in_time_strip_L2_S1_x", "(in time) hit counts vs stripID in L2, S1, x view", MAXNOFSTRIP_L2_x, 0, MAXNOFSTRIP_L2_x);
320 h_hit_in_time_strip_L2_S2_x = new TH1F("h_hit_in_time_strip_L2_S2_x", "(in time) hit counts vs stripID in L2, S2, x view", MAXNOFSTRIP_L2_x, 0, MAXNOFSTRIP_L2_x);
321 // v
322 h_hit_in_time_strip_L1_S1_v = new TH1F("h_hit_in_time_strip_L1_S1_v", "(in time) hit counts vs stripID in L1, S1, v view", MAXNOFSTRIP_L1_v, 0, MAXNOFSTRIP_L1_v);
323 h_hit_in_time_strip_L2_S1_v = new TH1F("h_hit_in_time_strip_L2_S1_v", "(in time) hit counts vs stripID in L2, S1, v view", MAXNOFSTRIP_L2_v, 0, MAXNOFSTRIP_L2_v);
324 h_hit_in_time_strip_L2_S2_v = new TH1F("h_hit_in_time_strip_L2_S2_v", "(in time) hit counts vs stripID in L2, S2, v view", MAXNOFSTRIP_L2_v, 0, MAXNOFSTRIP_L2_v);
325
326 // charge vs strip ID
327 // x
328 h_hit_in_time_charge_vs_strip_L1_S1_x = new TH2F("h_hit_in_time_charge_vs_strip_L1_S1_x", "(in time) hit charge (fC) vs stripID in L1, S1, x view", MAXNOFSTRIP_L1_x, 0, MAXNOFSTRIP_L1_x, hit_maxcharge, 0, hit_maxcharge);
329 h_hit_in_time_charge_vs_strip_L1_S2_x = new TH2F("h_hit_in_time_charge_vs_strip_L1_S2_x", "(in time) hit charge (fC) vs stripID in L1, S2, x view", MAXNOFSTRIP_L1_x, 0, MAXNOFSTRIP_L1_x, hit_maxcharge, 0, hit_maxcharge);
330 h_hit_in_time_charge_vs_strip_L2_S1_x = new TH2F("h_hit_in_time_charge_vs_strip_L2_S1_x", "(in time) hit charge (fC) vs stripID in L2, S1, x view", MAXNOFSTRIP_L2_x, 0, MAXNOFSTRIP_L2_x, hit_maxcharge, 0, hit_maxcharge);
331 h_hit_in_time_charge_vs_strip_L2_S2_x = new TH2F("h_hit_in_time_charge_vs_strip_L2_S2_x", "(in time) hit charge (fC) vs stripID in L2, S2, x view", MAXNOFSTRIP_L2_x, 0, MAXNOFSTRIP_L2_x, hit_maxcharge, 0, hit_maxcharge);
332 // v
333 h_hit_in_time_charge_vs_strip_L1_S1_v = new TH2F("h_hit_in_time_charge_vs_strip_L1_S1_v", "(in time) hit charge (fC) vs stripID in L1, S1, v view", MAXNOFSTRIP_L1_v, 0, MAXNOFSTRIP_L1_v, hit_maxcharge, 0, hit_maxcharge);
334 h_hit_in_time_charge_vs_strip_L2_S1_v = new TH2F("h_hit_in_time_charge_vs_strip_L2_S1_v", "(in time) hit charge (fC) vs stripID in L2, S1, v view", MAXNOFSTRIP_L2_v, 0, MAXNOFSTRIP_L2_v, hit_maxcharge, 0, hit_maxcharge);
335 h_hit_in_time_charge_vs_strip_L2_S2_v = new TH2F("h_hit_in_time_charge_vs_strip_L2_S2_v", "(in time) hit charge (fC) vs stripID in L2, S2, v view", MAXNOFSTRIP_L2_v, 0, MAXNOFSTRIP_L2_v, hit_maxcharge, 0, hit_maxcharge);
336
337 output->cd();
338
339}
340
341void CgemCosmicRayQA::define_cluster1d_histo() {
342 cout << "define_cluster1d_histo" << endl;
343
344 double cluster1d_max_charge = 250;
345 double cluster1d_max_size = 40;
346
347 TString dname = "cluster1d_histo";
348 output-> mkdir(dname, "histograms of all cluster 1D");
349 output->cd(dname);
350
351 // nof cluster 1d
352 // x
353 h_nofcluster1d_L1_S1_x = new TH1F("h_nofcluster1d_L1_S1_x", "number of cluster1d in L1, S1, x view", MAXNOFCLUSTERS, 0, MAXNOFCLUSTERS);
354 h_nofcluster1d_L1_S2_x = new TH1F("h_nofcluster1d_L1_S2_x", "number of cluster1d in L1, S2, x view", MAXNOFCLUSTERS, 0, MAXNOFCLUSTERS);
355 h_nofcluster1d_L2_S1_x = new TH1F("h_nofcluster1d_L2_S1_x", "number of cluster1d in L2, S1, x view", MAXNOFCLUSTERS, 0, MAXNOFCLUSTERS);
356 h_nofcluster1d_L2_S2_x = new TH1F("h_nofcluster1d_L2_S2_x", "number of cluster1d in L2, S2, x view", MAXNOFCLUSTERS, 0, MAXNOFCLUSTERS);
357 // v
358 h_nofcluster1d_L1_S1_v = new TH1F("h_nofcluster1d_L1_S1_v", "number of cluster1d in L1, S1, v view", MAXNOFCLUSTERS, 0, MAXNOFCLUSTERS);
359 h_nofcluster1d_L2_S1_v = new TH1F("h_nofcluster1d_L2_S1_v", "number of cluster1d in L2, S1, v view", MAXNOFCLUSTERS, 0, MAXNOFCLUSTERS);
360 h_nofcluster1d_L2_S2_v = new TH1F("h_nofcluster1d_L2_S2_v", "number of cluster1d in L2, S2, v view", MAXNOFCLUSTERS, 0, MAXNOFCLUSTERS);
361
362 // cl. size
363 // x
364 h_cluster1d_size_L1_S1_x = new TH1F("h_cluster1d_size_L1_S1_x", "cluster1d size in L1, S1, x view", cluster1d_max_size, 0, cluster1d_max_size);
365 h_cluster1d_size_L1_S2_x = new TH1F("h_cluster1d_size_L1_S2_x", "cluster1d size in L1, S2, x view", cluster1d_max_size, 0, cluster1d_max_size);
366 h_cluster1d_size_L2_S1_x = new TH1F("h_cluster1d_size_L2_S1_x", "cluster1d size in L2, S1, x view", cluster1d_max_size, 0, cluster1d_max_size);
367 h_cluster1d_size_L2_S2_x = new TH1F("h_cluster1d_size_L2_S2_x", "cluster1d size in L2, S2, x view", cluster1d_max_size, 0, cluster1d_max_size);
368 // v
369 h_cluster1d_size_L1_S1_v = new TH1F("h_cluster1d_size_L1_S1_v", "cluster1d size in L1, S1, v view", cluster1d_max_size, 0, cluster1d_max_size);
370 h_cluster1d_size_L2_S1_v = new TH1F("h_cluster1d_size_L2_S1_v", "cluster1d size in L2, S1, v view", cluster1d_max_size, 0, cluster1d_max_size);
371 h_cluster1d_size_L2_S2_v = new TH1F("h_cluster1d_size_L2_S2_v", "cluster1d size in L2, S2, v view", cluster1d_max_size, 0, cluster1d_max_size);
372
373 // charge vs phi (x view)
374 h_cluster1d_charge_vs_phi_L1_S1_x = new TH2F("h_cluster1d_charge_vs_phi_L1_S1_x", "cluster1d charge (fC) vs phi (deg) in L1, S1, x view", 180, -180, 0, cluster1d_max_charge, 0, cluster1d_max_charge);
375 h_cluster1d_charge_vs_phi_L1_S2_x = new TH2F("h_cluster1d_charge_vs_phi_L1_S2_x", "cluster1d charge (fC) vs phi (deg) in L1, S2, x view", 180, 0, 180, cluster1d_max_charge, 0, cluster1d_max_charge);
376 h_cluster1d_charge_vs_phi_L2_S1_x = new TH2F("h_cluster1d_charge_vs_phi_L2_S1_x", "cluster1d charge (fC) vs phi (deg) in L2, S1, x view", 180, -180, 0, cluster1d_max_charge, 0, cluster1d_max_charge);
377 h_cluster1d_charge_vs_phi_L2_S2_x = new TH2F("h_cluster1d_charge_vs_phi_L2_S2_x", "cluster1d charge (fC) vs phi (deg) in L2, S2, x view", 180, 0, 180, cluster1d_max_charge, 0, cluster1d_max_charge);
378
379 // charge vs v (v view)
380 h_cluster1d_charge_vs_v_L1_S1_v = new TH2F("h_cluster1d_charge_vs_v_L1_S1_v", "cluster1d charge (fC) vs v (mm) in L1, S1, v view", 0.7*MAXNOFSTRIP_L1_v, 0, 0.7*MAXNOFSTRIP_L1_v, cluster1d_max_charge, 0, cluster1d_max_charge);
381 h_cluster1d_charge_vs_v_L2_S1_v = new TH2F("h_cluster1d_charge_vs_v_L2_S1_v", "cluster1d charge (fC) vs v (mm) in L2, S1, v view", 0.7*MAXNOFSTRIP_L2_v, 0, 0.7*MAXNOFSTRIP_L2_v, cluster1d_max_charge, 0, cluster1d_max_charge);
382 h_cluster1d_charge_vs_v_L2_S2_v = new TH2F("h_cluster1d_charge_vs_v_L2_S2_v", "cluster1d charge (fC) vs v (mm) in L2, S2, v view", 0.7*MAXNOFSTRIP_L2_v, 0, 0.7*MAXNOFSTRIP_L2_v, cluster1d_max_charge, 0, cluster1d_max_charge);
383
384 output->cd();
385
386}
387
388void CgemCosmicRayQA::define_cluster_selected_histo() {
389 cout << "define_cluster_selected_histo" << endl;
390
391 double cluster1d_max_charge = 250;
392 double cluster1d_max_size = 40;
393
394 TString dname = "cluster_selected_histo";
395 output->mkdir(dname, "histograms of cluster 1D and 2D from selected clusters");
396 output->cd(dname);
397
398 // charge of the selected clusters
399 h_cluster1d_charge_selected_L1_S1_x = new TH1F("h_cluster1d_charge_selected_L1_S1_x", "cluster1d charge (fC) in L1, S1, x view (only selected clusters)", cluster1d_max_charge, 0, cluster1d_max_charge);
400 h_cluster1d_charge_selected_L1_S2_x = new TH1F("h_cluster1d_charge_selected_L1_S2_x", "cluster1d charge (fC) in L1, S2, x view (only selected clusters)", cluster1d_max_charge, 0, cluster1d_max_charge);
401 h_cluster1d_charge_selected_L2_S1_x = new TH1F("h_cluster1d_charge_selected_L2_S1_x", "cluster1d charge (fC) in L2, S1, x view (only selected clusters)", cluster1d_max_charge, 0, cluster1d_max_charge);
402 h_cluster1d_charge_selected_L2_S2_x = new TH1F("h_cluster1d_charge_selected_L2_S2_x", "cluster1d charge (fC) in L2, S2, x view (only selected clusters)", cluster1d_max_charge, 0, cluster1d_max_charge);
403
404 h_cluster1d_charge_selected_L1_S1_v = new TH1F("h_cluster1d_charge_selected_L1_S1_v", "cluster1d charge (fC) in L1, S1, v view (only selected clusters)", cluster1d_max_charge, 0, cluster1d_max_charge);
405 h_cluster1d_charge_selected_L1_S2_v = new TH1F("h_cluster1d_charge_selected_L1_S2_v", "cluster1d charge (fC) in L1, S2, v view (only selected clusters)", cluster1d_max_charge, 0, cluster1d_max_charge);
406 h_cluster1d_charge_selected_L2_S1_v = new TH1F("h_cluster1d_charge_selected_L2_S1_v", "cluster1d charge (fC) in L2, S1, v view (only selected clusters)", cluster1d_max_charge, 0, cluster1d_max_charge);
407 h_cluster1d_charge_selected_L2_S2_v = new TH1F("h_cluster1d_charge_selected_L2_S2_v", "cluster1d charge (fC) in L2, S2, v view (only selected clusters)", cluster1d_max_charge, 0, cluster1d_max_charge);
408
409 // cl. size
410 // x
411 h_cluster1d_size_selected_L1_S1_x = new TH1F("h_cluster1d_size_selected_L1_S1_x", "cluster1d size in L1, S1, x view (only selected clusters)", cluster1d_max_size, 0, cluster1d_max_size);
412 h_cluster1d_size_selected_L1_S2_x = new TH1F("h_cluster1d_size_selected_L1_S2_x", "cluster1d size in L1, S2, x view (only selected clusters)", cluster1d_max_size, 0, cluster1d_max_size);
413 h_cluster1d_size_selected_L2_S1_x = new TH1F("h_cluster1d_size_selected_L2_S1_x", "cluster1d size in L2, S1, x view (only selected clusters)", cluster1d_max_size, 0, cluster1d_max_size);
414 h_cluster1d_size_selected_L2_S2_x = new TH1F("h_cluster1d_size_selected_L2_S2_x", "cluster1d size in L2, S2, x view (only selected clusters)", cluster1d_max_size, 0, cluster1d_max_size);
415 // v
416 h_cluster1d_size_selected_L1_S1_v = new TH1F("h_cluster1d_size_selected_L1_S1_v", "cluster1d size in L1, S1, v view (only selected clusters)", cluster1d_max_size, 0, cluster1d_max_size);
417 h_cluster1d_size_selected_L1_S2_v = new TH1F("h_cluster1d_size_selected_L1_S2_v", "cluster1d size in L1, S2, v view (only selected clusters)", cluster1d_max_size, 0, cluster1d_max_size);
418 h_cluster1d_size_selected_L2_S1_v = new TH1F("h_cluster1d_size_selected_L2_S1_v", "cluster1d size in L2, S1, v view (only selected clusters)", cluster1d_max_size, 0, cluster1d_max_size);
419 h_cluster1d_size_selected_L2_S2_v = new TH1F("h_cluster1d_size_selected_L2_S2_v", "cluster1d size in L2, S2, v view (only selected clusters)", cluster1d_max_size, 0, cluster1d_max_size);
420
421 double cluster2d_max_charge = 500;
422 double cluster2d_max_size = 40;
423
424 // charge of the selected clusters 2D
425 h_cluster2d_charge_selected_L1_S1 = new TH1F("h_cluster2d_charge_selected_L1_S1", "cluster2d charge (fC) in L1, S1 (only selected clusters)", cluster2d_max_charge, 0, cluster2d_max_charge);
426 h_cluster2d_charge_selected_L1_S2 = new TH1F("h_cluster2d_charge_selected_L1_S2", "cluster2d charge (fC) in L1, S2 (only selected clusters)", cluster2d_max_charge, 0, cluster2d_max_charge);
427 h_cluster2d_charge_selected_L2_S1 = new TH1F("h_cluster2d_charge_selected_L2_S1", "cluster2d charge (fC) in L2, S1 (only selected clusters)", cluster2d_max_charge, 0, cluster2d_max_charge);
428 h_cluster2d_charge_selected_L2_S2 = new TH1F("h_cluster2d_charge_selected_L2_S2", "cluster2d charge (fC) in L2, S2 (only selected clusters)", cluster2d_max_charge, 0, cluster2d_max_charge);
429
430 output->cd();
431
432}
433
434void CgemCosmicRayQA::define_cluster1d_from_fit_histo() {
435 cout << "define_cluster1d_from_fit_histo" << endl;
436
437 double cluster1d_max_charge = 250;
438 double cluster1d_max_size = 40;
439
440 TString dname = "cluster1d_from_fit_histo";
441 output->mkdir(dname, "histograms of cluster 1D from fitted tracks");
442 output->cd(dname);
443
444 // charge of the used clusters used in the fit
445 h_cluster1d_charge_fitted_L1_S1_x = new TH1F("h_cluster1d_charge_fitted_L1_S1_x", "cluster1d charge (fC) in L1, S1, x view (only clusters in the fit)", cluster1d_max_charge, 0, cluster1d_max_charge);
446 h_cluster1d_charge_fitted_L1_S2_x = new TH1F("h_cluster1d_charge_fitted_L1_S2_x", "cluster1d charge (fC) in L1, S2, x view (only clusters in the fit)", cluster1d_max_charge, 0, cluster1d_max_charge);
447 h_cluster1d_charge_fitted_L2_S1_x = new TH1F("h_cluster1d_charge_fitted_L2_S1_x", "cluster1d charge (fC) in L2, S1, x view (only clusters in the fit)", cluster1d_max_charge, 0, cluster1d_max_charge);
448 h_cluster1d_charge_fitted_L2_S2_x = new TH1F("h_cluster1d_charge_fitted_L2_S2_x", "cluster1d charge (fC) in L2, S2, x view (only clusters in the fit)", cluster1d_max_charge, 0, cluster1d_max_charge);
449
450 h_cluster1d_charge_fitted_L1_S1_v = new TH1F("h_cluster1d_charge_fitted_L1_S1_v", "cluster1d charge (fC) in L1, S1, v view (only clusters in the fit)", cluster1d_max_charge, 0, cluster1d_max_charge);
451 h_cluster1d_charge_fitted_L1_S2_v = new TH1F("h_cluster1d_charge_fitted_L1_S2_v", "cluster1d charge (fC) in L1, S2, v view (only clusters in the fit)", cluster1d_max_charge, 0, cluster1d_max_charge);
452 h_cluster1d_charge_fitted_L2_S1_v = new TH1F("h_cluster1d_charge_fitted_L2_S1_v", "cluster1d charge (fC) in L2, S1, v view (only clusters in the fit)", cluster1d_max_charge, 0, cluster1d_max_charge);
453 h_cluster1d_charge_fitted_L2_S2_v = new TH1F("h_cluster1d_charge_fitted_L2_S2_v", "cluster1d charge (fC) in L2, S2, v view (only clusters in the fit)", cluster1d_max_charge, 0, cluster1d_max_charge);
454
455 output->cd();
456
457}
458
459void CgemCosmicRayQA::define_cluster2d_histo() {
460 cout << "define_cluster2d_histo" << endl;
461
462 double cluster2d_max_charge = 500;
463 double cluster2d_max_size = 40;
464
465 TString dname = "cluster2d_histo";
466 output-> mkdir(dname, "histograms of all cluster 2D");
467 output->cd(dname);
468
469 // nof cluster 2d
470 h_nofcluster2d_L1_S1 = new TH1F("h_nofcluster2d_L1_S1", "number of cluster2d in L1, S1", MAXNOFCLUSTERS, 0, MAXNOFCLUSTERS);
471 h_nofcluster2d_L1_S2 = new TH1F("h_nofcluster2d_L1_S2", "number of cluster2d in L1, S2", MAXNOFCLUSTERS, 0, MAXNOFCLUSTERS);
472 h_nofcluster2d_L2_S1 = new TH1F("h_nofcluster2d_L2_S1", "number of cluster2d in L2, S1", MAXNOFCLUSTERS, 0, MAXNOFCLUSTERS);
473 h_nofcluster2d_L2_S2 = new TH1F("h_nofcluster2d_L2_S2", "number of cluster2d in L2, S2", MAXNOFCLUSTERS, 0, MAXNOFCLUSTERS);
474
475 // cl. size CHECK
476 // h_cluster2d_size_L1_S1 = new TH1F("h_cluster2d_size_L1_S1", "cluster2d size in L1, S1", 40, 0, 40);
477 // h_cluster2d_size_L1_S2 = new TH1F("h_cluster2d_size_L1_S2", "cluster2d size in L1, S2", 40, 0, 40);
478 // h_cluster2d_size_L2_S1 = new TH1F("h_cluster2d_size_L2_S1", "cluster2d size in L2, S1", 40, 0, 40);
479 // h_cluster2d_size_L2_S2 = new TH1F("h_cluster2d_size_L2_S2", "cluster2d size in L2, S2", 40, 0, 40);
480
481 // charge vs phi
482 h_cluster2d_charge_vs_phi_L1_S1 = new TH2F("h_cluster2d_charge_vs_phi_L1_S1", "cluster2d charge (fC) vs phi (deg) in L1, S1", 180, -180, 0, cluster2d_max_charge, 0, cluster2d_max_charge);
483 h_cluster2d_charge_vs_phi_L1_S2 = new TH2F("h_cluster2d_charge_vs_phi_L1_S2", "cluster2d charge (fC) vs phi (deg) in L1, S2", 180, 0, 180, cluster2d_max_charge, 0, cluster2d_max_charge);
484 h_cluster2d_charge_vs_phi_L2_S1 = new TH2F("h_cluster2d_charge_vs_phi_L2_S1", "cluster2d charge (fC) vs phi (deg) in L2, S1", 180, -180, 0, cluster2d_max_charge, 0, cluster2d_max_charge);
485 h_cluster2d_charge_vs_phi_L2_S2 = new TH2F("h_cluster2d_charge_vs_phi_L2_S2", "cluster2d charge (fC) vs phi (deg) in L2, S2", 180, 0, 180, cluster2d_max_charge, 0, cluster2d_max_charge);
486
487 // charge vs z
488 h_cluster2d_charge_vs_z_L1_S1 = new TH2F("h_cluster2d_charge_vs_z_L1_S1", "cluster2d charge (fC) vs z (mm) in L1, S1", MAXLENGTH_L1_x, -0.5*MAXLENGTH_L1_x-20, 0.5*MAXLENGTH_L1_x+20, cluster2d_max_charge, 0, cluster2d_max_charge);
489 h_cluster2d_charge_vs_z_L1_S2 = new TH2F("h_cluster2d_charge_vs_z_L1_S2", "cluster2d charge (fC) vs z (mm) in L1, S2", MAXLENGTH_L1_x, -0.5*MAXLENGTH_L1_x-20, 0.5*MAXLENGTH_L1_x+20, cluster2d_max_charge, 0, cluster2d_max_charge);
490 h_cluster2d_charge_vs_z_L2_S1 = new TH2F("h_cluster2d_charge_vs_z_L2_S1", "cluster2d charge (fC) vs z (mm) in L2, S1", MAXLENGTH_L2_x, -0.5*MAXLENGTH_L2_x-20, 0.5*MAXLENGTH_L2_x+20, cluster2d_max_charge, 0, cluster2d_max_charge);
491 h_cluster2d_charge_vs_z_L2_S2 = new TH2F("h_cluster2d_charge_vs_z_L2_S2", "cluster2d charge (fC) vs z (mm) in L2, S2", MAXLENGTH_L2_x, -0.5*MAXLENGTH_L2_x-20, 0.5*MAXLENGTH_L2_x+20, cluster2d_max_charge, 0, cluster2d_max_charge);
492
493 // z vs phi
494 h_cluster2d_z_vs_phi_L1_S1 = new TH2F("h_cluster2d_z_vs_phi_L1_S1", "cluster2d z (mm) vs phi (deg) in L1, S1", 180, -180, 0, MAXLENGTH_L1_x, -0.5*MAXLENGTH_L1_x-20, 0.5*MAXLENGTH_L1_x+20);
495 h_cluster2d_z_vs_phi_L1_S2 = new TH2F("h_cluster2d_z_vs_phi_L1_S2", "cluster2d z (mm) vs phi (deg) in L1, S2", 180, 0, 180, MAXLENGTH_L1_x, -0.5*MAXLENGTH_L1_x-20, 0.5*MAXLENGTH_L1_x+20);
496 h_cluster2d_z_vs_phi_L2_S1 = new TH2F("h_cluster2d_z_vs_phi_L2_S1", "cluster2d z (mm) vs phi (deg) in L2, S1", 180, -180, 0, MAXLENGTH_L2_x, -0.5*MAXLENGTH_L2_x-20, 0.5*MAXLENGTH_L2_x+20);
497 h_cluster2d_z_vs_phi_L2_S2 = new TH2F("h_cluster2d_z_vs_phi_L2_S2", "cluster2d z (mm) vs phi (deg) in L2, S2", 180, 0, 180, MAXLENGTH_L2_x, -0.5*MAXLENGTH_L2_x-20, 0.5*MAXLENGTH_L2_x+20);
498
499 output->cd();
500
501}
502
503
504void CgemCosmicRayQA::define_track_histo() {
505 cout << "define_track_histo" << endl;
506
507 TString dname = "track_histo";
508 output-> mkdir(dname, "histograms of fitted tracks");
509 output->cd(dname);
510
511 // nof fitted tracks
512 h_noftrack = new TH1F("h_noftrack", "number of fitted tracks", MAXNOFTRACKS, 0, MAXNOFTRACKS);
513
514 // chi2 distro
515 h_track_chi2 = new TH1F("h_track_chi2", "fitted track chi2", 8000, 0, 8000);
516
517 // xpoca, ypoca, zpoca distro
518 h_track_xpoca = new TH1F("h_track_xpoca", "fitted track point of closest approach x (mm)", 100, -100, 100);
519 h_track_ypoca = new TH1F("h_track_ypoca", "fitted track point of closest approach y (mm)", 100, -100, 100);
520 h_track_zpoca = new TH1F("h_track_zpoca", "fitted track point of closest approach z (mm)", 100, -150, 150);
521
522 // TEST PLANE
523 // residual in rphi
524 h_test_residual_rphi = new TH1F("h_test_residual_rphi", "test plane: residual in R * phi (mm)", 100, -5, 5);
525 h_test_res_rphi_vs_chi2 = new TH2F("h_test_res_rphi_vs_chi2", "test plane: residual in R * phi (mm) vs chi2", 100, -5, 5, 1000, 0, 1000);
526
527 // TPC residual in rphi
528 h_test_residual_rphi_tpc = new TH1F("h_test_residual_rphi_tpc", "test plane: TPC residual in R * phi (mm)", 100, -5, 5);
529 h_test_res_rphi_vs_chi2_tpc = new TH2F("h_test_res_rphi_vs_chi2_tpc", "test plane: TPC residual in R * phi (mm) vs chi2", 100, -5, 5, 1000, 0, 1000);
530
531 // residual in z
532 h_test_residual_z = new TH1F("h_test_residual_z", "test plane: residual in z (mm)", 100, -5, 5);
533 h_test_res_z_vs_chi2 = new TH2F("h_test_res_z_vs_chi2", "test plane: residual in z (mm) vs chi2", 100, -5, 5, 1000, 0, 1000);
534
535 // TPC residual in z
536 h_test_residual_z_tpc = new TH1F("h_test_residual_z_tpc", "test plane: TPC residual in z (mm)", 100, -5, 5);
537 h_test_res_z_vs_chi2_tpc = new TH2F("h_test_res_z_vs_chi2_tpc", "test plane: TPC residual in z (mm) vs chi2", 100, -5, 5, 1000, 0, 1000);
538
539 // TRACKER
540 // residual in phi
541 h_residual_rphi_L1_S1 = new TH1F("h_residual_rphi_L1_S1", "residual in R * phi (mm) on L1, S1", 100, -5, 5);
542 h_residual_rphi_L1_S2 = new TH1F("h_residual_rphi_L1_S2", "residual in R * phi (mm) on L1, S2", 100, -5, 5);
543 h_residual_rphi_L2_S1 = new TH1F("h_residual_rphi_L2_S1", "residual in R * phi (mm) on L2, S1", 100, -5, 5);
544 h_residual_rphi_L2_S2 = new TH1F("h_residual_rphi_L2_S2", "residual in R * phi (mm) on L2, S2", 100, -5, 5);
545
546 h_res_rphi_vs_chi2_L1_S1 = new TH2F("h_res_rphi_vs_chi2_L1_S1", "residual in R * phi (mm) on L1, S1 vs chi2", 100, -5, 5, 1000, 0, 1000);
547 h_res_rphi_vs_chi2_L1_S2 = new TH2F("h_res_rphi_vs_chi2_L1_S2", "residual in R * phi (mm) on L1, S2 vs chi2", 100, -5, 5, 1000, 0, 1000);
548 h_res_rphi_vs_chi2_L2_S1 = new TH2F("h_res_rphi_vs_chi2_L2_S1", "residual in R * phi (mm) on L2, S1 vs chi2", 100, -5, 5, 1000, 0, 1000);
549 h_res_rphi_vs_chi2_L2_S2 = new TH2F("h_res_rphi_vs_chi2_L2_S2", "residual in R * phi (mm) on L2, S2 vs chi2", 100, -5, 5, 1000, 0, 1000);
550
551 // TPC residual in phi
552 h_residual_rphi_L1_S1_tpc = new TH1F("h_residual_rphi_L1_S1_tpc", "TPC residual in R * phi (mm) on L1, S1", 100, -5, 5);
553 h_residual_rphi_L1_S2_tpc = new TH1F("h_residual_rphi_L1_S2_tpc", "TPC residual in R * phi (mm) on L1, S2", 100, -5, 5);
554 h_residual_rphi_L2_S1_tpc = new TH1F("h_residual_rphi_L2_S1_tpc", "TPC residual in R * phi (mm) on L2, S1", 100, -5, 5);
555 h_residual_rphi_L2_S2_tpc = new TH1F("h_residual_rphi_L2_S2_tpc", "TPC residual in R * phi (mm) on L2, S2", 100, -5, 5);
556
557 h_res_rphi_vs_chi2_L1_S1_tpc = new TH2F("h_res_rphi_vs_chi2_L1_S1_tpc", "TPC residual in R * phi (mm) on L1, S1 vs chi2", 100, -5, 5, 1000, 0, 1000);
558 h_res_rphi_vs_chi2_L1_S2_tpc = new TH2F("h_res_rphi_vs_chi2_L1_S2_tpc", "TPC residual in R * phi (mm) on L1, S2 vs chi2", 100, -5, 5, 1000, 0, 1000);
559 h_res_rphi_vs_chi2_L2_S1_tpc = new TH2F("h_res_rphi_vs_chi2_L2_S1_tpc", "TPC residual in R * phi (mm) on L2, S1 vs chi2", 100, -5, 5, 1000, 0, 1000);
560 h_res_rphi_vs_chi2_L2_S2_tpc = new TH2F("h_res_rphi_vs_chi2_L2_S2_tpc", "TPC residual in R * phi (mm) on L2, S2 vs chi2", 100, -5, 5, 1000, 0, 1000);
561
562 // residual in z
563 h_residual_z_L1_S1 = new TH1F("h_residual_z_L1_S1", "residual in z (mm) on L1, S1", 100, -5, 5);
564 h_residual_z_L1_S2 = new TH1F("h_residual_z_L1_S2", "residual in z (mm) on L1, S2", 100, -5, 5);
565 h_residual_z_L2_S1 = new TH1F("h_residual_z_L2_S1", "residual in z (mm) on L2, S1", 100, -5, 5);
566 h_residual_z_L2_S2 = new TH1F("h_residual_z_L2_S2", "residual in z (mm) on L2, S2", 100, -5, 5);
567
568 h_res_z_vs_chi2_L1_S1 = new TH2F("h_res_z_vs_chi2_L1_S1", "residual in z (mm) on L1, S1 vs chi2", 100, -5, 5, 1000, 0, 1000);
569 h_res_z_vs_chi2_L1_S2 = new TH2F("h_res_z_vs_chi2_L1_S2", "residual in z (mm) on L1, S2 vs chi2", 100, -5, 5, 1000, 0, 1000);
570 h_res_z_vs_chi2_L2_S1 = new TH2F("h_res_z_vs_chi2_L2_S1", "residual in z (mm) on L2, S1 vs chi2", 100, -5, 5, 1000, 0, 1000);
571 h_res_z_vs_chi2_L2_S2 = new TH2F("h_res_z_vs_chi2_L2_S2", "residual in z (mm) on L2, S2 vs chi2", 100, -5, 5, 1000, 0, 1000);
572
573 // TPC residual in z
574 h_residual_z_L1_S1_tpc = new TH1F("h_residual_z_L1_S1_tpc", "TPC residual in z (mm) on L1, S1", 100, -5, 5);
575 h_residual_z_L1_S2_tpc = new TH1F("h_residual_z_L1_S2_tpc", "TPC residual in z (mm) on L1, S2", 100, -5, 5);
576 h_residual_z_L2_S1_tpc = new TH1F("h_residual_z_L2_S1_tpc", "TPC residual in z (mm) on L2, S1", 100, -5, 5);
577 h_residual_z_L2_S2_tpc = new TH1F("h_residual_z_L2_S2_tpc", "TPC residual in z (mm) on L2, S2", 100, -5, 5);
578
579 h_res_z_vs_chi2_L1_S1_tpc = new TH2F("h_res_z_vs_chi2_L1_S1_tpc", "TPC residual in z (mm) on L1, S1 vs chi2", 100, -5, 5, 1000, 0, 1000);
580 h_res_z_vs_chi2_L1_S2_tpc = new TH2F("h_res_z_vs_chi2_L1_S2_tpc", "TPC residual in z (mm) on L1, S2 vs chi2", 100, -5, 5, 1000, 0, 1000);
581 h_res_z_vs_chi2_L2_S1_tpc = new TH2F("h_res_z_vs_chi2_L2_S1_tpc", "TPC residual in z (mm) on L2, S1 vs chi2", 100, -5, 5, 1000, 0, 1000);
582 h_res_z_vs_chi2_L2_S2_tpc = new TH2F("h_res_z_vs_chi2_L2_S2_tpc", "TPC residual in z (mm) on L2, S2 vs chi2", 100, -5, 5, 1000, 0, 1000);
583
584 // signal
585 h_signal_charge_2d = new TH1F("h_signal_charge_2d", "total charge of the signal (fC)", 500, 0, 500);
586 h_signal_charge_1d_x = new TH1F("h_signal_charge_1d_x", "charge x of the signal (fC)", 500, 0, 500);
587 h_signal_charge_1d_v = new TH1F("h_signal_charge_1d_v", "charge v of the signal (fC)", 500, 0, 500);
588 h_signal_size_1d_x = new TH1F("h_signal_size_1d_x", "cl.size x of the signal", 50, 0, 50);
589 h_signal_size_1d_v = new TH1F("h_signal_size_1d_v", "cl.size v of the signal", 50, 0, 50);
590
591 h_signal_charge_2d_chi2cut = new TH1F("h_signal_charge_2d_chi2cut", "total charge of the signal (fC) with a chi2 cut", 500, 0, 500);
592 h_signal_charge_1d_x_chi2cut = new TH1F("h_signal_charge_1d_x_chi2cut", "charge x of the signal (fC) with a chi2 cut", 500, 0, 500);
593 h_signal_charge_1d_v_chi2cut = new TH1F("h_signal_charge_1d_v_chi2cut", "charge v of the signal (fC) with a chi2 cut", 500, 0, 500);
594 h_signal_size_1d_x_chi2cut = new TH1F("h_signal_size_1d_x_chi2cut", "cl.size x of the signal with a chi2 cut", 50, 0, 50);
595 h_signal_size_1d_v_chi2cut = new TH1F("h_signal_size_1d_v_chi2cut", "cl.size v of the signal with a chi2 cut", 50, 0, 50);
596
597 // background
598 h_background_charge_2d = new TH1F("h_background_charge_2d", "total charge of the background (fC)", 500, 0, 500);
599 h_background_charge_1d_x = new TH1F("h_background_charge_1d_x", "charge x of the background (fC)", 500, 0, 500);
600 h_background_charge_1d_v = new TH1F("h_background_charge_1d_v", "charge v of the background (fC)", 500, 0, 500);
601 h_background_size_1d_x = new TH1F("h_background_size_1d_x", "cl.size x of the background", 50, 0, 50);
602 h_background_size_1d_v = new TH1F("h_background_size_1d_v", "cl.size v of the background", 50, 0, 50);
603
604 // efficiency vs phi
605 h_eff_vs_phi = new TH1F("h_eff_vs_phi", "Track efficiency vs phi (deg)", MAX_PHI_BIN, MIN_PHI, MAX_PHI);
606 h_eff_vs_phi_vs_chi2 = new TH2F("h_eff_vs_phi_vs_chi2", "Track efficiency vs phi (deg) vs chi2", MAX_PHI_BIN, MIN_PHI, MAX_PHI, CHI_BIN, CHI_MIN, CHI_MAX);
607
608 // efficiency vs z
609 h_eff_vs_z = new TH1F("h_eff_vs_z", "Track efficiency vs z (mm)", MAX_Z_BIN, MIN_Z, MAX_Z);
610 h_eff_vs_z_vs_chi2 = new TH2F("h_eff_vs_z_vs_chi2", "Track efficiency vs z (mm) vs chi2", MAX_Z_BIN, MIN_Z, MAX_Z, CHI_BIN, CHI_MIN, CHI_MAX);
611
612 // efficiency vs phi vs z
613 h_eff_vs_phi_vs_z = new TH2F("h_eff_vs_phi_vs_z", "Track efficiency vs phi (deg) vs z (mm)", MAX_PHI_BIN, MIN_PHI, MAX_PHI, MAX_Z_BIN, MIN_Z, MAX_Z);
614
615 // cc resolution in rphi vs chi2
616 h_resolution_rphi_vs_chi2 = new TH1F("h_resolution_rphi_vs_chi2" , "test plane: cc resolution in R * phi (mm) vs chi2" , CHI_BIN, CHI_MIN, CHI_MAX);
617 h_resolution_rphi_tpc_vs_chi2 = new TH1F("h_resolution_rphi_tpc_vs_chi2", "test plane: tpc resolution in R * phi (mm) vs chi2", CHI_BIN, CHI_MIN, CHI_MAX);
618
619 h_resolution_z_vs_chi2 = new TH1F("h_resolution_z_vs_chi2" , "test plane: cc resolution in z (mm) vs chi2" , CHI_BIN, CHI_MIN, CHI_MAX);
620 h_resolution_z_tpc_vs_chi2 = new TH1F("h_resolution_z_tpc_vs_chi2", "test plane: tpc resolution in z (mm) vs chi2", CHI_BIN, CHI_MIN, CHI_MAX);
621
622 // Incident Angles
623
624 // L1 residual
625 h_test_residual_rphi_vs_angxy_L1_cc = new TH2F("h_test_residual_rphi_vs_angxy_L1_cc" , "test plane: CC residual in R * phi (mm) vs ang_{xy}_{L1}", 100, -5, 5, ANG_BIN, ANG_MIN, ANG_MAX);
626 h_test_residual_rphi_vs_angxy_L1_tpc = new TH2F("h_test_residual_rphi_vs_angxy_L1_tpc", "test plane: TPC residual in R * phi (mm) vs ang_{xy}_{L1}", 100, -5, 5, ANG_BIN, ANG_MIN, ANG_MAX);
627
628 // L2 residual
629 h_test_residual_rphi_vs_angxy_L2_cc = new TH2F("h_test_residual_rphi_vs_angxy_L2_cc" , "test plane: CC residual in R * phi (mm) vs ang_{xy}_{L2}", 100, -5, 5, ANG_BIN, ANG_MIN, ANG_MAX);
630 h_test_residual_rphi_vs_angxy_L2_tpc = new TH2F("h_test_residual_rphi_vs_angxy_L2_tpc", "test plane: TPC residual in R * phi (mm) vs ang_{xy}_{L2}", 100, -5, 5, ANG_BIN, ANG_MIN, ANG_MAX);
631
632 // cc resolution
633 h_resolution_vs_ang_xy_L1 = new TH1F("h_resolution_vs_ang_xy_L1", "cc resolution in R * phi (mm) vs L1 ang_{xy}", ANG_BIN, ANG_MIN, ANG_MAX);
634 h_resolution_vs_ang_xy_L2 = new TH1F("h_resolution_vs_ang_xy_L2", "cc resolution in R * phi (mm) vs L2 ang_{xy}", ANG_BIN, ANG_MIN, ANG_MAX);
635
636 // tpc resolution
637 h_resolution_tpc_vs_ang_xy_L1 = new TH1F("h_resolution_tpc_vs_ang_xy_L1", "tpc resolution in R * phi (mm) vs L1 ang_{xy}", ANG_BIN, ANG_MIN, ANG_MAX);
638 h_resolution_tpc_vs_ang_xy_L2 = new TH1F("h_resolution_tpc_vs_ang_xy_L2", "tpc resolution in R * phi (mm) vs L2 ang_{xy}", ANG_BIN, ANG_MIN, ANG_MAX);
639
640 output->cd();
641}
642
643
644StatusCode CgemCosmicRayQA::initialize(){
645 MsgStream log(msgSvc(), name());
646 log << MSG::INFO << "CgemCosmicRayQA initialize()" << endreq;
647 cout << "QA info" << endl;
648 cout << "filename " << filename << endl;
649 cout << "alignment flag " << align_flag << endl;
650 cout << "cut_chi2 " << cut_chi2 << endl;
651
652 // SmartIF<IProperty> propMgr( IID_IProperty, iface );
653 // SmartIF<IAppMgrUI> appMgr( IID_IAppMgrUI, iface );
654
655 // CgemGeomSvc
656 StatusCode sc;
657 sc = service("CgemGeomSvc", m_SvcCgem);
658 if(sc != StatusCode::SUCCESS) {
659 log << MSG::ERROR << "can not use CgemGeomSvc" << endreq;
660 return StatusCode::FAILURE;
661 }
662 middleplane = m_SvcCgem->getMidDriftPtr() ;
663 alignment = m_SvcCgem->getAlignPtr() ;
664
665 // cout << "ALLINEAMENTO L1 " << alignment->getDx(0) << " " << alignment->getDy(0) << " " << alignment->getDz(0) << endl;
666 // cout << "ALLINEAMENTO L2 " << alignment->getDx(1) << " " << alignment->getDy(1) << " " << alignment->getDz(1) << endl;
667
668 anode = m_SvcCgem->getReadoutPlane(0, 0);
669 anode_mid_gap_radius[0] = anode->getMidRAtGap();
670 anode = m_SvcCgem->getReadoutPlane(1, 0);
671 anode_mid_gap_radius[1] = anode->getMidRAtGap();
672 anode = m_SvcCgem->getReadoutPlane(2, 0);
673 anode_mid_gap_radius[2] = anode->getMidRAtGap();
674
675
676 cout << "CgemLineFit alignment " << align_flag << endl;
677 if(align_flag) {
678 for(int ilay=0; ilay<3; ilay++) {
679 cout << "LAYER " << ilay+1 << endl;
680 cout << "midplane radius " << middleplane->getR(ilay) << endl;
681 for(int ishe=0; ishe<2; ishe++) {
682 cout << "shift dx " << alignment->getDx(ilay,ishe) << " dy " << alignment->getDy(ilay,ishe) << " dz " << alignment->getDz(ilay,ishe) << endl;
683 cout << "rotation Rx " << alignment->getRx(ilay,ishe) << " Ry " << alignment->getRy(ilay,ishe) << " Rz " << alignment->getRz(ilay,ishe) << endl;
684 }
685 }
686 }
687
688 // LUT
689 // LUT
690 // cout << "CREATING LUT" << endl;
691 // lutreader = new CgemLUTReader(lutfile);
692 // lutreader->ReadLUT();
693 // lutreader->PrintMap(1,0,0);
694 // lutreader->PrintMap(1,1,0);
695
696 event=0;
697 define_output();
698 // all
699 define_hit_histo();
700 define_cluster1d_histo();
701 define_cluster2d_histo();
702 define_track_histo();
703 // selection
704 define_hit_in_time_histo();
705 define_cluster_selected_histo();
706 define_cluster1d_from_fit_histo();
707
708 nchi2 = 0;
709 nstopped = 0;
710 nfittedtrack = 0;
711 nvalidtrack = 0;
712 ninside = 0;
713 noutside = 0;
714
715 test_l = 0;
716 test_ml = 0;
717 test_mh = 0;
718 test_h = 0;
719 test_p = 0;
720
721 // For efficiency vs phi
722 for(int i = 0; i < MAX_PHI_BIN; i++) {
723 n_validtrack_phi_arr[i] = 0;
724 n_insidetrack_phi_arr[i] = 0;
725 }
726
727 // For efficiency vs z
728 for(int i = 0; i < MAX_Z_BIN; i++) {
729 n_validtrack_z_arr[i] = 0;
730 n_insidetrack_z_arr[i] = 0;
731 }
732
733 for(int i = 0; i < MAX_PHI_BIN; i++) {
734 for(int j = 0; j < MAX_Z_BIN; j++) {
735 n_validtrack_phi_z_mat[i][j] = 0;
736 n_insidetrack_phi_z_mat[i][j] = 0;
737 }
738 }
739
740 for(int j = 0; j < CHI_BIN; j++) {
741 for(int i = 0; i < MAX_PHI_BIN; i++) {
742 n_validtrack_phi_chi2_mat[j][i] = 0;
743 n_insidetrack_phi_chi2_mat[j][i] = 0;
744 }
745 for(int i = 0; i < MAX_Z_BIN; i++) {
746 n_validtrack_z_chi2_mat[j][i] = 0;
747 n_insidetrack_z_chi2_mat[j][i] = 0;
748 }
749 }
750
751 int chi_dummy_var = 0;
752 int j_dum = 0;
753
754 for(int i = 0; i < CHI_BIN; i++) {
755
756 chi_dummy_var = j_dum*j_dum + 1;
757 if(i/3 == 1) chi_dummy_var = chi_dummy_var*10;
758 if(i/3 == 2) chi_dummy_var = chi_dummy_var*100;
759 if(i/3 == 3) chi_dummy_var = chi_dummy_var*1000;
760
761 chi2_cut_arr[i] = chi_dummy_var;
762 if(j_dum < 2) j_dum++;
763 else j_dum = 0;
764 }
765
766 read_file(filename);
767 cout << "entries " << tree->GetEntries() << endl;
768
769 return StatusCode::SUCCESS;
770
771}
772
773StatusCode CgemCosmicRayQA::execute(){
774
775 MsgStream log(msgSvc(), name());
776
777 if(event%1000==0) cout << "--------->CgemCosmicRayQA::execute " << event << endl;
778
779 /**
780 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),"/Event/EventHeader");
781 if (!eventHeader) {
782 log << MSG::FATAL << "Could not find Event Header" << endreq;
783 return StatusCode::FAILURE;
784 }
785 cout << "EVENTO " << eventHeader->eventNumber() << " RUN " << eventHeader->runNumber() << endl;
786 **/
787
788
789
790 //interface to event data service
791 ISvcLocator* svcLocator = Gaudi::svcLocator();
792 StatusCode sc=svcLocator->service("EventDataSvc", m_evtSvc);
793 if (sc.isFailure())
794 cout<<"Could not accesss EventDataSvc!"<<endl;
795
796 // for(int ievt = 0; ievt < tree->GetEntries(); ievt++) {
797 tree->GetEntry(event);
798
799 // ..................................................................................................
800 // loop on hits
801 fill_hit_histo();
802
803 // ..................................................................................................
804 // loop on hits in time window
805 fill_hit_in_time_histo();
806
807 // ..................................................................................................
808 // loop on cluster1d
809 fill_cluster1d_histo();
810
811 // ..................................................................................................
812 // loop on cluster2d
813 fill_cluster2d_histo();
814
815 // ..................................................................................................
816 // loop on cluster selected
817 fill_cluster_selected_histo();
818
819 // ..................................................................................................
820 // loop on track
821 fill_track_histo();
822
823 event++;
824 if(event==tree->GetEntries()) return StatusCode::FAILURE; // CHECK HACK THIS
825
826 return StatusCode::SUCCESS;
827}
828
829StatusCode CgemCosmicRayQA::finalize(){
830 MsgStream log(msgSvc(),name());
831 log << MSG::INFO << "CgemCosmicRayQA finalize()" << endreq;
832 cout << "finalize" << endl;
833
834 // ..................................................................................................
835 // analisys
836 analyze();
837 histo_cosmetics();
838
839 output->cd();
840 double efficiency = (double) ninside/nvalidtrack;
841 double background = (double) noutside/nvalidtrack;
842 int nentries = tree->GetEntries();
843 int nof_fitted_track = nfittedtrack;
844 TTree *res_tree = new TTree("qa_results","QA results");
845 res_tree->Branch("ntrigger", &nentries, "ntrigger/I");
846 res_tree->Branch("ntrack_fitted", &nfittedtrack, "ntrack_fitted/I");
847 res_tree->Branch("ntrack_not_passing_chi2_cut", &nchi2, "ntrack_not_passing_chi2_cut/I");
848 res_tree->Branch("ntrack_not_passing_ene_cut", &nstopped, "ntrack_not_passing_ene_cut/I");
849 res_tree->Branch("ntrack_valid", &nvalidtrack, "ntrack_valid/I");
850 res_tree->Branch("efficiency", &efficiency, "efficiency/D");
851 res_tree->Branch("background", &background,"background/D");
852 res_tree->Branch("nentries_residual_rphi", &nentries_residual_rphi, "nentries_residual_rphi/I");
853 res_tree->Branch("mean_residual_rphi", &mean_residual_rphi, "mean_residual_rphi/D");
854 res_tree->Branch("sigma_residual_rphi", &sigma_residual_rphi, "sigma_residual_rphi/D");
855 res_tree->Branch("nentries_residual_z", &nentries_residual_z, "nentries_residual_z/I");
856 res_tree->Branch("mean_residual_z", &mean_residual_z, "mean_residual_z/D");
857 res_tree->Branch("sigma_residual_z", &sigma_residual_z, "sigma_residual_z/D");
858 res_tree->Fill();
859 // res_tree->Write();
860
861 output->Write();
862 output->Close();
863
864
865 return StatusCode::SUCCESS;
866}
867
868
869void CgemCosmicRayQA::fill_hit_histo() {
870 // cout << "event " << event << " nhit " << nhit << endl;
871 // ..................................................................................................
872 // loop on hits
873 int nhit_L1_S1_x = 0; int nhit_L1_S2_x = 0; int nhit_L2_S1_x = 0; int nhit_L2_S2_x = 0;
874 int nhit_L1_S1_v = 0; int nhit_L2_S1_v = 0; int nhit_L2_S2_v = 0;
875 for(int ihit = 0; ihit < nhit; ihit++) {
876
877 // L1
878 if(hit_layer[ihit] == 0) { // ................................................. LAYER 1
879
880 if(hit_view[ihit] == 0) {
881
882 // x view
883 if(hit_sheet[ihit] == 0) { // ................................... S1 (CHECK aggiungi ed usa phi)
884 h_hit_charge_L1_S1_x->Fill(hit_q[ihit]);
885 h_hit_time_L1_S1_x->Fill(hit_t[ihit]);
886 h_hit_charge_vs_strip_L1_S1_x->Fill(hit_strip[ihit], hit_q[ihit]);
887 h_hit_time_vs_strip_L1_S1_x->Fill(hit_strip[ihit], hit_t[ihit]);
888 h_hit_charge_vs_time_L1_S1_x->Fill(hit_t[ihit], hit_q[ihit]);
889 nhit_L1_S1_x++;
890 }
891 else { // ................................... S2
892 h_hit_charge_L1_S2_x->Fill(hit_q[ihit]);
893 h_hit_time_L1_S2_x->Fill(hit_t[ihit]);
894 h_hit_charge_vs_strip_L1_S2_x->Fill(hit_strip[ihit], hit_q[ihit]);
895 h_hit_time_vs_strip_L1_S2_x->Fill(hit_strip[ihit], hit_t[ihit]);
896 h_hit_charge_vs_time_L1_S2_x->Fill(hit_t[ihit], hit_q[ihit]);
897 nhit_L1_S2_x++;
898 }
899 }
900 else {
901
902 // v view
903 h_hit_charge_L1_S1_v->Fill(hit_q[ihit]);
904 h_hit_time_L1_S1_v->Fill(hit_t[ihit]);
905 h_hit_charge_vs_strip_L1_S1_v->Fill(hit_strip[ihit], hit_q[ihit]);
906 h_hit_time_vs_strip_L1_S1_v->Fill(hit_strip[ihit], hit_t[ihit]);
907 h_hit_charge_vs_time_L1_S1_v->Fill(hit_t[ihit], hit_q[ihit]);
908 h_hit_charge_vs_length_L1_S1_v->Fill(hit_length[ihit], hit_q[ihit]);
909 nhit_L1_S1_v++;
910 }
911 }
912 else if(hit_layer[ihit] == 1) { // ................................................. LAYER 2
913
914 if(hit_view[ihit] == 0) {
915
916 // x view
917 if(hit_sheet[ihit] == 0) { // ................................... S1
918 h_hit_charge_L2_S1_x->Fill(hit_q[ihit]);
919 h_hit_time_L2_S1_x->Fill(hit_t[ihit]);
920 h_hit_charge_vs_strip_L2_S1_x->Fill(hit_strip[ihit], hit_q[ihit]);
921 h_hit_time_vs_strip_L2_S1_x->Fill(hit_strip[ihit], hit_t[ihit]);
922 h_hit_charge_vs_time_L2_S1_x->Fill(hit_t[ihit], hit_q[ihit]);
923 nhit_L2_S1_x++;
924 }
925 else { // ................................... S2
926 h_hit_charge_L2_S2_x->Fill(hit_q[ihit]);
927 h_hit_time_L2_S2_x->Fill(hit_t[ihit]);
928 h_hit_charge_vs_strip_L2_S2_x->Fill(hit_strip[ihit], hit_q[ihit]);
929 h_hit_time_vs_strip_L2_S2_x->Fill(hit_strip[ihit], hit_t[ihit]);
930 h_hit_charge_vs_time_L2_S2_x->Fill(hit_t[ihit], hit_q[ihit]);
931 nhit_L2_S2_x++;
932 }
933 }
934 else {
935
936 // v view
937 if(hit_sheet[ihit] == 0) { // ................................... S1
938 h_hit_charge_L2_S1_v->Fill(hit_q[ihit]);
939 h_hit_time_L2_S1_v->Fill(hit_t[ihit]);
940 h_hit_charge_vs_strip_L2_S1_v->Fill(hit_strip[ihit], hit_q[ihit]);
941 h_hit_time_vs_strip_L2_S1_v->Fill(hit_strip[ihit], hit_t[ihit]);
942 h_hit_charge_vs_time_L2_S1_v->Fill(hit_t[ihit], hit_q[ihit]);
943 h_hit_charge_vs_length_L2_S1_v->Fill(hit_length[ihit], hit_q[ihit]);
944 nhit_L2_S1_v++;
945 }
946 else { // ................................... S2
947 h_hit_charge_L2_S2_v->Fill(hit_q[ihit]);
948 h_hit_time_L2_S2_v->Fill(hit_t[ihit]);
949 h_hit_charge_vs_strip_L2_S2_v->Fill(hit_strip[ihit], hit_q[ihit]);
950 h_hit_time_vs_strip_L2_S2_v->Fill(hit_strip[ihit], hit_t[ihit]);
951 h_hit_charge_vs_time_L2_S2_v->Fill(hit_t[ihit], hit_q[ihit]);
952 h_hit_charge_vs_length_L2_S2_v->Fill(hit_length[ihit], hit_q[ihit]);
953 nhit_L2_S2_v++;
954 }
955 }
956 }
957 }
958
959 h_nofhit_L1_S1_x->Fill(nhit_L1_S1_x);
960 h_nofhit_L1_S2_x->Fill(nhit_L1_S2_x);
961 h_nofhit_L2_S1_x->Fill(nhit_L2_S1_x);
962 h_nofhit_L2_S2_x->Fill(nhit_L2_S2_x);
963 h_nofhit_L1_S1_v->Fill(nhit_L1_S1_v);
964 h_nofhit_L2_S1_v->Fill(nhit_L2_S1_v);
965 h_nofhit_L2_S2_v->Fill(nhit_L2_S2_v);
966
967}
968
969void CgemCosmicRayQA::fill_hit_in_time_histo() {
970 // cout << "event " << event << " nhit " << nhit << endl;
971 // ..................................................................................................
972 // loop on hits
973 int nhit_L1_S1_x = 0; int nhit_L1_S2_x = 0; int nhit_L2_S1_x = 0; int nhit_L2_S2_x = 0;
974 int nhit_L1_S1_v = 0; int nhit_L2_S1_v = 0; int nhit_L2_S2_v = 0;
975 for(int ihit = 0; ihit < nhit; ihit++) {
976
977 // if hit outside time window -> continue
978 if(hit_t[ihit] < minDigiTime || hit_t[ihit] > maxDigiTime) continue;
979
980 // L1
981 if(hit_layer[ihit] == 0) { // ................................................. LAYER 1
982
983 if(hit_view[ihit] == 0) {
984
985 // x view
986 if(hit_sheet[ihit] == 0) { // ................................... S1 (CHECK aggiungi ed usa phi)
987 h_hit_in_time_charge_L1_S1_x->Fill(hit_q[ihit]);
988 h_hit_in_time_strip_L1_S1_x->Fill(hit_strip[ihit]);
989 h_hit_in_time_charge_vs_strip_L1_S1_x->Fill(hit_strip[ihit], hit_q[ihit]);
990 }
991 else { // ................................... S2
992 h_hit_in_time_charge_L1_S2_x->Fill(hit_q[ihit]);
993 h_hit_in_time_strip_L1_S2_x->Fill(hit_strip[ihit]);
994 h_hit_in_time_charge_vs_strip_L1_S2_x->Fill(hit_strip[ihit], hit_q[ihit]);
995 }
996 }
997 else {
998
999 // v view
1000
1001 h_hit_in_time_charge_L1_S1_v->Fill(hit_q[ihit]);
1002 h_hit_in_time_strip_L1_S1_v->Fill(hit_strip[ihit]);
1003 h_hit_in_time_charge_vs_strip_L1_S1_v->Fill(hit_strip[ihit], hit_q[ihit]);
1004 }
1005 }
1006 else if(hit_layer[ihit] == 1) { // ................................................. LAYER 2
1007
1008 if(hit_view[ihit] == 0) {
1009
1010 // x view
1011 if(hit_sheet[ihit] == 0) { // ................................... S1
1012 h_hit_in_time_charge_L2_S1_x->Fill(hit_q[ihit]);
1013 h_hit_in_time_strip_L2_S1_x->Fill(hit_strip[ihit]);
1014 h_hit_in_time_charge_vs_strip_L2_S1_x->Fill(hit_strip[ihit], hit_q[ihit]);
1015 }
1016 else { // ................................... S2
1017 h_hit_in_time_charge_L2_S2_x->Fill(hit_q[ihit]);
1018 h_hit_in_time_strip_L2_S2_x->Fill(hit_strip[ihit]);
1019 h_hit_in_time_charge_vs_strip_L2_S2_x->Fill(hit_strip[ihit], hit_q[ihit]);
1020 }
1021 }
1022 else {
1023
1024 // v view
1025 if(hit_sheet[ihit] == 0) { // ................................... S1
1026 h_hit_in_time_charge_L2_S1_v->Fill(hit_q[ihit]);
1027 h_hit_in_time_strip_L2_S1_v->Fill(hit_strip[ihit]);
1028 h_hit_in_time_charge_vs_strip_L2_S1_v->Fill(hit_strip[ihit], hit_q[ihit]);
1029 }
1030 else { // ................................... S2
1031 h_hit_in_time_charge_L2_S2_v->Fill(hit_q[ihit]);
1032 h_hit_in_time_strip_L2_S2_v->Fill(hit_strip[ihit]);
1033 h_hit_in_time_charge_vs_strip_L2_S2_v->Fill(hit_strip[ihit], hit_q[ihit]);
1034 }
1035 }
1036 }
1037 }
1038}
1039
1040void CgemCosmicRayQA::fill_cluster1d_histo() {
1041
1042 // loop on cluster 1d
1043 int counter=0;
1044 for(int iclu = 0; iclu < ncluster; iclu++) {
1045
1046 if(cluster_1d_view[iclu] == -1) continue;
1047 if(cluster_1d_view[iclu] != 0 && cluster_1d_view[iclu] != 1) cout << "ERROR " << cluster_1d_view[iclu] << endl;
1048 counter++;
1049 // L1
1050 if(cluster_1d_layerid[iclu] == 0) { // ................................................. LAYER 1
1051
1052 if(cluster_1d_view[iclu] == 0) {
1053
1054 // x view
1055 if(cluster_1d_phi[iclu] < 0) { // ................................... S1
1056 h_cluster1d_size_L1_S1_x->Fill(cluster_1d_size[iclu]);
1057 h_cluster1d_charge_vs_phi_L1_S1_x->Fill(TMath::RadToDeg()*cluster_1d_phi[iclu], cluster_1d_q[iclu]);
1058 }
1059 else { // ................................... S2
1060 h_cluster1d_size_L1_S2_x->Fill(cluster_1d_size[iclu]);
1061 h_cluster1d_charge_vs_phi_L1_S2_x->Fill(TMath::RadToDeg()*cluster_1d_phi[iclu], cluster_1d_q[iclu]);
1062 }
1063 }
1064 else {
1065 // v view
1066 h_cluster1d_size_L1_S1_v->Fill(cluster_1d_size[iclu]);
1067 h_cluster1d_charge_vs_v_L1_S1_v->Fill(cluster_1d_v[iclu], cluster_1d_q[iclu]);
1068 }
1069 }
1070 else if(cluster_1d_layerid[iclu] == 1) { // ................................................. LAYER 2
1071
1072 if(cluster_1d_view[iclu] == 0) {
1073
1074 // x view
1075 if(cluster_1d_sheetid[iclu] == 0) { // ................................... S1
1076 h_cluster1d_size_L2_S1_x->Fill(cluster_1d_size[iclu]);
1077 h_cluster1d_charge_vs_phi_L2_S1_x->Fill(TMath::RadToDeg()*cluster_1d_phi[iclu], cluster_1d_q[iclu]);
1078 }
1079 else { // ................................... S2
1080 h_cluster1d_size_L2_S2_x->Fill(cluster_1d_size[iclu]);
1081 h_cluster1d_charge_vs_phi_L2_S2_x->Fill(TMath::RadToDeg()*cluster_1d_phi[iclu], cluster_1d_q[iclu]);
1082 }
1083 }
1084 else {
1085
1086 // v view
1087 if(cluster_1d_sheetid[iclu] == 0) { // ................................... S1
1088 h_cluster1d_size_L2_S1_v->Fill(cluster_1d_size[iclu]);
1089 h_cluster1d_charge_vs_v_L2_S1_v->Fill(cluster_1d_v[iclu], cluster_1d_q[iclu]);
1090 }
1091 else { // ................................... S2
1092 h_cluster1d_size_L2_S2_v->Fill(cluster_1d_size[iclu]);
1093 h_cluster1d_charge_vs_v_L2_S2_v->Fill(cluster_1d_v[iclu], cluster_1d_q[iclu]);
1094 }
1095 }
1096 }
1097 }
1098
1099 h_nofcluster1d_L1_S1_x->Fill(ncluster_1d_L1_S1_x);
1100 // h_nofcluster1d_L1_S2_x->Fill(ncluster_1d_L1_S2_x); // CHECK
1101 h_nofcluster1d_L2_S1_x->Fill(ncluster_1d_L2_S1_x);
1102 h_nofcluster1d_L2_S2_x->Fill(ncluster_1d_L2_S2_x);
1103 h_nofcluster1d_L1_S1_v->Fill(ncluster_1d_L1_S1_v);
1104 h_nofcluster1d_L2_S1_v->Fill(ncluster_1d_L2_S1_v);
1105 h_nofcluster1d_L2_S2_v->Fill(ncluster_1d_L2_S2_v);
1106
1107 if(ncluster_1d != counter) cout << "ERROR in fill_cluster1d_histo: ncluster_1d " << ncluster_1d << ", counter " << counter <<endl;
1108
1109}
1110
1111/**
1112* SELECTION 1 - require all firing planes, i.e. at least one cluster 2d for each plane
1113* get cluster 2d with max charge on each plane
1114* SELECTION 2 - require all charges 1d above cuts
1115* SELECTION 3 - require cl.size 1d > 1 strip
1116**/
1117void CgemCosmicRayQA::fill_cluster_selected_histo() {
1118
1119 double tmp_charge_L1bot = 0;
1120 double tmp_charge_L1top = 0;
1121 double tmp_charge_L2bot = 0;
1122 double tmp_charge_L2top = 0;
1123 int tmp_cluster_L1bot = -1;
1124 int tmp_cluster_L1top = -1;
1125 int tmp_cluster_L2bot = -1;
1126 int tmp_cluster_L2top = -1;
1127 // loop on cluster 2d
1128 for(int iclu = 0; iclu < ncluster; iclu++) {
1129
1130 if(cluster_2d_view[iclu] == -1) continue;
1131 if(cluster_2d_view[iclu] != 2) cout << "ERROR in fill_cluster2d_histo " << cluster_2d_view[iclu] << endl;
1132
1133 if(cluster_2d_layerid[iclu] == 0) { // ................................................. LAYER 1
1134
1135 if(cluster_2d_phi[iclu] < 0) { // ................................... S1
1136 if(cluster_2d_q[iclu] > tmp_charge_L1bot) {
1137 tmp_charge_L1bot = cluster_2d_q[iclu];
1138 tmp_cluster_L1bot = iclu;
1139 }
1140 }
1141 else { // ................................... S2
1142 if(cluster_2d_q[iclu] > tmp_charge_L1top) {
1143 tmp_charge_L1top = cluster_2d_q[iclu];
1144 tmp_cluster_L1top = iclu;
1145 }
1146 }
1147 }
1148 else if(cluster_2d_layerid[iclu] == 1) { // ............................................ LAYER 2
1149
1150 if(cluster_2d_sheetid[iclu] == 0) { // ............................... S1
1151 if(cluster_2d_q[iclu] > tmp_charge_L2bot) {
1152 tmp_charge_L2bot = cluster_2d_q[iclu];
1153 tmp_cluster_L2bot = iclu;
1154 }
1155 }
1156 else { // ................................... S2
1157 if(cluster_2d_q[iclu] > tmp_charge_L2top) {
1158 tmp_charge_L2top = cluster_2d_q[iclu];
1159 tmp_cluster_L2top = iclu;
1160 }
1161 }
1162 }
1163 }
1164
1165 // SELECTION 1 - four firing
1166 if(tmp_cluster_L1bot == -1 || tmp_cluster_L1top == -1 || tmp_cluster_L2bot == -1 || tmp_cluster_L2top == -1) return;
1167
1168 // check it is not using the same v strip on L1
1169 if(cluster_2d_idv[tmp_cluster_L1bot] == cluster_2d_idv[tmp_cluster_L1top]) return;
1170
1171 // SELECTION 2 - energy
1172 if(cluster_1d_q[cluster_2d_idx[tmp_cluster_L1bot]] < cut_ene_L1_x || cluster_1d_q[cluster_2d_idv[tmp_cluster_L1bot]] < cut_ene_L1_v ||
1173 cluster_1d_q[cluster_2d_idx[tmp_cluster_L1top]] < cut_ene_L1_x || cluster_1d_q[cluster_2d_idv[tmp_cluster_L1top]] < cut_ene_L1_v ||
1174 cluster_1d_q[cluster_2d_idx[tmp_cluster_L2bot]] < cut_ene_L2_x || cluster_1d_q[cluster_2d_idv[tmp_cluster_L2bot]] < cut_ene_L2_v ||
1175 cluster_1d_q[cluster_2d_idx[tmp_cluster_L2top]] < cut_ene_L2_x || cluster_1d_q[cluster_2d_idv[tmp_cluster_L2top]] < cut_ene_L2_v) return;
1176
1177 // SELECTION 3 - size
1178 if(cluster_1d_size[cluster_2d_idx[tmp_cluster_L1bot]] <= select_size || cluster_1d_size[cluster_2d_idv[tmp_cluster_L1bot]] <= select_size ||
1179 cluster_1d_size[cluster_2d_idx[tmp_cluster_L1top]] <= select_size || cluster_1d_size[cluster_2d_idv[tmp_cluster_L1top]] <= select_size ||
1180 cluster_1d_size[cluster_2d_idx[tmp_cluster_L2bot]] <= select_size || cluster_1d_size[cluster_2d_idv[tmp_cluster_L2bot]] <= select_size ||
1181 cluster_1d_size[cluster_2d_idx[tmp_cluster_L2top]] <= select_size || cluster_1d_size[cluster_2d_idv[tmp_cluster_L2top]] <= select_size) return;
1182
1183 // fill charge 2d
1184 h_cluster2d_charge_selected_L1_S1->Fill(tmp_charge_L1bot);
1185 h_cluster2d_charge_selected_L1_S2->Fill(tmp_charge_L1top);
1186 h_cluster2d_charge_selected_L2_S1->Fill(tmp_charge_L2bot);
1187 h_cluster2d_charge_selected_L2_S2->Fill(tmp_charge_L2top);
1188
1189 // fill charge 1d x
1190 h_cluster1d_charge_selected_L1_S1_x->Fill(cluster_1d_q[cluster_2d_idx[tmp_cluster_L1bot]]);
1191 h_cluster1d_charge_selected_L1_S2_x->Fill(cluster_1d_q[cluster_2d_idx[tmp_cluster_L1top]]);
1192 h_cluster1d_charge_selected_L2_S1_x->Fill(cluster_1d_q[cluster_2d_idx[tmp_cluster_L2bot]]);
1193 h_cluster1d_charge_selected_L2_S2_x->Fill(cluster_1d_q[cluster_2d_idx[tmp_cluster_L2top]]);
1194
1195 // fill charge 1d v
1196 h_cluster1d_charge_selected_L1_S1_v->Fill(cluster_1d_q[cluster_2d_idv[tmp_cluster_L1bot]]);
1197 h_cluster1d_charge_selected_L1_S2_v->Fill(cluster_1d_q[cluster_2d_idv[tmp_cluster_L1top]]);
1198 h_cluster1d_charge_selected_L2_S1_v->Fill(cluster_1d_q[cluster_2d_idv[tmp_cluster_L2bot]]);
1199 h_cluster1d_charge_selected_L2_S2_v->Fill(cluster_1d_q[cluster_2d_idv[tmp_cluster_L2top]]);
1200
1201 // fill size 1d x
1202 h_cluster1d_size_selected_L1_S1_x->Fill(cluster_1d_size[cluster_2d_idx[tmp_cluster_L1bot]]);
1203 h_cluster1d_size_selected_L1_S2_x->Fill(cluster_1d_size[cluster_2d_idx[tmp_cluster_L1top]]);
1204 h_cluster1d_size_selected_L2_S1_x->Fill(cluster_1d_size[cluster_2d_idx[tmp_cluster_L2bot]]);
1205 h_cluster1d_size_selected_L2_S2_x->Fill(cluster_1d_size[cluster_2d_idx[tmp_cluster_L2top]]);
1206
1207 // fill size 1d v
1208 h_cluster1d_size_selected_L1_S1_v->Fill(cluster_1d_size[cluster_2d_idv[tmp_cluster_L1bot]]);
1209 h_cluster1d_size_selected_L1_S2_v->Fill(cluster_1d_size[cluster_2d_idv[tmp_cluster_L1top]]);
1210 h_cluster1d_size_selected_L2_S1_v->Fill(cluster_1d_size[cluster_2d_idv[tmp_cluster_L2bot]]);
1211 h_cluster1d_size_selected_L2_S2_v->Fill(cluster_1d_size[cluster_2d_idv[tmp_cluster_L2top]]);
1212
1213
1214}
1215
1216void CgemCosmicRayQA::fill_cluster2d_histo() {
1217
1218 // loop on cluster 2d
1219 int ncluster_2d_L1_S1 = 0; int ncluster_2d_L1_S2 = 0; int ncluster_2d_L2_S1 = 0; int ncluster_2d_L2_S2 = 0;
1220 int counter=0;
1221 for(int iclu = 0; iclu < ncluster; iclu++) {
1222
1223 if(cluster_2d_view[iclu] == -1) continue;
1224 if(cluster_2d_view[iclu] != 2) cout << "ERROR in fill_cluster2d_histo " << cluster_2d_view[iclu] << endl;
1225 counter++;
1226 // L1
1227 if(cluster_2d_layerid[iclu] == 0) { // ................................................. LAYER 1
1228
1229 if(cluster_2d_phi[iclu] < 0) { // ................................... S1
1230 // h_cluster2d_size_L1_S1->Fill(cluster_2d_size[iclu]);
1231 h_cluster2d_charge_vs_phi_L1_S1->Fill(TMath::RadToDeg()*cluster_2d_phi[iclu], cluster_2d_q[iclu]);
1232 h_cluster2d_charge_vs_z_L1_S1->Fill(cluster_2d_z[iclu], cluster_2d_q[iclu]);
1233 h_cluster2d_z_vs_phi_L1_S1->Fill(TMath::RadToDeg()*cluster_2d_phi[iclu], cluster_2d_z[iclu]);
1234 ncluster_2d_L1_S1++;
1235 }
1236 else { // ................................... S2
1237 // h_cluster2d_size_L1_S2->Fill(cluster_2d_size[iclu]);
1238 h_cluster2d_charge_vs_phi_L1_S2->Fill(TMath::RadToDeg()*cluster_2d_phi[iclu], cluster_2d_q[iclu]);
1239 h_cluster2d_charge_vs_z_L1_S2->Fill(cluster_2d_z[iclu], cluster_2d_q[iclu]);
1240 h_cluster2d_z_vs_phi_L1_S2->Fill(TMath::RadToDeg()*cluster_2d_phi[iclu], cluster_2d_z[iclu]);
1241 ncluster_2d_L1_S2++;
1242 }
1243 }
1244 else if(cluster_2d_layerid[iclu] == 1) { // ................................................. LAYER 2
1245
1246 if(cluster_2d_sheetid[iclu] == 0) { // ................................... S1
1247 // h_cluster2d_size_L2_S1->Fill(cluster_2d_size[iclu]);
1248 h_cluster2d_charge_vs_phi_L2_S1->Fill(TMath::RadToDeg()*cluster_2d_phi[iclu], cluster_2d_q[iclu]);
1249 h_cluster2d_charge_vs_z_L2_S1->Fill(cluster_2d_z[iclu], cluster_2d_q[iclu]);
1250 h_cluster2d_z_vs_phi_L2_S1->Fill(TMath::RadToDeg()*cluster_2d_phi[iclu], cluster_2d_z[iclu]);
1251 ncluster_2d_L2_S1++;
1252 }
1253 else { // ................................... S2
1254 // h_cluster2d_size_L2_S2->Fill(cluster_2d_size[iclu]);
1255 h_cluster2d_charge_vs_phi_L2_S2->Fill(TMath::RadToDeg()*cluster_2d_phi[iclu], cluster_2d_q[iclu]);
1256 h_cluster2d_charge_vs_z_L2_S2->Fill(cluster_2d_z[iclu], cluster_2d_q[iclu]);
1257 h_cluster2d_z_vs_phi_L2_S2->Fill(TMath::RadToDeg()*cluster_2d_phi[iclu], cluster_2d_z[iclu]);
1258 ncluster_2d_L2_S2++;
1259 }
1260 }
1261 }
1262
1263
1264 h_nofcluster2d_L1_S1->Fill(ncluster_2d_L1_S1);
1265 h_nofcluster2d_L1_S2->Fill(ncluster_2d_L1_S2); // CHECK
1266 h_nofcluster2d_L2_S1->Fill(ncluster_2d_L2_S1);
1267 h_nofcluster2d_L2_S2->Fill(ncluster_2d_L2_S2);
1268 if(ncluster_2d != counter) cout << "ERROR in fill_cluster2d_histo: ncluster_2d " << ncluster_2d << ", counter " << counter << endl;
1269}
1270
1271void CgemCosmicRayQA::fill_track_histo() {
1272
1273 // is fitted?
1274 if(track_chi2 == 9999) return;
1275 nfittedtrack++;
1276
1277 h_track_chi2->Fill(track_chi2);
1278
1279 bool isgood = apply_cuts(); // CHECK
1280 if(isgood == false) return;
1281 nvalidtrack++;
1282
1283 h_track_xpoca->Fill(track_xpoca_glo);
1284 h_track_ypoca->Fill(track_ypoca_glo);
1285 h_track_zpoca->Fill(track_zpoca_glo);
1286
1287 // -----------------------------
1288 for(int iclu = 0; iclu < track_nfitpoint; iclu++) {
1289 int ilay = track_layerid[iclu];
1290 int ishe = track_sheetid[iclu];
1291 int cluster2d_id = track_clusterid[iclu];
1292 int idx = cluster_2d_idx[cluster2d_id];
1293 int idv = cluster_2d_idv[cluster2d_id];
1294 if(ilay==0 && ishe==0) {
1295 h_cluster1d_charge_fitted_L1_S1_x->Fill(cluster_1d_q[idx]);
1296 h_cluster1d_charge_fitted_L1_S1_v->Fill(cluster_1d_q[idv]);
1297 }
1298 else if(ilay==0 && ishe==1) {
1299 h_cluster1d_charge_fitted_L1_S2_x->Fill(cluster_1d_q[idx]);
1300 h_cluster1d_charge_fitted_L1_S2_v->Fill(cluster_1d_q[idv]);
1301 }
1302 else if(ilay==1 && ishe==0) {
1303 h_cluster1d_charge_fitted_L2_S1_x->Fill(cluster_1d_q[idx]);
1304 h_cluster1d_charge_fitted_L2_S1_v->Fill(cluster_1d_q[idv]);
1305 }
1306 else if(ilay==1 && ishe==1) {
1307 h_cluster1d_charge_fitted_L2_S2_x->Fill(cluster_1d_q[idx]);
1308 h_cluster1d_charge_fitted_L2_S2_v->Fill(cluster_1d_q[idv]);
1309 }
1310 }
1311 // fitted ..............................
1312 double fit_phi[MAXNOFLAYER][MAXNOFSHEET];
1313 double fit_z[MAXNOFLAYER][MAXNOFSHEET];
1314
1315 // phi
1316 fit_phi[0][0] = track_phi1bot_loc;
1317 fit_phi[0][1] = track_phi1top_loc;
1318 fit_phi[1][0] = track_phi2bot_loc;
1319 fit_phi[1][1] = track_phi2top_loc;
1320
1321 // z
1322 if(align_flag==true) {
1323
1324 HepPoint3D fit_L1_S1_glo(track_x1bot_glo, track_y1bot_glo, track_z1bot_glo);
1325 HepPoint3D fit_L1_S1_loc = alignment->point_transform(0, 0, fit_L1_S1_glo);
1326 HepPoint3D fit_L1_S2_glo(track_x1top_glo, track_y1top_glo, track_z1top_glo);
1327 HepPoint3D fit_L1_S2_loc = alignment->point_transform(0, 1, fit_L1_S2_glo);
1328
1329 HepPoint3D fit_L2_S1_glo(track_x2bot_glo, track_y2bot_glo, track_z2bot_glo);
1330 HepPoint3D fit_L2_S1_loc = alignment->point_transform(1, 0, fit_L2_S1_glo);
1331 HepPoint3D fit_L2_S2_glo(track_x2top_glo, track_y2top_glo, track_z2top_glo);
1332 HepPoint3D fit_L2_S2_loc = alignment->point_transform(1, 1, fit_L2_S2_glo);
1333
1334 fit_z[0][0] = fit_L1_S1_loc.z();
1335 fit_z[0][1] = fit_L1_S2_loc.z();
1336 fit_z[1][0] = fit_L2_S1_loc.z();
1337 fit_z[1][1] = fit_L2_S2_loc.z();
1338 }
1339 else {
1340
1341 fit_z[0][0] = track_z1bot_glo;
1342 fit_z[0][1] = track_z1top_glo;
1343 fit_z[1][0] = track_z2bot_glo;
1344 fit_z[1][1] = track_z2top_glo;
1345 }
1346 // ------------------------
1347
1348
1349 // experimental ........................
1350 double exp_phi[MAXNOFLAYER][MAXNOFSHEET];
1351 double exp_z[MAXNOFLAYER][MAXNOFSHEET];
1352
1353 double exp_phi_tpc[MAXNOFLAYER][MAXNOFSHEET];
1354 double exp_z_tpc[MAXNOFLAYER][MAXNOFSHEET];
1355
1356 bool stop_here=false;
1357 // trackers
1358 for(int iclu = 0; iclu < track_nfitpoint; iclu++) {
1359 int exp_lay = track_layerid[iclu];
1360 int exp_she = track_sheetid[iclu];
1361
1362 exp_phi[exp_lay][exp_she] = cluster_2d_phi[track_clusterid[iclu]];
1363 exp_z[exp_lay][exp_she] = cluster_2d_z[track_clusterid[iclu]];
1364
1365 exp_phi_tpc[exp_lay][exp_she] = cluster_2d_phi_tpc[track_clusterid[iclu]];
1366 exp_z_tpc[exp_lay][exp_she] = cluster_2d_z_tpc[track_clusterid[iclu]];
1367
1368 }
1369
1370 // residuals
1371 for(int ilay=0; ilay<MAXNOFLAYER; ilay++) {
1372 for(int ishe=0; ishe<MAXNOFLAYER; ishe++) {
1373
1374 double residual_rphi = anode_mid_gap_radius[ilay] * (exp_phi[ilay][ishe] - fit_phi[ilay][ishe]);
1375 double residual_z = exp_z[ilay][ishe] - fit_z[ilay][ishe];
1376
1377 double residual_rphi_tpc = anode_mid_gap_radius[ilay] * (exp_phi_tpc[ilay][ishe] - fit_phi[ilay][ishe]);
1378 double residual_z_tpc = exp_z_tpc[ilay][ishe] - fit_z[ilay][ishe];
1379
1380 if(ilay==0 && ishe==0) {
1381 h_residual_rphi_L1_S1->Fill(residual_rphi);
1382 h_residual_z_L1_S1->Fill(residual_z);
1383
1384 h_res_rphi_vs_chi2_L1_S1->Fill(residual_rphi, track_chi2);
1385 h_res_z_vs_chi2_L1_S1->Fill(residual_z, track_chi2);
1386
1387 h_residual_rphi_L1_S1_tpc->Fill(residual_rphi_tpc);
1388 h_residual_z_L1_S1_tpc->Fill(residual_z_tpc);
1389
1390 h_res_rphi_vs_chi2_L1_S1_tpc->Fill(residual_rphi_tpc, track_chi2);
1391 h_res_z_vs_chi2_L1_S1_tpc->Fill(residual_z_tpc, track_chi2);
1392 }
1393 else if(ilay==0 && ishe==1) {
1394 h_residual_rphi_L1_S2->Fill(residual_rphi);
1395 h_residual_z_L1_S2->Fill(residual_z);
1396
1397 h_res_rphi_vs_chi2_L1_S2->Fill(residual_rphi, track_chi2);
1398 h_res_z_vs_chi2_L1_S2->Fill(residual_z, track_chi2);
1399
1400 h_residual_rphi_L1_S2_tpc->Fill(residual_rphi_tpc);
1401 h_residual_z_L1_S2_tpc->Fill(residual_z_tpc);
1402
1403 h_res_rphi_vs_chi2_L1_S2_tpc->Fill(residual_rphi_tpc, track_chi2);
1404 h_res_z_vs_chi2_L1_S2_tpc->Fill(residual_z_tpc, track_chi2);
1405 }
1406 else if(ilay==1 && ishe==0) {
1407 h_residual_rphi_L2_S1->Fill(residual_rphi);
1408 h_residual_z_L2_S1->Fill(residual_z);
1409
1410 h_res_rphi_vs_chi2_L2_S1->Fill(residual_rphi, track_chi2);
1411 h_res_z_vs_chi2_L2_S1->Fill(residual_z, track_chi2);
1412
1413 h_residual_rphi_L2_S1_tpc->Fill(residual_rphi_tpc);
1414 h_residual_z_L2_S1_tpc->Fill(residual_z_tpc);
1415
1416 h_res_rphi_vs_chi2_L2_S1_tpc->Fill(residual_rphi_tpc, track_chi2);
1417 h_res_z_vs_chi2_L2_S1_tpc->Fill(residual_z_tpc, track_chi2);
1418 }
1419 else if(ilay==1 && ishe==1) {
1420 h_residual_rphi_L2_S2->Fill(residual_rphi);
1421 h_residual_z_L2_S2->Fill(residual_z);
1422
1423 h_res_rphi_vs_chi2_L2_S2->Fill(residual_rphi, track_chi2);
1424 h_res_z_vs_chi2_L2_S2->Fill(residual_z, track_chi2);
1425
1426 h_residual_rphi_L2_S2_tpc->Fill(residual_rphi_tpc);
1427 h_residual_z_L2_S2_tpc->Fill(residual_z_tpc);
1428
1429 h_res_rphi_vs_chi2_L2_S2_tpc->Fill(residual_rphi_tpc, track_chi2);
1430 h_res_z_vs_chi2_L2_S2_tpc->Fill(residual_z_tpc, track_chi2);
1431 }
1432 }
1433 }
1434
1435 // --------------------------- IF WE ARE TESTING --------------------
1436 // test plane
1437 int exp_lay = track_test_layerid;
1438 int exp_she = track_test_sheetid;
1439 if(exp_lay != -1 && exp_she != -1) {
1440 int closest_cluster2d_id = find_closest_exp_cluster2d(fit_phi[exp_lay][exp_she], fit_z[exp_lay][exp_she]);
1441 if(closest_cluster2d_id == -1) return;
1442
1443 exp_phi[exp_lay][exp_she] = cluster_2d_phi[closest_cluster2d_id];
1444 exp_z[exp_lay][exp_she] = cluster_2d_z[closest_cluster2d_id];
1445
1446 exp_phi_tpc[exp_lay][exp_she] = cluster_2d_phi_tpc[closest_cluster2d_id];
1447 exp_z_tpc[exp_lay][exp_she] = cluster_2d_z_tpc[closest_cluster2d_id];
1448
1449 // test plane (redundant)
1450 double residual_rphi = anode_mid_gap_radius[track_test_layerid] * (exp_phi[track_test_layerid][track_test_sheetid] - fit_phi[track_test_layerid][track_test_sheetid]);
1451 double residual_z = exp_z[track_test_layerid][track_test_sheetid] - fit_z[track_test_layerid][track_test_sheetid];
1452
1453 double exp_phi_test_layer = TMath::RadToDeg()*exp_phi[track_test_layerid][track_test_sheetid];
1454 vector_test_phi.push_back(exp_phi_test_layer);
1455
1456 double exp_z_test_layer = exp_z[track_test_layerid][track_test_sheetid];
1457 vector_test_z.push_back(exp_z_test_layer);
1458
1459 vector_chi2.push_back(track_chi2);
1460
1461 if(track_chi2 <= cut_chi2) {
1462
1463 for(int jj = 0; jj < MAX_PHI_BIN; jj++) {
1464 if(exp_phi_test_layer >= MIN_PHI+(PHI_STEP*jj) && exp_phi_test_layer < MIN_PHI+(PHI_STEP*(jj+1))) n_validtrack_phi_arr[jj]++;
1465 }
1466
1467 for(int jj = 0; jj < MAX_Z_BIN; jj++) {
1468 if(exp_z_test_layer >= MIN_Z+(Z_STEP*jj) && exp_z_test_layer < MIN_Z+(Z_STEP*(jj+1))) n_validtrack_z_arr[jj]++;
1469 }
1470
1471 for(int ii = 0; ii < MAX_PHI_BIN; ii++) {
1472 for(int jj = 0; jj < MAX_Z_BIN; jj++) {
1473 if( (exp_phi_test_layer >= MIN_PHI+(PHI_STEP*ii) && exp_phi_test_layer < MIN_PHI+(PHI_STEP*(ii+1))) && (exp_z_test_layer >= MIN_Z+(Z_STEP*jj) && exp_z_test_layer < MIN_Z+(Z_STEP*(jj+1))) ) n_validtrack_phi_z_mat[ii][jj]++;
1474 }
1475 }
1476 }
1477
1478 for(int chi_i = 0; chi_i < CHI_BIN; chi_i++) {
1479 for(int jj = 0; jj < MAX_PHI_BIN; jj++) {
1480 if( track_chi2 < chi2_cut_arr[chi_i] && (exp_phi_test_layer >= MIN_PHI+(PHI_STEP*jj) && exp_phi_test_layer < MIN_PHI+(PHI_STEP*(jj+1)))) n_validtrack_phi_chi2_mat[chi_i][jj]++;
1481 }
1482 for(int jj = 0; jj < MAX_Z_BIN; jj++) {
1483 if( track_chi2 < chi2_cut_arr[chi_i] && (exp_z_test_layer >= MIN_Z+(Z_STEP*jj) && exp_z_test_layer < MIN_Z+(Z_STEP*(jj+1)))) n_validtrack_z_chi2_mat[chi_i][jj]++;
1484 }
1485 }
1486
1487 double residual_rphi_tpc = anode_mid_gap_radius[track_test_layerid] * (exp_phi_tpc[track_test_layerid][track_test_sheetid] - fit_phi[track_test_layerid][track_test_sheetid]);
1488 double residual_z_tpc = exp_z_tpc[track_test_layerid][track_test_sheetid] - fit_z[track_test_layerid][track_test_sheetid];
1489
1490 h_test_residual_rphi->Fill(residual_rphi);
1491 h_test_residual_z->Fill(residual_z);
1492
1493 h_test_res_rphi_vs_chi2->Fill(residual_rphi, track_chi2);
1494 h_test_res_z_vs_chi2->Fill(residual_z, track_chi2);
1495
1496 h_test_residual_rphi_tpc->Fill(residual_rphi_tpc);
1497 h_test_residual_z_tpc->Fill(residual_z_tpc);
1498
1499 h_test_res_rphi_vs_chi2_tpc->Fill(residual_rphi_tpc, track_chi2);
1500 h_test_res_z_vs_chi2_tpc->Fill(residual_z_tpc, track_chi2);
1501
1502 h_test_residual_rphi_vs_angxy_L1_cc ->Fill(residual_rphi , ang_xy_L1*TMath::RadToDeg());
1503 h_test_residual_rphi_vs_angxy_L1_tpc->Fill(residual_rphi_tpc, ang_xy_L1*TMath::RadToDeg());
1504 h_test_residual_rphi_vs_angxy_L2_cc ->Fill(residual_rphi , ang_xy_L2*TMath::RadToDeg());
1505 h_test_residual_rphi_vs_angxy_L2_tpc->Fill(residual_rphi_tpc, ang_xy_L2*TMath::RadToDeg());
1506
1507 // save info for analysis
1508 double q = cluster_2d_q[closest_cluster2d_id];
1509 // 1d components
1510 int idx = cluster_2d_idx[closest_cluster2d_id];
1511 int idv = cluster_2d_idv[closest_cluster2d_id];
1512 double qx = cluster_1d_q[idx];
1513 double qv = cluster_1d_q[idv];
1514 double sizex = cluster_1d_size[idx];
1515 double sizev = cluster_1d_size[idv];
1516
1517 vector_test_residual_rphi.push_back(residual_rphi);
1518 vector_test_residual_z.push_back(residual_z);
1519 vector_test_charge.push_back(q);
1520 vector_test_charge_x.push_back(qx);
1521 vector_test_charge_v.push_back(qv);
1522 vector_test_size_x.push_back(sizex);
1523 vector_test_size_v.push_back(sizev);
1524
1525 vector_test_clusterid.push_back(closest_cluster2d_id);
1526 vector_test_idx.push_back(idx);
1527 vector_test_idv.push_back(idv);
1528 }
1529}
1530
1531bool CgemCosmicRayQA::apply_cuts() {
1532
1533 if(track_chi2 > cut_chi2) {
1534 nchi2++;
1535 return false;
1536 }
1537
1538 // ------------------------------
1539 // trackers
1540 bool stop_here=false;
1541 for(int iclu = 0; iclu < track_nfitpoint; iclu++) {
1542 int exp_lay = track_layerid[iclu];
1543 int exp_she = track_sheetid[iclu];
1544
1545 int idx = cluster_2d_idx[track_clusterid[iclu]];
1546 int idv = cluster_2d_idv[track_clusterid[iclu]];
1547 double qx = cluster_1d_q[idx];
1548 double qv = cluster_1d_q[idv];
1549
1550 if(exp_lay==0 && qx < cut_ene_L1_x) stop_here=true;
1551 if(exp_lay==0 && qv < cut_ene_L1_v) stop_here=true;
1552 if(exp_lay==1 && qx < cut_ene_L2_x) stop_here=true;
1553 if(exp_lay==1 && qv < cut_ene_L2_v) stop_here=true;
1554 }
1555 if(stop_here==true) {
1556 nstopped++;
1557 return false;
1558 }
1559 return true;
1560
1561}
1562// ------------------- ANALYSIS --------------------
1563void CgemCosmicRayQA::analyze() {
1564
1565 // hnoftrack->Fill(ntrack);
1566
1567 nentries_residual_rphi = h_test_residual_rphi->GetEntries();
1568 nentries_residual_z = h_test_residual_z->GetEntries();
1569
1570 // fit R * phi residual
1571 if(nentries_residual_rphi > 0) {
1572 TF1* f1 = new TF1("f1", "gaus", -5, 5);
1573 f1->SetParameters(1, 0, 1);
1574 h_test_residual_rphi->Fit("f1", "WRN");
1575 double param1[3];
1576 f1->GetParameters(&param1[0]);
1577 mean_residual_rphi = param1[1];
1578 sigma_residual_rphi = param1[2];
1579 delete f1;
1580 }
1581
1582 // fit z residual
1583 if(nentries_residual_z > 0) {
1584 TF1 *f2 = new TF1("f2", "gaus", -5, 5);
1585 f2->SetParameters(1, 0, 1);
1586 h_test_residual_z->Fit("f2", "WRN");
1587 double param2[3];
1588 f2->GetParameters(&param2[0]);
1589 mean_residual_z = param2[1];
1590 sigma_residual_z = param2[2];
1591 delete f2;
1592 }
1593
1594 // fit residual vs chi2
1595 TF1* f3 = new TF1("f3", "gaus", -5, 5); TF1* f4 = new TF1("f4", "gaus", -5, 5);
1596 TF1* f5 = new TF1("f5", "gaus", -5, 5); TF1* f6 = new TF1("f6", "gaus", -5, 5);
1597
1598 TH1D *h1_dummy = new TH1D("h1_dummy", "h1_dummy", 100, -5, 5); TH1D *h2_dummy = new TH1D("h2_dummy", "h2_dummy", 100, -5, 5);
1599 TH1D *h3_dummy = new TH1D("h3_dummy", "h3_dummy", 100, -5, 5); TH1D *h4_dummy = new TH1D("h4_dummy", "h4_dummy", 100, -5, 5);
1600
1601 int i_bin = 0;
1602 int j_dumbin = 0;
1603
1604 double param3[3]; double param4[3];
1605 double param5[3]; double param6[3];
1606
1607 for(int i_chi = 0; i_chi < CHI_BIN; i_chi++) {
1608 f3->SetParameters(1, 0, 1); f4->SetParameters(1, 0, 1);
1609 f5->SetParameters(1, 0, 1); f6->SetParameters(1, 0, 1);
1610
1611 i_bin = j_dumbin*j_dumbin + 1;
1612 if(i_chi/3 == 1) i_bin = i_bin*10 ;
1613 if(i_chi/3 == 2) i_bin = i_bin*100 ;
1614 if(i_chi/3 == 3) i_bin = i_bin*1000;
1615
1616 if(j_dumbin < 2) j_dumbin++;
1617 else j_dumbin = 0;
1618
1619 h1_dummy = h_test_res_rphi_vs_chi2->ProjectionX("h1_dummy", 1, i_bin);
1620 h2_dummy = h_test_res_z_vs_chi2 ->ProjectionX("h2_dummy", 1, i_bin);
1621
1622 h3_dummy = h_test_res_rphi_vs_chi2_tpc->ProjectionX("h3_dummy", 1, i_bin);
1623 h4_dummy = h_test_res_z_vs_chi2_tpc ->ProjectionX("h4_dummy", 1, i_bin);
1624
1625 if(h1_dummy->GetEntries() != 0) {
1626 h1_dummy->Fit("f3", "WRN");
1627 f3->GetParameters(&param3[0]);
1628 h_resolution_rphi_vs_chi2->SetBinContent(i_chi+1, param3[2]);
1629 }
1630 else h_resolution_rphi_vs_chi2->SetBinContent(i_chi+1, -1);
1631
1632 if(h2_dummy->GetEntries() != 0) {
1633 h2_dummy->Fit("f4", "WRN");
1634 f4->GetParameters(&param4[0]);
1635 h_resolution_z_vs_chi2 ->SetBinContent(i_chi+1, param4[2]);
1636 }
1637 else h_resolution_z_vs_chi2 ->SetBinContent(i_chi+1, -1);
1638
1639 if(h3_dummy->GetEntries() != 0) {
1640 h3_dummy->Fit("f5", "WRN");
1641 f5->GetParameters(&param5[0]);
1642 h_resolution_rphi_tpc_vs_chi2->SetBinContent(i_chi+1, param5[2]);
1643 }
1644 else h_resolution_rphi_tpc_vs_chi2->SetBinContent(i_chi+1, -1);
1645
1646 if(h4_dummy->GetEntries() != 0) {
1647 h4_dummy->Fit("f6", "WRN");
1648 f6->GetParameters(&param6[0]);
1649 h_resolution_z_tpc_vs_chi2 ->SetBinContent(i_chi+1, param6[2]);
1650 }
1651 else h_resolution_z_tpc_vs_chi2 ->SetBinContent(i_chi+1, -1);
1652
1653 h1_dummy->Reset(); h2_dummy->Reset();
1654 h3_dummy->Reset(); h4_dummy->Reset();
1655 }
1656
1657 for(int i_clean = 0; i_clean < 3; i_clean++) {
1658 param3[i_clean] = 0; param4[i_clean] = 0;
1659 param5[i_clean] = 0; param6[i_clean] = 0;
1660 }
1661 h1_dummy->Reset(); h2_dummy->Reset();
1662 h3_dummy->Reset(); h4_dummy->Reset();
1663
1664 // fit residual vs ang_xy
1665 for(int i_ang = 0; i_ang < 10; i_ang++) {
1666
1667 f3->SetParameters(1, 0, 1); f4->SetParameters(1, 0, 1);
1668 f5->SetParameters(1, 0, 1); f6->SetParameters(1, 0, 1);
1669
1670 h1_dummy = h_test_residual_rphi_vs_angxy_L1_cc->ProjectionX("h1_dummy", i_ang+1, i_ang+1);
1671 h2_dummy = h_test_residual_rphi_vs_angxy_L2_cc->ProjectionX("h2_dummy", i_ang+1, i_ang+1);
1672
1673 h3_dummy = h_test_residual_rphi_vs_angxy_L1_tpc->ProjectionX("h3_dummy", i_ang+1, i_ang+1);
1674 h4_dummy = h_test_residual_rphi_vs_angxy_L2_tpc->ProjectionX("h4_dummy", i_ang+1, i_ang+1);
1675
1676 if(h1_dummy->GetEntries() != 0) {
1677 h1_dummy->Fit("f3", "WRN");
1678 f3->GetParameters(&param3[0]);
1679 h_resolution_vs_ang_xy_L1->SetBinContent(i_ang+1, param3[2]);
1680 }
1681 else h_resolution_vs_ang_xy_L1->SetBinContent(i_ang+1, -1);
1682
1683 if(h2_dummy->GetEntries() != 0) {
1684 h2_dummy->Fit("f4", "WRN");
1685 f4->GetParameters(&param4[0]);
1686 h_resolution_vs_ang_xy_L2->SetBinContent(i_ang+1, param4[2]);
1687 }
1688 else h_resolution_vs_ang_xy_L2->SetBinContent(i_ang+1, -1);
1689
1690 if(h3_dummy->GetEntries() != 0) {
1691 h3_dummy->Fit("f5", "WRN");
1692 f5->GetParameters(&param5[0]);
1693 h_resolution_tpc_vs_ang_xy_L1->SetBinContent(i_ang+1, param5[2]);
1694 }
1695 else h_resolution_tpc_vs_ang_xy_L1->SetBinContent(i_ang+1, -1);
1696
1697 if(h4_dummy->GetEntries() != 0) {
1698 h4_dummy->Fit("f6", "WRN");
1699 f6->GetParameters(&param6[0]);
1700 h_resolution_tpc_vs_ang_xy_L2->SetBinContent(i_ang+1, param6[2]);
1701 }
1702 else h_resolution_tpc_vs_ang_xy_L2->SetBinContent(i_ang+1, -1);
1703
1704 h1_dummy->Reset(); h2_dummy->Reset();
1705 h3_dummy->Reset(); h4_dummy->Reset();
1706 }
1707 delete f3; delete f4;
1708 delete f5; delete f6;
1709
1710 delete h1_dummy; delete h2_dummy;
1711 delete h3_dummy; delete h4_dummy;
1712
1713 // inside and outside n-sigma on test plane
1714 int const nsigma_inside = 5; // CHECK cout to ttree
1715 int const nsigma_outside = 10; // CHECK cout to ttree
1716
1717 for(int ivec=0; ivec < vector_test_charge.size(); ivec++) {
1718
1719 double residual_rphi = vector_test_residual_rphi.at(ivec);
1720 double residual_z = vector_test_residual_z.at(ivec);
1721
1722 double q = vector_test_charge.at(ivec);
1723 double qx = vector_test_charge_x.at(ivec);
1724 double qv = vector_test_charge_v.at(ivec);
1725 double sizex = vector_test_size_x.at(ivec);
1726 double sizev = vector_test_size_v.at(ivec);
1727 double phi_t = vector_test_phi.at(ivec);
1728 double z_t = vector_test_z.at(ivec);
1729
1730 double t_chi2 = vector_chi2.at(ivec);
1731
1732 int clusterid = vector_test_clusterid.at(ivec);
1733 int idx = vector_test_idx.at(ivec);
1734 int idv = vector_test_idv.at(ivec);
1735
1736
1737 // signal
1738 if((residual_rphi > mean_residual_rphi - nsigma_inside * sigma_residual_rphi && residual_rphi < mean_residual_rphi + nsigma_inside * sigma_residual_rphi)
1739 && (residual_z > mean_residual_z - nsigma_inside * sigma_residual_z && residual_z < mean_residual_z + nsigma_inside * sigma_residual_z)) {
1740
1741 h_signal_charge_2d->Fill(q);
1742 h_signal_charge_1d_x->Fill(qx);
1743 h_signal_charge_1d_v->Fill(qv);
1744 h_signal_size_1d_x->Fill(sizex);
1745 h_signal_size_1d_v->Fill(sizev);
1746
1747 if(t_chi2 <= cut_chi2) {
1748 h_signal_charge_2d_chi2cut->Fill(q);
1749 h_signal_charge_1d_x_chi2cut->Fill(qx);
1750 h_signal_charge_1d_v_chi2cut->Fill(qv);
1751 h_signal_size_1d_x_chi2cut->Fill(sizex);
1752 h_signal_size_1d_v_chi2cut->Fill(sizev);
1753
1754 for(int jj = 0; jj < MAX_PHI_BIN; jj++) {
1755 if(phi_t >= MIN_PHI+(PHI_STEP*jj) && phi_t < MIN_PHI+(PHI_STEP*(jj+1))) n_insidetrack_phi_arr[jj]++;
1756 }
1757 for(int jj = 0; jj < MAX_Z_BIN; jj++) {
1758 if(z_t >= MIN_Z+(Z_STEP*jj) && z_t < MIN_Z+(Z_STEP*(jj+1))) n_insidetrack_z_arr[jj]++;
1759 }
1760
1761 for(int ii = 0; ii < MAX_PHI_BIN; ii++) {
1762 for(int jj = 0; jj < MAX_Z_BIN; jj++) {
1763 if((phi_t >= MIN_PHI+(PHI_STEP*ii) && phi_t < MIN_PHI+(PHI_STEP*(ii+1))) && (z_t >= MIN_Z+(Z_STEP*jj) && z_t < MIN_Z+(Z_STEP*(jj+1)))) n_insidetrack_phi_z_mat[ii][jj]++;
1764 }
1765 }
1766
1767 }
1768
1769 for(int chi_i = 0; chi_i < CHI_BIN; chi_i++) {
1770 for(int jj = 0; jj < MAX_PHI_BIN; jj++) {
1771 if(t_chi2 < chi2_cut_arr[chi_i] && (phi_t >= MIN_PHI+(PHI_STEP*jj) && phi_t < MIN_PHI+(PHI_STEP*(jj+1)))) n_insidetrack_phi_chi2_mat[chi_i][jj]++;
1772 }
1773 for(int jj = 0; jj < MAX_Z_BIN; jj++) {
1774 if(t_chi2 < chi2_cut_arr[chi_i] && (z_t >= MIN_Z+(Z_STEP*jj) && z_t < MIN_Z+(Z_STEP*(jj+1)))) n_insidetrack_z_chi2_mat[chi_i][jj]++;
1775 }
1776 }
1777
1778 ninside++;
1779 }
1780
1781 // outside
1782 if((residual_rphi < mean_residual_rphi - nsigma_outside * sigma_residual_rphi || residual_rphi > mean_residual_rphi + nsigma_outside * sigma_residual_rphi)
1783 && (residual_z < mean_residual_z - nsigma_outside * sigma_residual_z || residual_z > mean_residual_z + nsigma_outside * sigma_residual_z)) {
1784
1785 h_background_charge_2d->Fill(q);
1786 h_background_charge_1d_x->Fill(qx);
1787 h_background_charge_1d_v->Fill(qv);
1788 h_background_size_1d_x->Fill(sizex);
1789 h_background_size_1d_v->Fill(sizev);
1790
1791 noutside++;
1792 }
1793 }
1794
1795
1796 for(int i_bin_phi = 0; i_bin_phi < MAX_PHI_BIN; i_bin_phi++) {
1797 if(n_validtrack_phi_arr[i_bin_phi] != 0) h_eff_vs_phi->SetBinContent(i_bin_phi+1, ((double)n_insidetrack_phi_arr[i_bin_phi])/((double)n_validtrack_phi_arr[i_bin_phi]));
1798 else h_eff_vs_phi->SetBinContent(i_bin_phi+1, 0);
1799
1800 // cout << "Eff_n vs phi " << i_bin_phi << " " << n_insidetrack_phi_arr[i_bin_phi] << "/" << n_validtrack_phi_arr[i_bin_phi] << endl;
1801 // cout << "Eff vs phi " << i_bin_phi << " " << h_eff_vs_phi->GetBinContent(i_bin_phi+1) << endl;
1802 // cout << "Eff vs phi " << i_bin_phi << " " << ((double)n_insidetrack_phi_arr[i_bin_phi])/((double)n_validtrack_phi_arr[i_bin_phi]) << endl;
1803 }
1804
1805 for(int i_bin_z = 0; i_bin_z < MAX_Z_BIN; i_bin_z++) {
1806 if(n_validtrack_z_arr[i_bin_z] != 0) h_eff_vs_z->SetBinContent(i_bin_z+1, ((double)n_insidetrack_z_arr[i_bin_z])/((double)n_validtrack_z_arr[i_bin_z]));
1807 else h_eff_vs_z->SetBinContent(i_bin_z+1, 0);
1808
1809 // cout << "Eff_n vs z " << i_bin_z << " " << n_insidetrack_z_arr[i_bin_z] << "/" << n_validtrack_z_arr[i_bin_z] << endl;
1810 // cout << "Eff vs z " << i_bin_z << " " << h_eff_vs_z->GetBinContent(i_bin_z+1) << endl;
1811 // cout << "Eff vs z " << i_bin_z << " " << ((double)n_insidetrack_z_arr[i_bin_z])/((double)n_validtrack_z_arr[i_bin_z]) << endl;
1812 }
1813
1814 for(int i_bin_phi = 0; i_bin_phi < MAX_PHI_BIN; i_bin_phi++) {
1815 for(int i_bin_z = 0; i_bin_z < MAX_Z_BIN; i_bin_z++) {
1816 if(n_validtrack_phi_z_mat[i_bin_phi][i_bin_z] != 0) h_eff_vs_phi_vs_z->SetBinContent(i_bin_phi+1, i_bin_z+1, ((double)n_insidetrack_phi_z_mat[i_bin_phi][i_bin_z])/((double)n_validtrack_phi_z_mat[i_bin_phi][i_bin_z]));
1817 else h_eff_vs_phi_vs_z->SetBinContent(i_bin_phi+1, i_bin_z+1, 0);
1818
1819 // cout << "Eff_n vs phi_z " << i_bin_z << " " << i_bin_phi << " " << n_insidetrack_phi_z_mat[i_bin_phi][i_bin_z] << "/" << n_validtrack_phi_z_mat[i_bin_phi][i_bin_z] << endl;
1820 // cout << "Eff vs phi_z " << i_bin_z << " " << i_bin_phi << " " << ((double)n_insidetrack_phi_z_mat[i_bin_phi][i_bin_z])/((double)n_validtrack_phi_z_mat[i_bin_phi][i_bin_z]) << endl;
1821 }
1822 }
1823
1824 for(int i_bin_chi = 0; i_bin_chi < CHI_BIN; i_bin_chi++) {
1825 for(int i_bin_phi = 0; i_bin_phi < MAX_PHI_BIN; i_bin_phi++) {
1826 if(n_validtrack_phi_chi2_mat[i_bin_chi][i_bin_phi] != 0) h_eff_vs_phi_vs_chi2->SetBinContent(i_bin_phi+1, i_bin_chi+1, ((double)n_insidetrack_phi_chi2_mat[i_bin_chi][i_bin_phi])/((double)n_validtrack_phi_chi2_mat[i_bin_chi][i_bin_phi]));
1827 else h_eff_vs_phi_vs_chi2->SetBinContent(i_bin_phi+1, i_bin_chi+1, 0);
1828
1829 // cout << "Eff_n vs phi_chi " << i_bin_chi << " " << i_bin_phi << " " << n_insidetrack_phi_chi2_mat[i_bin_chi][i_bin_phi] << "/" << n_validtrack_phi_chi2_mat[i_bin_chi][i_bin_phi] << endl;
1830 // cout << "Eff vs phi_chi " << i_bin_chi << " " << i_bin_phi << " " << ((double)n_insidetrack_phi_chi2_mat[i_bin_chi][i_bin_phi])/((double)n_validtrack_phi_chi2_mat[i_bin_chi][i_bin_phi]) << endl;
1831 }
1832 for(int i_bin_z = 0; i_bin_z < MAX_Z_BIN; i_bin_z++) {
1833 if(n_validtrack_z_chi2_mat[i_bin_chi][i_bin_z] != 0) h_eff_vs_z_vs_chi2->SetBinContent(i_bin_z+1, i_bin_chi+1, ((double)n_insidetrack_z_chi2_mat[i_bin_chi][i_bin_z])/((double)n_validtrack_z_chi2_mat[i_bin_chi][i_bin_z]));
1834 else h_eff_vs_z_vs_chi2->SetBinContent(i_bin_z+1, i_bin_chi+1, 0);
1835
1836 // cout << "Eff_n vs z_chi " << i_bin_chi << " " << i_bin_z << " " << n_insidetrack_z_chi2_mat[i_bin_chi][i_bin_z] << "/" << n_validtrack_z_chi2_mat[i_bin_chi][i_bin_z] << endl;
1837 // cout << "Eff vs z_chi " << i_bin_chi << " " << i_bin_z << " " << ((double)n_insidetrack_z_chi2_mat[i_bin_chi][i_bin_z])/((double)n_validtrack_z_chi2_mat[i_bin_chi][i_bin_z]) << endl;
1838 }
1839 }
1840
1841}
1842
1843
1844void CgemCosmicRayQA::histo_cosmetics() {
1845
1846 // ---------------------------------------
1847 // reset range
1848 int lastbin = h_nofhit_L1_S1_x->FindLastBinAbove(0);
1849 h_nofhit_L1_S1_x->GetXaxis()->SetRangeUser(0, lastbin);
1850 lastbin = h_nofhit_L1_S2_x->FindLastBinAbove(0);
1851 h_nofhit_L1_S2_x->GetXaxis()->SetRangeUser(0, lastbin);
1852 lastbin = h_nofhit_L2_S1_x->FindLastBinAbove(0);
1853 h_nofhit_L2_S1_x->GetXaxis()->SetRangeUser(0, lastbin);
1854 lastbin = h_nofhit_L2_S2_x->FindLastBinAbove(0);
1855 h_nofhit_L2_S2_x->GetXaxis()->SetRangeUser(0, lastbin);
1856
1857 lastbin = h_nofhit_L1_S1_v->FindLastBinAbove(0);
1858 h_nofhit_L1_S1_v->GetXaxis()->SetRangeUser(0, lastbin);
1859 lastbin = h_nofhit_L2_S1_v->FindLastBinAbove(0);
1860 h_nofhit_L2_S1_v->GetXaxis()->SetRangeUser(0, lastbin);
1861 lastbin = h_nofhit_L2_S2_v->FindLastBinAbove(0);
1862 h_nofhit_L2_S2_v->GetXaxis()->SetRangeUser(0, lastbin);
1863
1864 lastbin = h_nofcluster1d_L1_S1_x->FindLastBinAbove(0);
1865 h_nofcluster1d_L1_S1_x->GetXaxis()->SetRangeUser(0, lastbin);
1866 lastbin = h_nofcluster1d_L1_S2_x->FindLastBinAbove(0);
1867 h_nofcluster1d_L1_S2_x->GetXaxis()->SetRangeUser(0, lastbin);
1868 lastbin = h_nofcluster1d_L2_S1_x->FindLastBinAbove(0);
1869 h_nofcluster1d_L2_S1_x->GetXaxis()->SetRangeUser(0, lastbin);
1870 lastbin = h_nofcluster1d_L2_S2_x->FindLastBinAbove(0);
1871 h_nofcluster1d_L2_S2_x->GetXaxis()->SetRangeUser(0, lastbin);
1872
1873 lastbin = h_nofcluster1d_L1_S1_v->FindLastBinAbove(0);
1874 h_nofcluster1d_L1_S1_v->GetXaxis()->SetRangeUser(0, lastbin);
1875 lastbin = h_nofcluster1d_L2_S1_v->FindLastBinAbove(0);
1876 h_nofcluster1d_L2_S1_v->GetXaxis()->SetRangeUser(0, lastbin);
1877 lastbin = h_nofcluster1d_L2_S2_v->FindLastBinAbove(0);
1878 h_nofcluster1d_L2_S2_v->GetXaxis()->SetRangeUser(0, lastbin);
1879
1880}
1881
1882int CgemCosmicRayQA::find_closest_exp_cluster2d(double fphi, double fz) {
1883
1884 int thisclusterid = -1;
1885 double thisdistance = 1000;
1886 for(int iclu = 0; iclu < ncluster; iclu++) {
1887
1888 if(cluster_2d_view[iclu] == -1) continue;
1889 if(cluster_2d_view[iclu] != 2) cout << "ERROR in find_closest_exp_cluster2d " << cluster_2d_view[iclu] << endl;
1890 if(cluster_2d_layerid[iclu] != track_test_layerid) continue;
1891 if(cluster_2d_sheetid[iclu] != track_test_sheetid) continue;
1892
1893 double closest_phi = cluster_2d_phi[iclu];
1894 double closest_z = cluster_2d_z[iclu];
1895
1896 double distance_rphi = fabs(closest_phi - fphi) *anode_mid_gap_radius[track_test_layerid];
1897 double distance_z = fabs(closest_z - fz);
1898 double distance = sqrt(distance_rphi*distance_rphi + distance_z*distance_z);
1899
1900 if(distance < thisdistance) {
1901 thisdistance = distance;
1902 thisclusterid = iclu;
1903 }
1904 }
1905
1906 return thisclusterid;
1907}
MAX_Z
#define MAX_Z
Definition CgemCosmicRayQA.h:66
MAXLENGTH_L2_x
#define MAXLENGTH_L2_x
Definition CgemCosmicRayQA.h:56
MAXNOFSHEET
#define MAXNOFSHEET
Definition CgemCosmicRayQA.h:48
CHI_BIN
#define CHI_BIN
Definition CgemCosmicRayQA.h:72
MAX_PHI
#define MAX_PHI
Definition CgemCosmicRayQA.h:65
MAXNOFSTRIP_L2_v
#define MAXNOFSTRIP_L2_v
Definition CgemCosmicRayQA.h:53
MAXNOFTRACKS
#define MAXNOFTRACKS
Definition CgemCosmicRayQA.h:42
MAXNOFSTRIP_L1_x
#define MAXNOFSTRIP_L1_x
Definition CgemCosmicRayQA.h:50
Z_STEP
#define Z_STEP
Definition CgemCosmicRayQA.h:76
MAXNOFCLUSTERS
#define MAXNOFCLUSTERS
Definition CgemCosmicRayQA.h:43
MAX_PHI_BIN
#define MAX_PHI_BIN
Definition CgemCosmicRayQA.h:70
CHI_MAX
#define CHI_MAX
Definition CgemCosmicRayQA.h:67
ANG_MAX
#define ANG_MAX
Definition CgemCosmicRayQA.h:68
MAXNOFSTRIP_L1_v
#define MAXNOFSTRIP_L1_v
Definition CgemCosmicRayQA.h:51
MIN_Z
#define MIN_Z
Definition CgemCosmicRayQA.h:61
MAX_Z_BIN
#define MAX_Z_BIN
Definition CgemCosmicRayQA.h:71
MAXLENGTH_L1_x
#define MAXLENGTH_L1_x
Definition CgemCosmicRayQA.h:55
MAXNOFHITS
#define MAXNOFHITS
Definition CgemCosmicRayQA.h:45
MIN_PHI
#define MIN_PHI
Definition CgemCosmicRayQA.h:60
MAXLENGTH_L2_v
#define MAXLENGTH_L2_v
Definition CgemCosmicRayQA.h:58
ANG_MIN
#define ANG_MIN
Definition CgemCosmicRayQA.h:63
PHI_STEP
#define PHI_STEP
Definition CgemCosmicRayQA.h:75
MAXNOFLAYER
#define MAXNOFLAYER
Definition CgemCosmicRayQA.h:47
MAXNOFSTRIP_L2_x
#define MAXNOFSTRIP_L2_x
Definition CgemCosmicRayQA.h:52
CHI_MIN
#define CHI_MIN
Definition CgemCosmicRayQA.h:62
ANG_BIN
#define ANG_BIN
Definition CgemCosmicRayQA.h:73
MAXLENGTH_L1_v
#define MAXLENGTH_L1_v
Definition CgemCosmicRayQA.h:57
f1
TFile * f1
Definition DataBase/tau_mode.c:4
nentries
Int_t nentries
Definition DataBase/tau_mode.c:8
q
****INTEGER imax DOUBLE PRECISION m_pi *DOUBLE PRECISION m_amfin DOUBLE PRECISION m_Chfin DOUBLE PRECISION m_Xenph DOUBLE PRECISION m_sinw2 DOUBLE PRECISION m_GFermi DOUBLE PRECISION m_MfinMin DOUBLE PRECISION m_ta2 INTEGER m_out INTEGER m_KeyFSR INTEGER m_KeyQCD *COMMON c_Semalib $ !copy of input $ !CMS energy $ !beam mass $ !final mass $ !beam charge $ !final charge $ !smallest final mass $ !Z mass $ !Z width $ !EW mixing angle $ !Gmu Fermi $ alphaQED at q
Definition KKsem.h:33
msgSvc
IMessageSvc * msgSvc()
Definition ServiceAccessor.h:13
CgemCosmicRayQA::CgemCosmicRayQA
CgemCosmicRayQA(const std::string &name, ISvcLocator *pSvcLocator)
Definition CgemCosmicRayQA.cxx:46
CgemCosmicRayQA::define_cluster_selected_histo
void define_cluster_selected_histo()
Definition CgemCosmicRayQA.cxx:388
CgemCosmicRayQA::find_closest_exp_cluster2d
int find_closest_exp_cluster2d(double fphi, double fz)
Definition CgemCosmicRayQA.cxx:1882
CgemCosmicRayQA::initialize
StatusCode initialize()
Definition CgemCosmicRayQA.cxx:644
CgemCosmicRayQA::define_cluster1d_from_fit_histo
void define_cluster1d_from_fit_histo()
Definition CgemCosmicRayQA.cxx:434
CgemCosmicRayQA::finalize
StatusCode finalize()
Definition CgemCosmicRayQA.cxx:829
CgemCosmicRayQA::read_file
bool read_file(TString name)
Definition CgemCosmicRayQA.cxx:71
CgemGeoAlign::getDx
double getDx(int layer_vir)
Definition CgemGeoAlign.h:69
CgemGeoAlign::getRx
double getRx(int layer_vir)
Definition CgemGeoAlign.h:72
CgemGeoAlign::getRy
double getRy(int layer_vir)
Definition CgemGeoAlign.h:73
CgemGeoAlign::getRz
double getRz(int layer_vir)
Definition CgemGeoAlign.h:74
CgemGeoAlign::getDz
double getDz(int layer_vir)
Definition CgemGeoAlign.h:71
CgemGeoAlign::point_transform
HepPoint3D point_transform(int layer_vir, HepPoint3D pos)
Definition CgemGeoAlign.cxx:241
CgemGeoAlign::getDy
double getDy(int layer_vir)
Definition CgemGeoAlign.h:70
CgemMidDriftPlane::getR
double getR(int layer)
Definition CgemMidDriftPlane.h:19
ICgemGeomSvc::getReadoutPlane
virtual CgemGeoReadoutPlane * getReadoutPlane(int iLayer, int iSheet) const =0
ICgemGeomSvc::getAlignPtr
virtual CgemGeoAlign * getAlignPtr() const =0
ICgemGeomSvc::getMidDriftPtr
virtual CgemMidDriftPlane * getMidDriftPtr() const =0