CGEM BOSS 6.6.5.f
BESIII Offline Software System
Loading...
Searching...
No Matches
bin_peak3D.cxx
Go to the documentation of this file.
1#include <TChain.h>
2#include <TFile.h>
3#include <vector>
4#include <math.h>
5#include "TH1.h"
6#include "TH2.h"
7#include "TGraph.h"
8#include "TCanvas.h"
9#include <iostream>
10#include <fstream>
11#include <string>
12#include <sstream>
13#include <typeinfo.h>
14using namespace std;
15
16
17
18int bin_peak3D()
19{
20 double M_PI = 3.1415926;
21///*
22//////////////////////////////////////////////////////////
23// This file has been automatically generated
24// (Thu Feb 1 20:29:58 2018 by ROOT version5.34/09)
25// from TTree hit/hit
26// found on file: hough_hist_test.root
27//////////////////////////////////////////////////////////
28
29
30//Reset ROOT and connect tree file
31 //gROOT->Reset();
32 //TFile *f = (TFile*)gROOT->GetListOfFiles()->FindObject("hough_hist_test.root");
33 //if (!f) {
34 //f = new TFile("hough_hist_test.root");
35 //}
36 //f->GetObject("hit",tree);
37 TChain * hit = new TChain("hit");
38 hit->Add("hough_hist_test.root");
39
40//Declaration of leaves types
41 Int_t hit_run;
42 Int_t hit_evt;
43 Int_t hit_nhit;
44 Int_t hit_hitid[10000];
45 Int_t hit_layer[10000];
46 Int_t hit_wire[10000];
47 Double_t hit_x[10000];
48 Double_t hit_y[10000];
49 Double_t hit_z[10000];
50 Double_t hit_driftdist[10000];
51 Double_t hit_drifttime[10000];
52 Int_t hit_flag[10000];
53 Double_t hit_truth_x[10000];
54 Double_t hit_truth_y[10000];
55 Double_t hit_truth_z[10000];
56 Double_t hit_truth_drift[10000];
57 Int_t hit_truth_ambig[10000];
58
59 // Set branch addresses.
60 hit->SetBranchAddress("hit_run",&hit_run);
61 hit->SetBranchAddress("hit_evt",&hit_evt);
62 hit->SetBranchAddress("hit_nhit",&hit_nhit);
63 hit->SetBranchAddress("hit_hitid",hit_hitid);
64 hit->SetBranchAddress("hit_layer",hit_layer);
65 hit->SetBranchAddress("hit_wire",hit_wire);
66 hit->SetBranchAddress("hit_x",hit_x);
67 hit->SetBranchAddress("hit_y",hit_y);
68 hit->SetBranchAddress("hit_z",hit_z);
69 hit->SetBranchAddress("hit_driftdist",hit_driftdist);
70 hit->SetBranchAddress("hit_drifttime",hit_drifttime);
71 hit->SetBranchAddress("hit_flag",hit_flag);
72 hit->SetBranchAddress("hit_truth_x",hit_truth_x);
73 hit->SetBranchAddress("hit_truth_y",hit_truth_y);
74 hit->SetBranchAddress("hit_truth_z",hit_truth_z);
75 hit->SetBranchAddress("hit_truth_drift",hit_truth_drift);
76 hit->SetBranchAddress("hit_truth_ambig",hit_truth_ambig);
77
78// This is the loop skeleton
79// To read only selected branches, Insert statements like:
80// hit->SetBranchStatus("*",0); // disable all branches
81// TTreePlayer->SetBranchStatus("branchname",1); // activate branchname
82//*/
83
84//////////////////////////////////////////////////////////
85// This file has been automatically generated
86// (Thu Feb 1 20:29:46 2018 by ROOT version5.34/09)
87// from TTree hot/hot
88// found on file: hough_hist_test.root
89//////////////////////////////////////////////////////////
90
91
92//Reset ROOT and connect tree file
93 //gROOT->Reset();
94 //TFile *f = (TFile*)gROOT->GetListOfFiles()->FindObject("hough_hist_test.root");
95 //if (!f) {
96 //f = new TFile("hough_hist_test.root");
97 //}
98 //f->GetObject("hot",tree);
99 TChain * hot = new TChain("hot");
100 hot->Add("hough_hist_test.root");
101
102//Declaration of leaves types
103 Int_t hot_run;
104 Int_t hot_evt;
105 Int_t hot_trk;
106 Int_t hot_nhot;
107 Int_t hot_hitid[1000];
108 Int_t hot_layer[1000];
109 Int_t hot_wire[1000];
110 Double_t hot_x[1000];
111 Double_t hot_y[1000];
112 Double_t hot_z[1000];
113 Double_t hot_x0[1000];
114 Double_t hot_y0[1000];
115 Double_t hot_z0[1000];
116 Double_t hot_s0[1000];
117 Double_t hot_x1[1000];
118 Double_t hot_y1[1000];
119 Double_t hot_z1[1000];
120 Double_t hot_s1[1000];
121 Double_t hot_drift[1000];
122 Int_t hot_flag[1000];
123 Double_t hot_deltaD[1000];
124 Double_t hot_truth_x[1000];
125 Double_t hot_truth_y[1000];
126 Double_t hot_truth_z[1000];
127 Double_t hot_truth_drift[1000];
128 Int_t hot_truth_ambig[1000];
129
130 // Set branch addresses.
131 hot->SetBranchAddress("hot_run",&hot_run);
132 hot->SetBranchAddress("hot_evt",&hot_evt);
133 hot->SetBranchAddress("hot_trk",&hot_trk);
134 hot->SetBranchAddress("hot_nhot",&hot_nhot);
135 hot->SetBranchAddress("hot_hitid",hot_hitid);
136 hot->SetBranchAddress("hot_layer",hot_layer);
137 hot->SetBranchAddress("hot_wire",hot_wire);
138 hot->SetBranchAddress("hot_x",hot_x);
139 hot->SetBranchAddress("hot_y",hot_y);
140 hot->SetBranchAddress("hot_z",hot_z);
141 hot->SetBranchAddress("hot_x0",hot_x0);
142 hot->SetBranchAddress("hot_y0",hot_y0);
143 hot->SetBranchAddress("hot_z0",hot_z0);
144 hot->SetBranchAddress("hot_s0",hot_s0);
145 hot->SetBranchAddress("hot_x1",hot_x1);
146 hot->SetBranchAddress("hot_y1",hot_y1);
147 hot->SetBranchAddress("hot_z1",hot_z1);
148 hot->SetBranchAddress("hot_s1",hot_s1);
149 hot->SetBranchAddress("hot_drift",hot_drift);
150 hot->SetBranchAddress("hot_flag",hot_flag);
151 hot->SetBranchAddress("hot_deltaD",hot_deltaD);
152 hot->SetBranchAddress("hot_truth_x",hot_truth_x);
153 hot->SetBranchAddress("hot_truth_y",hot_truth_y);
154 hot->SetBranchAddress("hot_truth_z",hot_truth_z);
155 hot->SetBranchAddress("hot_truth_drift",hot_truth_drift);
156 hot->SetBranchAddress("hot_truth_ambig",hot_truth_ambig);
157
158// This is the loop skeleton
159// To read only selected branches, Insert statements like:
160// hot->SetBranchStatus("*",0); // disable all branches
161// TTreePlayer->SetBranchStatus("branchname",1); // activate branchname
162
163
164//////////////////////////////////////////////////////////
165// This file has been automatically generated
166// (Thu Feb 1 19:24:40 2018 by ROOT version5.34/09)
167// from TTree trk/trk
168// found on file: hough_hist_test.root
169//////////////////////////////////////////////////////////
170
171///*
172//Reset ROOT and connect tree file
173 //gROOT->Reset();
174 //TFile *f = (TFile*)gROOT->GetListOfFiles()->FindObject("hough_hist_test.root");
175 //if (!f) {
176 //f = new TFile("hough_hist_test.root");
177 //}
178 //TTree* tree;
179 //f->GetObject("trk",tree);
180 TChain * trk = new TChain("trk");
181 trk->Add("hough_hist_test.root");
182
183//Declaration of leaves types
184 Int_t trk_run;
185 Int_t trk_evt;
186 Int_t trk_ntrk;
187 Int_t trk_trackId[100];
188 Int_t trk_charge[100];
189 Double_t trk_dr[100];
190 Double_t trk_phi0[100];
191 Double_t trk_kappa[100];
192 Double_t trk_dz[100];
193 Double_t trk_tanl[100];
194 Double_t trk_pxy[100];
195 Double_t trk_px[100];
196 Double_t trk_py[100];
197 Double_t trk_pz[100];
198 Double_t trk_p[100];
199 Double_t trk_theta[100];
200 Double_t trk_phi[100];
201 Double_t trk_x[100];
202 Double_t trk_y[100];
203 Double_t trk_z[100];
204 Double_t trk_r[100];
205 Double_t trk_chi2[100];
206 Double_t trk_fiTerm[100];
207 Double_t trk_matchChi2[100];
208 Int_t trk_nhit[100];
209 Int_t trk_ncluster[100];
210 Int_t trk_stat[100];
211 Int_t trk_ndof[100];
212 Int_t trk_nster[100];
213 Int_t trk_nlayer[100];
214 Int_t trk_firstLayer[100];
215 Int_t trk_lastLayer[100];
216 Int_t trk_nCgemXClusters[100];
217 Int_t trk_nCgemVClusters[100];
218 Int_t trk_nhop[100];
219 Int_t trk_nhot[100];
220 Double_t trk_Xc[100];
221 Double_t trk_Yc[100];
222 Double_t trk_R[100];
223 Int_t trk_truth_charge[100];
224 Double_t trk_truth_dr[100];
225 Double_t trk_truth_phi0[100];
226 Double_t trk_truth_kappa[100];
227 Double_t trk_truth_dz[100];
228 Double_t trk_truth_tanl[100];
229 Double_t trk_truth_pxy[100];
230 Double_t trk_truth_px[100];
231 Double_t trk_truth_py[100];
232 Double_t trk_truth_pz[100];
233 Double_t trk_truth_p[100];
234 Double_t trk_truth_theta[100];
235 Double_t trk_truth_phi[100];
236 Double_t trk_truth_x[100];
237 Double_t trk_truth_y[100];
238 Double_t trk_truth_z[100];
239 Double_t trk_truth_r[100];
240 Double_t trk_truth_cosTheta[100];
241 Double_t trk_truth_Xc[100];
242 Double_t trk_truth_Yc[100];
243 Double_t trk_truth_R[100];
244
245 // Set branch addresses.
246 trk->SetBranchAddress("trk_run",&trk_run);
247 trk->SetBranchAddress("trk_evt",&trk_evt);
248 trk->SetBranchAddress("trk_ntrk",&trk_ntrk);
249 trk->SetBranchAddress("trk_trackId",trk_trackId);
250 trk->SetBranchAddress("trk_charge",trk_charge);
251 trk->SetBranchAddress("trk_dr",trk_dr);
252 trk->SetBranchAddress("trk_phi0",trk_phi0);
253 trk->SetBranchAddress("trk_kappa",trk_kappa);
254 trk->SetBranchAddress("trk_dz",trk_dz);
255 trk->SetBranchAddress("trk_tanl",trk_tanl);
256 trk->SetBranchAddress("trk_pxy",trk_pxy);
257 trk->SetBranchAddress("trk_px",trk_px);
258 trk->SetBranchAddress("trk_py",trk_py);
259 trk->SetBranchAddress("trk_pz",trk_pz);
260 trk->SetBranchAddress("trk_p",trk_p);
261 trk->SetBranchAddress("trk_theta",trk_theta);
262 trk->SetBranchAddress("trk_phi",trk_phi);
263 trk->SetBranchAddress("trk_x",trk_x);
264 trk->SetBranchAddress("trk_y",trk_y);
265 trk->SetBranchAddress("trk_z",trk_z);
266 trk->SetBranchAddress("trk_r",trk_r);
267 trk->SetBranchAddress("trk_chi2",trk_chi2);
268 trk->SetBranchAddress("trk_fiTerm",trk_fiTerm);
269 trk->SetBranchAddress("trk_matchChi2",trk_matchChi2);
270 trk->SetBranchAddress("trk_nhit",trk_nhit);
271 trk->SetBranchAddress("trk_ncluster",trk_ncluster);
272 trk->SetBranchAddress("trk_stat",trk_stat);
273 trk->SetBranchAddress("trk_ndof",trk_ndof);
274 trk->SetBranchAddress("trk_nster",trk_nster);
275 trk->SetBranchAddress("trk_nlayer",trk_nlayer);
276 trk->SetBranchAddress("trk_firstLayer",trk_firstLayer);
277 trk->SetBranchAddress("trk_lastLayer",trk_lastLayer);
278 trk->SetBranchAddress("trk_nCgemXClusters",trk_nCgemXClusters);
279 trk->SetBranchAddress("trk_nCgemVClusters",trk_nCgemVClusters);
280 trk->SetBranchAddress("trk_nhop",trk_nhop);
281 trk->SetBranchAddress("trk_nhot",trk_nhot);
282 trk->SetBranchAddress("trk_Xc",trk_Xc);
283 trk->SetBranchAddress("trk_Yc",trk_Yc);
284 trk->SetBranchAddress("trk_R",trk_R);
285 trk->SetBranchAddress("trk_truth_charge",trk_truth_charge);
286 trk->SetBranchAddress("trk_truth_dr",trk_truth_dr);
287 trk->SetBranchAddress("trk_truth_phi0",trk_truth_phi0);
288 trk->SetBranchAddress("trk_truth_kappa",trk_truth_kappa);
289 trk->SetBranchAddress("trk_truth_dz",trk_truth_dz);
290 trk->SetBranchAddress("trk_truth_tanl",trk_truth_tanl);
291 trk->SetBranchAddress("trk_truth_pxy",trk_truth_pxy);
292 trk->SetBranchAddress("trk_truth_px",trk_truth_px);
293 trk->SetBranchAddress("trk_truth_py",trk_truth_py);
294 trk->SetBranchAddress("trk_truth_pz",trk_truth_pz);
295 trk->SetBranchAddress("trk_truth_p",trk_truth_p);
296 trk->SetBranchAddress("trk_truth_theta",trk_truth_theta);
297 trk->SetBranchAddress("trk_truth_phi",trk_truth_phi);
298 trk->SetBranchAddress("trk_truth_x",trk_truth_x);
299 trk->SetBranchAddress("trk_truth_y",trk_truth_y);
300 trk->SetBranchAddress("trk_truth_z",trk_truth_z);
301 trk->SetBranchAddress("trk_truth_r",trk_truth_r);
302 trk->SetBranchAddress("trk_truth_cosTheta",trk_truth_cosTheta);
303 trk->SetBranchAddress("trk_truth_Xc",trk_truth_Xc);
304 trk->SetBranchAddress("trk_truth_Yc",trk_truth_Yc);
305 trk->SetBranchAddress("trk_truth_R",trk_truth_R);
306
307// This is the loop skeleton
308// To read only selected branches, Insert statements like:
309// trk->SetBranchStatus("*",0); // disable all branches
310// TTreePlayer->SetBranchStatus("branchname",1); // activate branchname
311//*/
312 TAxis *x_axis,*y_axis;
313 stringstream ssname;
314 string sname;
315 const char * name;
316 int point =10 ;
317 int yMax;
318 double peak_rate_mean_kb = 0;
319 double peak_rate_mean_tr = 0;
320 double peak_rate_RMS_kb = 0;
321 double peak_rate_RMS_tr = 0;
322 TGraphErrors *gr_peak_rate_mean_kb = new TGraphErrors();
323 TGraphErrors *gr_peak_rate_mean_tr = new TGraphErrors();
324 int ibin = 0;
325 int nbin = 100;
326 while(nbin <= 2000){
327 int x_bin = nbin;
328 int y_bin = nbin;
329 double x_max = 0.93/sqrt(1-0.93*0.93);
330 double y_max = 50;
331 TH1D *h_peak_rate_kb = new TH1D("peak_rate","peak_rate", 50,0,1.1);
332 TH1D *h_peak_rate_tr = new TH1D("peak_rate_tr","peak_rate_tr", 50,0,1.1);
333 TH1D *h_layersOnPeak_kb = new TH1D("layersOnPeak_","layersOnPeak_",50,0,50);
334 TH1D *h_layersOnPeak_tr = new TH1D("layersOnPeak_tr","layersOnPeak_tr",50,0,50);
335 TH1D *h_layersOfEvent = new TH1D("layersOfEvent","layersOfEvent",50,0,50);
336
337 Long64_t nentries = hot->GetEntries();
338 int event =11;
339 Long64_t nbytes = 0;
340 for (Long64_t i=0; i<nentries;i++)
341 //for (Long64_t i=0; i<1000;i++)
342 //for (Long64_t i=event; i<event+1;i++)
343 {
344 nbytes += hot->GetEntry(i);
345///*
346 TH2D* hough = new TH2D("hough","hough",x_bin,-x_max,x_max,y_bin,-y_max,y_max);
347 TH2D *houghhit[1000];
348 //TH2D *houghhit0[1000];
349 //TH2D *houghhit1[1000];
350 int layersOnPeak[1000];
351 int nhit(0),ihit(0);
352 for(int j=0;j<hot_nhot;j++){
353 if(hot_flag[j] ==0)continue;
354 if(hot_layer[j] <3 && hot_flag[j] ==1)continue;
355 nhit++;
356 //cout<<hot_layer[j]<<endl;
357 //fill_histogram(hot_x[j],hot_y[j],hot_drift[j],x_bin*2,hough);
358 layersOnPeak[ihit] = hot_layer[j];
359 h_layersOfEvent->Fill(hot_layer[j]);
360 ssname.str("");
361 ssname <<"layer:"<<hot_layer[j]<<" ,wire"<<hot_wire[j];
362 sname = ssname.str();
363 name = sname.c_str();
364 houghhit[ihit] = new TH2D(name,name,x_bin,-x_max,x_max,y_bin,-y_max,y_max);
365 fill_histogram(hot_s0[j],hot_z0[j],x_max,nbin*point,houghhit[ihit]);
366 //houghhit0[ihit] = new TH2D(name,name,x_bin,-x_max,x_max,y_bin,-y_max,y_max);
367 //fill_histogram(hot_s0[j],hot_z0[j],x_max,nbin,houghhit0[ihit]);
368 //hough->Add(houghhit0[ihit]);
369 if(hot_layer[j] >3)fill_histogram(hot_s1[j],hot_z1[j],x_max,nbin*point,houghhit[ihit]);
370 hough->Add(houghhit[ihit]);
371 //houghhit1[ihit] = new TH2D(name,name,x_bin,-x_max,x_max,y_bin,-y_max,y_max);
372 //fill_histogram(hot_s1[j],hot_z1[j],x_max,nbin,houghhit1[ihit]);
373 //hough->Add(houghhit1[ihit]);
374 ihit++;
375 }
376
377 //TCanvas *c1 = new TCanvas("houghspace","houghspace",1000,500 );
378 //c1->Divide(2,1);
379 //c1->cd(1);
380 //hough->Draw();
381 //c1->cd(2);
382 //hough->Draw("LEGO2Z");
383 //c1->SaveAs("houghspace.pdf");
384
385 int binx,biny,binz;
386 hough->GetMaximumBin(binx,biny,binz);
387 int height = hough->GetBinContent(binx,biny);
388 //double x = hough->GetXaxis()->GetBinCenter(binx);
389 //double y = hough->GetYaxis()->GetBinCenter(biny);
390 //cout<<x<<" "<<y<<endl;
391
392///*
393 //TCanvas *c;
394 int nhot(0);
395 for(int ihit=0;ihit<nhit;ihit++){
396 //ssname.str("");
397 //ssname <<ihit<<".pdf";
398 //sname = ssname.str();
399 //name = sname.c_str();
400 //c = new TCanvas(name,name,700,500);
401 //houghhit[ihit]->Draw();
402 //cout<<houghhit[ihit]->GetBinContent(binx,biny)<<endl;
403 //cout<<layersOnPeak[ihit]<<endl;
404 if(houghhit[ihit]->GetBinContent(binx,biny)>0 ){
405 nhot++;
406 h_layersOnPeak_kb->Fill(layersOnPeak[ihit]);
407 //cout<<layersOnPeak[ihit]<<endl;
408 }
409 //cout<<endl;
410 //cout<<layersOnPeak[ihit]<<" "<<houghhit[ihit]->GetBinContent(binx,biny)<<endl;
411 //c->SaveAs(name);
412 delete houghhit[ihit];
413 delete houghhit[ihit];
414 }
415 //cout<<nhit<<" "<<height<<" "<<nhot<<" "<<(double)height/(double)nhit<<" "<<(double)nhot/(double)nhit<<endl;
416 double peak_rate = (double)nhot/(double)nhit;
417 double peak_height_rate = (double)height/(double)nhit;
418 h_peak_rate_kb->Fill(peak_rate);
419 //h_peak_rate_tr->Fill(peak_height_rate);
420 //cout<<nhot<<" "<<height<<" "<<nhit<<endl;
421
422 delete hough;
423//*/
424
425///*
426 TH2D* hough2 = new TH2D("hough2","hough2",x_bin,-x_max,x_max,y_bin,-y_max,y_max);
427 TH2D *houghhit2[1000];
428 int nhit2(0),ihit2(0);
429 for(int j=0;j<hot_nhot;j++){
430 if(hot_flag[j] ==0)continue;
431 if(hot_layer[j] <3 && hot_flag[j] ==1)continue;
432 nhit2++;
433 ssname.str("");
434 ssname <<"2layer:"<<hot_layer[j]<<" ,wire"<<hot_wire[j];
435 sname = ssname.str();
436 name = sname.c_str();
437 houghhit2[ihit2] = new TH2D(name,name,x_bin,-x_max,x_max,y_bin,-y_max,y_max);
438 fill_histogram2(hot_s0[j],hot_z0[j],x_max,nbin*point,houghhit2[ihit2]);
439 if(hot_layer[j] >3)fill_histogram2(hot_s1[j],hot_z1[j],x_max,nbin*point,houghhit2[ihit2]);
440 hough2->Add(houghhit2[ihit2]);
441 ihit2++;
442 }
443
444///*
445
446 int binx2,biny2,binz2;
447 hough2->GetMaximumBin(binx2,biny2,binz2);
448 int nhot2(0);
449 for(int ihit2=0;ihit2<nhit2;ihit2++){
450 //ssname.str("");
451 //ssname <<ihit<<".pdf";
452 //sname = ssname.str();
453 //name = sname.c_str();
454 //c = new TCanvas(name,name,700,500);
455 //houghhit[ihit]->Draw();
456 //cout<<houghhit2[ihit2]->GetBinContent(binx,biny)<<endl;
457 //cout<<layersOnPeak[ihit]<<endl;
458 if(houghhit2[ihit2]->GetBinContent(binx2,biny2)>0 ){
459 nhot2++;
460 h_layersOnPeak_tr->Fill(layersOnPeak[ihit2]);
461 }
462 delete houghhit2[ihit2];
463 delete houghhit2[ihit2];
464 }
465 double peak_rate2 = (double)nhot2/(double)nhit2;
466 h_peak_rate_tr->Fill(peak_rate2);
467
468//*/
469 delete hough2;
470 }
471///*
472 ssname.str("");
473 ssname <<x_bin<<"bins_layersOnPeak";
474 sname = ssname.str();
475 name = sname.c_str();
476 TCanvas *c2 = new TCanvas(name,name,700,500 );
477 yMax = h_layersOnPeak_kb->GetMaximum() > h_layersOnPeak_tr->GetMaximum()? h_layersOnPeak_kb->GetMaximum():h_layersOnPeak_tr->GetMaximum();
478 yMax = yMax > h_layersOfEvent->GetMaximum()? yMax:h_layersOfEvent->GetMaximum();
479 h_layersOnPeak_kb->GetYaxis()->SetRangeUser(0,yMax*1.2);
480 h_layersOnPeak_kb->Draw();
481 h_layersOnPeak_kb->SetLineColor(4);
482 x_axis = h_layersOnPeak_kb->GetXaxis();
483 y_axis = h_layersOnPeak_kb->GetYaxis();
484 x_axis->SetTitle("Layer");
485 x_axis->CenterTitle();
486 y_axis->SetTitle("Entries");
487 y_axis->CenterTitle();
488 h_layersOnPeak_tr->Draw("same");
489 h_layersOnPeak_tr->SetLineColor(2);
490 h_layersOfEvent->Draw("same");
491 h_layersOfEvent->SetLineColor(1);
492 TLegend* leg2= new TLegend(0.15 ,0.77,0.40,0.88);
493 leg2->AddEntry(h_layersOfEvent,"all","l");
494 leg2->AddEntry(h_layersOnPeak_kb,"k - b","l");
495 leg2->AddEntry(h_layersOnPeak_tr,"#theta - #rho","l");
496 leg2->SetFillColor(kWhite);
497 leg2->SetLineColor(kWhite);
498 leg2->Draw();
499 ssname.str("");
500 ssname <<x_bin<<"bins_layersOnPeak_sz"<<".pdf";
501 sname = ssname.str();
502 name = sname.c_str();
503 c2->SaveAs(name);
504 ssname.str("");
505 ssname <<x_bin<<"bins_layersOnPeak_sz"<<".eps";
506 sname = ssname.str();
507 name = sname.c_str();
508 c2->SaveAs(name);
509 ssname.str("");
510 ssname <<x_bin<<"bins_layersOnPeak_sz"<<".png";
511 sname = ssname.str();
512 name = sname.c_str();
513 c2->SaveAs(name);
514
515
516
517
518
519 ssname.str("");
520 ssname <<x_bin<<"bins_rate";
521 sname = ssname.str();
522 name = sname.c_str();
523 TCanvas *c3 = new TCanvas(name,name,700,500 );
524 h_peak_rate_kb->Draw();
525 h_peak_rate_tr->Draw("same");
526 h_peak_rate_kb->SetLineColor(4);
527 h_peak_rate_tr->SetLineColor(2);
528 x_axis = h_peak_rate_kb->GetXaxis();
529 y_axis = h_peak_rate_kb->GetYaxis();
530 //x_axis->SetLimits(0,80 );
531 x_axis->SetTitle("Rate");
532 x_axis->CenterTitle();
533 //x_axis->SetTitleSize(0.05);
534 //x_axis->SetLabelSize(0.04);
535 //x_axis->SetTitleOffset(1.05);
536 y_axis->SetTitle("Entries");
537 y_axis->CenterTitle();
538 //y_axis->SetTitleSize(0.05);
539 //y_axis->SetLabelSize(0.04);
540 //y_axis->SetTitleOffset(1.25);
541 yMax = h_peak_rate_kb->GetMaximum() > h_peak_rate_tr->GetMaximum()? h_peak_rate_kb->GetMaximum():h_peak_rate_tr->GetMaximum();
542 h_peak_rate_kb->GetYaxis()->SetRangeUser(0,yMax*1.2);
543 TLegend* leg3= new TLegend(0.12 ,0.80,0.33,0.89);
544 leg3->AddEntry(h_peak_rate_kb,"k - b","l");
545 leg3->AddEntry(h_peak_rate_tr,"#theta - #rho","l");
546 leg3->SetFillColor(kWhite);
547 leg3->SetLineColor(kWhite);
548 leg3->Draw();
549 ssname.str("");
550 ssname <<x_bin<<"bins_rate_sz"<<".pdf";
551 sname = ssname.str();
552 name = sname.c_str();
553 c3->SaveAs(name);
554 ssname.str("");
555 ssname <<x_bin<<"bins_rate_sz"<<".eps";
556 sname = ssname.str();
557 name = sname.c_str();
558 c3->SaveAs(name);
559 ssname.str("");
560 ssname <<x_bin<<"bins_rate_sz"<<".png";
561 sname = ssname.str();
562 name = sname.c_str();
563 c3->SaveAs(name);
564 peak_rate_mean_kb = h_peak_rate_kb->GetMean();
565 peak_rate_mean_tr = h_peak_rate_tr->GetMean();
566 peak_rate_RMS_kb = h_peak_rate_kb->GetRMS();
567 peak_rate_RMS_tr = h_peak_rate_tr->GetRMS();
568 cout<<x_bin<<" "<<peak_rate_mean_kb<<" "<<peak_rate_RMS_kb<<" "<<peak_rate_mean_tr<<" "<<peak_rate_RMS_tr<<endl;
569 gr_peak_rate_mean_kb->SetPoint(ibin,x_bin,peak_rate_mean_kb);
570 gr_peak_rate_mean_kb->SetPointError(ibin,0,peak_rate_RMS_kb);
571 gr_peak_rate_mean_tr->SetPoint(ibin,x_bin,peak_rate_mean_tr);
572 gr_peak_rate_mean_tr->SetPointError(ibin,0,peak_rate_RMS_tr);
573
574 int bin;
575 if(nbin<500) bin = 50;
576 else if(nbin<1000) bin =100;
577 else if(nbin<2000)bin =200;
578 else if(nbin<5000)bin =500;
579 else bin = 1000;
580 nbin += bin;
581 ibin++;
582
583 delete h_layersOnPeak_kb;
584 delete h_layersOnPeak_tr;
585 delete h_layersOfEvent;
586 delete h_peak_rate_kb;
587 delete h_peak_rate_tr;
588//*/
589 }
590///*
591 TCanvas *c4 = new TCanvas("peak_rate_mean_kb","peak_rate_mean_kb",700,500 );
592 gr_peak_rate_mean_kb->Draw("AP");
593 x_axis = gr_peak_rate_mean_kb->GetXaxis();
594 y_axis = gr_peak_rate_mean_kb->GetYaxis();
595 //x_axis->SetLimits(0,80 );
596 x_axis->SetTitle("nbin");
597 x_axis->CenterTitle();
598 //x_axis->SetTitleSize(0.05);
599 //x_axis->SetLabelSize(0.04);
600 //x_axis->SetTitleOffset(1.05);
601 y_axis->SetTitle("mean rate");
602 y_axis->CenterTitle();
603 //y_axis->SetTitleSize(0.05);
604 //y_axis->SetLabelSize(0.04);
605 //y_axis->SetTitleOffset(1.25);
606 gr_peak_rate_mean_kb->SetMinimum(0.0);
607 gr_peak_rate_mean_kb->SetMaximum(1.2);
608 gr_peak_rate_mean_kb->SetMarkerSize(1.0);
609 gr_peak_rate_mean_kb->SetMarkerStyle(24);
610 gr_peak_rate_mean_kb->SetMarkerColor(4);
611 gr_peak_rate_mean_kb->SetLineColor(4);
612 gr_peak_rate_mean_tr->Draw("Psame");
613 gr_peak_rate_mean_tr->SetMarkerSize(1.0);
614 gr_peak_rate_mean_tr->SetMarkerStyle(24);
615 gr_peak_rate_mean_tr->SetMarkerColor(2);
616 gr_peak_rate_mean_tr->SetLineColor(2);
617 TLegend* leg4= new TLegend(0.70 ,0.75,0.85,0.85);
618 leg4->SetFillColor(kWhite);
619 leg4->SetLineColor(kWhite);
620 leg4->AddEntry(gr_peak_rate_mean_kb,"k - b","p");
621 leg4->AddEntry(gr_peak_rate_mean_tr,"#theta - #rho","p");
622 leg4->Draw();
623
624 c4->SaveAs("peak_rate_mean_sz.eps");
625 c4->SaveAs("peak_rate_mean_sz.png");
626 c4->SaveAs("peak_rate_mean_sz.pdf");
627//*/
628}
629///*
630void fill_histogram(double s, double z, double x_max, int nbin, TH2D* hough)
631{
632
633 double dx = 2*x_max/nbin;
634 double x = -x_max - 0.5*dx;
635 double y = -1*s*x+z;
636 double dy = -1*s*dx;
637 for(int i =0; i<nbin; i++){
638 x += dx;
639 y += dy;
640 hough->Fill(x,y);
641 }
642}
643//*/
644///*
645void fill_histogram2(double s, double z, double x_max, int nbin, TH2D* hough)
646{
647
648 double step = 2*x_max/nbin;
649 double theta = -x_max - 0.5*step;
650 double rho = 0;
651 for(int i =0; i<nbin; i++){
652 theta += step;
653 rho = s*cos(theta) + z*sin(theta);
654 hough->Fill(theta,rho);
655 }
656}
657//*/
658/*
659void fill_histogram(double x, double y, double r, int nbin, TH2D* hough)
660{
661 double M_PI = 3.1415926;
662 double U = 2*x/(x*x+y*y-r*r);
663 double V = 2*y/(x*x+y*y-r*r);
664 double R = 2*r/(x*x+y*y-r*r);
665 double delta_phi = 2.*M_PI/nbin;
666 double phi = -delta_phi/2;
667 for(int i =0; i<nbin; i++){
668 phi += delta_phi;
669 double u = U + R*cos(phi);
670 double v = V + R*sin(phi);
671 double normal = (v-V)/(u-U);
672 double k = -1./normal;
673 double b = v - k*u;
674 double x_cross = -b/(k+1/k);
675 double y_cross = b/(1+k*k);
676 // std::cout<<"x y centerX centerY "<<x<<" "<<y<<" "<<centerX<<" "<<centerY<<std::endl;
677 // std::cout<<"k b "<<k<<" "<<b<<std::endl;
678 // std::cout<<"xcross ycross "<<x_cross<<" "<<y_cross<<std::endl;
679
680 double rho=sqrt(x_cross*x_cross+y_cross*y_cross);
681 double theta=atan2(y_cross,x_cross);
682 if(theta<0) {
683 theta=theta+M_PI;
684 rho=-rho;
685 }
686 if( normal ==0 && x>0) {
687 rho = fabs(x);
688 theta = 0;
689 }
690 if( normal ==0 && x<0) {
691 rho = -fabs(x);
692 theta = M_PI;
693 }
694 hough->Fill(theta,rho);
695 }
696 //double slant = _y*cos(_theta)-_x*sin(_theta);
697 //_slant = slant;
698// std::cout<<"THETA RHO "<<theta_temp<<" "<<rho_temp<<std::endl;
699// std::cout<<std::endl;
700}
701*/
x
Double_t x[10]
Definition: DataBase/tau_mode.c:57
nentries
Int_t nentries
Definition: DataBase/tau_mode.c:8
bin
*******INTEGER m_nBinMax INTEGER m_NdiMax !No of bins in histogram for cell exploration division $ !Last vertex $ !Last active cell $ !Last cell in buffer $ !No of sampling when dividing cell $ !No of function total $ !Flag for random ceel for $ !Flag for type of for WtMax $ !Flag which decides whether vertices are included in the sampling $ entire domain is hyp !Maximum effective eevents per bin
Definition: FoamA.h:85
s
XmlRpcServer s
Definition: HelloServer.cpp:11
sin
double sin(const BesAngle a)
Definition: InstallArea/include/MdcGeom/MdcGeom/BesAngle.h:210
cos
double cos(const BesAngle a)
Definition: InstallArea/include/MdcGeom/MdcGeom/BesAngle.h:213
M_PI
#define M_PI
Definition: TConstant.h:4
fill_histogram2
void fill_histogram2(double s, double z, double x_max, int nbin, TH2D *hough)
Definition: bin_peak3D.cxx:645
bin_peak3D
int bin_peak3D()
Definition: bin_peak3D.cxx:18
fill_histogram
void fill_histogram(double s, double z, double x_max, int nbin, TH2D *hough)
Definition: bin_peak3D.cxx:630