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,""),
84 // Input slot #0 works with an Ntuple
85 DefineInput(0, TChain::Class());
86 // Output slot #0 writes into a TH1 container
87 DefineOutput(0, TObjArray::Class()) ;
90 //______________________________________________________________________________
91 AliPMDQATask::~AliPMDQATask()
94 fOutputContainer->Clear() ;
95 delete fOutputContainer ;
138 //______________________________________________________________________________
139 void AliPMDQATask::ConnectInputData(const Option_t*)
141 // Initialisation of branch container and histograms
143 AliInfo(Form("*** Initialization of %s", GetName())) ;
146 fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
148 AliError(Form("Input 0 for %s not found\n", GetName()));
152 // One should first check if the branch address was taken by some other task
153 char ** address = (char **)GetBranchAddress(0, "ESD");
155 fESD = (AliESD*)(*address);
158 SetBranchAddress(0, "ESD", &fESD);
162 //________________________________________________________________________
163 void AliPMDQATask::CreateOutputObjects()
169 fhPMDP1 = new TH2F("fhPMDP1","XY of Clusters",100,-100.,100.,100,-100.,100.);
170 fhPMDC2 = new TH1F("fhPMDC2","CPV PHI",200,-1,9);
171 fhPMDP2 = new TH1F("fhPMDP2","PRE PHI",200,-1,9);
172 fhPMDC3 = new TH1F("fhPMDC3","CPV Clus",30,0.,500.);
173 fhPMDP3 = new TH1F("fhPMDP3","PRE N-gammalike",20,0.,500.);
174 fhPMDP4 = new TH1F("fhPMDP4","PRE EDEP",30,0.,1000.);
175 fhPMDC5 = new TH1F("fhPMDC5","CPV n-cell",20,0.,100.);
176 fhPMDP5 = new TH1F("fhPMDP5","PMD n-cell",20,0.,100.);
177 fhPMDCP0 = new TH2F("fhPMDCP0","PRE CLUS Quad.1 vs 2",150,0.,300.,150,0.,300.);
178 fhPMDCP1 = new TH2F("fhPMDCP1","PRE CLUS Quad.3 vs 4",150,0.,300.,150,0.,300.);
179 fhPMDCP2 = new TH2F("fhPMDCP2","PRE EDEP Quad.3 vs 4",50,0.,300.,50,0.,300.);
180 fhPMDCP3 = new TH2F("fhPMDCP3","PRE EDEP vs Tot Clus ",10,0.,1000.,10,0.,300.);
181 fhPMDCP4 = new TH2F("fhPMDCP4","PRE Clus vs CPV Clus ",150,0.,200.,150,0.,200.);
183 fhPMDSM1 = new TH2F("fhPMDSM1","PRE Cluster XY",200,-100,100,200,-100,100);
184 fhPMDSM2 = new TH2F("fhPMDSM2","",999,-100.0,100.0,999,-100.0,100.0);
185 fhPMDSM3 = new TH2F("fhPMDSM3","",999,-100.0,100.0,999,-100.0,100.0);
186 fhPMDSM4 = new TH2F("fhPMDSM4","",999,-100.0,100.0,999,-100.0,100.0);
187 fhPMDSM5 = new TH2F("fhPMDSM5","",999,-100.0,100.0,999,-100.0,100.0);
188 fhPMDSM6 = new TH2F("fhPMDSM6","",999,-100.0,100.0,999,-100.0,100.0);
189 fhPMDSM7 = new TH2F("fhPMDSM7","",999,-100.0,100.0,999,-100.0,100.0);
190 fhPMDSM8 = new TH2F("fhPMDSM8","",999,-100.0,100.0,999,-100.0,100.0);
191 fhPMDSM9 = new TH2F("fhPMDSM9","",999,-100.0,100.0,999,-100.0,100.0);
192 fhPMDSM10 = new TH2F("fhPMDSM10","",999,-100.0,100.0,999,-100.0,100.0);
193 fhPMDSM11 = new TH2F("fhPMDSM11","",999,-100.0,100.0,999,-100.0,100.0);
194 fhPMDSM12 = new TH2F("fhPMDSM12","",999,-100.0,100.0,999,-100.0,100.0);
195 fhPMDSM13 = new TH2F("fhPMDSM13","",999,-100.0,100.0,999,-100.0,100.0);
196 fhPMDSM14 = new TH2F("fhPMDSM14","",999,-100.0,100.0,999,-100.0,100.0);
197 fhPMDSM15 = new TH2F("fhPMDSM15","",999,-100.0,100.0,999,-100.0,100.0);
198 fhPMDSM16 = new TH2F("fhPMDSM16","",999,-100.0,100.0,999,-100.0,100.0);
199 fhPMDSM17 = new TH2F("fhPMDSM17","",999,-100.0,100.0,999,-100.0,100.0);
200 fhPMDSM18 = new TH2F("fhPMDSM18","",999,-100.0,100.0,999,-100.0,100.0);
201 fhPMDSM19 = new TH2F("fhPMDSM19","",999,-100.0,100.0,999,-100.0,100.0);
202 fhPMDSM20 = new TH2F("fhPMDSM20","",999,-100.0,100.0,999,-100.0,100.0);
203 fhPMDSM21 = new TH2F("fhPMDSM21","",999,-100.0,100.0,999,-100.0,100.0);
204 fhPMDSM22 = new TH2F("fhPMDSM22","",999,-100.0,100.0,999,-100.0,100.0);
205 fhPMDSM23 = new TH2F("fhPMDSM23","",999,-100.0,100.0,999,-100.0,100.0);
206 fhPMDSM24 = new TH2F("fhPMDSM24","",999,-100.0,100.0,999,-100.0,100.0);
207 fhPMDSM = new TH1F("fhPMDSM","Plot of all 24 Super Modules",24,0,24);
209 // create output container
211 fOutputContainer = new TObjArray(38) ;
212 fOutputContainer->SetName("PMD") ;
214 fOutputContainer->AddAt(fhPMDP1, 0 );
215 fOutputContainer->AddAt(fhPMDC2, 1 );
216 fOutputContainer->AddAt(fhPMDP2, 2 );
217 fOutputContainer->AddAt(fhPMDC3, 3 );
218 fOutputContainer->AddAt(fhPMDP3, 4 );
219 fOutputContainer->AddAt(fhPMDP4, 5 );
220 fOutputContainer->AddAt(fhPMDC5, 6 );
221 fOutputContainer->AddAt(fhPMDP5, 7 );
222 fOutputContainer->AddAt(fhPMDCP0, 8 );
223 fOutputContainer->AddAt(fhPMDCP1, 9);
224 fOutputContainer->AddAt(fhPMDCP2, 10 );
225 fOutputContainer->AddAt(fhPMDCP3, 11 );
226 fOutputContainer->AddAt(fhPMDCP4, 12 );
228 fOutputContainer->AddAt(fhPMDSM1, 13 );
229 fOutputContainer->AddAt(fhPMDSM2, 14 );
230 fOutputContainer->AddAt(fhPMDSM3, 15 );
231 fOutputContainer->AddAt(fhPMDSM4, 16 );
232 fOutputContainer->AddAt(fhPMDSM5, 17 );
233 fOutputContainer->AddAt(fhPMDSM6, 18 );
234 fOutputContainer->AddAt(fhPMDSM7, 19 );
235 fOutputContainer->AddAt(fhPMDSM8, 20 );
236 fOutputContainer->AddAt(fhPMDSM9, 21 );
237 fOutputContainer->AddAt(fhPMDSM10, 22 );
238 fOutputContainer->AddAt(fhPMDSM11, 23 );
239 fOutputContainer->AddAt(fhPMDSM12, 24 );
240 fOutputContainer->AddAt(fhPMDSM13, 25 );
241 fOutputContainer->AddAt(fhPMDSM14, 26 );
242 fOutputContainer->AddAt(fhPMDSM15, 27 );
243 fOutputContainer->AddAt(fhPMDSM16, 28 );
244 fOutputContainer->AddAt(fhPMDSM17, 29 );
245 fOutputContainer->AddAt(fhPMDSM18, 30 );
246 fOutputContainer->AddAt(fhPMDSM19, 31 );
247 fOutputContainer->AddAt(fhPMDSM20, 32 );
248 fOutputContainer->AddAt(fhPMDSM21, 33 );
249 fOutputContainer->AddAt(fhPMDSM22, 34 );
250 fOutputContainer->AddAt(fhPMDSM23, 35 );
251 fOutputContainer->AddAt(fhPMDSM24, 36 );
252 fOutputContainer->AddAt(fhPMDSM, 37 );
255 //______________________________________________________________________________
256 void AliPMDQATask::Exec(Option_t *)
258 // Processing of one event
260 Long64_t entry = fChain->GetReadEntry() ;
263 AliError("fESD is not connected to the input!") ;
267 if ( !((entry-1)%100) )
268 AliInfo(Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(fChain))->GetFile()->GetName(), entry)) ;
270 // ************************ PMD
272 AliPMDUtility *cc = new AliPMDUtility();
281 Float_t preCluQUAD[4] ;
282 Float_t cpvCluQUAD[4] ;
283 Float_t preADCQUAD[4] ;
284 Float_t cpvADCQUAD[4] ;
285 Float_t preCelQUAD[4] ;
286 Float_t cpvCelQUAD[4] ;
288 Int_t npmdCl = fESD->GetNumberOfPmdTracks();
290 // ****** The loop over PMD clusters
292 for (Int_t kk = 0; kk < 4 ; kk++) {
293 cpvCluQUAD[kk] = 0.0 ;
294 preCluQUAD[kk] = 0.0 ;
295 cpvCelQUAD[kk] = 0.0 ;
296 preCelQUAD[kk] = 0.0 ;
297 preADCQUAD[kk] = 0.0 ;
302 AliESDPmdTrack * pmdtr = fESD->GetPmdTrack(npmdCl);
303 Int_t det = pmdtr->GetDetector();
304 Float_t clsX = pmdtr->GetClusterX();
305 Float_t clsY = pmdtr->GetClusterY();
306 Float_t clsZ = pmdtr->GetClusterZ();
307 Float_t ncell = pmdtr->GetClusterCells();
308 Float_t adc = pmdtr->GetClusterADC();
310 cc->SetXYZ(clsX,clsY,clsZ);
313 Float_t eta = cc->GetEta();
314 Float_t phi = cc->GetPhi();
316 // Calculating S.Module Number from Cluster .
318 CalculateSMN(clsX, clsY, smn);
321 if(smn >= 0 && smn <= 5) {
323 cpvADCQUAD[0] =+ adc ;
324 cpvCelQUAD[0] =+ ncell ;
326 if(smn >= 6 && smn <=11) {
328 cpvADCQUAD[1] =+ adc ;
329 cpvCelQUAD[1] =+ ncell ;
331 if(smn >=12 && smn <=17) {
333 cpvADCQUAD[2] =+ adc ;
334 cpvCelQUAD[2] =+ ncell ;
336 if(smn >=18 && smn <=23) {
338 cpvADCQUAD[3] =+ adc ;
339 cpvCelQUAD[3] =+ ncell ;
342 if(eta >= 2.3 && eta <= 3.5)
349 if(smn >= 0 && smn <= 5) {
351 preADCQUAD[0] =+ adc ;
352 preCelQUAD[0] =+ ncell ;
354 if(smn >= 6 && smn <=11) {
356 preADCQUAD[1] =+ adc ;
357 preCelQUAD[1] =+ ncell ;
359 if(smn >=12 && smn <=17) {
361 preADCQUAD[2] =+ adc ;
362 preCelQUAD[2] =+ ncell ;
364 if(smn >=18 && smn <=23) {
366 preADCQUAD[3] =+ adc ;
367 preCelQUAD[3] =+ ncell ;
371 if(smn == 0) fhPMDSM1->Fill(-clsX,clsY);
372 if(smn == 0) fhPMDSM1->Fill(-clsX,clsY);
373 if(smn == 1) fhPMDSM2->Fill(-clsX,clsY);
374 if(smn == 2) fhPMDSM3->Fill(-clsX,clsY);
375 if(smn == 3) fhPMDSM4->Fill(-clsX,clsY);
376 if(smn == 4) fhPMDSM5->Fill(-clsX,clsY);
377 if(smn == 5) fhPMDSM6->Fill(-clsX,clsY);
378 if(smn == 6) fhPMDSM7->Fill(-clsX,clsY);
379 if(smn == 7) fhPMDSM8->Fill(-clsX,clsY);
380 if(smn == 8) fhPMDSM9->Fill(-clsX,clsY);
381 if(smn == 9) fhPMDSM10->Fill(-clsX,clsY);
382 if(smn ==10) fhPMDSM11->Fill(-clsX,clsY);
383 if(smn ==11) fhPMDSM12->Fill(-clsX,clsY);
384 if(smn ==12) fhPMDSM13->Fill(-clsX,clsY);
385 if(smn ==13) fhPMDSM14->Fill(-clsX,clsY);
386 if(smn ==14) fhPMDSM15->Fill(-clsX,clsY);
387 if(smn ==15) fhPMDSM16->Fill(-clsX,clsY);
388 if(smn ==16) fhPMDSM17->Fill(-clsX,clsY);
389 if(smn ==17) fhPMDSM18->Fill(-clsX,clsY);
390 if(smn ==18) fhPMDSM19->Fill(-clsX,clsY);
391 if(smn ==19) fhPMDSM20->Fill(-clsX,clsY);
392 if(smn ==20) fhPMDSM21->Fill(-clsX,clsY);
393 if(smn ==21) fhPMDSM22->Fill(-clsX,clsY);
394 if(smn ==22) fhPMDSM23->Fill(-clsX,clsY);
395 if(smn ==23) fhPMDSM24->Fill(-clsX,clsY);
397 if(eta >= 2.3 && eta <= 3.5)
401 fhPMDP1->Fill(clsX,clsY);
404 for (Int_t k = 0 ; k < 4 ; k++) {
405 totCPVClus =+ cpvCluQUAD [k] ;
406 totPREClus =+ preCluQUAD [k] ;
407 totCPVCell =+ cpvCelQUAD [k] ;
408 totPRECell =+ preCelQUAD [k] ;
409 totPREEdep =+ preADCQUAD [k] ;
411 Float_t totCPVpreClus = totPREClus + totCPVClus ;
413 // if(eta >= 2.3 && eta <= 3.5) {
414 fhPMDC3->Fill(totCPVClus);
415 fhPMDP3->Fill(totPREClus);
416 fhPMDP4->Fill(totPREEdep);
417 fhPMDP5->Fill(totPRECell);
418 fhPMDCP0->Fill(preCluQUAD[0],preCluQUAD[1]);
419 fhPMDCP1->Fill(preCluQUAD[2],preCluQUAD[3]);
420 fhPMDCP2->Fill(preADCQUAD[2],preADCQUAD[3]);
421 fhPMDCP3->Fill(totPREEdep,totCPVpreClus);
422 fhPMDCP4->Fill(totPREClus,totCPVClus);
430 PostData(0, fOutputContainer);
433 //______________________________________________________________________________
434 void AliPMDQATask::Terminate(Option_t *)
436 // Processing when the event loop is ended
437 fOutputContainer = (TObjArray*)GetOutputData(0);
439 fhPMDP1 = (TH2F*)fOutputContainer->At(0);
440 fhPMDC2 = (TH1F*)fOutputContainer->At(1);
441 fhPMDP2 = (TH1F*)fOutputContainer->At(2);
442 fhPMDC3 = (TH1F*)fOutputContainer->At(3);
443 fhPMDP3 = (TH1F*)fOutputContainer->At(4);
444 fhPMDP4 = (TH1F*)fOutputContainer->At(5);
445 fhPMDC5 = (TH1F*)fOutputContainer->At(6);
446 fhPMDP5 = (TH1F*)fOutputContainer->At(7);
447 fhPMDCP0 = (TH2F*)fOutputContainer->At(8);
448 fhPMDCP1 = (TH2F*)fOutputContainer->At(9);
449 fhPMDCP2 = (TH2F*)fOutputContainer->At(10);
450 fhPMDCP3 = (TH2F*)fOutputContainer->At(11);
451 fhPMDCP4 = (TH2F*)fOutputContainer->At(12);
453 fhPMDSM1 = (TH2F*)fOutputContainer->At(13);
454 fhPMDSM2 = (TH2F*)fOutputContainer->At(14);
455 fhPMDSM3 = (TH2F*)fOutputContainer->At(15);
456 fhPMDSM4 = (TH2F*)fOutputContainer->At(16);
457 fhPMDSM5 = (TH2F*)fOutputContainer->At(17);
458 fhPMDSM6 = (TH2F*)fOutputContainer->At(18);
459 fhPMDSM7 = (TH2F*)fOutputContainer->At(19);
460 fhPMDSM8 = (TH2F*)fOutputContainer->At(20);
461 fhPMDSM9 = (TH2F*)fOutputContainer->At(21);
462 fhPMDSM10 = (TH2F*)fOutputContainer->At(22);
463 fhPMDSM11 = (TH2F*)fOutputContainer->At(23);
464 fhPMDSM12 = (TH2F*)fOutputContainer->At(24);
465 fhPMDSM13 = (TH2F*)fOutputContainer->At(25);
466 fhPMDSM14 = (TH2F*)fOutputContainer->At(26);
467 fhPMDSM15 = (TH2F*)fOutputContainer->At(27);
468 fhPMDSM16 = (TH2F*)fOutputContainer->At(28);
469 fhPMDSM17 = (TH2F*)fOutputContainer->At(29);
470 fhPMDSM18 = (TH2F*)fOutputContainer->At(30);
471 fhPMDSM19 = (TH2F*)fOutputContainer->At(31);
472 fhPMDSM20 = (TH2F*)fOutputContainer->At(32);
473 fhPMDSM21 = (TH2F*)fOutputContainer->At(33);
474 fhPMDSM22 = (TH2F*)fOutputContainer->At(34);
475 fhPMDSM23 = (TH2F*)fOutputContainer->At(35);
476 fhPMDSM24 = (TH2F*)fOutputContainer->At(36);
477 fhPMDSM = (TH1F*)fOutputContainer->At(37);
479 Bool_t problem = kFALSE ;
480 AliInfo(Form(" *** %s Report:", GetName())) ;
482 gStyle->SetOptStat(110000);
483 gStyle->SetOptFit(1);
485 TCanvas *cPMD0 = new TCanvas("cPMD0","PMD ESD Test #1", 10,10, 600, 600);
486 cPMD0->Range(-100, -100,100 ,100 );
487 fhPMDSM1->SetMarkerColor(2);
489 fhPMDSM1->GetXaxis()->SetTitle("Cluster X");
490 fhPMDSM1->GetYaxis()->SetTitle("Cluster Y");
491 fhPMDSM2->SetMarkerColor(2);
492 fhPMDSM2->Draw("same");
493 fhPMDSM3->SetMarkerColor(2);
494 fhPMDSM3->Draw("same");
495 fhPMDSM4->SetMarkerColor(2);
496 fhPMDSM4->Draw("same");
497 fhPMDSM5->SetMarkerColor(2);
498 fhPMDSM5->Draw("same");
499 fhPMDSM6->SetMarkerColor(2);
500 fhPMDSM6->Draw("same");
501 fhPMDSM7->SetMarkerColor(4);
502 fhPMDSM7->Draw("same");
503 fhPMDSM8->SetMarkerColor(4);
504 fhPMDSM8->Draw("same");
505 fhPMDSM9->SetMarkerColor(4);
506 fhPMDSM9->Draw("same");
507 fhPMDSM10->SetMarkerColor(4);
508 fhPMDSM10->Draw("same");
509 fhPMDSM11->SetMarkerColor(4);
510 fhPMDSM11->Draw("same");
511 fhPMDSM12->SetMarkerColor(4);
512 fhPMDSM12->Draw("same");
513 fhPMDSM13->SetMarkerColor(6);
514 fhPMDSM13->Draw("same");
515 fhPMDSM14->SetMarkerColor(6);
516 fhPMDSM14->Draw("same");
517 fhPMDSM15->SetMarkerColor(6);
518 fhPMDSM15->Draw("same");
519 fhPMDSM16->SetMarkerColor(6);
520 fhPMDSM16->Draw("same");
521 fhPMDSM17->SetMarkerColor(6);
522 fhPMDSM17->Draw("same");
523 fhPMDSM18->SetMarkerColor(6);
524 fhPMDSM18->Draw("same");
525 fhPMDSM19->SetMarkerColor(8);
526 fhPMDSM19->Draw("same");
527 fhPMDSM20->SetMarkerColor(8);
528 fhPMDSM20->Draw("same");
529 fhPMDSM21->SetMarkerColor(8);
530 fhPMDSM21->Draw("same");
531 fhPMDSM22->SetMarkerColor(8);
532 fhPMDSM22->Draw("same");
533 fhPMDSM23->SetMarkerColor(8);
534 fhPMDSM23->Draw("same");
535 fhPMDSM24->SetMarkerColor(8);
536 fhPMDSM24->Draw("same");
539 DrawPMDBoundarySM1();
540 DrawPMDBoundarySM2();
541 DrawPMDBoundarySM3();
542 DrawPMDBoundarySM4();
543 cPMD0->Print("ClusterXY.eps");
545 TCanvas *cPMD1 = new TCanvas("cPMD1"," PMD ESD Test #2",10, 10, 600,600);
548 cPMD1->SetFillColor(0);
549 fhPMDC2->SetLineColor(4);
552 fhPMDP2->SetLineColor(2);
554 cPMD1->Print("CPVPREphi.eps");
556 TCanvas *cPMD2 = new TCanvas("cPMD2","PMD ESD test #3",10, 10, 600, 600);
558 fhPMDSM->SetFillColor(2);
560 cPMD2->Print("AllSMN.eps");
562 TCanvas *cPMD3 = new TCanvas("cPMD3", "PMD ESD test #4",10, 10, 600, 600);
565 fhPMDCP0->SetMarkerColor(9);
568 fhPMDCP1->SetMarkerColor(6);
571 fhPMDP3->SetLineColor(2);
574 fhPMDCP4->SetMarkerColor(3);
576 cPMD3->Print("CPVPREClus.eps");
578 TCanvas *cPMD4 = new TCanvas("cPMD4","PMD ESD test #5", 10, 10, 600, 600);
581 fhPMDC3->SetLineColor(4);
584 fhPMDP4->SetLineColor(2);
586 cPMD4->Print("CPVPREAdc.eps");
589 sprintf(line, ".!tar -zcf %s.tar.gz *.eps", GetName()) ;
590 gROOT->ProcessLine(line);
592 AliInfo(Form("!!! All the eps files are in %s.tar.gz !!!", GetName())) ;
596 report="Problems found, please check!!!";
600 AliInfo(Form("*** %s Summary Report: %s \n",GetName(), report)) ;
603 //______________________________________________________________________________
604 void AliPMDQATask::CalculateSMN( Float_t clsX, Float_t clsY, Int_t & smn)
606 Double_t xcon[96] = {75.133, 54.204, 53.254, 32.326, 31.376,10.447,
607 75.133, 54.204, 53.254, 32.326, 31.376,10.447,
608 75.133, 54.204, 53.254, 32.326, 31.376,10.447,
609 75.133, 54.204, 53.254, 32.326, 31.376,10.447,
610 -75.133, -54.204, -53.254, -32.326, -31.376,-10.447,
611 -75.133, -54.204, -53.254, -32.326, -31.376,-10.447,
612 -75.133, -54.204, -53.254, -32.326, -31.376,-10.447,
613 -75.133, -54.204, -53.254, -32.326, -31.376,-10.447,
614 9.167, -32.543, -33.493, -75.133,
615 9.167, -32.543, -33.493, -75.133,
616 9.167, -32.543, -33.493, -75.133,
617 9.167, -32.543, -33.493, -75.133,
618 9.167, -32.543, -33.493, -75.133,
619 9.167, -32.543, -33.493, -75.133,
620 -9.167, 32.543, 33.493, 75.133,
621 -9.167, 32.543, 33.493, 75.133,
622 -9.167, 32.543, 33.493, 75.133,
623 -9.167, 32.543, 33.493, 75.133,
624 -9.167, 32.543, 33.493, 75.133,
625 -9.167, 32.543, 33.493, 75.133};
627 Double_t ycon[96] = {86.475, 86.475, 86.475, 86.475, 86.475, 86.475,
628 38.225, 38.225, 38.225, 38.225, 38.225, 38.225,
629 37.325, 37.325, 37.325, 37.325, 37.325, 37.325,
630 -10.925, -10.925, -10.925, -10.925, -10.925, -10.925,
631 -86.475, -86.475, -86.475, -86.475, -86.475, -86.475,
632 -38.225, -38.225, -38.225, -38.225, -38.225, -38.225,
633 -37.325, -37.325, -37.325, -37.325, -37.325, -37.325,
634 10.925, 10.925, 10.925, 10.925, 10.925, 10.925,
635 86.475, 86.475, 86.475, 86.475,
636 62.225, 62.225, 62.225, 62.225,
637 61.325, 61.325, 61.325, 61.325,
638 37.075, 37.075, 37.075, 37.075,
639 36.175, 36.175, 36.175, 36.175,
640 11.925, 11.925, 11.925 , 11.925,
641 -86.475, -86.475, -86.475, -86.475,
642 -62.225, -62.225, -62.225, -62.225,
643 -61.325, -61.325, -61.325, -61.325,
644 -37.075, -37.075, -37.075, -37.075,
645 -36.175, -36.175, -36.175, -36.175,
646 -11.925, -11.925, -11.925 , -11.925 };
648 if((clsX <= xcon[0]) && (clsX >= xcon[1]) &&
649 (clsY <= ycon[0]) && (clsY >= ycon[6])) smn = 0 ;
651 else if((clsX <=xcon[2]) && (clsX >= xcon[3]) &&
652 (clsY <= ycon[1]) && (clsY >= ycon[7]))smn = 1 ;
654 else if((clsX <=xcon[4]) && (clsX >= xcon[5]) &&
655 (clsY <= ycon[3]) && (clsY >= ycon[8]))smn = 2 ;
657 else if((clsX <= xcon[0]) && (clsX >= xcon[1]) &&
658 (clsY <= ycon[12]) && (clsY >= ycon[18])) smn = 3 ;
660 else if((clsX <=xcon[2]) && (clsX >= xcon[3]) &&
661 (clsY <= ycon[12]) && (clsY >= ycon[18]))smn = 4 ;
663 else if((clsX <=xcon[4]) && (clsX >= xcon[5]) &&
664 (clsY <= ycon[12]) && (clsY >= ycon[18]))smn = 5 ;
665 //------------------------------------------------------------------
666 else if((clsX >= xcon[24]) && (clsX <= xcon[25]) &&
667 (clsY >= ycon[24]) && (clsY <= ycon[30])) smn = 6 ;
669 else if((clsX >=xcon[26]) && (clsX <= xcon[27]) &&
670 (clsY >= ycon[25]) && (clsY <= ycon[31]))smn = 7 ;
672 else if((clsX >=xcon[28]) && (clsX <= xcon[29]) &&
673 (clsY >= ycon[26]) && (clsY <= ycon[32]))smn = 8 ;
675 else if((clsX >= xcon[24]) && (clsX <= xcon[25]) &&
676 (clsY >= ycon[36]) && (clsY <= ycon[42])) smn = 9 ;
678 else if((clsX >=xcon[26]) && (clsX <= xcon[27]) &&
679 (clsY >= ycon[36]) && (clsY <= ycon[42]))smn = 10;
681 else if((clsX >=xcon[28]) && (clsX <= xcon[29]) &&
682 (clsY >= ycon[36]) && (clsY <= ycon[42]))smn = 11;
683 //------------------------------------------------------------------
684 else if((clsX <= xcon[48]) && (clsX >= xcon[49]) &&
685 (clsY <= ycon[48]) && (clsY >= ycon[52])) smn = 12 ;
687 else if((clsX <=xcon[50]) && (clsX >= xcon[51]) &&
688 (clsY <= ycon[48]) && (clsY >= ycon[52]))smn = 13 ;
690 else if((clsX <=xcon[48]) && (clsX >= xcon[49]) &&
691 (clsY <= ycon[56]) && (clsY >= ycon[60]))smn = 14 ;
693 else if((clsX <=xcon[50]) && (clsX >= xcon[51]) &&
694 (clsY <= ycon[56]) && (clsY >= ycon[60]))smn = 15 ;
696 else if((clsX <=xcon[48]) && (clsX >= xcon[49]) &&
697 (clsY <= ycon[64]) && (clsY >= ycon[68]))smn = 16 ;
699 else if((clsX <=xcon[50]) && (clsX >= xcon[51]) &&
700 (clsY <= ycon[64]) && (clsY >= ycon[68]))smn = 17 ;
701 //--------------------------------------------------------------
702 else if((clsX >= xcon[72]) && (clsX <= xcon[73]) &&
703 (clsY >= ycon[72]) && (clsY <= ycon[76])) smn = 18 ;
705 else if((clsX >=xcon[74]) && (clsX <= xcon[75]) &&
706 (clsY >= ycon[72]) && (clsY <= ycon[76]))smn = 19 ;
708 else if((clsX >=xcon[72]) && (clsX <= xcon[73]) &&
709 (clsY >= ycon[80]) && (clsY <= ycon[84]))smn = 20 ;
711 else if((clsX >=xcon[74]) && (clsX <= xcon[75]) &&
712 (clsY >= ycon[80]) && (clsY <= ycon[84]))smn = 21;
714 else if((clsX >= xcon[72]) && (clsX <= xcon[73]) &&
715 (clsY >= ycon[88]) && (clsY <= ycon[92])) smn = 22 ;
717 else if((clsX >=xcon[74]) && (clsX <= xcon[75]) &&
718 (clsY >= ycon[88]) && (clsY <= ycon[92]))smn = 23 ;
722 //______________________________________________________________________________
723 void AliPMDQATask::DrawPMDBoundary()
725 // Draw PMD boundaries
727 gStyle->SetLineWidth(2);
728 gStyle->SetLineColor(2);
730 l = new TLine(75.1333, 86.475, -75.1333, 86.475); l->Draw("same");
731 l = new TLine(-75.1333, 86.470,-75.1333, -86.475); l->Draw("same");
732 l = new TLine(-75.1333, -86.475,75.1333, -86.475); l->Draw("same");
733 l = new TLine(75.1333, -86.475,75.1333, 86.475); l->Draw("same");
736 //______________________________________________________________________________
737 void AliPMDQATask::DrawPMDBoundarySM1()
739 // Draw boundaries of Super Module 1
741 gStyle->SetLineWidth(1);
742 gStyle->SetLineColor(4);
744 l = new TLine(-75.1333, 86.475, -10.447, 86.475); l->Draw("same");
745 l = new TLine(-10.447, 86.475, -10.446, -10.925); l->Draw("same");
746 l = new TLine(-10.446, -10.925, -75.1333,-10.925); l->Draw("same");
747 l = new TLine(-75.1333,-10.925, -75.1333, 86.475); l->Draw("same");
750 //______________________________________________________________________________
751 void AliPMDQATask::DrawPMDBoundarySM2()
753 // Draw boundaries of Super Module 2
755 gStyle->SetLineWidth(1);
756 gStyle->SetLineColor(4);
758 l = new TLine(75.1333, -86.475, 10.446, -86.475); l->Draw("same");
759 l = new TLine(10.446, -86.475, 10.446, 10.925); l->Draw("same");
760 l = new TLine(10.446, 10.925, 75.1333, 10.925); l->Draw("same");
761 l = new TLine(75.1333, 10.925, 75.1333, -86.475); l->Draw("same");
765 //______________________________________________________________________________
766 void AliPMDQATask::DrawPMDBoundarySM3()
768 // Draw boundaries of Super Module 3
770 gStyle->SetLineWidth(1);
771 gStyle->SetLineColor(1);
773 l = new TLine( -9.167, 86.475, 75.1333, 86.475); l->Draw("same");
774 l = new TLine(75.1333,86.475, 75.1333, 11.925); l->Draw("same");
775 l = new TLine(75.1333,11.925, -9.167, 11.925); l->Draw("same");
776 l = new TLine( -9.167, 11.925, -9.167, 86.475); l->Draw("same");
779 //______________________________________________________________________________
780 void AliPMDQATask::DrawPMDBoundarySM4()
782 // Draw boundaries of Super Module 4
784 gStyle->SetLineWidth(1);
785 gStyle->SetLineColor(1);
787 l = new TLine(9.167, -86.475, -75.1333,-86.475); l->Draw("same");
788 l = new TLine(-75.1333,-86.475, -75.1333,-11.925); l->Draw("same");
789 l = new TLine(-75.1333,-11.925, 9.167, -11.925); l->Draw("same");
790 l = new TLine(9.167, -11.925, 9.167, -86.475); l->Draw("same");