2 The class for calculating the global (not detector specific) quality assurance.
3 It reuses the following TLists from its base class
4 AliQADataMaker::fRecPointsQAList (for keeping the track residuals)
5 AliQADataMaker::fESDsQAList (for keeping global ESD QA data)
12 #include "AliCDBPath.h"
13 #include "AliCDBEntry.h"
14 #include "AliCDBManager.h"
15 #include "AliDetectorRecoParam.h"
16 #include "AliQAChecker.h"
17 #include "AliGlobalQADataMaker.h"
18 #include "AliGeomManager.h"
19 #include "AliESDEvent.h"
21 #include "AliRawReader.h"
22 #include "AliESDVZERO.h"
23 #include "AliMultiplicity.h"
25 ClassImp(AliGlobalQADataMaker)
27 //____________________________________________________________________________
28 void AliGlobalQADataMaker::EndOfDetectorCycle(AliQAv1::TASKINDEX_t task, TObjArray ** list)
30 //Detector specific actions at end of cycle
32 AliQAChecker::Instance()->Run(AliQAv1::kGLOBAL, task, list) ;
35 //____________________________________________________________________________
36 void AliGlobalQADataMaker::InitRaws()
38 // create Raws histograms in Raws subdir
39 ClonePerTrigClass(AliQAv1::kRAWS); // this should be the last line
42 //____________________________________________________________________________
43 void AliGlobalQADataMaker::InitRecoParams()
45 // Get the recoparam form the OCDB
48 AliDebug(AliQAv1::GetQADebugLevel(), Form("Loading reconstruction parameter objects for detector %s", name.Data()));
49 AliCDBPath path(name.Data(),"Calib","RecoParam");
50 AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath());
53 AliDebug(AliQAv1::GetQADebugLevel(), Form("Couldn't find RecoParam entry in OCDB for detector %s",name.Data()));
56 TObject * recoParamObj = entry->GetObject() ;
57 if ( strcmp(recoParamObj->ClassName(), "TObjArray") == 0 ) {
58 // The detector has only one set of reco parameters
59 AliDebug(AliQAv1::GetQADebugLevel(), Form("Array of reconstruction parameters found for detector %s",name.Data()));
60 TObjArray *recoParamArray = static_cast<TObjArray*>(recoParamObj) ;
61 for (Int_t iRP=0; iRP<recoParamArray->GetEntriesFast(); iRP++) {
62 fRecoParam = static_cast<AliDetectorRecoParam*>(recoParamArray->At(iRP)) ;
65 else if (fRecoParam->IsDefault())
69 else if (recoParamObj->InheritsFrom("AliDetectorRecoParam")) {
70 // The detector has only one set of reco parameters
71 // Registering it in AliRecoParam
72 AliDebug(AliQAv1::GetQADebugLevel(), Form("Single set of reconstruction parameters found for detector %s",name.Data()));
73 fRecoParam = static_cast<AliDetectorRecoParam*>(recoParamObj) ;
74 static_cast<AliDetectorRecoParam*>(recoParamObj)->SetAsDefault();
76 AliError(Form("No valid RecoParam object found in the OCDB for detector %s",name.Data()));
82 //____________________________________________________________________________
83 void AliGlobalQADataMaker::InitRecPoints() {
84 //------------------------------------------------------
85 // This function books the histograms of *track*residuals*
86 // as a part of global QA
87 //------------------------------------------------------
88 static Bool_t first = kTRUE ;
89 if ( ! first ) return;
90 const Char_t *name[]={
91 "hGlobalSPD1ResidualsY","SPD1ResidualsZ",
92 "hGlobalSPD2ResidualsY","SPD2ResidualsZ",
93 "hGlobalSDD1ResidualsY","SDD1ResidualsZ",
94 "hGlobalSDD2ResidualsY","SDD2ResidualsZ",
95 "hGlobalSSD1ResidualsY","SSD1ResidualsZ",
96 "hGlobalSSD2ResidualsY","SSD2ResidualsZ",
98 "hGlobalTPC1ResidualsY","TPC1ResidualsZ",
99 "hGlobalTPC2ResidualsY","TPC2ResidualsZ",
101 "hGlobalTRD1ResidualsY","TRD1ResidualsZ",
102 "hGlobalTRD2ResidualsY","TRD2ResidualsZ",
103 "hGlobalTRD3ResidualsY","TRD3ResidualsZ",
104 "hGlobalTRD4ResidualsY","TRD4ResidualsZ",
105 "hGlobalTRD5ResidualsY","TRD5ResidualsZ",
106 "hGlobalTRD6ResidualsY","TRD6ResidualsZ",
108 "hGlobalTOFResidualsY","TOFResidualsZ",
110 "hGlobalPHOS1ResidualsY","PHOS1ResidualsZ",
111 "hGlobalPHOS2ResidualsY","PHOS2ResidualsZ",
113 "hGlobalHMPIDResidualsY","HMPIDResidualsZ",
115 "hGlobalMUONResidualsY","MUONResidualsZ",
117 "hGlobalEMCALResidualsY","EMCALResidualsZ"
119 const Char_t *title[]={
120 "SPD1 residuals Y","SPD1 residuals Z",
121 "SPD2 residuals Y","SPD2 residuals Z",
122 "SDD1 residuals Y","SDD1 residuals Z",
123 "SDD2 residuals Y","SDD2 residuals Z",
124 "SSD1 residuals Y","SSD1 residuals Z",
125 "SSD2 residuals Y","SSD2 residuals Z",
127 "TPC1 residuals Y","TPC1 residuals Z",
128 "TPC2 residuals Y","TPC2 residuals Z",
130 "TRD1 residuals Y","TRD1 residuals Z",
131 "TRD2 residuals Y","TRD2 residuals Z",
132 "TRD3 residuals Y","TRD3 residuals Z",
133 "TRD4 residuals Y","TRD4 residuals Z",
134 "TRD5 residuals Y","TRD5 residuals Z",
135 "TRD6 residuals Y","TRD6 residuals Z",
137 "TOF residuals Y","TOF residuals Z",
139 "PHOS1 residuals Y","PHOS1 residuals Z",
140 "PHOS2 residuals Y","PHOS2 residuals Z",
142 "HMPID residuals Y","HMPID residuals Z",
144 "MUON residuals Y","MUON residuals Z",
146 "EMCAL residuals Y","EMCAL residuals Z"
149 for (Int_t m=1; m<AliGeomManager::kLastLayer; m++) {
151 TH1F *h=new TH1F(name[i],title[i],100,-5.,5.);
152 Add2RecPointsList(h,i);
153 h=new TH1F(name[i+1],title[i+1],100,-5.,5.);
154 Add2RecPointsList(h,i+1);
158 new TH1F("hGlobalSSD1AbsoluteResidualsYNegZ",
159 "SSD1 Absolute Residuals Y Neg Z",100,-2.,2.),40);
161 new TH1F("hGlobalSSD1AbsoluteResidualsZNegZ",
162 "SSD1 Absolute Residuals Z Neg Z",100,-2.,2.),41);
164 new TH1F("hGlobalSSD1AbsoluteResidualsYPosZ",
165 "SSD1 Absolute Residuals Y Pos Z",100,-2.,2.),42);
167 new TH1F("hGlobalSSD1AbsoluteResidualsZPosZ",
168 "SSD1 Absolute Residuals Z Pos Z",100,-2.,2.),43);
172 new TH1F("hGlobalSSD2AbsoluteResidualsYNegZ",
173 "SSD2 Absolute Residuals Y Neg Z",100,-3.,3.),44);
175 new TH1F("hGlobalSSD2AbsoluteResidualsZNegZ",
176 "SSD2 Absolute Residuals Z Neg Z",100,-3.,3.),45);
178 new TH1F("hGlobalSSD2AbsoluteResidualsYPosZ",
179 "SSD2 Absolute Residuals Y Pos Z",100,-3.,3.),46);
181 new TH1F("hGlobalSSD2AbsoluteResidualsZPosZ",
182 "SSD2Absolute Residuals Z Pos Z",100,-3.,3.),47);
186 ClonePerTrigClass(AliQAv1::kRECPOINTS); // this should be the last line
189 //____________________________________________________________________________
190 void AliGlobalQADataMaker::InitESDs() {
191 //------------------------------------------------------
192 // This function books the ESD QA histograms
193 // as a part of global QA
194 //------------------------------------------------------
196 const Bool_t expert = kTRUE ;
197 const Bool_t image = kTRUE ;
200 const Char_t *name[]={
201 "hGlobalPrimaryVertex"
203 const Char_t *title[]={
204 "Z-distribution of the primary vertex"
206 Add2ESDsList(new TH1F(name[0],title[0],100,-20.,20.),kEvt0,!expert,image);
209 {// Cluster related QA
210 const Char_t *name[]={
211 "hGlobalFractionAssignedClustersITS",
212 "hGlobalFractionAssignedClustersTPC",
213 "hGlobalFractionAssignedClustersTRD",
214 "hGlobalClustersPerITSModule"
216 const Char_t *title[]={
217 "Fraction of the assigned clusters in ITS",
218 "Fraction of the assigned clusters in TPC",
219 "Fraction of the assigned clusters in TRD",
220 "Number of clusters per an ITS module"
222 Add2ESDsList(new TH1F(name[0],title[0],100,0.,2.),kClr0, !expert, image);
223 Add2ESDsList(new TH1F(name[1],title[1],100,0.,2.),kClr1, !expert, image);
224 Add2ESDsList(new TH1F(name[2],title[2],100,0.,2.),kClr2, !expert, image);
225 Add2ESDsList(new TH1F(name[3],title[3],2201,-0.5,2200.5),kClr3, !expert, image);
229 const Char_t *name[]={
230 "hGlobalTrackAzimuthe", // kTrk0
231 "hGlobalTrackEta", // kTrk1
232 "hGlobalTPCTrackpT", // kTrk2
233 "hGlobalTPCITSMatchedpT", // kTrk3
234 "hGlobalTPCTOFMatchedpT", // kTrk4
235 "hGlobalTPCITSMatchingProbability", // kTrk5
236 "hGlobalTPCTOFMatchingProbability", // kTrk6
237 "hGlobalTPCsideAposDCA", // kTrk7
238 "hGlobalTPCsideAnegDCA", // kTrk8
239 "hGlobalTPCsideCposDCA", // kTrk9
240 "hGlobalTPCsideCnegDCA" // kTrk10
242 const Char_t *title[]={
243 "Track azimuthal distribution (rad)", // kTrk0
244 "Track pseudo-rapidity distribution", // kTrk1
245 "TPC: track momentum distribution (GeV)", // kTrk2
246 "TPC-ITS matched: track momentum distribution (GeV)", // kTrk3
247 "TPC-TOF matched: track momentum distribution (GeV)", // kTrk4
248 "TPC-ITS track-matching probability", // kTrk5
249 "TPC-TOF track-matching probability", // kTrk6
250 "TPC side A: DCA for the positive tracks (mm)", // kTrk7
251 "TPC side A: DCA for the negative tracks (mm)", // kTrk8
252 "TPC side C: DCA for the positive tracks (mm)", // kTrk9
253 "TPC side C: DCA for the negative tracks (mm)" // kTrk10
255 Add2ESDsList(new TH1F(name[0],title[0],100, 0.,TMath::TwoPi()),kTrk0, !expert, image);
256 Add2ESDsList(new TH1F(name[1],title[1],100,-2.00,2.00),kTrk1, !expert, image);
257 Add2ESDsList(new TH1F(name[2],title[2],50, 0.20,5.00),kTrk2, !expert, image);
258 Add2ESDsList(new TH1F(name[3],title[3],50, 0.20,5.00),kTrk3, !expert, image);
259 Add2ESDsList(new TH1F(name[4],title[4],50, 0.20,5.00),kTrk4, !expert, image);
260 Add2ESDsList(new TH1F(name[5],title[5],50, 0.20,5.00),kTrk5, !expert, image);
261 Add2ESDsList(new TH1F(name[6],title[6],50, 0.20,5.00),kTrk6, !expert, image);
262 Add2ESDsList(new TH1F(name[7],title[7],50, -25.0,25.0),kTrk7, !expert, image);
263 Add2ESDsList(new TH1F(name[8],title[8],50, -25.0,25.0),kTrk8, !expert, image);
264 Add2ESDsList(new TH1F(name[9],title[9],50, -25.0,25.0),kTrk9, !expert, image);
265 Add2ESDsList(new TH1F(name[10],title[10],50, -25.0,25.0),kTrk10, !expert, image);
269 const Char_t *name[]={
270 "hGlobalPromptK0sMass",
271 "hGlobalOfflineK0sMass",
272 "hGlobalPromptLambda0Lambda0BarMass",
273 "hGlobalOfflineLambda0Lambda0BarMass"
275 const Char_t *title[]={
276 "On-the-fly K0s mass (GeV)",
277 "Offline K0s mass (GeV)",
278 "On-the-fly Lambda0 + Lambda0Bar mass (GeV)",
279 "Offline Lambda0 + Lambda0Bar mass (GeV)"
281 Add2ESDsList(new TH1F(name[0],title[0],50, 0.4477,0.5477),kK0on, !expert, image);
282 Add2ESDsList(new TH1F(name[1],title[1],50, 0.4477,0.5477),kK0off, !expert, image);
283 Add2ESDsList(new TH1F(name[2],title[2],50, 1.0657,1.1657),kL0on, !expert, image);
284 Add2ESDsList(new TH1F(name[3],title[3],50, 1.0657,1.1657),kL0off, !expert, image);
288 const Char_t *name[]={
291 "hGlobalTOFTrackingvsMeasured",
292 "hGlobalTPCdEdxvsMomentum"
294 const Char_t *title[]={
295 "ITS: dEdx (ADC) for particles with momentum 0.4 - 0.5 (GeV)",
296 "TPC: dEdx (ADC) for particles with momentum 0.4 - 0.5 (GeV)",
297 "TOF: tracking - measured (ps)",
298 "TPC: dEdx (A.U.) vs momentum (GeV)"
300 Add2ESDsList(new TH1F(name[0],title[0],50,0.00,200.),kPid0, !expert, image);
301 Add2ESDsList(new TH1F(name[1],title[1],50,0.00,100.),kPid1, !expert, image);
302 Add2ESDsList(new TH1F(name[2],title[2],50,-3500.,3500.),kPid2, !expert, image);
303 Add2ESDsList(new TH2F(name[3],title[3],1500,0.05,15.,700,0.,700.),kPid3,!expert,image);
305 {// Multiplicity related QA
306 const Char_t *name[]={
310 const Char_t *title[]={
311 "Multiplicity: V0A vs ITS",
312 "Multiplicity: V0C vs ITS"
314 TH2F *h0=new TH2F(name[0],title[0],41,-0.5,40.5, 33,-0.5,32.5);
315 Add2ESDsList(h0,kMlt0, !expert, image);
316 TH2F *h1=new TH2F(name[1],title[1],41,-0.5,40.5, 33,-0.5,32.5);
317 Add2ESDsList(h1,kMlt1, !expert, image);
320 ClonePerTrigClass(AliQAv1::kESDS); // this should be the last line
323 //____________________________________________________________________________
324 void AliGlobalQADataMaker::MakeRaws(AliRawReader* rawReader)
326 //Fill prepared histograms with Raw digit properties
328 IncEvCountCycleRaws();
329 IncEvCountTotalRaws();
333 //____________________________________________________________________________
334 void AliGlobalQADataMaker::MakeESDs(AliESDEvent * event) {
335 //-----------------------------------------------------------
336 // This function fills the ESD QA histograms
337 // as a part of global QA
338 //-----------------------------------------------------------
340 const AliESDEvent *esd=event;
343 const AliESDVertex *vtx=esd->GetPrimaryVertex();
344 if (!vtx->GetStatus()) return;
346 Double_t xv=vtx->GetXv();
347 Double_t yv=vtx->GetYv();
348 Double_t zv=vtx->GetZv();
349 FillESDsData(kEvt0,zv);
352 Int_t ntrk=esd->GetNumberOfTracks() ;
353 for (Int_t i=0; i<ntrk; i++) {
354 const AliESDtrack *track=esd->GetTrack(i);
356 // Cluster related QA
357 if (track->IsOn(AliESDtrack::kITSrefit)) {
358 Int_t n=track->GetITSclusters(0);
359 FillESDsData(kClr0,Float_t(n)/6.); //6 is the number of ITS layers
362 for (Int_t j=0; j<6; ++j) {
365 if (!track->GetITSModuleIndexInfo(j,idet,sts,xloc,zloc)) continue;
368 if ((sts==1)||(sts==2)||(sts==4)) FillESDsData(kClr3,idet);
371 if (track->IsOn(AliESDtrack::kTPCrefit)) {
372 Int_t n =track->GetTPCNcls();
373 Int_t nf=track->GetTPCNclsF(); // number of crossed TPC pad rows
375 Double_t val = n*1.0/nf;
376 FillESDsData(kClr1,val);
380 if (track->IsOn(AliESDtrack::kTRDrefit)) {
381 Int_t n=track->GetTRDclusters(0);
382 FillESDsData(kClr2,Float_t(n)/(6*24));//(6*24) is the number of TRD time bins
385 Double_t p=track->GetP();
388 if (track->IsOn(AliESDtrack::kTPCrefit)) {
390 track->GetDZ(xv,yv,zv,esd->GetMagneticField(),dz);
391 if ((TMath::Abs(dz[0])<3.) && (TMath::Abs(dz[1])<3.)) { // beam pipe
392 Double_t phi=track->Phi();
393 FillESDsData(kTrk0,phi);
394 Double_t y=track->Eta();
395 FillESDsData(kTrk1,y);
397 if (TMath::Abs(y)<0.9) {
398 FillESDsData(kTrk2,p);
399 if (track->IsOn(AliESDtrack::kITSrefit)) FillESDsData(kTrk3,p);
400 //if (track->IsOn(AliESDtrack::kTOFout)) FillESDsData(kTrk4,p);
401 if (track->GetTOFsignal()>0) FillESDsData(kTrk4,p);
405 const AliExternalTrackParam *tpcTrack=track->GetTPCInnerParam();
406 const AliExternalTrackParam *innTrack=track->GetInnerParam();
410 tpcTrack->GetDZ(xv,yv,zv,esd->GetMagneticField(),dz);
412 if (innTrack->GetZ() > 0)
413 if (innTrack->GetTgl()> 0) { // TPC side A
414 if (tpcTrack->GetSign() > 0) FillESDsData(kTrk7,dz[0]);
415 else FillESDsData(kTrk8,dz[0]);
417 if (innTrack->GetZ() < 0)
418 if (innTrack->GetTgl()< 0) { // TPC side C
419 if (tpcTrack->GetSign() > 0) FillESDsData(kTrk9,dz[0]);
420 else FillESDsData(kTrk10,dz[0]);
425 if ((p>0.4) && (p<0.5)) {
426 if (track->IsOn(AliESDtrack::kITSpid)) {
427 Double_t dedx=track->GetITSsignal();
428 FillESDsData(kPid0,dedx);
430 if (track->IsOn(AliESDtrack::kTPCpid)) {
431 Double_t dedx=track->GetTPCsignal();
432 FillESDsData(kPid1,dedx);
436 if (track->IsOn(AliESDtrack::kITSrefit))
437 if (track->IsOn(AliESDtrack::kTPCrefit))
438 if (track->IsOn(AliESDtrack::kTOFout)) {
440 track->GetDZ(xv,yv,zv,esd->GetMagneticField(),dz);
443 track->GetIntegratedTimes(times);
444 Double_t tof=track->GetTOFsignal()/*-847055 -1771207*/;
445 FillESDsData(kPid2,times[2]-tof);
449 const AliExternalTrackParam *par=track->GetInnerParam();
451 Double_t pp=par->GetP();
452 Double_t dedx=track->GetTPCsignal();
453 FillESDsData(kPid3,pp,dedx);
458 // Multiplicity related QA
459 AliESDVZERO *mltV0 =esd->GetVZEROData();
460 const AliMultiplicity *mltITS=esd->GetMultiplicity();
463 Short_t nv0a=mltV0->GetNbPMV0A();
464 Short_t nv0c=mltV0->GetNbPMV0C();
465 Int_t nits=mltITS->GetNumberOfTracklets();
466 FillESDsData(kMlt0,nits,nv0a);
467 FillESDsData(kMlt1,nits,nv0c);
471 for (int itr = -1; itr<GetNEventTrigClasses(); itr++) {
472 TH1 *tpc = GetMatchingESDsHisto(kTrk2,itr);
473 TH1 *its = GetMatchingESDsHisto(kTrk3,itr);
474 TH1 *tof = GetMatchingESDsHisto(kTrk4,itr);
475 TH1* h5 = GetMatchingESDsHisto(kTrk5,itr);
476 TH1* h6 = GetMatchingESDsHisto(kTrk6,itr);
477 if (h5 && h6 && tpc && its && tof) {
481 h5->Divide(its,tpc,1,1.,"b");
482 h6->Divide(tof,tpc,1,1.,"b");
487 Int_t nV0=esd->GetNumberOfV0s();
488 for (Int_t i=0; i<nV0; i++) {
490 AliESDv0 v0(*esd->GetV0(i));
492 Int_t nidx=TMath::Abs(v0.GetNindex());
493 AliESDtrack *ntrack1=esd->GetTrack(nidx);
494 if (!ntrack1->IsOn(AliESDtrack::kTPCrefit)) continue;
496 Int_t pidx=TMath::Abs(v0.GetPindex());
497 AliESDtrack *ptrack1=esd->GetTrack(pidx);
498 if (!ptrack1->IsOn(AliESDtrack::kTPCrefit)) continue;
500 v0.ChangeMassHypothesis(kK0Short);
501 mass=v0.GetEffMass();
502 if (v0.GetOnFlyStatus())
503 FillESDsData(kK0on,mass);
505 FillESDsData(kK0off,mass);
507 v0.ChangeMassHypothesis(kLambda0);
508 mass=v0.GetEffMass();
509 if (v0.GetOnFlyStatus())
510 FillESDsData(kL0on,mass);
512 FillESDsData(kL0off,mass);
514 v0.ChangeMassHypothesis(kLambda0Bar);
515 mass=v0.GetEffMass();
516 if (v0.GetOnFlyStatus())
517 FillESDsData(kL0on,mass);
519 FillESDsData(kL0off,mass);
522 IncEvCountCycleESDs();
523 IncEvCountTotalESDs();