BOSS 7.0.7
BESIII Offline Software System
Loading...
Searching...
No Matches
Public Member Functions | List of all members
DedxCalibEAng Class Reference

#include <DedxCalibEAng.h>

+ Inheritance diagram for DedxCalibEAng:

Public Member Functions

 DedxCalibEAng (const std::string &name, ISvcLocator *pSvcLocator)
 
 ~DedxCalibEAng ()
 
void initializing ()
 
void BookHists ()
 
void genNtuple ()
 
void FillHists ()
 
void AnalyseHists ()
 
void WriteHists ()
 
- Public Member Functions inherited from DedxCalib
 DedxCalib (const std::string &name, ISvcLocator *pSvcLocator)
 
 ~DedxCalib ()
 
StatusCode initialize ()
 
StatusCode execute ()
 
StatusCode finalize ()
 
virtual void genNtuple ()=0
 
virtual void initializing ()=0
 
virtual void BookHists ()=0
 
virtual void FillHists ()=0
 
virtual void AnalyseHists ()=0
 
virtual void WriteHists ()=0
 

Additional Inherited Members

- Protected Member Functions inherited from DedxCalib
double cal_dedx_bitrunc (float truncate, std::vector< double > phlist)
 
double cal_dedx (float truncate, std::vector< double > phlist)
 
void getCurvePar ()
 
void set_dEdx (int landau, float dEdx, int trkalg, int runflag, int dedxhit_use, float ptrk, float theta, float pl_rp, int vflag[3], double t0)
 
void ReadRecFileList ()
 
- Protected Attributes inherited from DedxCalib
IMdcGeomSvc * geosvc
 
IDedxCorrecSvc * exsvc
 
float truncate
 
vector< double > Curve_Para
 
vector< double > Sigma_Para
 
int vFlag [3]
 
double m_dedx_exp [5]
 
double m_ex_sigma [5]
 
double m_pid_prob [5]
 
double m_chi_ex [5]
 
vector< string > m_recFileVector
 
int ParticleFlag
 
int m_calibflag
 
int m_phShapeFlag
 
std::string m_eventType
 
std::string m_recFileList
 
std::string m_rootfile
 
std::string m_curvefile
 

Detailed Description

Definition at line 10 of file DedxCalibEAng.h.

Constructor & Destructor Documentation

◆ DedxCalibEAng()

DedxCalibEAng::DedxCalibEAng ( const std::string &  name,
ISvcLocator *  pSvcLocator 
)

Definition at line 18 of file DedxCalibEAng.cxx.

18 : DedxCalib(name, pSvcLocator){
19 declareProperty("Truncate", m_truncate = 0.7);
20 declareProperty("PMin", pMin=0.2);
21 declareProperty("PMax", pMax=2.3);
22}

◆ ~DedxCalibEAng()

DedxCalibEAng::~DedxCalibEAng ( )
inline

Definition at line 14 of file DedxCalibEAng.h.

14{};

Member Function Documentation

◆ AnalyseHists()

void DedxCalibEAng::AnalyseHists ( )
virtual

Implements DedxCalib.

Definition at line 151 of file DedxCalibEAng.cxx.

152{
153 MsgStream log(msgSvc(), name());
154 log<<MSG::INFO<<"DedxCalibEAng::AnalyseHists()"<<endreq;
155
156 TF1* func;
157 Double_t entry=0,mean=0,rms=0;
158 Double_t binmax=0;
159 Double_t lp[3]={0};
160 gStyle->SetOptFit(1111);
161 for(Int_t ip=0;ip<NumSlices;ip++)
162 {
163 entry = m_eangle[ip]->Integral();
164 if(entry<10) continue;
165 mean = m_eangle[ip]->GetMean();
166 rms = m_eangle[ip]->GetRMS();
167 binmax = m_eangle[ip]->GetMaximumBin();
168 lp[1] = m_eangle[ip]->GetBinCenter(binmax);
169 lp[2] = 200;
170 lp[0] = (m_eangle[ip]->GetBinContent(binmax)+m_eangle[ip]->GetBinContent(binmax-1)+m_eangle[ip]->GetBinContent(binmax+1))/3;
171
172 if(m_phShapeFlag==0)
173 {
174 func = new TF1("mylan",mylan,10,3000,4);
175 func->SetParLimits(0, 0, entry);
176 func->SetParLimits(1, 0, mean);
177 func->SetParLimits(2, 10, rms);
178 func->SetParLimits(3, -1, 0);
179 func->SetParameters(entry/36., mean-200, rms/2.3, -0.5);
180 }
181 else if(m_phShapeFlag==1)
182 {
183 func = new TF1("landaun",landaun,10,3000,3);
184 func->SetParameters(lp[0], lp[1], lp[2] );
185 }
186 else if(m_phShapeFlag==2)
187 {
188 func = new TF1("Landau",Landau,10,3000,2);
189 func->SetParameters(0.7*mean, rms );
190 }
191 else if(m_phShapeFlag==3)
192 {
193 func = new TF1("Vavilov",Vavilov,10,3000,2);
194 func->SetParameters(0.05, 0.6);
195 }
196 else if(m_phShapeFlag==4)
197 {
198 func = new TF1("AsymGauss",AsymGauss,10,3000,4);
199 func->SetParameters(entry, mean, rms, 1.0 );
200 }
201 func->SetLineWidth(2.1);
202 func->SetLineColor(2);
203
204 m_eangle[ip]->Fit(func, "rq" );// , "", mean-rms, mean+rms/2);
205 m_pos_eangle[ip]->Fit(func, "rq");// , "", mean-rms, mean+rms/2);
206 m_neg_eangle[ip]->Fit(func, "rq");// , "", mean-rms, mean+rms/2);
207 }
208}
Landau
double Landau(double *x, double *par)
Definition: DedxCalibParameters.cxx:34
Vavilov
double Vavilov(double *x, double *par)
Definition: DedxCalibParameters.cxx:41
NumSlices
#define NumSlices
Definition: DedxCalibParameters.h:61
landaun
double landaun(double *x, double *par)
Definition: DedxCalibParameters.cxx:25
AsymGauss
double AsymGauss(double *x, double *par)
Definition: DedxCalibParameters.cxx:50
mylan
double mylan(double *x, double *par)
Definition: DedxCalibParameters.cxx:17
msgSvc
IMessageSvc * msgSvc()
Definition: ServiceAccessor.h:13
DedxCalib::m_phShapeFlag
int m_phShapeFlag
Definition: DedxCalib.h:52

◆ BookHists()

void DedxCalibEAng::BookHists ( )
virtual

Implements DedxCalib.

Definition at line 24 of file DedxCalibEAng.cxx.

25{
26 MsgStream log(msgSvc(), name());
27 log<<MSG::INFO<<"DedxCalibEAng::BookHists()"<<endreq;
28
29 ReadRecFileList();
30
31 m_eangle = new TH1F*[NumSlices];
32 m_trunc_eangle = new TH1F*[NumSlices];
33 m_pos_eangle = new TH1F*[NumSlices];
34 m_neg_eangle = new TH1F*[NumSlices];
35 stringstream hlabel;
36 for(int i=0;i<NumSlices;i++)
37 {
38 hlabel.str("");
39 hlabel<<"dEdx"<<"_eang"<<i;
40 m_eangle[i] = new TH1F(hlabel.str().c_str(), hlabel.str().c_str(), NumHistBins, MinHistValue, MaxHistValue);
41 hlabel.str("");
42 hlabel<<"trunc_dEdx"<<"_eang"<<i;
43 m_trunc_eangle[i] = new TH1F(hlabel.str().c_str(), hlabel.str().c_str(), NumHistBins, MinHistValue, MaxHistValue);
44 hlabel.str("");
45 hlabel<<"pos_dEdx"<<"_eang"<<i;
46 m_pos_eangle[i] = new TH1F(hlabel.str().c_str(), hlabel.str().c_str(), NumHistBins, MinHistValue, MaxHistValue);
47 hlabel.str("");
48 hlabel<<"neg_dEdx"<<"_eang"<<i;
49 m_neg_eangle[i] = new TH1F(hlabel.str().c_str(), hlabel.str().c_str(), NumHistBins, MinHistValue, MaxHistValue);
50 }
51 hlabel.str("");
52 hlabel<<"dEdxVsEAng";
53 m_EAngAverdE = new TH1F(hlabel.str().c_str(),"dEdx & EAng",NumSlices,PhiMin,PhiMax);
54 hlabel.str("");
55 hlabel<<"trunc_dEdxVsEAng";
56 m_trunc_EAngAverdE = new TH1F(hlabel.str().c_str(),"trunc_dEdx & EAng",NumSlices,PhiMin,PhiMax);
57 hlabel.str("");
58 hlabel<<"pos_dEdxVsEAng";
59 m_pos_EAngAverdE = new TH1F(hlabel.str().c_str(),"pos_dEdx & EAng",NumSlices,PhiMin,PhiMax);
60 hlabel.str("");
61 hlabel<<"neg_dEdxVsEAng";
62 m_neg_EAngAverdE = new TH1F(hlabel.str().c_str(),"neg_dEdx & EAng",NumSlices,PhiMin,PhiMax);
63
64}
MinHistValue
#define MinHistValue
Definition: DedxCalibParameters.h:46
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
DedxCalib::ReadRecFileList
void ReadRecFileList()
Definition: DedxCalib.cxx:89

◆ FillHists()

void DedxCalibEAng::FillHists ( )
virtual

Implements DedxCalib.

Definition at line 66 of file DedxCalibEAng.cxx.

67{
68 MsgStream log(msgSvc(), name());
69 log<<MSG::INFO<<"DedxCalibEAng::FillHists()"<<endreq;
70
71 TFile* f;
72 TTree* n102;
73 string runlist;
74
75 double dedx=0;
76 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;
77 cout<<"m_recFileVector.size()= "<<m_recFileVector.size()<<endl;
78
79 vector<double> m_phlist(0);
80 vector<double> m_phlist_test(0);
81 vector<int> m_iEAng_list(0);
82 double old_costheta(0), old_ptrk(0);
83
84 for(unsigned int i=0; i<m_recFileVector.size(); i++)
85 {
86 runlist = m_recFileVector[i];
87 cout<<"runlist: "<<runlist.c_str()<<endl;
88 f = new TFile(runlist.c_str());
89 n102 = (TTree*)f->Get("n102");
90 n102-> SetBranchAddress("adc_raw",&dedx);
91 n102-> SetBranchAddress("path_rphi",&pathlength);
92 n102-> SetBranchAddress("wire",&wid);
93 n102-> SetBranchAddress("layer",&layid);
94 n102-> SetBranchAddress("doca_in",&dd_in);
95 n102-> SetBranchAddress("driftdist",&driftdist);
96 n102-> SetBranchAddress("eangle",&eangle);
97 n102-> SetBranchAddress("zhit",&zhit);
98 n102-> SetBranchAddress("runNO",&runNO);
99 n102-> SetBranchAddress("evtNO",&evtNO);
100 n102-> SetBranchAddress("runFlag",&runFlag);
101 n102-> SetBranchAddress("costheta1",&costheta);
102 n102-> SetBranchAddress("t01",&tes);
103 n102->SetBranchAddress("ptrk1",&ptrk);
104 n102->SetBranchAddress("charge",&charge);
105 for(int j=0;j<n102->GetEntries();j++)
106 {
107 n102->GetEntry(j);
108 // begin of save truncated hits information for EACH track
109 if((fabs(costheta-old_costheta)>0.001 || fabs(ptrk-old_ptrk>0.001)) && m_phlist.size()>=2){
110 sort(m_phlist_test.begin(), m_phlist_test.end());
111 unsigned int nsample = (unsigned int)(m_phlist_test.size()*m_truncate);
112 double m_trunc_value = m_phlist_test[nsample];
113 for(unsigned int i=0; i<m_phlist.size(); i++){
114 if(m_phlist[i]>m_trunc_value) continue;
115 m_trunc_eangle[m_iEAng_list[i]]->Fill(m_phlist[i]);
116 }
117 m_phlist.clear();
118 m_phlist_test.clear();
119 m_iEAng_list.clear();
120 }
121 old_costheta = costheta;
122 old_ptrk = ptrk;
123 // end of save truncated information
124 if(eangle >0.25) eangle = eangle -0.5;
125 else if(eangle <-0.25) eangle = eangle +0.5;
126 if(eangle>PhiMin && eangle<PhiMax)
127 {
128 if(ptrk>pMax || ptrk<pMin) continue;
129 if(tes>1400) continue;
130 // if (ptrk<0.1) continue;
131 if(layid <4) continue;
132 else if(layid <8)
133 {if(dedx<MinADCValuecut || dedx>MaxADCValuecut || pathlength>RPhi_PathMaxCut || pathlength<Iner_RPhi_PathMinCut || driftdist>Iner_DriftDistCut) continue;}
134 else
135 {if(dedx<MinADCValuecut || dedx>MaxADCValuecut || pathlength>RPhi_PathMaxCut || pathlength<Out_RPhi_PathMinCut || driftdist>Out_DriftDistCut) continue;}
136 int iEAng = 0;
137 iEAng = (int) ((eangle-PhiMin)/StepEAng);
138 dedx = exsvc->StandardHitCorrec(0,runFlag,2,runNO,evtNO,pathlength,wid,layid,dedx,dd_in,eangle,costheta);
139 dedx = exsvc->StandardTrackCorrec(0, runFlag, 2, runNO, evtNO, dedx, costheta, tes);
140 m_eangle[iEAng]->Fill(dedx);
141 m_phlist.push_back(dedx);
142 m_phlist_test.push_back(dedx);
143 m_iEAng_list.push_back(iEAng);
144 if(charge>0) m_pos_eangle[iEAng]->Fill(dedx);
145 if(charge<0) m_neg_eangle[iEAng]->Fill(dedx);
146 }
147 }
148 }
149}
StepEAng
const double StepEAng
Definition: DedxCalibEAng.cxx:16
MaxADCValuecut
#define MaxADCValuecut
Definition: DedxCalibParameters.h:54
Iner_DriftDistCut
#define Iner_DriftDistCut
Definition: DedxCalibParameters.h:58
RPhi_PathMaxCut
#define RPhi_PathMaxCut
Definition: DedxCalibParameters.h:55
Out_DriftDistCut
#define Out_DriftDistCut
Definition: DedxCalibParameters.h:59
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
f
TFile f("ana_bhabha660a_dqa_mcPat_zy_old.root")
SetBranchAddress
curve SetBranchAddress("CurveSize",&CurveSize)
costheta
float costheta
Definition: prepare_data.cxx:28
charge
float charge
Definition: prepare_data.cxx:28
ptrk
float ptrk
Definition: prepare_data.cxx:28

◆ genNtuple()

void DedxCalibEAng::genNtuple ( )
inlinevirtual

Implements DedxCalib.

Definition at line 17 of file DedxCalibEAng.h.

17{}

◆ initializing()

void DedxCalibEAng::initializing ( )
inlinevirtual

Implements DedxCalib.

Definition at line 15 of file DedxCalibEAng.h.

15{}

◆ WriteHists()

void DedxCalibEAng::WriteHists ( )
virtual

Implements DedxCalib.

Definition at line 210 of file DedxCalibEAng.cxx.

211{
212 MsgStream log(msgSvc(), name());
213 log<<MSG::INFO<<"DedxCalibEAng::WriteHists()"<<endreq;
214
215 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};
216 double trunc_gain[NumSlices]={0};
217 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};
218 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};
219 double Ip_eangle[NumSlices]={0},eangle[NumSlices]={0};
220
221 double Norm(0), trunc_Norm(0);
222 int count(0);
223 for(Int_t ip=0;ip<NumSlices;ip++)
224 {
225 Ip_eangle[ip] = ip;
226 eangle[ip] = (ip+0.5)*StepEAng + PhiMin;
227
228 entryNo[ip]= m_eangle[ip]->Integral();
229 trunc_mean[ip] = m_trunc_eangle[ip]->GetMean();
230 trunc_sigma[ip] = m_trunc_eangle[ip]->GetRMS();
231 pos_entryNo[ip]= m_pos_eangle[ip]->Integral();
232 pos_mean[ip] = m_pos_eangle[ip]->GetMean();
233 pos_sigma[ip] = m_pos_eangle[ip]->GetRMS();
234 neg_entryNo[ip]= m_neg_eangle[ip]->Integral();
235 neg_mean[ip] = m_neg_eangle[ip]->GetMean();
236 neg_sigma[ip] = m_neg_eangle[ip]->GetRMS();
237 if(entryNo[ip]<10 || pos_entryNo[ip]<10 || neg_entryNo[ip]<10) continue;
238
239 if(m_phShapeFlag==0)
240 {
241 fitmean[ip] = m_eangle[ip]->GetFunction("mylan")->GetParameter(1);
242 fitmeanerr[ip] = m_eangle[ip]->GetFunction("mylan")->GetParError(1);
243 fitsigma[ip] = m_eangle[ip]->GetFunction("mylan")->GetParameter(2);
244 chisq[ip] = (m_eangle[ip]->GetFunction("mylan")-> GetChisquare())/(m_eangle[ip]->GetFunction("mylan")-> GetNDF());
245
246 pos_fitmean[ip] = m_pos_eangle[ip]->GetFunction("mylan")->GetParameter(1);
247 pos_fitmeanerr[ip] = m_pos_eangle[ip]->GetFunction("mylan")->GetParError(1);
248 pos_fitsigma[ip] = m_pos_eangle[ip]->GetFunction("mylan")->GetParameter(2);
249 pos_chisq[ip] = (m_pos_eangle[ip]->GetFunction("mylan")-> GetChisquare())/(m_pos_eangle[ip]->GetFunction("mylan")-> GetNDF());
250
251 neg_fitmean[ip] = m_neg_eangle[ip]->GetFunction("mylan")->GetParameter(1);
252 neg_fitmeanerr[ip] = m_neg_eangle[ip]->GetFunction("mylan")->GetParError(1);
253 neg_fitsigma[ip] = m_neg_eangle[ip]->GetFunction("mylan")->GetParameter(2);
254 neg_chisq[ip] = (m_neg_eangle[ip]->GetFunction("mylan")-> GetChisquare())/(m_neg_eangle[ip]->GetFunction("mylan")-> GetNDF());
255 }
256 if(m_phShapeFlag==1)
257 {
258 fitmean[ip] = m_eangle[ip]->GetFunction("landaun")->GetParameter(1);
259 fitmeanerr[ip] = m_eangle[ip]->GetFunction("landaun")->GetParError(1);
260 fitsigma[ip] = m_eangle[ip]->GetFunction("landaun")->GetParameter(2);
261 chisq[ip] = (m_eangle[ip]->GetFunction("landaun")-> GetChisquare())/(m_eangle[ip]->GetFunction("mylan")-> GetNDF());
262
263 pos_fitmean[ip] = m_pos_eangle[ip]->GetFunction("landaun")->GetParameter(1);
264 pos_fitmeanerr[ip] = m_pos_eangle[ip]->GetFunction("landaun")->GetParError(1);
265 pos_fitsigma[ip] = m_pos_eangle[ip]->GetFunction("landaun")->GetParameter(2);
266 pos_chisq[ip] = (m_pos_eangle[ip]->GetFunction("landaun")-> GetChisquare())/(m_pos_eangle[ip]->GetFunction("landaun")-> GetNDF());
267
268 neg_fitmean[ip] = m_neg_eangle[ip]->GetFunction("landaun")->GetParameter(1);
269 neg_fitmeanerr[ip] = m_neg_eangle[ip]->GetFunction("landaun")->GetParError(1);
270 neg_fitsigma[ip] = m_neg_eangle[ip]->GetFunction("landaun")->GetParameter(2);
271 neg_chisq[ip] = (m_neg_eangle[ip]->GetFunction("landaun")-> GetChisquare())/(m_neg_eangle[ip]->GetFunction("landaun")-> GetNDF());
272 }
273 if(m_phShapeFlag==2)
274 {
275 fitmean[ip] = m_eangle[ip]->GetFunction("Landau")->GetParameter(1);
276 fitmeanerr[ip] = m_eangle[ip]->GetFunction("Landau")->GetParError(1);
277 fitsigma[ip] = m_eangle[ip]->GetFunction("Landau")->GetParameter(2);
278 chisq[ip] = (m_eangle[ip]->GetFunction("Landau")-> GetChisquare())/(m_eangle[ip]->GetFunction("Landau")-> GetNDF());
279
280 pos_fitmean[ip] = m_pos_eangle[ip]->GetFunction("Landau")->GetParameter(1);
281 pos_fitmeanerr[ip] = m_pos_eangle[ip]->GetFunction("Landau")->GetParError(1);
282 pos_fitsigma[ip] = m_pos_eangle[ip]->GetFunction("Landau")->GetParameter(2);
283 pos_chisq[ip] = (m_pos_eangle[ip]->GetFunction("Landau")-> GetChisquare())/(m_pos_eangle[ip]->GetFunction("Landau")-> GetNDF());
284
285 neg_fitmean[ip] = m_neg_eangle[ip]->GetFunction("Landau")->GetParameter(1);
286 neg_fitmeanerr[ip] = m_neg_eangle[ip]->GetFunction("Landau")->GetParError(1);
287 neg_fitsigma[ip] = m_neg_eangle[ip]->GetFunction("Landau")->GetParameter(2);
288 neg_chisq[ip] = (m_neg_eangle[ip]->GetFunction("Landau")-> GetChisquare())/(m_neg_eangle[ip]->GetFunction("Landau")-> GetNDF());
289 }
290 if(m_phShapeFlag==3)
291 {
292 fitmean[ip] = m_eangle[ip]->GetFunction("Vavilov")->GetParameter(1);
293 fitmeanerr[ip] = m_eangle[ip]->GetFunction("Vavilov")->GetParError(1);
294 fitsigma[ip] = m_eangle[ip]->GetFunction("Vavilov")->GetParameter(2);
295 chisq[ip] = (m_eangle[ip]->GetFunction("Vavilov")-> GetChisquare())/(m_eangle[ip]->GetFunction("Vavilov")-> GetNDF());
296
297 pos_fitmean[ip] = m_pos_eangle[ip]->GetFunction("Vavilov")->GetParameter(1);
298 pos_fitmeanerr[ip] = m_pos_eangle[ip]->GetFunction("Vavilov")->GetParError(1);
299 pos_fitsigma[ip] = m_pos_eangle[ip]->GetFunction("Vavilov")->GetParameter(2);
300 pos_chisq[ip] = (m_pos_eangle[ip]->GetFunction("Vavilov")-> GetChisquare())/(m_pos_eangle[ip]->GetFunction("Vavilov")-> GetNDF());
301
302 neg_fitmean[ip] = m_neg_eangle[ip]->GetFunction("Vavilov")->GetParameter(1);
303 neg_fitmeanerr[ip] = m_neg_eangle[ip]->GetFunction("Vavilov")->GetParError(1);
304 neg_fitsigma[ip] = m_neg_eangle[ip]->GetFunction("Vavilov")->GetParameter(2);
305 neg_chisq[ip] = (m_neg_eangle[ip]->GetFunction("Vavilov")-> GetChisquare())/(m_neg_eangle[ip]->GetFunction("Vavilov")-> GetNDF());
306 }
307 if(m_phShapeFlag==4)
308 {
309 fitmean[ip] = m_eangle[ip]->GetFunction("AsymGauss")->GetParameter(1);
310 fitmeanerr[ip] = m_eangle[ip]->GetFunction("AsymGauss")->GetParError(1);
311 fitsigma[ip] = m_eangle[ip]->GetFunction("AsymGauss")->GetParameter(2);
312 chisq[ip] = (m_eangle[ip]->GetFunction("AsymGauss")-> GetChisquare())/(m_eangle[ip]->GetFunction("AsymGauss")-> GetNDF());
313
314 pos_fitmean[ip] = m_pos_eangle[ip]->GetFunction("AsymGauss")->GetParameter(1);
315 pos_fitmeanerr[ip] = m_pos_eangle[ip]->GetFunction("AsymGauss")->GetParError(1);
316 pos_fitsigma[ip] = m_pos_eangle[ip]->GetFunction("AsymGauss")->GetParameter(2);
317 pos_chisq[ip] = (m_pos_eangle[ip]->GetFunction("AsymGauss")-> GetChisquare())/(m_pos_eangle[ip]->GetFunction("AsymGauss")-> GetNDF());
318
319 neg_fitmean[ip] = m_neg_eangle[ip]->GetFunction("AsymGauss")->GetParameter(1);
320 neg_fitmeanerr[ip] = m_neg_eangle[ip]->GetFunction("AsymGauss")->GetParError(1);
321 neg_fitsigma[ip] = m_neg_eangle[ip]->GetFunction("AsymGauss")->GetParameter(2);
322 neg_chisq[ip] = (m_neg_eangle[ip]->GetFunction("AsymGauss")-> GetChisquare())/(m_neg_eangle[ip]->GetFunction("AsymGauss")-> GetNDF());
323 }
324 if(entryNo[ip]>1500 && fitmean[ip]>0 && fitmeanerr[ip]>0 && fitmeanerr[ip]<10 && fitsigma[ip]>0 && fitsigma[ip]<200){
325 Norm += fitmean[ip];
326 trunc_Norm += trunc_mean[ip];
327 count++;
328 }
329 }
330 Norm = Norm/count;
331 trunc_Norm = trunc_Norm/count;
332 // Norm = 550; // use an absolute normalization
333 cout<<"count= "<<count<<" average value: "<< Norm << " truncted: " << trunc_Norm << endl;
334 for(Int_t i=0;i<NumSlices;i++)
335 { // the difference between gain and mean is the normalization
336 if(entryNo[i]>10 && fitmean[i]>0 && fitsigma[i]>0) gain[i]=fitmean[i]/Norm;
337 else gain[i] = 1;
338 if(entryNo[i]>10) trunc_gain[i]=trunc_mean[i]/trunc_Norm;
339 else trunc_mean[i] = 1;
340 if(pos_entryNo[i]>10 && pos_fitmean[i]>0 && pos_fitsigma[i]>0) pos_gain[i]=pos_fitmean[i]/Norm;
341 else pos_gain[i] = 1;
342 if(neg_entryNo[i]>10 && neg_fitmean[i]>0 && neg_fitsigma[i]>0) neg_gain[i]=neg_fitmean[i]/Norm;
343 else neg_gain[i] = 1;
344 if(gain[i] != 1){
345 m_EAngAverdE -> SetBinContent(i+1, gain[i]);
346 m_EAngAverdE->SetBinError(i+1, fitmeanerr[i]/Norm);
347 }
348 if(trunc_gain[i] != 1){
349 m_trunc_EAngAverdE -> SetBinContent(i+1, trunc_gain[i]);
350 m_trunc_EAngAverdE->SetBinError(i+1, trunc_sigma[i]/trunc_Norm);
351 }
352 if(pos_gain[i] != 1){
353 m_pos_EAngAverdE -> SetBinContent(i+1, pos_gain[i]);
354 m_pos_EAngAverdE->SetBinError(i+1, pos_fitmeanerr[i]/Norm);
355 }
356 if(neg_gain[i] != 1){
357 m_neg_EAngAverdE -> SetBinContent(i+1, neg_gain[i]);
358 m_neg_EAngAverdE->SetBinError(i+1, neg_fitmeanerr[i]/Norm);
359 }
360 }
361
362 log<<MSG::INFO<<"begin getting calibration constant!!! "<<endreq;
363 double denangle[NumSlices]={0};
364 int denangle_entry[1] = {100};
365 m_EAngAverdE->Fit("pol1","r","",-0.25,0.25);
366 double par0 = m_EAngAverdE->GetFunction("pol1")->GetParameter(0);
367 double par1 = m_EAngAverdE->GetFunction("pol1")->GetParameter(1);
368 // double par2 = m_EAngAverdE->GetFunction("pol2")->GetParameter(2);
369 // par1 = -par1; // a trick is used here, because we don't really understand
370 // for(int i=0;i<NumSlices;i++) denangle[i] = (par0+par1*eangle[i])/par0;//+par2*pow(eangle[i],2))/par0;
371 for(int i=0;i<NumSlices;i++) denangle[i] = gain[i]; // use raw values instead of fit
372
373 log<<MSG::INFO<<"begin generating root file!!! "<<endreq;
374 TFile* f = new TFile(m_rootfile.c_str(), "recreate");
375 for(Int_t ip=0;ip<NumSlices;ip++){
376 m_eangle[ip]->Write();
377 m_trunc_eangle[ip]->Write();
378 m_pos_eangle[ip]->Write();
379 m_neg_eangle[ip]->Write();
380 }
381 m_EAngAverdE->Write();
382 m_trunc_EAngAverdE->Write();
383 m_pos_EAngAverdE->Write();
384 m_neg_EAngAverdE->Write();
385
386 TTree* entracalib = new TTree("entracalib","entracalib");
387 entracalib->Branch("1denangle", denangle, "denangle[100]/D");
388 entracalib->Branch("1denangle_entry", denangle_entry, "denangle_entry[1]/I");
389 entracalib->Fill();
390 entracalib-> Write();
391
392 TTree* ddgcalib = new TTree("ddgcalib", "ddgcalib");
393 ddgcalib -> Branch("hits", entryNo, "entryNo[100]/D");
394 ddgcalib -> Branch("truncmean", trunc_mean, "trunc_mean[100]/D");
395 ddgcalib -> Branch("truncsigma", trunc_sigma, "trunc_sigma[100]/D");
396 ddgcalib -> Branch("fitmean", fitmean, "fitmean[100]/D");
397 ddgcalib -> Branch("fitmeanerr", fitmeanerr, "fitmeanerr[100]/D");
398 ddgcalib -> Branch("chi", fitsigma, "fitsigma[100]/D");
399 ddgcalib -> Branch("gain", gain, "gain[100]/D");
400 ddgcalib -> Branch("truncgain", trunc_gain, "trunc_gain[100]/D");
401 ddgcalib -> Branch("chisq", chisq, "chisq[100]/D");
402 ddgcalib -> Branch("pos_hits", pos_entryNo, "pos_entryNo[100]/D");
403 ddgcalib -> Branch("pos_mean", pos_mean, "pos_mean[100]/D");
404 ddgcalib -> Branch("pos_sigma", pos_sigma, "pos_sigma[100]/D");
405 ddgcalib -> Branch("pos_fitmean", pos_fitmean, "pos_fitmean[100]/D");
406 ddgcalib -> Branch("pos_fitmeanerr", pos_fitmeanerr, "pos_fitmeanerr[100]/D");
407 ddgcalib -> Branch("pos_chi", pos_fitsigma, "pos_fitsigma[100]/D");
408 ddgcalib -> Branch("pos_gain", pos_gain, "pos_gain[100]/D");
409 ddgcalib -> Branch("pos_chisq", pos_chisq, "pos_chisq[100]/D");
410 ddgcalib -> Branch("neg_hits", neg_entryNo, "neg_entryNo[100]/D");
411 ddgcalib -> Branch("neg_mean", neg_mean, "neg_mean[100]/D");
412 ddgcalib -> Branch("neg_sigma", neg_sigma, "neg_sigma[100]/D");
413 ddgcalib -> Branch("neg_fitmean", neg_fitmean, "neg_fitmean[100]/D");
414 ddgcalib -> Branch("neg_fitmeanerr", neg_fitmeanerr, "neg_fitmeanerr[100]/D");
415 ddgcalib -> Branch("neg_chi", neg_fitsigma, "neg_fitsigma[100]/D");
416 ddgcalib -> Branch("neg_gain", neg_gain, "neg_gain[100]/D");
417 ddgcalib -> Branch("neg_chisq", neg_chisq, "neg_chisq[100]/D");
418 ddgcalib -> Branch("Ip_eangle", Ip_eangle, "Ip_eangle[100]/D");
419 ddgcalib -> Branch("eangle", eangle, "eangle[100]/D");
420 ddgcalib -> Fill();
421 ddgcalib -> Write();
422 f->Close();
423
424 TCanvas c1("c1", "canvas", 500, 400);
425 c1.Print((m_rootfile+".ps[").c_str());
426 m_EAngAverdE->Draw();
427 c1.Print((m_rootfile+".ps").c_str());
428 m_trunc_EAngAverdE->Draw();
429 c1.Print((m_rootfile+".ps").c_str());
430 m_pos_EAngAverdE->Draw();
431 c1.Print((m_rootfile+".ps").c_str());
432 m_neg_EAngAverdE->Draw();
433 c1.Print((m_rootfile+".ps").c_str());
434 for(Int_t ip=0;ip<NumSlices;ip++)
435 {
436 m_eangle[ip]->Draw();
437 c1.Print((m_rootfile+".ps").c_str());
438 }
439 for(Int_t ip=0;ip<NumSlices;ip++)
440 {
441 m_trunc_eangle[ip]->Draw();
442 c1.Print((m_rootfile+".ps").c_str());
443 }
444 for(Int_t ip=0;ip<NumSlices;ip++)
445 {
446 m_pos_eangle[ip]->Draw();
447 c1.Print((m_rootfile+".ps").c_str());
448 }
449 for(Int_t ip=0;ip<NumSlices;ip++)
450 {
451 m_neg_eangle[ip]->Draw();
452 c1.Print((m_rootfile+".ps").c_str());
453 }
454 c1.Print((m_rootfile+".ps]").c_str());
455
456 for(Int_t ip=0;ip<NumSlices;ip++)
457 {
458 delete m_eangle[ip];
459 delete m_trunc_eangle[ip];
460 delete m_pos_eangle[ip];
461 delete m_neg_eangle[ip];
462 }
463 delete m_EAngAverdE;
464 delete m_trunc_EAngAverdE;
465 delete m_pos_EAngAverdE;
466 delete m_neg_EAngAverdE;
467
468}
Branch
curve Branch("CurveSize",&CurveSize,"CurveSize/I")
Fill
curve Fill()
Write
curve Write()
count
DOUBLE_PRECISION count[3]
Definition: EvtPhokharaDef.hh:64
DedxCalib::m_rootfile
std::string m_rootfile
Definition: DedxCalib.h:55
EvtCyclic3::c_str
char * c_str(Index i)
Definition: EvtCyclic3.cc:252
The documentation for this class was generated from the following files: