BOSS 7.1.0
BESIII Offline Software System
Loading...
Searching...
No Matches
DedxCalibEAng.cxx
Go to the documentation of this file.
1#include "GaudiKernel/MsgStream.h"
2
3#include <sstream>
4#include <string>
5#include "TTree.h"
6#include "TH1F.h"
7#include "TF1.h"
8#include "TFile.h"
9#include "TStyle.h"
10#include "TCanvas.h"
11
12#include "DedxCalibAlg/DedxCalibEAng.h"
13DECLARE_COMPONENT(DedxCalibEAng)
14
15using namespace std;
16
17const double StepEAng = (PhiMax-PhiMin)/NumSlices;
18
19DedxCalibEAng::DedxCalibEAng(const std::string& name, ISvcLocator* pSvcLocator): DedxCalib(name, pSvcLocator){
20 declareProperty("Truncate", m_truncate = 0.7);
21 declareProperty("PMin", pMin=0.2);
22 declareProperty("PMax", pMax=2.3);
23}
24
25void DedxCalibEAng::BookHists()
26{
27 MsgStream log(msgSvc(), name());
28 log<<MSG::INFO<<"DedxCalibEAng::BookHists()"<<endreq;
29
30 ReadRecFileList();
31
32 m_eangle = new TH1F*[NumSlices];
33 m_trunc_eangle = new TH1F*[NumSlices];
34 m_pos_eangle = new TH1F*[NumSlices];
35 m_neg_eangle = new TH1F*[NumSlices];
36 stringstream hlabel;
37 for(int i=0;i<NumSlices;i++)
38 {
39 hlabel.str("");
40 hlabel<<"dEdx"<<"_eang"<<i;
41 m_eangle[i] = new TH1F(hlabel.str().c_str(), hlabel.str().c_str(), NumHistBins, MinHistValue, MaxHistValue);
42 hlabel.str("");
43 hlabel<<"trunc_dEdx"<<"_eang"<<i;
44 m_trunc_eangle[i] = new TH1F(hlabel.str().c_str(), hlabel.str().c_str(), NumHistBins, MinHistValue, MaxHistValue);
45 hlabel.str("");
46 hlabel<<"pos_dEdx"<<"_eang"<<i;
47 m_pos_eangle[i] = new TH1F(hlabel.str().c_str(), hlabel.str().c_str(), NumHistBins, MinHistValue, MaxHistValue);
48 hlabel.str("");
49 hlabel<<"neg_dEdx"<<"_eang"<<i;
50 m_neg_eangle[i] = new TH1F(hlabel.str().c_str(), hlabel.str().c_str(), NumHistBins, MinHistValue, MaxHistValue);
51 }
52 hlabel.str("");
53 hlabel<<"dEdxVsEAng";
54 m_EAngAverdE = new TH1F(hlabel.str().c_str(),"dEdx & EAng",NumSlices,PhiMin,PhiMax);
55 hlabel.str("");
56 hlabel<<"trunc_dEdxVsEAng";
57 m_trunc_EAngAverdE = new TH1F(hlabel.str().c_str(),"trunc_dEdx & EAng",NumSlices,PhiMin,PhiMax);
58 hlabel.str("");
59 hlabel<<"pos_dEdxVsEAng";
60 m_pos_EAngAverdE = new TH1F(hlabel.str().c_str(),"pos_dEdx & EAng",NumSlices,PhiMin,PhiMax);
61 hlabel.str("");
62 hlabel<<"neg_dEdxVsEAng";
63 m_neg_EAngAverdE = new TH1F(hlabel.str().c_str(),"neg_dEdx & EAng",NumSlices,PhiMin,PhiMax);
64
65}
66
67void DedxCalibEAng::FillHists()
68{
69 MsgStream log(msgSvc(), name());
70 log<<MSG::INFO<<"DedxCalibEAng::FillHists()"<<endreq;
71
72 TFile* f;
73 TTree* n102;
74 string runlist;
75
76 double dedx=0;
77 float runNO=0,evtNO=0,runFlag=0,pathlength=0,wid=0,layid=0,dd_in=0,driftdist=0,eangle=0,zhit=0,costheta=0,tes=0,ptrk=0,charge=0;
78 cout<<"m_recFileVector.size()= "<<m_recFileVector.size()<<endl;
79
80 vector<double> m_phlist(0);
81 vector<double> m_phlist_test(0);
82 vector<int> m_iEAng_list(0);
83 double old_costheta(0), old_ptrk(0);
84
85 for(unsigned int i=0; i<m_recFileVector.size(); i++)
86 {
87 runlist = m_recFileVector[i];
88 cout<<"runlist: "<<runlist.c_str()<<endl;
89 f = new TFile(runlist.c_str());
90 n102 = (TTree*)f->Get("n102");
91 n102-> SetBranchAddress("adc_raw",&dedx);
92 n102-> SetBranchAddress("path_rphi",&pathlength);
93 n102-> SetBranchAddress("wire",&wid);
94 n102-> SetBranchAddress("layer",&layid);
95 n102-> SetBranchAddress("doca_in",&dd_in);
96 n102-> SetBranchAddress("driftdist",&driftdist);
97 n102-> SetBranchAddress("eangle",&eangle);
98 n102-> SetBranchAddress("zhit",&zhit);
99 n102-> SetBranchAddress("runNO",&runNO);
100 n102-> SetBranchAddress("evtNO",&evtNO);
101 n102-> SetBranchAddress("runFlag",&runFlag);
102 n102-> SetBranchAddress("costheta1",&costheta);
103 n102-> SetBranchAddress("t01",&tes);
104 n102->SetBranchAddress("ptrk1",&ptrk);
105 n102->SetBranchAddress("charge",&charge);
106 for(int j=0;j<n102->GetEntries();j++)
107 {
108 n102->GetEntry(j);
109 // begin of save truncated hits information for EACH track
110 if((fabs(costheta-old_costheta)>0.001 || fabs(ptrk-old_ptrk>0.001)) && m_phlist.size()>=2){
111 sort(m_phlist_test.begin(), m_phlist_test.end());
112 unsigned int nsample = (unsigned int)(m_phlist_test.size()*m_truncate);
113 double m_trunc_value = m_phlist_test[nsample];
114 for(unsigned int i=0; i<m_phlist.size(); i++){
115 if(m_phlist[i]>m_trunc_value) continue;
116 m_trunc_eangle[m_iEAng_list[i]]->Fill(m_phlist[i]);
117 }
118 m_phlist.clear();
119 m_phlist_test.clear();
120 m_iEAng_list.clear();
121 }
122 old_costheta = costheta;
123 old_ptrk = ptrk;
124 // end of save truncated information
125 if(eangle >0.25) eangle = eangle -0.5;
126 else if(eangle <-0.25) eangle = eangle +0.5;
127 if(eangle>PhiMin && eangle<PhiMax)
128 {
129 if(ptrk>pMax || ptrk<pMin) continue;
130 if(tes>1400) continue;
131 // if (ptrk<0.1) continue;
132 if(layid <4) continue;
133 else if(layid <8)
134 {if(dedx<MinADCValuecut || dedx>MaxADCValuecut || pathlength>RPhi_PathMaxCut || pathlength<Iner_RPhi_PathMinCut || driftdist>Iner_DriftDistCut) continue;}
135 else
136 {if(dedx<MinADCValuecut || dedx>MaxADCValuecut || pathlength>RPhi_PathMaxCut || pathlength<Out_RPhi_PathMinCut || driftdist>Out_DriftDistCut) continue;}
137 int iEAng = 0;
138 iEAng = (int) ((eangle-PhiMin)/StepEAng);
139 dedx = exsvc->StandardHitCorrec(0,runFlag,2,runNO,evtNO,pathlength,wid,layid,dedx,dd_in,eangle,costheta);
140 dedx = exsvc->StandardTrackCorrec(0, runFlag, 2, runNO, evtNO, dedx, costheta, tes);
141 m_eangle[iEAng]->Fill(dedx);
142 m_phlist.push_back(dedx);
143 m_phlist_test.push_back(dedx);
144 m_iEAng_list.push_back(iEAng);
145 if(charge>0) m_pos_eangle[iEAng]->Fill(dedx);
146 if(charge<0) m_neg_eangle[iEAng]->Fill(dedx);
147 }
148 }
149 }
150}
151
152void DedxCalibEAng::AnalyseHists() // fit histograms
153{
154 MsgStream log(msgSvc(), name());
155 log<<MSG::INFO<<"DedxCalibEAng::AnalyseHists()"<<endreq;
156
157 TF1* func;
158 Double_t entry=0,mean=0,rms=0;
159 Double_t binmax=0;
160 Double_t lp[3]={0};
161 gStyle->SetOptFit(1111);
162 for(Int_t ip=0;ip<NumSlices;ip++)
163 {
164 entry = m_eangle[ip]->Integral();
165 if(entry<10) continue;
166 mean = m_eangle[ip]->GetMean();
167 rms = m_eangle[ip]->GetRMS();
168 binmax = m_eangle[ip]->GetMaximumBin();
169 lp[1] = m_eangle[ip]->GetBinCenter(binmax);
170 lp[2] = 200;
171 lp[0] = (m_eangle[ip]->GetBinContent(binmax)+m_eangle[ip]->GetBinContent(binmax-1)+m_eangle[ip]->GetBinContent(binmax+1))/3;
172
173 if(m_phShapeFlag==0)
174 {
175 func = new TF1("mylan",mylan,10,3000,4);
176 func->SetParLimits(0, 0, entry);
177 func->SetParLimits(1, 0, mean);
178 func->SetParLimits(2, 10, rms);
179 func->SetParLimits(3, -1, 0);
180 func->SetParameters(entry/36., mean-200, rms/2.3, -0.5);
181 }
182 else if(m_phShapeFlag==1)
183 {
184 func = new TF1("landaun",landaun,10,3000,3);
185 func->SetParameters(lp[0], lp[1], lp[2] );
186 }
187 else if(m_phShapeFlag==2)
188 {
189 func = new TF1("Landau",Landau,10,3000,2);
190 func->SetParameters(0.7*mean, rms );
191 }
192 else if(m_phShapeFlag==3)
193 {
194 func = new TF1("Vavilov",Vavilov,10,3000,2);
195 func->SetParameters(0.05, 0.6);
196 }
197 else if(m_phShapeFlag==4)
198 {
199 func = new TF1("AsymGauss",AsymGauss,10,3000,4);
200 func->SetParameters(entry, mean, rms, 1.0 );
201 }
202 func->SetLineWidth(2.1);
203 func->SetLineColor(2);
204
205 m_eangle[ip]->Fit(func, "rq" );// , "", mean-rms, mean+rms/2);
206 m_pos_eangle[ip]->Fit(func, "rq");// , "", mean-rms, mean+rms/2);
207 m_neg_eangle[ip]->Fit(func, "rq");// , "", mean-rms, mean+rms/2);
208 }
209}
210
211void DedxCalibEAng::WriteHists()
212{
213 MsgStream log(msgSvc(), name());
214 log<<MSG::INFO<<"DedxCalibEAng::WriteHists()"<<endreq;
215
216 double entryNo[NumSlices]={0},trunc_mean[NumSlices]={0},trunc_sigma[NumSlices]={0},fitmean[NumSlices]={0},fitmeanerr[NumSlices]={0},fitsigma[NumSlices]={0},gain[NumSlices]={0},chisq[NumSlices]={0};
217 double trunc_gain[NumSlices]={0};
218 double pos_entryNo[NumSlices]={0},pos_mean[NumSlices]={0},pos_sigma[NumSlices]={0},pos_fitmean[NumSlices]={0},pos_fitmeanerr[NumSlices]={0},pos_fitsigma[NumSlices]={0},pos_gain[NumSlices]={0},pos_chisq[NumSlices]={0};
219 double neg_entryNo[NumSlices]={0},neg_mean[NumSlices]={0},neg_sigma[NumSlices]={0},neg_fitmean[NumSlices]={0},neg_fitmeanerr[NumSlices]={0},neg_fitsigma[NumSlices]={0},neg_gain[NumSlices]={0},neg_chisq[NumSlices]={0};
220 double Ip_eangle[NumSlices]={0},eangle[NumSlices]={0};
221
222 double Norm(0), trunc_Norm(0);
223 int count(0);
224 for(Int_t ip=0;ip<NumSlices;ip++)
225 {
226 Ip_eangle[ip] = ip;
227 eangle[ip] = (ip+0.5)*StepEAng + PhiMin;
228
229 entryNo[ip]= m_eangle[ip]->Integral();
230 trunc_mean[ip] = m_trunc_eangle[ip]->GetMean();
231 trunc_sigma[ip] = m_trunc_eangle[ip]->GetRMS();
232 pos_entryNo[ip]= m_pos_eangle[ip]->Integral();
233 pos_mean[ip] = m_pos_eangle[ip]->GetMean();
234 pos_sigma[ip] = m_pos_eangle[ip]->GetRMS();
235 neg_entryNo[ip]= m_neg_eangle[ip]->Integral();
236 neg_mean[ip] = m_neg_eangle[ip]->GetMean();
237 neg_sigma[ip] = m_neg_eangle[ip]->GetRMS();
238 if(entryNo[ip]<10 || pos_entryNo[ip]<10 || neg_entryNo[ip]<10) continue;
239
240 if(m_phShapeFlag==0)
241 {
242 fitmean[ip] = m_eangle[ip]->GetFunction("mylan")->GetParameter(1);
243 fitmeanerr[ip] = m_eangle[ip]->GetFunction("mylan")->GetParError(1);
244 fitsigma[ip] = m_eangle[ip]->GetFunction("mylan")->GetParameter(2);
245 chisq[ip] = (m_eangle[ip]->GetFunction("mylan")-> GetChisquare())/(m_eangle[ip]->GetFunction("mylan")-> GetNDF());
246
247 pos_fitmean[ip] = m_pos_eangle[ip]->GetFunction("mylan")->GetParameter(1);
248 pos_fitmeanerr[ip] = m_pos_eangle[ip]->GetFunction("mylan")->GetParError(1);
249 pos_fitsigma[ip] = m_pos_eangle[ip]->GetFunction("mylan")->GetParameter(2);
250 pos_chisq[ip] = (m_pos_eangle[ip]->GetFunction("mylan")-> GetChisquare())/(m_pos_eangle[ip]->GetFunction("mylan")-> GetNDF());
251
252 neg_fitmean[ip] = m_neg_eangle[ip]->GetFunction("mylan")->GetParameter(1);
253 neg_fitmeanerr[ip] = m_neg_eangle[ip]->GetFunction("mylan")->GetParError(1);
254 neg_fitsigma[ip] = m_neg_eangle[ip]->GetFunction("mylan")->GetParameter(2);
255 neg_chisq[ip] = (m_neg_eangle[ip]->GetFunction("mylan")-> GetChisquare())/(m_neg_eangle[ip]->GetFunction("mylan")-> GetNDF());
256 }
257 if(m_phShapeFlag==1)
258 {
259 fitmean[ip] = m_eangle[ip]->GetFunction("landaun")->GetParameter(1);
260 fitmeanerr[ip] = m_eangle[ip]->GetFunction("landaun")->GetParError(1);
261 fitsigma[ip] = m_eangle[ip]->GetFunction("landaun")->GetParameter(2);
262 chisq[ip] = (m_eangle[ip]->GetFunction("landaun")-> GetChisquare())/(m_eangle[ip]->GetFunction("mylan")-> GetNDF());
263
264 pos_fitmean[ip] = m_pos_eangle[ip]->GetFunction("landaun")->GetParameter(1);
265 pos_fitmeanerr[ip] = m_pos_eangle[ip]->GetFunction("landaun")->GetParError(1);
266 pos_fitsigma[ip] = m_pos_eangle[ip]->GetFunction("landaun")->GetParameter(2);
267 pos_chisq[ip] = (m_pos_eangle[ip]->GetFunction("landaun")-> GetChisquare())/(m_pos_eangle[ip]->GetFunction("landaun")-> GetNDF());
268
269 neg_fitmean[ip] = m_neg_eangle[ip]->GetFunction("landaun")->GetParameter(1);
270 neg_fitmeanerr[ip] = m_neg_eangle[ip]->GetFunction("landaun")->GetParError(1);
271 neg_fitsigma[ip] = m_neg_eangle[ip]->GetFunction("landaun")->GetParameter(2);
272 neg_chisq[ip] = (m_neg_eangle[ip]->GetFunction("landaun")-> GetChisquare())/(m_neg_eangle[ip]->GetFunction("landaun")-> GetNDF());
273 }
274 if(m_phShapeFlag==2)
275 {
276 fitmean[ip] = m_eangle[ip]->GetFunction("Landau")->GetParameter(1);
277 fitmeanerr[ip] = m_eangle[ip]->GetFunction("Landau")->GetParError(1);
278 fitsigma[ip] = m_eangle[ip]->GetFunction("Landau")->GetParameter(2);
279 chisq[ip] = (m_eangle[ip]->GetFunction("Landau")-> GetChisquare())/(m_eangle[ip]->GetFunction("Landau")-> GetNDF());
280
281 pos_fitmean[ip] = m_pos_eangle[ip]->GetFunction("Landau")->GetParameter(1);
282 pos_fitmeanerr[ip] = m_pos_eangle[ip]->GetFunction("Landau")->GetParError(1);
283 pos_fitsigma[ip] = m_pos_eangle[ip]->GetFunction("Landau")->GetParameter(2);
284 pos_chisq[ip] = (m_pos_eangle[ip]->GetFunction("Landau")-> GetChisquare())/(m_pos_eangle[ip]->GetFunction("Landau")-> GetNDF());
285
286 neg_fitmean[ip] = m_neg_eangle[ip]->GetFunction("Landau")->GetParameter(1);
287 neg_fitmeanerr[ip] = m_neg_eangle[ip]->GetFunction("Landau")->GetParError(1);
288 neg_fitsigma[ip] = m_neg_eangle[ip]->GetFunction("Landau")->GetParameter(2);
289 neg_chisq[ip] = (m_neg_eangle[ip]->GetFunction("Landau")-> GetChisquare())/(m_neg_eangle[ip]->GetFunction("Landau")-> GetNDF());
290 }
291 if(m_phShapeFlag==3)
292 {
293 fitmean[ip] = m_eangle[ip]->GetFunction("Vavilov")->GetParameter(1);
294 fitmeanerr[ip] = m_eangle[ip]->GetFunction("Vavilov")->GetParError(1);
295 fitsigma[ip] = m_eangle[ip]->GetFunction("Vavilov")->GetParameter(2);
296 chisq[ip] = (m_eangle[ip]->GetFunction("Vavilov")-> GetChisquare())/(m_eangle[ip]->GetFunction("Vavilov")-> GetNDF());
297
298 pos_fitmean[ip] = m_pos_eangle[ip]->GetFunction("Vavilov")->GetParameter(1);
299 pos_fitmeanerr[ip] = m_pos_eangle[ip]->GetFunction("Vavilov")->GetParError(1);
300 pos_fitsigma[ip] = m_pos_eangle[ip]->GetFunction("Vavilov")->GetParameter(2);
301 pos_chisq[ip] = (m_pos_eangle[ip]->GetFunction("Vavilov")-> GetChisquare())/(m_pos_eangle[ip]->GetFunction("Vavilov")-> GetNDF());
302
303 neg_fitmean[ip] = m_neg_eangle[ip]->GetFunction("Vavilov")->GetParameter(1);
304 neg_fitmeanerr[ip] = m_neg_eangle[ip]->GetFunction("Vavilov")->GetParError(1);
305 neg_fitsigma[ip] = m_neg_eangle[ip]->GetFunction("Vavilov")->GetParameter(2);
306 neg_chisq[ip] = (m_neg_eangle[ip]->GetFunction("Vavilov")-> GetChisquare())/(m_neg_eangle[ip]->GetFunction("Vavilov")-> GetNDF());
307 }
308 if(m_phShapeFlag==4)
309 {
310 fitmean[ip] = m_eangle[ip]->GetFunction("AsymGauss")->GetParameter(1);
311 fitmeanerr[ip] = m_eangle[ip]->GetFunction("AsymGauss")->GetParError(1);
312 fitsigma[ip] = m_eangle[ip]->GetFunction("AsymGauss")->GetParameter(2);
313 chisq[ip] = (m_eangle[ip]->GetFunction("AsymGauss")-> GetChisquare())/(m_eangle[ip]->GetFunction("AsymGauss")-> GetNDF());
314
315 pos_fitmean[ip] = m_pos_eangle[ip]->GetFunction("AsymGauss")->GetParameter(1);
316 pos_fitmeanerr[ip] = m_pos_eangle[ip]->GetFunction("AsymGauss")->GetParError(1);
317 pos_fitsigma[ip] = m_pos_eangle[ip]->GetFunction("AsymGauss")->GetParameter(2);
318 pos_chisq[ip] = (m_pos_eangle[ip]->GetFunction("AsymGauss")-> GetChisquare())/(m_pos_eangle[ip]->GetFunction("AsymGauss")-> GetNDF());
319
320 neg_fitmean[ip] = m_neg_eangle[ip]->GetFunction("AsymGauss")->GetParameter(1);
321 neg_fitmeanerr[ip] = m_neg_eangle[ip]->GetFunction("AsymGauss")->GetParError(1);
322 neg_fitsigma[ip] = m_neg_eangle[ip]->GetFunction("AsymGauss")->GetParameter(2);
323 neg_chisq[ip] = (m_neg_eangle[ip]->GetFunction("AsymGauss")-> GetChisquare())/(m_neg_eangle[ip]->GetFunction("AsymGauss")-> GetNDF());
324 }
325 if(entryNo[ip]>1500 && fitmean[ip]>0 && fitmeanerr[ip]>0 && fitmeanerr[ip]<10 && fitsigma[ip]>0 && fitsigma[ip]<200){
326 Norm += fitmean[ip];
327 trunc_Norm += trunc_mean[ip];
328 count++;
329 }
330 }
331 Norm = Norm/count;
332 trunc_Norm = trunc_Norm/count;
333 // Norm = 550; // use an absolute normalization
334 cout<<"count= "<<count<<" average value: "<< Norm << " truncted: " << trunc_Norm << endl;
335 for(Int_t i=0;i<NumSlices;i++)
336 { // the difference between gain and mean is the normalization
337 if(entryNo[i]>10 && fitmean[i]>0 && fitsigma[i]>0) gain[i]=fitmean[i]/Norm;
338 else gain[i] = 1;
339 if(entryNo[i]>10) trunc_gain[i]=trunc_mean[i]/trunc_Norm;
340 else trunc_mean[i] = 1;
341 if(pos_entryNo[i]>10 && pos_fitmean[i]>0 && pos_fitsigma[i]>0) pos_gain[i]=pos_fitmean[i]/Norm;
342 else pos_gain[i] = 1;
343 if(neg_entryNo[i]>10 && neg_fitmean[i]>0 && neg_fitsigma[i]>0) neg_gain[i]=neg_fitmean[i]/Norm;
344 else neg_gain[i] = 1;
345 if(gain[i] != 1){
346 m_EAngAverdE -> SetBinContent(i+1, gain[i]);
347 m_EAngAverdE->SetBinError(i+1, fitmeanerr[i]/Norm);
348 }
349 if(trunc_gain[i] != 1){
350 m_trunc_EAngAverdE -> SetBinContent(i+1, trunc_gain[i]);
351 m_trunc_EAngAverdE->SetBinError(i+1, trunc_sigma[i]/trunc_Norm);
352 }
353 if(pos_gain[i] != 1){
354 m_pos_EAngAverdE -> SetBinContent(i+1, pos_gain[i]);
355 m_pos_EAngAverdE->SetBinError(i+1, pos_fitmeanerr[i]/Norm);
356 }
357 if(neg_gain[i] != 1){
358 m_neg_EAngAverdE -> SetBinContent(i+1, neg_gain[i]);
359 m_neg_EAngAverdE->SetBinError(i+1, neg_fitmeanerr[i]/Norm);
360 }
361 }
362
363 log<<MSG::INFO<<"begin getting calibration constant!!! "<<endreq;
364 double denangle[NumSlices]={0};
365 int denangle_entry[1] = {100};
366 m_EAngAverdE->Fit("pol1","r","",-0.25,0.25);
367 double par0 = m_EAngAverdE->GetFunction("pol1")->GetParameter(0);
368 double par1 = m_EAngAverdE->GetFunction("pol1")->GetParameter(1);
369 // double par2 = m_EAngAverdE->GetFunction("pol2")->GetParameter(2);
370 // par1 = -par1; // a trick is used here, because we don't really understand
371 // for(int i=0;i<NumSlices;i++) denangle[i] = (par0+par1*eangle[i])/par0;//+par2*pow(eangle[i],2))/par0;
372 for(int i=0;i<NumSlices;i++) denangle[i] = gain[i]; // use raw values instead of fit
373
374 log<<MSG::INFO<<"begin generating root file!!! "<<endreq;
375 TFile* f = new TFile(m_rootfile.c_str(), "recreate");
376 for(Int_t ip=0;ip<NumSlices;ip++){
377 m_eangle[ip]->Write();
378 m_trunc_eangle[ip]->Write();
379 m_pos_eangle[ip]->Write();
380 m_neg_eangle[ip]->Write();
381 }
382 m_EAngAverdE->Write();
383 m_trunc_EAngAverdE->Write();
384 m_pos_EAngAverdE->Write();
385 m_neg_EAngAverdE->Write();
386
387 TTree* entracalib = new TTree("entracalib","entracalib");
388 entracalib->Branch("1denangle", denangle, "denangle[100]/D");
389 entracalib->Branch("1denangle_entry", denangle_entry, "denangle_entry[1]/I");
390 entracalib->Fill();
391 entracalib-> Write();
392
393 TTree* ddgcalib = new TTree("ddgcalib", "ddgcalib");
394 ddgcalib -> Branch("hits", entryNo, "entryNo[100]/D");
395 ddgcalib -> Branch("truncmean", trunc_mean, "trunc_mean[100]/D");
396 ddgcalib -> Branch("truncsigma", trunc_sigma, "trunc_sigma[100]/D");
397 ddgcalib -> Branch("fitmean", fitmean, "fitmean[100]/D");
398 ddgcalib -> Branch("fitmeanerr", fitmeanerr, "fitmeanerr[100]/D");
399 ddgcalib -> Branch("chi", fitsigma, "fitsigma[100]/D");
400 ddgcalib -> Branch("gain", gain, "gain[100]/D");
401 ddgcalib -> Branch("truncgain", trunc_gain, "trunc_gain[100]/D");
402 ddgcalib -> Branch("chisq", chisq, "chisq[100]/D");
403 ddgcalib -> Branch("pos_hits", pos_entryNo, "pos_entryNo[100]/D");
404 ddgcalib -> Branch("pos_mean", pos_mean, "pos_mean[100]/D");
405 ddgcalib -> Branch("pos_sigma", pos_sigma, "pos_sigma[100]/D");
406 ddgcalib -> Branch("pos_fitmean", pos_fitmean, "pos_fitmean[100]/D");
407 ddgcalib -> Branch("pos_fitmeanerr", pos_fitmeanerr, "pos_fitmeanerr[100]/D");
408 ddgcalib -> Branch("pos_chi", pos_fitsigma, "pos_fitsigma[100]/D");
409 ddgcalib -> Branch("pos_gain", pos_gain, "pos_gain[100]/D");
410 ddgcalib -> Branch("pos_chisq", pos_chisq, "pos_chisq[100]/D");
411 ddgcalib -> Branch("neg_hits", neg_entryNo, "neg_entryNo[100]/D");
412 ddgcalib -> Branch("neg_mean", neg_mean, "neg_mean[100]/D");
413 ddgcalib -> Branch("neg_sigma", neg_sigma, "neg_sigma[100]/D");
414 ddgcalib -> Branch("neg_fitmean", neg_fitmean, "neg_fitmean[100]/D");
415 ddgcalib -> Branch("neg_fitmeanerr", neg_fitmeanerr, "neg_fitmeanerr[100]/D");
416 ddgcalib -> Branch("neg_chi", neg_fitsigma, "neg_fitsigma[100]/D");
417 ddgcalib -> Branch("neg_gain", neg_gain, "neg_gain[100]/D");
418 ddgcalib -> Branch("neg_chisq", neg_chisq, "neg_chisq[100]/D");
419 ddgcalib -> Branch("Ip_eangle", Ip_eangle, "Ip_eangle[100]/D");
420 ddgcalib -> Branch("eangle", eangle, "eangle[100]/D");
421 ddgcalib -> Fill();
422 ddgcalib -> Write();
423 f->Close();
424
425 TCanvas c1("c1", "canvas", 500, 400);
426 c1.Print((m_rootfile+".ps[").c_str());
427 m_EAngAverdE->Draw();
428 c1.Print((m_rootfile+".ps").c_str());
429 m_trunc_EAngAverdE->Draw();
430 c1.Print((m_rootfile+".ps").c_str());
431 m_pos_EAngAverdE->Draw();
432 c1.Print((m_rootfile+".ps").c_str());
433 m_neg_EAngAverdE->Draw();
434 c1.Print((m_rootfile+".ps").c_str());
435 for(Int_t ip=0;ip<NumSlices;ip++)
436 {
437 m_eangle[ip]->Draw();
438 c1.Print((m_rootfile+".ps").c_str());
439 }
440 for(Int_t ip=0;ip<NumSlices;ip++)
441 {
442 m_trunc_eangle[ip]->Draw();
443 c1.Print((m_rootfile+".ps").c_str());
444 }
445 for(Int_t ip=0;ip<NumSlices;ip++)
446 {
447 m_pos_eangle[ip]->Draw();
448 c1.Print((m_rootfile+".ps").c_str());
449 }
450 for(Int_t ip=0;ip<NumSlices;ip++)
451 {
452 m_neg_eangle[ip]->Draw();
453 c1.Print((m_rootfile+".ps").c_str());
454 }
455 c1.Print((m_rootfile+".ps]").c_str());
456
457 for(Int_t ip=0;ip<NumSlices;ip++)
458 {
459 delete m_eangle[ip];
460 delete m_trunc_eangle[ip];
461 delete m_pos_eangle[ip];
462 delete m_neg_eangle[ip];
463 }
464 delete m_EAngAverdE;
465 delete m_trunc_EAngAverdE;
466 delete m_pos_EAngAverdE;
467 delete m_neg_EAngAverdE;
468
469}
Branch
curve Branch("CurveSize",&CurveSize,"CurveSize/I")
Fill
curve Fill()
Write
curve Write()
f
TFile f("ana_bhabha660a_dqa_mcPat_zy_old.root")
SetBranchAddress
data SetBranchAddress("time",&time)
StepEAng
const double StepEAng
Definition: DedxCalibEAng.cxx:17
Landau
double Landau(double *x, double *par)
Definition: DedxCalibParameters.cxx:34
MaxADCValuecut
#define MaxADCValuecut
Definition: DedxCalibParameters.h:54
Vavilov
double Vavilov(double *x, double *par)
Definition: DedxCalibParameters.cxx:41
NumSlices
#define NumSlices
Definition: DedxCalibParameters.h:61
Iner_DriftDistCut
#define Iner_DriftDistCut
Definition: DedxCalibParameters.h:58
RPhi_PathMaxCut
#define RPhi_PathMaxCut
Definition: DedxCalibParameters.h:55
MinHistValue
#define MinHistValue
Definition: DedxCalibParameters.h:46
landaun
double landaun(double *x, double *par)
Definition: DedxCalibParameters.cxx:25
Out_DriftDistCut
#define Out_DriftDistCut
Definition: DedxCalibParameters.h:59
AsymGauss
double AsymGauss(double *x, double *par)
Definition: DedxCalibParameters.cxx:50
PhiMax
#define PhiMax
Definition: DedxCalibParameters.h:67
PhiMin
#define PhiMin
Definition: DedxCalibParameters.h:66
NumHistBins
#define NumHistBins
Definition: DedxCalibParameters.h:45
MaxHistValue
#define MaxHistValue
Definition: DedxCalibParameters.h:47
mylan
double mylan(double *x, double *par)
Definition: DedxCalibParameters.cxx:17
count
DOUBLE_PRECISION count[3]
Definition: EvtPhokharaDef.hh:67
msgSvc
IMessageSvc * msgSvc()
Definition: ServiceAccessor.h:13
DedxCalibEAng::DedxCalibEAng
DedxCalibEAng(const std::string &name, ISvcLocator *pSvcLocator)
Definition: DedxCalibEAng.cxx:19
DedxCalib::ReadRecFileList
void ReadRecFileList()
Definition: DedxCalib.cxx:89
DedxCalib::m_phShapeFlag
int m_phShapeFlag
Definition: DedxCalib.h:52
DedxCalib::m_rootfile
std::string m_rootfile
Definition: DedxCalib.h:55
DedxCalib::exsvc
IDedxCorrecSvc * exsvc
Definition: DedxCalib.h:30
DedxCalib::m_recFileVector
vector< string > m_recFileVector
Definition: DedxCalib.h:48
IDedxCorrecSvc::StandardHitCorrec
virtual double StandardHitCorrec(int calib_rec_Flag, int runFlag, int ntpFlag, int runNO, int evtNO, double pathl, int wid, int layid, double adc, double dd, double eangle, double costheta) const =0
IDedxCorrecSvc::StandardTrackCorrec
virtual double StandardTrackCorrec(int calib_rec_Flag, int typFlag, int ntpFlag, int runNO, int evtNO, double ex, double costheta, double t0) const =0
costheta
float costheta
Definition: prepare_data.cxx:28
charge
float charge
Definition: prepare_data.cxx:28
ptrk
float ptrk
Definition: prepare_data.cxx:28