TENDER becomes Tender
[u/mrichter/AliRoot.git] / PWGPP / TPC / macros / MakeReportTPC.C
CommitLineData
fa3caf32 1/*
2bfe5463 2 //gSystem->AddIncludePath("-I$ALICE_ROOT/PWGPP/TPC");
3 //.L $ALICE_ROOT/PWGPP/TPC/macros/MakeReportTPC.C+
fa3caf32 4 //
5 //MakeReportTPC();
6
7 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
8 gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
9 AliXRDPROOFtoolkit tool;
10 TChain * chain = tool.MakeChain("summaryTPCQA.txt","tpcQA",0,200000);
fc567e9f 11
12 TTree *tree = chain->CopyTree("1");
fa3caf32 13*/
14
15
16#if !defined(__CINT__) || defined(__MAKECINT__)
17#include "TSystem.h"
fc567e9f 18#include "TMath.h"
19#include "TVectorD.h"
fa3caf32 20#include "TFile.h"
21#include "TGrid.h"
22#include "TF1.h"
23#include "TH1.h"
24#include "TH2.h"
25#include "TProfile.h"
26#include "THnSparse.h"
27#include "TTreeStream.h"
28#include "AliPerformanceTPC.h"
fc567e9f 29#include "AliPerformanceDEdx.h"
fa3caf32 30#endif
31
32TObject *pTPCObject=0;
fc567e9f 33TObject *pTPCObjectGain=0;
fa3caf32 34TTreeSRedirector *pcstream=0;
35void Init();
36void ReadObjects(const char * path=0);
37void MakeReportTPC(Int_t run, const char *path=0 );
38void AnalyzeNCL();
39void AnalyzeDrift();
fc567e9f 40void AnalyzeGain();
fa3caf32 41
42void Init(){
43 //
44 //
45 //
46 gSystem->Load("libANALYSIS.so");
47 gSystem->Load("libANALYSISalice.so");
af472fff 48 gSystem->Load("libTender.so");
fa3caf32 49 gSystem->Load("libCORRFW.so");
50 gSystem->Load("libPWG0base.so");
51 gSystem->Load("libPWG0dep.so");
52 gSystem->Load("libPWG0selectors.so");
2bfe5463 53 gSystem->Load("libPWGPP.so");
fa3caf32 54 gSystem->Load("libTPCcalib");
fc567e9f 55// gSystem->Setenv("alien_CLOSE_SE","ALICE::GSI::SE");
56// TGrid::Connect("alien://",0,0,"t");
57// gSystem->Setenv("alien_CLOSE_SE","ALICE::GSI::SE");
fa3caf32 58}
59
60
61void ReadObjects(const char * path){
62 //
63 //
64 //
65 TFile *f =0;
66 if (path) f=TFile::Open(Form("%s/TPC.PerformanceQAQA.root",path));
67 if (!path) f=TFile::Open("TPC.Performance.root");
68 if (!f){
69 printf("File %s not available\n", path);
70 return;
71 }
72 TList *list = (TList*)f->Get("TPC");
73 if (!list){
74 printf("QA %s not available\n", path);
75 return;
76 }
77 pTPCObject= (AliPerformanceTPC*)list->FindObject("AliPerformanceTPC");
fc567e9f 78 pTPCObjectGain = (AliPerformanceDEdx*)list->FindObject("AliPerformanceDEdxTPCInner");
fa3caf32 79}
80
81
82void MakeReportTPC(Int_t run, const char *path ){
83 //
84 // make a tpcQA report
85 // typical variables exported to the tree
86 //
87 Init();
88 ReadObjects(path);
89
90 pcstream= new TTreeSRedirector("TPC.PerformanceSummary.root");
91 (*pcstream)<<"tpcQA"<<"run="<<run;
92 AnalyzeNCL();
93 AnalyzeDrift();
fc567e9f 94 AnalyzeGain();
fa3caf32 95 (*pcstream)<<"tpcQA"<<"\n";
96 delete pcstream;
97}
98
99
100void AnalyzeNCL(){
101 //
102 // NCL statistic
103 //
104 // variables:
105 static Double_t meanTPCnclF=0;
106 static Double_t rmsTPCnclF=0;
107 static Double_t slopeATPCnclF=0;
108 static Double_t slopeCTPCnclF=0;
fc567e9f 109 static Double_t slopeATPCnclFErr=0;
110 static Double_t slopeCTPCnclFErr=0;
fa3caf32 111 static Double_t meanTPCncl=0;
112 static Double_t rmsTPCncl=0;
113 static Double_t slopeATPCncl=0;
114 static Double_t slopeCTPCncl=0;
fc567e9f 115 static Double_t slopeATPCnclErr=0;
116 static Double_t slopeCTPCnclErr=0;
fa3caf32 117 AliPerformanceTPC * pTPC= (AliPerformanceTPC *)pTPCObject;
118
119 TH1* his1D=0;
120 TProfile* hprof=0;
121 static TF1 *fpol1 = new TF1("fpol1","pol1");
122 //
123 // all clusters
124 // eta cut - +-1
125 // pt cut - +-0.250 GeV
126 pTPC->GetTPCTrackHisto()->GetAxis(5)->SetRangeUser(-1.,1.);
127 pTPC->GetTPCTrackHisto()->GetAxis(7)->SetRangeUser(0.25,10);
128 his1D = pTPC->GetTPCTrackHisto()->Projection(0);
129 meanTPCncl= his1D->GetMean();
130 rmsTPCncl= his1D->GetRMS();
131 delete his1D;
132 hprof = pTPC->GetTPCTrackHisto()->Projection(0,5)->ProfileX();
133 hprof->Fit(fpol1,"QNR","QNR",0,0.8);
134 slopeATPCncl= fpol1->GetParameter(1);
fc567e9f 135 slopeATPCnclErr= fpol1->GetParError(1);
fa3caf32 136 hprof->Fit(fpol1,"QNR","QNR",-0.8,0.0);
137 slopeCTPCncl= fpol1->GetParameter(1);
fc567e9f 138 slopeCTPCnclErr= fpol1->GetParameter(1);
fa3caf32 139 delete hprof;
140 //
141 // findable clusters
142 //
143 pTPC->GetTPCTrackHisto()->GetAxis(2)->SetRangeUser(0.4,1.1);
144 his1D = pTPC->GetTPCTrackHisto()->Projection(2);
145 meanTPCnclF= his1D->GetMean();
146 rmsTPCnclF= his1D->GetRMS();
fc567e9f 147 delete his1D;
fa3caf32 148 his1D = pTPC->GetTPCTrackHisto()->Projection(2,5)->ProfileX();
149 his1D->Fit(fpol1,"QNR","QNR",0,0.8);
150 slopeATPCnclF= fpol1->GetParameter(1);
fc567e9f 151 slopeATPCnclFErr= fpol1->GetParError(1);
fa3caf32 152 his1D->Fit(fpol1,"QNR","QNR",-0.8,0.0);
153 slopeCTPCnclF= fpol1->GetParameter(1);
fc567e9f 154 slopeCTPCnclFErr= fpol1->GetParameter(1);
fa3caf32 155 delete his1D;
156 printf("Cluster QA report\n");
157 printf("meanTPCnclF=\t%f\n",meanTPCnclF);
158 printf("rmsTPCnclF=\t%f\n",rmsTPCnclF);
159 printf("slopeATPCnclF=\t%f\n",slopeATPCnclF);
160 printf("slopeCTPCnclF=\t%f\n",slopeCTPCnclF);
161 printf("meanTPCncl=\t%f\n",meanTPCncl);
162 printf("rmsTPCncl=\t%f\n",rmsTPCncl);
163 printf("slopeATPCncl=\t%f\n",slopeATPCncl);
164 printf("slopeCTPCncl=\t%f\n",slopeCTPCncl);
165 //
166 // dump results to the tree
167 //
168 (*pcstream)<<"tpcQA"<<
fc567e9f 169 "meanTPCnclF="<<meanTPCnclF <<
fa3caf32 170 "rmsTPCnclF="<<rmsTPCnclF <<
171 "slopeATPCnclF="<< slopeATPCnclF<<
172 "slopeCTPCnclF="<< slopeCTPCnclF<<
fc567e9f 173 "slopeATPCnclFErr="<< slopeATPCnclFErr<<
174 "slopeCTPCnclFErr="<< slopeCTPCnclFErr<<
fa3caf32 175 "meanTPCncl="<<meanTPCncl <<
176 "rmsTPCncl="<< rmsTPCncl<<
177 "slopeATPCncl="<< slopeATPCncl<<
fc567e9f 178 "slopeCTPCncl="<< slopeCTPCncl<<
179 "slopeATPCnclErr="<< slopeATPCnclErr<<
180 "slopeCTPCnclErr="<< slopeCTPCnclErr;
fa3caf32 181
182
183}
184
185void AnalyzeDrift(){
186 //
187 // Analyze drift velocity imperfections
188 //
189 // variables:
190 static Double_t offsetdZA=0;
191 static Double_t slopedZA=0;
192 static Double_t offsetdZC=0;
193 static Double_t slopedZC=0;
fc567e9f 194 static Double_t offsetdZAErr=0;
195 static Double_t slopedZAErr=0;
196 static Double_t offsetdZCErr=0;
197 static Double_t slopedZCErr=0;
198 static Double_t offsetdZAchi2=0;
199 static Double_t slopedZAchi2=0;
358ebec5 200 static Double_t offsetdZCchi2=0;
fc567e9f 201 static Double_t slopedZCchi2=0;
fa3caf32 202 AliPerformanceTPC * pTPC= (AliPerformanceTPC *)pTPCObject;
203
204 TH1* his1D=0;
205 TH2* his2D=0;
206 static TF1 *fpol1 = new TF1("fpol1","pol1");
207 TObjArray arrayFit;
208 his2D = pTPC->GetTPCTrackHisto()->Projection(4,5);
209 his2D->FitSlicesY(0,0,-1,10,"QNR",&arrayFit);
210 delete his2D;
211 his1D = (TH1*) arrayFit.At(1);
fc567e9f 212 his1D->Fit(fpol1,"QNRROB=0.8","QNR",-0.8,-0.1);
fa3caf32 213 offsetdZC=fpol1->GetParameter(0);
214 slopedZC=fpol1->GetParameter(1);
358ebec5 215 //offsetdZCchi2=fpol1->GetChisquare();
fc567e9f 216 slopedZCchi2=fpol1->GetChisquare();
217 //
218 his1D->Fit(fpol1,"QNRROB=0.8","QNR",0.1,0.8);
fa3caf32 219 offsetdZA=fpol1->GetParameter(0);
220 slopedZA=fpol1->GetParameter(1);
fc567e9f 221 offsetdZAErr=fpol1->GetParError(0);
222 slopedZAErr=fpol1->GetParError(1);
223 offsetdZAchi2=fpol1->GetChisquare();
224 slopedZAchi2=fpol1->GetChisquare();
fa3caf32 225 //
226 printf("Drift velocity QA report\n");
227 printf("offsetdZA\t%f\n",offsetdZA);
228 printf("slopedZA\t%f\n",slopedZA);
229 printf("offsetdZC\t%f\n",offsetdZC);
230 printf("slopedZC\t%f\n",slopedZC);
231 //
232 // dump drift QA values
233 //
234 (*pcstream)<<"tpcQA"<<
235 "offsetdZA="<< offsetdZA<<
236 "slopedZA="<< slopedZA<<
237 "offsetdZC="<< offsetdZC<<
fc567e9f 238 "slopedZC="<<slopedZC<<
239 //
240 "offsetdZAErr="<< offsetdZAErr<<
241 "slopedZAErr="<< slopedZAErr<<
242 "offsetdZCErr="<< offsetdZCErr<<
243 "slopedZCErr="<<slopedZCErr<<
244 //
245 "offsetdZAchi2="<< offsetdZAchi2<<
246 "slopedZAchi2="<< slopedZAchi2<<
247 "offsetdZCchi2="<< offsetdZCchi2<<
248 "slopedZCchi2="<<slopedZCchi2;
fa3caf32 249
250}
fc567e9f 251
252
253void AnalyzeGain(){
254 //
255 // gain statistic
256 //
257 AliPerformanceDEdx * pTPCgain = (AliPerformanceDEdx *)pTPCObjectGain;
258 static TVectorD MIPvsSector(36);
259 static TVectorD sector(36);
260 static Float_t MIPmean = 0;
261 static Float_t MIPresolution = 0;
262 static Float_t attachSlopeC = 0;
263 static Float_t attachSlopeA = 0;
264
265 MIPvsSector.Zero();
266 //
267 // select MIP particles
268 //
358ebec5 269 if (pTPCgain){
270 pTPCgain->GetDeDxHisto()->GetAxis(7)->SetRangeUser(0.4,0.55);
271 pTPCgain->GetDeDxHisto()->GetAxis(0)->SetRangeUser(35,60);
272 pTPCgain->GetDeDxHisto()->GetAxis(6)->SetRangeUser(80,160);
273 pTPCgain->GetDeDxHisto()->GetAxis(5)->SetRangeUser(-1,1);
274 //
275 // MIP position and resolution
276 //
277 TF1 gausFit("gausFit","gaus");
278 TH1 * hisProj1D = pTPCgain->GetDeDxHisto()->Projection(0);
279 hisProj1D->Fit(&gausFit,"QN","QN");
280 delete hisProj1D;
281 MIPmean = gausFit.GetParameter(1);
282 MIPresolution = 0;
283 if (MIPmean!=0) MIPresolution = gausFit.GetParameter(2)/MIPmean;
284 //
285 // MIP position vs. dip angle (attachment)
286 //
287 pTPCgain->GetDeDxHisto()->GetAxis(5)->SetRangeUser(-3,0);
288 TH2* his2D = pTPCgain->GetDeDxHisto()->Projection(0,5);
289 TF1 * fpol = new TF1("fpol","pol1");
290 TObjArray arrayFit;
291 his2D->FitSlicesY(0,0,-1,10,"QN",&arrayFit);
292 delete his2D;
293 TH1 * his1D = (TH1*) arrayFit.At(1);
294 his1D->Fit(fpol,"QNROB=0.8","QNR",-1,0);
295 attachSlopeC = fpol->GetParameter(1);
296 //
297 pTPCgain->GetDeDxHisto()->GetAxis(5)->SetRangeUser(0,3);
298 TH2* his2DA = pTPCgain->GetDeDxHisto()->Projection(0,5);
299 TF1 * fpolA = new TF1("fpolA","pol1");
300 TObjArray arrayFitA;
301 his2DA->FitSlicesY(0,0,-1,10,"QN",&arrayFit);
302 delete his2DA;
303 TH1 * his1DA = (TH1*) arrayFit.At(1);
304 his1DA->Fit(fpolA,"QNROB=0.8","QN",0,1);
305 attachSlopeA = fpolA->GetParameter(1);
306 //
307 // MIP position vs. sector
308 //
309 pTPCgain->GetDeDxHisto()->GetAxis(5)->SetRangeUser(-3,0); // C side
310 for(Int_t i = 0; i < 18; i++) { // loop over sectors; correct mapping to be checked!
311 TH1* his1D=0;
312 Float_t phiLow = -TMath::Pi() + i*(20./360.)*(2*TMath::Pi());
313 Float_t phiUp = -TMath::Pi() + (i+1)*(20./360.)*(2*TMath::Pi());
314 pTPCgain->GetDeDxHisto()->GetAxis(1)->SetRangeUser(phiLow,phiUp);
315 his1D = pTPCgain->GetDeDxHisto()->Projection(0);
316 TF1 gausFunc("gausFunc","gaus");
317 his1D->Fit(&gausFunc, "QN");
318 MIPvsSector(i) = gausFunc.GetParameter(1);
319 sector(i)=i;
320 delete his1D;
321 }
322 //
323 pTPCgain->GetDeDxHisto()->GetAxis(5)->SetRangeUser(0,3); // A side
324 for(Int_t i = 0; i < 18; i++) { // loop over sectors; correct mapping to be checked!
325 TH1* his1D=0;
326 Float_t phiLow = -TMath::Pi() + i*(20./360.)*(2*TMath::Pi());
327 Float_t phiUp = -TMath::Pi() + (i+1)*(20./360.)*(2*TMath::Pi());
328 pTPCgain->GetDeDxHisto()->GetAxis(1)->SetRangeUser(phiLow,phiUp);
329 his1D = pTPCgain->GetDeDxHisto()->Projection(0);
330 TF1 gausFunc("gausFunc","gaus");
331 his1D->Fit(&gausFunc, "QN");
332 MIPvsSector(i+18) = gausFunc.GetParameter(1);
333 sector(i+18)=i+18;
334 delete his1D;
335 }
fc567e9f 336 }
337 //
338 printf("Gain QA report\n");
339 printf("MIP mean\t%f\n",MIPmean);
340 printf("MIP resolution\t%f\n",MIPresolution);
341 printf("MIPslopeA\t%f\n",attachSlopeA);
342 printf("MIPslopeC\t%f\n",attachSlopeC);
343 //
344 (*pcstream)<<"tpcQA"<<
345 "MIPattachSlopeC="<<attachSlopeC<<
346 "MIPattachSlopeA="<<attachSlopeA<<
347 "MIPresolution="<<MIPresolution<<
348 "MIPvsSector.="<<&MIPvsSector<<
349 "sector.="<<&sector<<
350 "MIPmean="<<MIPmean;
351
352}