7 #include <TGraphErrors.h>
8 #include <TMultiGraph.h>
13 int ModulesStat(
TString pathG=
"/home/karavdina/soft/pandaroot/macro/lmd/testPixelAlignSIM/mom_15/")
19 double v_tr_sc[nDt]={100};
32 double v_rt_sc[nDa]={0};
35 const int nSdstep = 10;
36 const int nuSd = nSdGL/nSdstep;
40 for(
int iSc=0;iSc<nSc;iSc++){
41 for(
int iDa=0;iDa<nDa;iDa++){
45 double ev_mod_b[10][nSd][nuSd+1][nDt];
46 double ev_mod_a[10][nSd][nuSd+1][nDt];
48 double resx_mod_b[10][nSd][nuSd+1][nDt];
49 double resx_mod_a[10][nSd][nuSd+1][nDt];
50 double resy_mod_b[10][nSd][nuSd+1][nDt];
51 double resy_mod_a[10][nSd][nuSd+1][nDt];
53 TString resname= pathG+
"/AlignmentLimits_EvModStat";
54 resname+=v_rt_sc[iDa];
57 resname+=
"constsSample";
58 TString resname_pdf = resname+
".pdf";
59 TString resname_pdf_o = resname_pdf+
"(";
60 TString resname_pdf_c = resname_pdf+
")";
62 for(
int idt=0;idt<nDt;idt++){
63 for(
int imod=0;imod<nMod;imod++){
65 for(
int ijsum=1;ijsum<nuSd+1;ijsum++){
66 int isum = 1+nSdstep*(ijsum-1);
68 for(
int ipart=0;ipart<nSd;ipart++){
70 double tr_sc = v_tr_sc[idt];
71 double rt_sc = v_rt_sc[iDa];
81 path+=
"Samples/Sample";
86 TString name1 = path +
"/hitsRes_before_alignment_sector";
89 TString name2 = path +
"/hitsRes_after_alignment_sector";
92 TFile *
f1 =
new TFile(name1,
"READ");
93 TNtuple *nb = (TNtuple*)f1->Get(
"nhits");
95 TH1I *ntrksb = (TH1I*)f1->Get(
"trks");
96 ev_mod_b[imod][ipart][ijsum][idt] = 1e-3*double(ntrksb->GetEntries());
97 TH1D *hresx_b =
new TH1D(
"hresx_b",
"",1e3,-1.,1.);
98 nb->Project(
"hresx_b",
"xrec-xmc");
99 resx_mod_b[imod][ipart][ijsum][idt] = 1e4*(hresx_b->GetMean());
100 TH1D *hresy_b =
new TH1D(
"hresy_b",
"",1e3,-1.,1.);
101 nb->Project(
"hresy_b",
"yrec-ymc");
102 resy_mod_b[imod][ipart][ijsum][idt] = 1e4*(hresy_b->GetMean());
103 TFile *
f2 =
new TFile(name2,
"READ");
104 TNtuple *na = (TNtuple*)f2->Get(
"nhits");
106 TH1I *ntrksa = (TH1I*)f2->Get(
"trks");
107 ev_mod_a[imod][ipart][ijsum][idt] = 1e-3*double(ntrksa->GetEntries());
108 TH1D *hresx_a =
new TH1D(
"hresx_a",
"",1e3,-1.,1.);
109 na->Project(
"hresx_a",
"xrec-xmc");
110 resx_mod_a[imod][ipart][ijsum][idt] = 1e4*(hresx_a->GetMean());
111 TH1D *hresy_a =
new TH1D(
"hresy_a",
"",1e3,-1.,1.);
112 na->Project(
"hresy_a",
"yrec-ymc");
113 resy_mod_a[imod][ipart][ijsum][idt] = 1e4*(hresy_a->GetMean());
114 cout<<
"imod:"<<imod<<
" dt:"<<idt<<
" isum:"<<isum<<
" ipart:"<<ipart<<
" nEv:"<<ev_mod_b[imod][ipart][ijsum][idt]<<
" "<<ev_mod_a[imod][ipart][ijsum][idt]<<
" "<<endl;
123 double av_evMod_b[10][nuSd+1][nDt];
124 double av_evMod_a[10][nuSd+1][nDt];
125 double err_av_evMod_b[10][nuSd+1][nDt];
126 double err_av_evMod_a[10][nuSd+1][nDt];
127 double av_resxMod_b[10][nuSd+1][nDt];
128 double av_resxMod_a[10][nuSd+1][nDt];
129 double err_av_resxMod_b[10][nuSd+1][nDt];
130 double err_av_resxMod_a[10][nuSd+1][nDt];
131 double av_resyMod_b[10][nuSd+1][nDt];
132 double av_resyMod_a[10][nuSd+1][nDt];
133 double err_av_resyMod_b[10][nuSd+1][nDt];
134 double err_av_resyMod_a[10][nuSd+1][nDt];
137 for(
int ijsum=1;ijsum<nuSd+1;ijsum++){
139 int isum = 1+nSdstep*(ijsum-1);
140 for(
int idt=0;idt<nDt;idt++){
141 for(
int imod=0;imod<nMod;imod++){
142 av_evMod_b[imod][ijsum][idt]=0;
143 av_evMod_a[imod][ijsum][idt]=0;
144 err_av_evMod_b[imod][ijsum][idt]=0;
145 err_av_evMod_a[imod][ijsum][idt]=0;
147 av_resxMod_b[imod][ijsum][idt]=0;
148 av_resxMod_a[imod][ijsum][idt]=0;
149 err_av_resxMod_b[imod][ijsum][idt]=0;
150 err_av_resxMod_a[imod][ijsum][idt]=0;
151 av_resyMod_b[imod][ijsum][idt]=0;
152 av_resyMod_a[imod][ijsum][idt]=0;
153 err_av_resyMod_b[imod][ijsum][idt]=0;
154 err_av_resyMod_a[imod][ijsum][idt]=0;
159 double i_nS = 1./nSd;
160 double i_nS_r = 1./(nSd-1);
164 for(
int ijsum=1;ijsum<nuSd+1;ijsum++){
166 int isum = 1+nSdstep*(ijsum-1);
167 for(
int idt=0;idt<nDt;idt++){
168 for(
int imod=0;imod<nMod;imod++){
169 for(
int ipart=0;ipart<isum;ipart++){
170 av_evMod_b[imod][ijsum][idt]+= ev_mod_b[imod][ipart][ijsum][idt];
171 av_evMod_a[imod][ijsum][idt]+= ev_mod_a[imod][ipart][ijsum][idt];
173 av_resxMod_b[imod][ijsum][idt]+= resx_mod_b[imod][ipart][ijsum][idt];
174 av_resxMod_a[imod][ijsum][idt]+= resx_mod_a[imod][ipart][ijsum][idt];
175 av_resyMod_b[imod][ijsum][idt]+= resy_mod_b[imod][ipart][ijsum][idt];
176 av_resyMod_a[imod][ijsum][idt]+= resy_mod_a[imod][ipart][ijsum][idt];
179 av_resxMod_b[imod][ijsum][idt] *= 1./isum;
180 av_resxMod_a[imod][ijsum][idt] *= 1./isum;
181 av_resyMod_b[imod][ijsum][idt] *= 1./isum;
182 av_resyMod_a[imod][ijsum][idt] *= 1./isum;
192 for(
int ijsum=1;ijsum<nuSd+1;ijsum++){
194 int isum = 1+nSdstep*(ijsum-1);
195 for(
int idt=0;idt<nDt;idt++){
196 for(
int imod=0;imod<nMod;imod++){
197 for(
int ipart=0;ipart<isum;ipart++){
198 err_av_resxMod_b[imod][ijsum][idt] += TMath::Power((resx_mod_b[imod][ipart][ijsum][idt]-av_resxMod_b[imod][ijsum][idt]),2);
199 err_av_resxMod_a[imod][ijsum][idt] += TMath::Power((resx_mod_a[imod][ipart][ijsum][idt]-av_resxMod_a[imod][ijsum][idt]),2);
200 err_av_resyMod_b[imod][ijsum][idt] += TMath::Power((resy_mod_b[imod][ipart][ijsum][idt]-av_resyMod_b[imod][ijsum][idt]),2);
201 err_av_resyMod_a[imod][ijsum][idt] += TMath::Power((resy_mod_a[imod][ipart][ijsum][idt]-av_resyMod_a[imod][ijsum][idt]),2);
203 err_av_resyMod_b[imod][ijsum][idt] *=1./isum;
204 err_av_resyMod_b[imod][ijsum][idt] =
sqrt(err_av_resyMod_b[imod][ijsum][idt]);
205 err_av_resyMod_a[imod][ijsum][idt] *=1./isum;
206 err_av_resyMod_a[imod][ijsum][idt] =
sqrt(err_av_resyMod_a[imod][ijsum][idt]);
226 TFile *
f =
new TFile(out,
"RECREATE");
230 for(
int imod=0;imod<nMod;imod++){
231 TGraphErrors *gr_av_b[nuSd+1];
232 TGraphErrors *gr_av_a[nuSd+1];
233 TMultiGraph *mgr_av_b =
new TMultiGraph();
234 TMultiGraph *mgr_av_a =
new TMultiGraph();
235 TGraphErrors *gr_avResx_b[nuSd+1];
236 TGraphErrors *gr_avResx_a[nuSd+1];
237 TMultiGraph *mgr_avResx_b =
new TMultiGraph();
238 TMultiGraph *mgr_avResx_a =
new TMultiGraph();
239 TGraphErrors *gr_avResy_b[nuSd+1];
240 TGraphErrors *gr_avResy_a[nuSd+1];
241 TMultiGraph *mgr_avResy_b =
new TMultiGraph();
242 TMultiGraph *mgr_avResy_a =
new TMultiGraph();
243 TLegend *leg_b =
new TLegend(0.91,0.65,0.99,0.98);
244 leg_b->SetHeader(
"BEFORE");
245 TLegend *leg_a =
new TLegend(0.91,0.65,0.99,0.98);
246 leg_a->SetHeader(
"AFTER");
249 for(
int ijsum=nuSd;ijsum>0;ijsum--){
251 int isum = 1+nSdstep*(ijsum-1);
255 gr_av_b[ijsum] =
new TGraphErrors(nDt,v_tr_sc, av_evMod_b[imod][ijsum],0,err_av_evMod_b[imod][ijsum]);
256 gr_av_a[ijsum] =
new TGraphErrors(nDt,v_tr_sc, av_evMod_a[imod][ijsum],0,err_av_evMod_a[imod][ijsum]);
257 gr_avResx_b[ijsum] =
new TGraphErrors(nDt,v_tr_sc, av_resxMod_b[imod][ijsum],0,err_av_resxMod_b[imod][ijsum]);
258 gr_avResx_a[ijsum] =
new TGraphErrors(nDt,v_tr_sc, av_resxMod_a[imod][ijsum],0,err_av_resxMod_a[imod][ijsum]);
259 gr_avResx_b[ijsum] =
new TGraphErrors(nDt,v_tr_sc, av_resxMod_b[imod][ijsum],0,err_av_resxMod_b[imod][ijsum]);
260 gr_avResx_a[ijsum] =
new TGraphErrors(nDt,v_tr_sc, av_resxMod_a[imod][ijsum],0,err_av_resxMod_a[imod][ijsum]);
261 gr_avResy_b[ijsum] =
new TGraphErrors(nDt,v_tr_sc, av_resyMod_b[imod][ijsum],0,err_av_resyMod_b[imod][ijsum]);
262 gr_avResy_a[ijsum] =
new TGraphErrors(nDt,v_tr_sc, av_resyMod_a[imod][ijsum],0,err_av_resyMod_a[imod][ijsum]);
263 gr_avResy_b[ijsum] =
new TGraphErrors(nDt,v_tr_sc, av_resyMod_b[imod][ijsum],0,err_av_resyMod_b[imod][ijsum]);
264 gr_avResy_a[ijsum] =
new TGraphErrors(nDt,v_tr_sc, av_resyMod_a[imod][ijsum],0,err_av_resyMod_a[imod][ijsum]);
269 gr_av_b[ijsum]->SetName(name_ev+
"_b");
270 gr_av_b[ijsum]->SetMarkerStyle(20+ijsum);
271 gr_av_b[ijsum]->SetMarkerColor(kGreen-3);
272 gr_av_b[ijsum]->SetMarkerSize(1.5);
273 gr_av_a[ijsum]->SetName(name_ev+
"_a");
274 gr_av_a[ijsum]->SetMarkerStyle(20+ijsum);
275 gr_av_a[ijsum]->SetMarkerColor(kOrange+7);
276 gr_av_a[ijsum]->SetMarkerSize(1.5);
277 mgr_av_b->Add(gr_av_b[ijsum]);
278 mgr_av_a->Add(gr_av_a[ijsum]);
279 gr_av_b[ijsum]->Write();
280 gr_av_a[ijsum]->Write();
281 leg_b->AddEntry(gr_av_b[ijsum],sname,
"P");
282 leg_a->AddEntry(gr_av_a[ijsum],sname,
"P");
284 gr_avResx_b[ijsum]->SetName(name_ev+
"_b");
285 gr_avResx_b[ijsum]->SetMarkerStyle(20+ijsum);
286 gr_avResx_b[ijsum]->SetMarkerColor(kGreen-3);
287 gr_avResx_b[ijsum]->SetMarkerSize(1.5);
288 gr_avResx_a[ijsum]->SetName(name_ev+
"_a");
289 gr_avResx_a[ijsum]->SetMarkerStyle(20+ijsum);
290 gr_avResx_a[ijsum]->SetMarkerColor(kOrange+7);
291 gr_avResx_a[ijsum]->SetMarkerSize(1.5);
292 mgr_avResx_b->Add(gr_avResx_b[ijsum]);
293 mgr_avResx_a->Add(gr_avResx_a[ijsum]);
294 gr_avResx_b[ijsum]->Write();
295 gr_avResx_a[ijsum]->Write();
296 gr_avResy_b[ijsum]->SetName(name_ev+
"_b");
297 gr_avResy_b[ijsum]->SetMarkerStyle(20+ijsum);
298 gr_avResy_b[ijsum]->SetMarkerColor(kGreen-3);
299 gr_avResy_b[ijsum]->SetMarkerSize(1.5);
300 gr_avResy_a[ijsum]->SetName(name_ev+
"_a");
301 gr_avResy_a[ijsum]->SetMarkerStyle(20+ijsum);
302 gr_avResy_a[ijsum]->SetMarkerColor(kOrange+7);
303 gr_avResy_a[ijsum]->SetMarkerSize(1.5);
304 mgr_avResy_b->Add(gr_avResy_b[ijsum]);
305 mgr_avResy_a->Add(gr_avResy_a[ijsum]);
306 gr_avResy_b[ijsum]->Write();
307 gr_avResy_a[ijsum]->Write();
313 mgr_av_b->SetTitle(smname+
"_b");
314 mgr_av_a->SetTitle(smname+
"_a");
315 leg_a->SetFillColor(0);
316 leg_b->SetFillColor(0);
317 mgr_av_b->SetName(smname+
"_b");
318 mgr_av_a->SetName(smname+
"_a");
324 mgr_avResx_b->SetName(smnamex+
"_b");
325 mgr_avResx_a->SetName(smnamex+
"_a");
326 mgr_avResy_b->SetName(smnamey+
"_b");
327 mgr_avResy_a->SetName(smnamey+
"_a");
328 mgr_av_b->Draw(
"AP");
329 mgr_av_b->GetXaxis()->SetTitle(
"#delta_{t}, #mum");
330 mgr_av_b->GetYaxis()->SetTitle(
"events, 10^{3}");
332 cRes.Print(resname_pdf_o);
334 mgr_av_a->Draw(
"AP");
335 mgr_av_a->GetXaxis()->SetTitle(
"#delta_{t}, #mum");
336 mgr_av_a->GetYaxis()->SetTitle(
"events, 10^{3}");
338 cRes.Print(resname_pdf_o);
343 leg_b->SetName(
"LegBefore");
344 leg_a->SetName(
"LegAfter");
349 mgr_avResx_b->Draw(
"AP");
350 mgr_avResx_b->GetXaxis()->SetTitle(
"#delta_{t}, #mum");
351 mgr_avResx_b->GetYaxis()->SetTitle(
"res_{X}, mean [#mum]");
353 cRes.Print(resname_pdf_o);
355 mgr_avResx_a->Draw(
"AP");
356 mgr_avResx_a->GetXaxis()->SetTitle(
"#delta_{t}, #mum");
357 mgr_avResx_a->GetYaxis()->SetTitle(
"res_{X}, mean [#mum]");
359 cRes.Print(resname_pdf_o);
361 mgr_avResx_a->Write();
362 mgr_avResx_b->Write();
364 mgr_avResy_b->Draw(
"AP");
365 mgr_avResy_b->GetXaxis()->SetTitle(
"#delta_{t}, #mum");
366 mgr_avResy_b->GetYaxis()->SetTitle(
"res_{Y}, mean [#mum]");
368 cRes.Print(resname_pdf_o);
370 mgr_avResy_a->Draw(
"AP");
371 mgr_avResy_a->GetXaxis()->SetTitle(
"#delta_{t}, #mum");
372 mgr_avResy_a->GetYaxis()->SetTitle(
"res_{Y}, mean [#mum]");
374 cRes.Print(resname_pdf_o);
376 mgr_avResy_a->Write();
377 mgr_avResy_b->Write();
380 cRes.Print(resname_pdf_c);
int ModulesStat(TString pathG="/home/karavdina/soft/pandaroot/macro/lmd/testPixelAlignSIM/mom_15/")
friend F32vec4 sqrt(const F32vec4 &a)