2 /**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
19 //_________________________________________________________________________
20 // An analysis task to check the PMD data in simulated data
23 //////////////////////////////////////////////////////////////////////////////
34 #include "AliPMDQATask.h"
35 #include "AliPMDUtility.h"
39 //______________________________________________________________________________
40 AliPMDQATask::AliPMDQATask(const char *name) :
41 AliAnalysisTask(name,""),
85 // Input slot #0 works with an Ntuple
86 DefineInput(0, TChain::Class());
87 // Output slot #0 writes into a TH1 container
88 DefineOutput(0, TObjArray::Class()) ;
91 //______________________________________________________________________________
92 AliPMDQATask::AliPMDQATask(const AliPMDQATask& ta) :
93 AliAnalysisTask(ta.GetName(), ""),
96 fOutputContainer(ta.fOutputContainer),
105 fhPMDCP0(ta.fhPMDCP0),
106 fhPMDCP1(ta.fhPMDCP1),
107 fhPMDCP2(ta.fhPMDCP2),
108 fhPMDCP3(ta.fhPMDCP3),
109 fhPMDCP4(ta.fhPMDCP4),
110 fhPMDSM1(ta.fhPMDSM1),
111 fhPMDSM2(ta.fhPMDSM2),
112 fhPMDSM3(ta.fhPMDSM3),
113 fhPMDSM4(ta.fhPMDSM4),
114 fhPMDSM5(ta.fhPMDSM5),
115 fhPMDSM6(ta.fhPMDSM6),
116 fhPMDSM7(ta.fhPMDSM7),
117 fhPMDSM8(ta.fhPMDSM8),
118 fhPMDSM9(ta.fhPMDSM9),
119 fhPMDSM10(ta.fhPMDSM10),
120 fhPMDSM11(ta.fhPMDSM11),
121 fhPMDSM12(ta.fhPMDSM12),
122 fhPMDSM13(ta.fhPMDSM13),
123 fhPMDSM14(ta.fhPMDSM14),
124 fhPMDSM15(ta.fhPMDSM15),
125 fhPMDSM16(ta.fhPMDSM16),
126 fhPMDSM17(ta.fhPMDSM17),
127 fhPMDSM18(ta.fhPMDSM18),
128 fhPMDSM19(ta.fhPMDSM19),
129 fhPMDSM20(ta.fhPMDSM20),
130 fhPMDSM21(ta.fhPMDSM21),
131 fhPMDSM22(ta.fhPMDSM22),
132 fhPMDSM23(ta.fhPMDSM23),
133 fhPMDSM24(ta.fhPMDSM24),
139 //_____________________________________________________________________________
140 AliPMDQATask& AliPMDQATask::operator = (const AliPMDQATask& ap)
142 // assignment operator
144 this->~AliPMDQATask();
145 new(this) AliPMDQATask(ap);
149 //______________________________________________________________________________
150 AliPMDQATask::~AliPMDQATask()
153 fOutputContainer->Clear() ;
154 delete fOutputContainer ;
197 //______________________________________________________________________________
198 void AliPMDQATask::ConnectInputData(const Option_t*)
200 // Initialisation of branch container and histograms
202 AliInfo(Form("*** Initialization of %s", GetName())) ;
205 fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
207 AliError(Form("Input 0 for %s not found\n", GetName()));
211 // One should first check if the branch address was taken by some other task
212 char ** address = (char **)GetBranchAddress(0, "ESD");
214 fESD = (AliESD*)(*address);
217 SetBranchAddress(0, "ESD", &fESD);
221 //________________________________________________________________________
222 void AliPMDQATask::CreateOutputObjects()
228 fhPMDP1 = new TH2F("fhPMDP1","XY of Clusters",100,-100.,100.,100,-100.,100.);
229 fhPMDC2 = new TH1F("fhPMDC2","CPV PHI",200,-1,9);
230 fhPMDP2 = new TH1F("fhPMDP2","PRE PHI",200,-1,9);
231 fhPMDC3 = new TH1F("fhPMDC3","CPV Clus",30,0.,500.);
232 fhPMDP3 = new TH1F("fhPMDP3","PRE N-gammalike",20,0.,500.);
233 fhPMDP4 = new TH1F("fhPMDP4","PRE EDEP",30,0.,1000.);
234 fhPMDC5 = new TH1F("fhPMDC5","CPV n-cell",20,0.,100.);
235 fhPMDP5 = new TH1F("fhPMDP5","PMD n-cell",20,0.,100.);
236 fhPMDCP0 = new TH2F("fhPMDCP0","PRE CLUS Quad.1 vs 2",150,0.,300.,150,0.,300.);
237 fhPMDCP1 = new TH2F("fhPMDCP1","PRE CLUS Quad.3 vs 4",150,0.,300.,150,0.,300.);
238 fhPMDCP2 = new TH2F("fhPMDCP2","PRE EDEP Quad.3 vs 4",50,0.,300.,50,0.,300.);
239 fhPMDCP3 = new TH2F("fhPMDCP3","PRE EDEP vs Tot Clus ",10,0.,1000.,10,0.,300.);
240 fhPMDCP4 = new TH2F("fhPMDCP4","PRE Clus vs CPV Clus ",150,0.,200.,150,0.,200.);
242 fhPMDSM1 = new TH2F("fhPMDSM1","PRE Cluster XY",200,-100,100,200,-100,100);
243 fhPMDSM2 = new TH2F("fhPMDSM2","",999,-100.0,100.0,999,-100.0,100.0);
244 fhPMDSM3 = new TH2F("fhPMDSM3","",999,-100.0,100.0,999,-100.0,100.0);
245 fhPMDSM4 = new TH2F("fhPMDSM4","",999,-100.0,100.0,999,-100.0,100.0);
246 fhPMDSM5 = new TH2F("fhPMDSM5","",999,-100.0,100.0,999,-100.0,100.0);
247 fhPMDSM6 = new TH2F("fhPMDSM6","",999,-100.0,100.0,999,-100.0,100.0);
248 fhPMDSM7 = new TH2F("fhPMDSM7","",999,-100.0,100.0,999,-100.0,100.0);
249 fhPMDSM8 = new TH2F("fhPMDSM8","",999,-100.0,100.0,999,-100.0,100.0);
250 fhPMDSM9 = new TH2F("fhPMDSM9","",999,-100.0,100.0,999,-100.0,100.0);
251 fhPMDSM10 = new TH2F("fhPMDSM10","",999,-100.0,100.0,999,-100.0,100.0);
252 fhPMDSM11 = new TH2F("fhPMDSM11","",999,-100.0,100.0,999,-100.0,100.0);
253 fhPMDSM12 = new TH2F("fhPMDSM12","",999,-100.0,100.0,999,-100.0,100.0);
254 fhPMDSM13 = new TH2F("fhPMDSM13","",999,-100.0,100.0,999,-100.0,100.0);
255 fhPMDSM14 = new TH2F("fhPMDSM14","",999,-100.0,100.0,999,-100.0,100.0);
256 fhPMDSM15 = new TH2F("fhPMDSM15","",999,-100.0,100.0,999,-100.0,100.0);
257 fhPMDSM16 = new TH2F("fhPMDSM16","",999,-100.0,100.0,999,-100.0,100.0);
258 fhPMDSM17 = new TH2F("fhPMDSM17","",999,-100.0,100.0,999,-100.0,100.0);
259 fhPMDSM18 = new TH2F("fhPMDSM18","",999,-100.0,100.0,999,-100.0,100.0);
260 fhPMDSM19 = new TH2F("fhPMDSM19","",999,-100.0,100.0,999,-100.0,100.0);
261 fhPMDSM20 = new TH2F("fhPMDSM20","",999,-100.0,100.0,999,-100.0,100.0);
262 fhPMDSM21 = new TH2F("fhPMDSM21","",999,-100.0,100.0,999,-100.0,100.0);
263 fhPMDSM22 = new TH2F("fhPMDSM22","",999,-100.0,100.0,999,-100.0,100.0);
264 fhPMDSM23 = new TH2F("fhPMDSM23","",999,-100.0,100.0,999,-100.0,100.0);
265 fhPMDSM24 = new TH2F("fhPMDSM24","",999,-100.0,100.0,999,-100.0,100.0);
266 fhPMDSM = new TH1F("fhPMDSM","Plot of all 24 Super Modules",24,0,24);
268 // create output container
270 fOutputContainer = new TObjArray(38) ;
271 fOutputContainer->SetName("PMD") ;
273 fOutputContainer->AddAt(fhPMDP1, 0 );
274 fOutputContainer->AddAt(fhPMDC2, 1 );
275 fOutputContainer->AddAt(fhPMDP2, 2 );
276 fOutputContainer->AddAt(fhPMDC3, 3 );
277 fOutputContainer->AddAt(fhPMDP3, 4 );
278 fOutputContainer->AddAt(fhPMDP4, 5 );
279 fOutputContainer->AddAt(fhPMDC5, 6 );
280 fOutputContainer->AddAt(fhPMDP5, 7 );
281 fOutputContainer->AddAt(fhPMDCP0, 8 );
282 fOutputContainer->AddAt(fhPMDCP1, 9);
283 fOutputContainer->AddAt(fhPMDCP2, 10 );
284 fOutputContainer->AddAt(fhPMDCP3, 11 );
285 fOutputContainer->AddAt(fhPMDCP4, 12 );
287 fOutputContainer->AddAt(fhPMDSM1, 13 );
288 fOutputContainer->AddAt(fhPMDSM2, 14 );
289 fOutputContainer->AddAt(fhPMDSM3, 15 );
290 fOutputContainer->AddAt(fhPMDSM4, 16 );
291 fOutputContainer->AddAt(fhPMDSM5, 17 );
292 fOutputContainer->AddAt(fhPMDSM6, 18 );
293 fOutputContainer->AddAt(fhPMDSM7, 19 );
294 fOutputContainer->AddAt(fhPMDSM8, 20 );
295 fOutputContainer->AddAt(fhPMDSM9, 21 );
296 fOutputContainer->AddAt(fhPMDSM10, 22 );
297 fOutputContainer->AddAt(fhPMDSM11, 23 );
298 fOutputContainer->AddAt(fhPMDSM12, 24 );
299 fOutputContainer->AddAt(fhPMDSM13, 25 );
300 fOutputContainer->AddAt(fhPMDSM14, 26 );
301 fOutputContainer->AddAt(fhPMDSM15, 27 );
302 fOutputContainer->AddAt(fhPMDSM16, 28 );
303 fOutputContainer->AddAt(fhPMDSM17, 29 );
304 fOutputContainer->AddAt(fhPMDSM18, 30 );
305 fOutputContainer->AddAt(fhPMDSM19, 31 );
306 fOutputContainer->AddAt(fhPMDSM20, 32 );
307 fOutputContainer->AddAt(fhPMDSM21, 33 );
308 fOutputContainer->AddAt(fhPMDSM22, 34 );
309 fOutputContainer->AddAt(fhPMDSM23, 35 );
310 fOutputContainer->AddAt(fhPMDSM24, 36 );
311 fOutputContainer->AddAt(fhPMDSM, 37 );
314 //______________________________________________________________________________
315 void AliPMDQATask::Exec(Option_t *)
317 // Processing of one event
319 Long64_t entry = fChain->GetReadEntry() ;
322 AliError("fESD is not connected to the input!") ;
326 if ( !((entry-1)%100) )
327 AliInfo(Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(fChain))->GetFile()->GetName(), entry)) ;
329 // ************************ PMD
331 AliPMDUtility *cc = new AliPMDUtility();
340 Float_t preCluQUAD[4] ;
341 Float_t cpvCluQUAD[4] ;
342 Float_t preADCQUAD[4] ;
343 Float_t cpvADCQUAD[4] ;
344 Float_t preCelQUAD[4] ;
345 Float_t cpvCelQUAD[4] ;
347 Int_t npmdCl = fESD->GetNumberOfPmdTracks();
349 // ****** The loop over PMD clusters
351 for (Int_t kk = 0; kk < 4 ; kk++) {
352 cpvCluQUAD[kk] = 0.0 ;
353 preCluQUAD[kk] = 0.0 ;
354 cpvCelQUAD[kk] = 0.0 ;
355 preCelQUAD[kk] = 0.0 ;
356 preADCQUAD[kk] = 0.0 ;
361 AliESDPmdTrack * pmdtr = fESD->GetPmdTrack(npmdCl);
362 Int_t det = pmdtr->GetDetector();
363 Float_t clsX = pmdtr->GetClusterX();
364 Float_t clsY = pmdtr->GetClusterY();
365 Float_t clsZ = pmdtr->GetClusterZ();
366 Float_t ncell = pmdtr->GetClusterCells();
367 Float_t adc = pmdtr->GetClusterADC();
369 cc->SetXYZ(clsX,clsY,clsZ);
372 Float_t eta = cc->GetEta();
373 Float_t phi = cc->GetPhi();
375 // Calculating S.Module Number from Cluster .
377 CalculateSMN(clsX, clsY, smn);
380 if(smn >= 0 && smn <= 5) {
382 cpvADCQUAD[0] =+ adc ;
383 cpvCelQUAD[0] =+ ncell ;
385 if(smn >= 6 && smn <=11) {
387 cpvADCQUAD[1] =+ adc ;
388 cpvCelQUAD[1] =+ ncell ;
390 if(smn >=12 && smn <=17) {
392 cpvADCQUAD[2] =+ adc ;
393 cpvCelQUAD[2] =+ ncell ;
395 if(smn >=18 && smn <=23) {
397 cpvADCQUAD[3] =+ adc ;
398 cpvCelQUAD[3] =+ ncell ;
401 if(eta >= 2.3 && eta <= 3.5)
408 if(smn >= 0 && smn <= 5) {
410 preADCQUAD[0] =+ adc ;
411 preCelQUAD[0] =+ ncell ;
413 if(smn >= 6 && smn <=11) {
415 preADCQUAD[1] =+ adc ;
416 preCelQUAD[1] =+ ncell ;
418 if(smn >=12 && smn <=17) {
420 preADCQUAD[2] =+ adc ;
421 preCelQUAD[2] =+ ncell ;
423 if(smn >=18 && smn <=23) {
425 preADCQUAD[3] =+ adc ;
426 preCelQUAD[3] =+ ncell ;
430 if(smn == 0) fhPMDSM1->Fill(-clsX,clsY);
431 if(smn == 0) fhPMDSM1->Fill(-clsX,clsY);
432 if(smn == 1) fhPMDSM2->Fill(-clsX,clsY);
433 if(smn == 2) fhPMDSM3->Fill(-clsX,clsY);
434 if(smn == 3) fhPMDSM4->Fill(-clsX,clsY);
435 if(smn == 4) fhPMDSM5->Fill(-clsX,clsY);
436 if(smn == 5) fhPMDSM6->Fill(-clsX,clsY);
437 if(smn == 6) fhPMDSM7->Fill(-clsX,clsY);
438 if(smn == 7) fhPMDSM8->Fill(-clsX,clsY);
439 if(smn == 8) fhPMDSM9->Fill(-clsX,clsY);
440 if(smn == 9) fhPMDSM10->Fill(-clsX,clsY);
441 if(smn ==10) fhPMDSM11->Fill(-clsX,clsY);
442 if(smn ==11) fhPMDSM12->Fill(-clsX,clsY);
443 if(smn ==12) fhPMDSM13->Fill(-clsX,clsY);
444 if(smn ==13) fhPMDSM14->Fill(-clsX,clsY);
445 if(smn ==14) fhPMDSM15->Fill(-clsX,clsY);
446 if(smn ==15) fhPMDSM16->Fill(-clsX,clsY);
447 if(smn ==16) fhPMDSM17->Fill(-clsX,clsY);
448 if(smn ==17) fhPMDSM18->Fill(-clsX,clsY);
449 if(smn ==18) fhPMDSM19->Fill(-clsX,clsY);
450 if(smn ==19) fhPMDSM20->Fill(-clsX,clsY);
451 if(smn ==20) fhPMDSM21->Fill(-clsX,clsY);
452 if(smn ==21) fhPMDSM22->Fill(-clsX,clsY);
453 if(smn ==22) fhPMDSM23->Fill(-clsX,clsY);
454 if(smn ==23) fhPMDSM24->Fill(-clsX,clsY);
456 if(eta >= 2.3 && eta <= 3.5)
460 fhPMDP1->Fill(clsX,clsY);
463 for (Int_t k = 0 ; k < 4 ; k++) {
464 totCPVClus =+ cpvCluQUAD [k] ;
465 totPREClus =+ preCluQUAD [k] ;
466 totCPVCell =+ cpvCelQUAD [k] ;
467 totPRECell =+ preCelQUAD [k] ;
468 totPREEdep =+ preADCQUAD [k] ;
470 Float_t totCPVpreClus = totPREClus + totCPVClus ;
472 // if(eta >= 2.3 && eta <= 3.5) {
473 fhPMDC3->Fill(totCPVClus);
474 fhPMDP3->Fill(totPREClus);
475 fhPMDP4->Fill(totPREEdep);
476 fhPMDP5->Fill(totPRECell);
477 fhPMDCP0->Fill(preCluQUAD[0],preCluQUAD[1]);
478 fhPMDCP1->Fill(preCluQUAD[2],preCluQUAD[3]);
479 fhPMDCP2->Fill(preADCQUAD[2],preADCQUAD[3]);
480 fhPMDCP3->Fill(totPREEdep,totCPVpreClus);
481 fhPMDCP4->Fill(totPREClus,totCPVClus);
489 PostData(0, fOutputContainer);
492 //______________________________________________________________________________
493 void AliPMDQATask::Terminate(Option_t *)
495 // Processing when the event loop is ended
496 fOutputContainer = (TObjArray*)GetOutputData(0);
498 fhPMDP1 = (TH2F*)fOutputContainer->At(0);
499 fhPMDC2 = (TH1F*)fOutputContainer->At(1);
500 fhPMDP2 = (TH1F*)fOutputContainer->At(2);
501 fhPMDC3 = (TH1F*)fOutputContainer->At(3);
502 fhPMDP3 = (TH1F*)fOutputContainer->At(4);
503 fhPMDP4 = (TH1F*)fOutputContainer->At(5);
504 fhPMDC5 = (TH1F*)fOutputContainer->At(6);
505 fhPMDP5 = (TH1F*)fOutputContainer->At(7);
506 fhPMDCP0 = (TH2F*)fOutputContainer->At(8);
507 fhPMDCP1 = (TH2F*)fOutputContainer->At(9);
508 fhPMDCP2 = (TH2F*)fOutputContainer->At(10);
509 fhPMDCP3 = (TH2F*)fOutputContainer->At(11);
510 fhPMDCP4 = (TH2F*)fOutputContainer->At(12);
512 fhPMDSM1 = (TH2F*)fOutputContainer->At(13);
513 fhPMDSM2 = (TH2F*)fOutputContainer->At(14);
514 fhPMDSM3 = (TH2F*)fOutputContainer->At(15);
515 fhPMDSM4 = (TH2F*)fOutputContainer->At(16);
516 fhPMDSM5 = (TH2F*)fOutputContainer->At(17);
517 fhPMDSM6 = (TH2F*)fOutputContainer->At(18);
518 fhPMDSM7 = (TH2F*)fOutputContainer->At(19);
519 fhPMDSM8 = (TH2F*)fOutputContainer->At(20);
520 fhPMDSM9 = (TH2F*)fOutputContainer->At(21);
521 fhPMDSM10 = (TH2F*)fOutputContainer->At(22);
522 fhPMDSM11 = (TH2F*)fOutputContainer->At(23);
523 fhPMDSM12 = (TH2F*)fOutputContainer->At(24);
524 fhPMDSM13 = (TH2F*)fOutputContainer->At(25);
525 fhPMDSM14 = (TH2F*)fOutputContainer->At(26);
526 fhPMDSM15 = (TH2F*)fOutputContainer->At(27);
527 fhPMDSM16 = (TH2F*)fOutputContainer->At(28);
528 fhPMDSM17 = (TH2F*)fOutputContainer->At(29);
529 fhPMDSM18 = (TH2F*)fOutputContainer->At(30);
530 fhPMDSM19 = (TH2F*)fOutputContainer->At(31);
531 fhPMDSM20 = (TH2F*)fOutputContainer->At(32);
532 fhPMDSM21 = (TH2F*)fOutputContainer->At(33);
533 fhPMDSM22 = (TH2F*)fOutputContainer->At(34);
534 fhPMDSM23 = (TH2F*)fOutputContainer->At(35);
535 fhPMDSM24 = (TH2F*)fOutputContainer->At(36);
536 fhPMDSM = (TH1F*)fOutputContainer->At(37);
538 Bool_t problem = kFALSE ;
539 AliInfo(Form(" *** %s Report:", GetName())) ;
541 gStyle->SetOptStat(110000);
542 gStyle->SetOptFit(1);
544 TCanvas *cPMD0 = new TCanvas("cPMD0","PMD ESD Test #1", 10,10, 600, 600);
545 cPMD0->Range(-100, -100,100 ,100 );
546 fhPMDSM1->SetMarkerColor(2);
548 fhPMDSM1->GetXaxis()->SetTitle("Cluster X");
549 fhPMDSM1->GetYaxis()->SetTitle("Cluster Y");
550 fhPMDSM2->SetMarkerColor(2);
551 fhPMDSM2->Draw("same");
552 fhPMDSM3->SetMarkerColor(2);
553 fhPMDSM3->Draw("same");
554 fhPMDSM4->SetMarkerColor(2);
555 fhPMDSM4->Draw("same");
556 fhPMDSM5->SetMarkerColor(2);
557 fhPMDSM5->Draw("same");
558 fhPMDSM6->SetMarkerColor(2);
559 fhPMDSM6->Draw("same");
560 fhPMDSM7->SetMarkerColor(4);
561 fhPMDSM7->Draw("same");
562 fhPMDSM8->SetMarkerColor(4);
563 fhPMDSM8->Draw("same");
564 fhPMDSM9->SetMarkerColor(4);
565 fhPMDSM9->Draw("same");
566 fhPMDSM10->SetMarkerColor(4);
567 fhPMDSM10->Draw("same");
568 fhPMDSM11->SetMarkerColor(4);
569 fhPMDSM11->Draw("same");
570 fhPMDSM12->SetMarkerColor(4);
571 fhPMDSM12->Draw("same");
572 fhPMDSM13->SetMarkerColor(6);
573 fhPMDSM13->Draw("same");
574 fhPMDSM14->SetMarkerColor(6);
575 fhPMDSM14->Draw("same");
576 fhPMDSM15->SetMarkerColor(6);
577 fhPMDSM15->Draw("same");
578 fhPMDSM16->SetMarkerColor(6);
579 fhPMDSM16->Draw("same");
580 fhPMDSM17->SetMarkerColor(6);
581 fhPMDSM17->Draw("same");
582 fhPMDSM18->SetMarkerColor(6);
583 fhPMDSM18->Draw("same");
584 fhPMDSM19->SetMarkerColor(8);
585 fhPMDSM19->Draw("same");
586 fhPMDSM20->SetMarkerColor(8);
587 fhPMDSM20->Draw("same");
588 fhPMDSM21->SetMarkerColor(8);
589 fhPMDSM21->Draw("same");
590 fhPMDSM22->SetMarkerColor(8);
591 fhPMDSM22->Draw("same");
592 fhPMDSM23->SetMarkerColor(8);
593 fhPMDSM23->Draw("same");
594 fhPMDSM24->SetMarkerColor(8);
595 fhPMDSM24->Draw("same");
598 DrawPMDBoundarySM1();
599 DrawPMDBoundarySM2();
600 DrawPMDBoundarySM3();
601 DrawPMDBoundarySM4();
602 cPMD0->Print("ClusterXY.eps");
604 TCanvas *cPMD1 = new TCanvas("cPMD1"," PMD ESD Test #2",10, 10, 600,600);
607 cPMD1->SetFillColor(0);
608 fhPMDC2->SetLineColor(4);
611 fhPMDP2->SetLineColor(2);
613 cPMD1->Print("CPVPREphi.eps");
615 TCanvas *cPMD2 = new TCanvas("cPMD2","PMD ESD test #3",10, 10, 600, 600);
617 fhPMDSM->SetFillColor(2);
619 cPMD2->Print("AllSMN.eps");
621 TCanvas *cPMD3 = new TCanvas("cPMD3", "PMD ESD test #4",10, 10, 600, 600);
624 fhPMDCP0->SetMarkerColor(9);
627 fhPMDCP1->SetMarkerColor(6);
630 fhPMDP3->SetLineColor(2);
633 fhPMDCP4->SetMarkerColor(3);
635 cPMD3->Print("CPVPREClus.eps");
637 TCanvas *cPMD4 = new TCanvas("cPMD4","PMD ESD test #5", 10, 10, 600, 600);
640 fhPMDC3->SetLineColor(4);
643 fhPMDP4->SetLineColor(2);
645 cPMD4->Print("CPVPREAdc.eps");
648 sprintf(line, ".!tar -zcf %s.tar.gz *.eps", GetName()) ;
649 gROOT->ProcessLine(line);
651 AliInfo(Form("!!! All the eps files are in %s.tar.gz !!!", GetName())) ;
655 report="Problems found, please check!!!";
659 AliInfo(Form("*** %s Summary Report: %s \n",GetName(), report.Data())) ;
662 //______________________________________________________________________________
663 void AliPMDQATask::CalculateSMN( Float_t clsX, Float_t clsY, Int_t & smn) const
665 Double_t xcon[96] = {75.133, 54.204, 53.254, 32.326, 31.376,10.447,
666 75.133, 54.204, 53.254, 32.326, 31.376,10.447,
667 75.133, 54.204, 53.254, 32.326, 31.376,10.447,
668 75.133, 54.204, 53.254, 32.326, 31.376,10.447,
669 -75.133, -54.204, -53.254, -32.326, -31.376,-10.447,
670 -75.133, -54.204, -53.254, -32.326, -31.376,-10.447,
671 -75.133, -54.204, -53.254, -32.326, -31.376,-10.447,
672 -75.133, -54.204, -53.254, -32.326, -31.376,-10.447,
673 9.167, -32.543, -33.493, -75.133,
674 9.167, -32.543, -33.493, -75.133,
675 9.167, -32.543, -33.493, -75.133,
676 9.167, -32.543, -33.493, -75.133,
677 9.167, -32.543, -33.493, -75.133,
678 9.167, -32.543, -33.493, -75.133,
679 -9.167, 32.543, 33.493, 75.133,
680 -9.167, 32.543, 33.493, 75.133,
681 -9.167, 32.543, 33.493, 75.133,
682 -9.167, 32.543, 33.493, 75.133,
683 -9.167, 32.543, 33.493, 75.133,
684 -9.167, 32.543, 33.493, 75.133};
686 Double_t ycon[96] = {86.475, 86.475, 86.475, 86.475, 86.475, 86.475,
687 38.225, 38.225, 38.225, 38.225, 38.225, 38.225,
688 37.325, 37.325, 37.325, 37.325, 37.325, 37.325,
689 -10.925, -10.925, -10.925, -10.925, -10.925, -10.925,
690 -86.475, -86.475, -86.475, -86.475, -86.475, -86.475,
691 -38.225, -38.225, -38.225, -38.225, -38.225, -38.225,
692 -37.325, -37.325, -37.325, -37.325, -37.325, -37.325,
693 10.925, 10.925, 10.925, 10.925, 10.925, 10.925,
694 86.475, 86.475, 86.475, 86.475,
695 62.225, 62.225, 62.225, 62.225,
696 61.325, 61.325, 61.325, 61.325,
697 37.075, 37.075, 37.075, 37.075,
698 36.175, 36.175, 36.175, 36.175,
699 11.925, 11.925, 11.925 , 11.925,
700 -86.475, -86.475, -86.475, -86.475,
701 -62.225, -62.225, -62.225, -62.225,
702 -61.325, -61.325, -61.325, -61.325,
703 -37.075, -37.075, -37.075, -37.075,
704 -36.175, -36.175, -36.175, -36.175,
705 -11.925, -11.925, -11.925 , -11.925 };
707 if((clsX <= xcon[0]) && (clsX >= xcon[1]) &&
708 (clsY <= ycon[0]) && (clsY >= ycon[6])) smn = 0 ;
710 else if((clsX <=xcon[2]) && (clsX >= xcon[3]) &&
711 (clsY <= ycon[1]) && (clsY >= ycon[7]))smn = 1 ;
713 else if((clsX <=xcon[4]) && (clsX >= xcon[5]) &&
714 (clsY <= ycon[3]) && (clsY >= ycon[8]))smn = 2 ;
716 else if((clsX <= xcon[0]) && (clsX >= xcon[1]) &&
717 (clsY <= ycon[12]) && (clsY >= ycon[18])) smn = 3 ;
719 else if((clsX <=xcon[2]) && (clsX >= xcon[3]) &&
720 (clsY <= ycon[12]) && (clsY >= ycon[18]))smn = 4 ;
722 else if((clsX <=xcon[4]) && (clsX >= xcon[5]) &&
723 (clsY <= ycon[12]) && (clsY >= ycon[18]))smn = 5 ;
724 //------------------------------------------------------------------
725 else if((clsX >= xcon[24]) && (clsX <= xcon[25]) &&
726 (clsY >= ycon[24]) && (clsY <= ycon[30])) smn = 6 ;
728 else if((clsX >=xcon[26]) && (clsX <= xcon[27]) &&
729 (clsY >= ycon[25]) && (clsY <= ycon[31]))smn = 7 ;
731 else if((clsX >=xcon[28]) && (clsX <= xcon[29]) &&
732 (clsY >= ycon[26]) && (clsY <= ycon[32]))smn = 8 ;
734 else if((clsX >= xcon[24]) && (clsX <= xcon[25]) &&
735 (clsY >= ycon[36]) && (clsY <= ycon[42])) smn = 9 ;
737 else if((clsX >=xcon[26]) && (clsX <= xcon[27]) &&
738 (clsY >= ycon[36]) && (clsY <= ycon[42]))smn = 10;
740 else if((clsX >=xcon[28]) && (clsX <= xcon[29]) &&
741 (clsY >= ycon[36]) && (clsY <= ycon[42]))smn = 11;
742 //------------------------------------------------------------------
743 else if((clsX <= xcon[48]) && (clsX >= xcon[49]) &&
744 (clsY <= ycon[48]) && (clsY >= ycon[52])) smn = 12 ;
746 else if((clsX <=xcon[50]) && (clsX >= xcon[51]) &&
747 (clsY <= ycon[48]) && (clsY >= ycon[52]))smn = 13 ;
749 else if((clsX <=xcon[48]) && (clsX >= xcon[49]) &&
750 (clsY <= ycon[56]) && (clsY >= ycon[60]))smn = 14 ;
752 else if((clsX <=xcon[50]) && (clsX >= xcon[51]) &&
753 (clsY <= ycon[56]) && (clsY >= ycon[60]))smn = 15 ;
755 else if((clsX <=xcon[48]) && (clsX >= xcon[49]) &&
756 (clsY <= ycon[64]) && (clsY >= ycon[68]))smn = 16 ;
758 else if((clsX <=xcon[50]) && (clsX >= xcon[51]) &&
759 (clsY <= ycon[64]) && (clsY >= ycon[68]))smn = 17 ;
760 //--------------------------------------------------------------
761 else if((clsX >= xcon[72]) && (clsX <= xcon[73]) &&
762 (clsY >= ycon[72]) && (clsY <= ycon[76])) smn = 18 ;
764 else if((clsX >=xcon[74]) && (clsX <= xcon[75]) &&
765 (clsY >= ycon[72]) && (clsY <= ycon[76]))smn = 19 ;
767 else if((clsX >=xcon[72]) && (clsX <= xcon[73]) &&
768 (clsY >= ycon[80]) && (clsY <= ycon[84]))smn = 20 ;
770 else if((clsX >=xcon[74]) && (clsX <= xcon[75]) &&
771 (clsY >= ycon[80]) && (clsY <= ycon[84]))smn = 21;
773 else if((clsX >= xcon[72]) && (clsX <= xcon[73]) &&
774 (clsY >= ycon[88]) && (clsY <= ycon[92])) smn = 22 ;
776 else if((clsX >=xcon[74]) && (clsX <= xcon[75]) &&
777 (clsY >= ycon[88]) && (clsY <= ycon[92]))smn = 23 ;
781 //______________________________________________________________________________
782 void AliPMDQATask::DrawPMDBoundary() const
784 // Draw PMD boundaries
786 gStyle->SetLineWidth(2);
787 gStyle->SetLineColor(2);
789 l = new TLine(75.1333, 86.475, -75.1333, 86.475); l->Draw("same");
790 l = new TLine(-75.1333, 86.470,-75.1333, -86.475); l->Draw("same");
791 l = new TLine(-75.1333, -86.475,75.1333, -86.475); l->Draw("same");
792 l = new TLine(75.1333, -86.475,75.1333, 86.475); l->Draw("same");
795 //______________________________________________________________________________
796 void AliPMDQATask::DrawPMDBoundarySM1() const
798 // Draw boundaries of Super Module 1
800 gStyle->SetLineWidth(1);
801 gStyle->SetLineColor(4);
803 l = new TLine(-75.1333, 86.475, -10.447, 86.475); l->Draw("same");
804 l = new TLine(-10.447, 86.475, -10.446, -10.925); l->Draw("same");
805 l = new TLine(-10.446, -10.925, -75.1333,-10.925); l->Draw("same");
806 l = new TLine(-75.1333,-10.925, -75.1333, 86.475); l->Draw("same");
809 //______________________________________________________________________________
810 void AliPMDQATask::DrawPMDBoundarySM2() const
812 // Draw boundaries of Super Module 2
814 gStyle->SetLineWidth(1);
815 gStyle->SetLineColor(4);
817 l = new TLine(75.1333, -86.475, 10.446, -86.475); l->Draw("same");
818 l = new TLine(10.446, -86.475, 10.446, 10.925); l->Draw("same");
819 l = new TLine(10.446, 10.925, 75.1333, 10.925); l->Draw("same");
820 l = new TLine(75.1333, 10.925, 75.1333, -86.475); l->Draw("same");
824 //______________________________________________________________________________
825 void AliPMDQATask::DrawPMDBoundarySM3() const
827 // Draw boundaries of Super Module 3
829 gStyle->SetLineWidth(1);
830 gStyle->SetLineColor(1);
832 l = new TLine( -9.167, 86.475, 75.1333, 86.475); l->Draw("same");
833 l = new TLine(75.1333,86.475, 75.1333, 11.925); l->Draw("same");
834 l = new TLine(75.1333,11.925, -9.167, 11.925); l->Draw("same");
835 l = new TLine( -9.167, 11.925, -9.167, 86.475); l->Draw("same");
838 //______________________________________________________________________________
839 void AliPMDQATask::DrawPMDBoundarySM4() const
841 // Draw boundaries of Super Module 4
843 gStyle->SetLineWidth(1);
844 gStyle->SetLineColor(1);
846 l = new TLine(9.167, -86.475, -75.1333,-86.475); l->Draw("same");
847 l = new TLine(-75.1333,-86.475, -75.1333,-11.925); l->Draw("same");
848 l = new TLine(-75.1333,-11.925, 9.167, -11.925); l->Draw("same");
849 l = new TLine(9.167, -11.925, 9.167, -86.475); l->Draw("same");