1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
16 //_________________________________________________________________________
17 // Algorythm class to analyze PHOSv1 events:
18 // Construct histograms and displays them.
19 // IHEP CPV/PHOS reconstruction algorithm used.
20 // Use the macro EditorBar.C for best access to the functionnalities
22 //*-- Author: B. Polichtchouk (IHEP)
23 //////////////////////////////////////////////////////////////////////////////
25 // --- ROOT system ---
31 #include "TParticle.h"
32 #include "TClonesArray.h"
38 // --- Standard library ---
40 // --- AliRoot header files ---
42 #include "AliRunLoader.h"
43 #include "AliHeader.h"
45 // --- PHOS header files ---
47 #include "AliPHOSIhepAnalyze.h"
48 #include "AliPHOSDigit.h"
49 #include "AliPHOSRecParticle.h"
50 #include "AliPHOSLoader.h"
51 #include "AliPHOSHit.h"
52 #include "AliPHOSImpact.h"
53 #include "AliPHOSvImpacts.h"
54 #include "AliPHOSCpvRecPoint.h"
56 #include "AliPHOSGeometry.h"
57 #include "AliPHOSEvalRecPoint.h"
59 ClassImp(AliPHOSIhepAnalyze)
61 //____________________________________________________________________________
63 AliPHOSIhepAnalyze::AliPHOSIhepAnalyze()
68 //____________________________________________________________________________
70 AliPHOSIhepAnalyze::AliPHOSIhepAnalyze(Text_t * name) : fFileName(name) {
71 // Constructor: open a header file
72 fRunLoader = AliRunLoader::Open(fFileName);
73 if (fRunLoader == 0x0)
75 AliFatal(Form("Can not load event from file %s",name));
79 //____________________________________________________________________________
80 void AliPHOSIhepAnalyze::AnalyzeCPV1(Int_t Nevents)
83 // Analyzes CPV characteristics: resolutions, cluster multiplicity,
84 // cluster lengths along Z and Phi.
85 // Author: Yuri Kharlov
87 // Modified by Boris Polichtchouk, 3.07.2001
92 TH1F *hDx = new TH1F("hDx" ,"CPV x-resolution@reconstruction",100,-5. , 5.);
93 TH1F *hDz = new TH1F("hDz" ,"CPV z-resolution@reconstruction",100,-5. , 5.);
94 TH1F *hChi2 = new TH1F("hChi2" ,"Chi2/dof of one-gamma fit",30, 0. , 10.);
95 TH1S *hNrp = new TH1S("hNrp" ,"CPV rec.point multiplicity", 21,-0.5,20.5);
96 TH1S *hNrpX = new TH1S("hNrpX","CPV rec.point Phi-length" , 21,-0.5,20.5);
97 TH1S *hNrpZ = new TH1S("hNrpZ","CPV rec.point Z-length" , 21,-0.5,20.5);
99 TH1F *hEg = new TH1F("hEg" ,"Energies of impacts",30,0.,6.);
100 TH1F *hEr = new TH1F("hEr" ,"Amplitudes of rec. points",50,0.,20.);
102 TList * fCpvImpacts ;
103 TBranch * branchCPVimpacts;
107 AliPHOSLoader* please = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
110 AliError(Form("Could not obtain the Loader object !"));
114 const AliPHOSGeometry * fGeom = please->PHOSGeometry();
116 AliInfo(Form("Start CPV Analysis-1. Resolutions, cluster multiplicity and lengths")) ;
117 for ( Int_t ievent=0; ievent<Nevents; ievent++) {
120 Int_t nChargedGen = 0;
122 Int_t ntracks = gAlice->GetEvent(ievent);
123 AliInfo(Form(">>>>>>>Event %d .<<<<<<<", ievent)) ;
125 /******************************************************************/
126 TTree* treeH = please->TreeH();
129 AliError(Form("Can not get TreeH"));
132 /******************************************************************/
134 // Get branch of CPV impacts
135 if (! (branchCPVimpacts =treeH->GetBranch("PHOSCpvImpacts")) ) {
136 AliWarning(Form("Couldn't find branch PHOSCpvImpacts. Exit.")) ;
140 // Create and fill arrays of hits for each CPV module
142 Int_t nOfModules = fGeom->GetNModules();
143 TClonesArray **hitsPerModule = new TClonesArray *[nOfModules];
145 for (iModule=0; iModule < nOfModules; iModule++)
146 hitsPerModule[iModule] = new TClonesArray("AliPHOSImpact",100);
148 TClonesArray *impacts;
149 AliPHOSImpact *impact;
152 // First go through all primary tracks and fill the arrays
153 // of hits per each CPV module
155 for (Int_t itrack=0; itrack<ntracks; itrack++) {
156 branchCPVimpacts ->SetAddress(&fCpvImpacts);
157 branchCPVimpacts ->GetEntry(itrack,0);
159 for (Int_t iModule=0; iModule < nOfModules; iModule++) {
160 impacts = (TClonesArray *)fCpvImpacts->At(iModule);
161 // Do loop over impacts in the module
162 for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
163 impact=(AliPHOSImpact*)impacts->At(iImpact);
164 hEg->Fill(impact->GetMomentum().E());
165 TClonesArray &lhits = *(TClonesArray *)hitsPerModule[iModule];
166 if(IsCharged(impact->GetPid())) nChargedGen++;
167 new(lhits[hitsPerModule[iModule]->GetEntriesFast()]) AliPHOSImpact(*impact);
170 fCpvImpacts->Clear();
172 for (iModule=0; iModule < nOfModules; iModule++) {
173 Int_t nsum = hitsPerModule[iModule]->GetEntriesFast();
174 printf("CPV module %d has %d impacts\n",iModule,nsum);
178 // Then go through reconstructed points and for each find
180 // The distance from the rec.point to the closest hit
181 // gives the coordinate resolution of the CPV
183 fRunLoader->GetEvent(ievent);
184 TIter nextRP(please->CpvRecPoints()) ;
185 AliPHOSCpvRecPoint *cpvRecPoint ;
186 Float_t xgen, ygen, zgen;
187 while( ( cpvRecPoint = (AliPHOSCpvRecPoint *)nextRP() ) ) {
189 Float_t chi2dof = ((AliPHOSEvalRecPoint*)cpvRecPoint)->Chi2Dof();
190 hChi2->Fill(chi2dof);
191 hEr->Fill(cpvRecPoint->GetEnergy());
194 cpvRecPoint->GetLocalPosition(locpos);
195 Int_t phosModule = cpvRecPoint->GetPHOSMod();
196 Int_t rpMult = cpvRecPoint->GetMultiplicity();
197 Int_t rpMultX, rpMultZ;
198 cpvRecPoint->GetClusterLengths(rpMultX,rpMultZ);
199 Float_t xrec = locpos.X();
200 Float_t zrec = locpos.Z();
201 Float_t dxmin = 1.e+10;
202 Float_t dzmin = 1.e+10;
203 Float_t r2min = 1.e+10;
206 Int_t nCPVhits = (hitsPerModule[phosModule-1])->GetEntriesFast();
207 Float_t locImpX=1e10,locImpZ=1e10; // local coords of closest impact
208 Float_t gImpX=1e10, gImpZ=1e10,gImpY=1e10; // global coords of closest impact
209 for (Int_t ihit=0; ihit<nCPVhits; ihit++) {
210 impact = (AliPHOSImpact*)(hitsPerModule[phosModule-1])->UncheckedAt(ihit);
215 //Transform to the local ref.frame
216 const AliPHOSGeometry* geom = please->PHOSGeometry();
217 Float_t phig = geom->GetPHOSAngle(phosModule);
218 Float_t phi = TMath::Pi()/180*phig;
219 Float_t distanceIPtoCPV = geom->GetIPtoOuterCoverDistance() -
220 (geom->GetFTPosition(1)+
221 geom->GetFTPosition(2)+
222 geom->GetCPVTextoliteThickness()
224 Float_t xoL,yoL,zoL ;
225 // xoL = xgen*TMath::Cos(phig)+ygen*TMath::Sin(phig) ;
226 // yoL = -xgen*TMath::Sin(phig)+ygen*TMath::Cos(phig) + distanceIPtoCPV;
227 xoL = xgen*TMath::Cos(phi)-ygen*TMath::Sin(phi) ;
228 yoL = xgen*TMath::Sin(phi)+ygen*TMath::Cos(phi) + distanceIPtoCPV;
231 r2 = TMath::Power((xoL-xrec),2) + TMath::Power((zoL-zrec),2);
243 AliInfo(Form("Impact global (X,Z,Y) = %f %f %f", gImpX, gImpZ, gImpY));
244 AliInfo(Form("Impact local (X,Z) = %f %f", locImpX, locImpZ));
245 AliInfo(Form("Reconstructed (X,Z) = %f %f", xrec, zrec));
246 AliInfo(Form("dxmin %f dzmin %f", dxmin, dzmin));
249 // hDr ->Fill(TMath::Sqrt(r2min));
251 hNrpX->Fill(rpMultX);
252 hNrpZ->Fill(rpMultZ);
254 delete [] hitsPerModule;
256 AliInfo(Form("++++ Event %d : total %d impacts, %d charged impacts and %d rec. points.",
257 ievent, nTotalGen, nChargedGen, please->CpvRecPoints()->GetEntries())) ;
261 Text_t outputname[80] ;
262 sprintf(outputname,"%s.analyzed",GetFileName().Data());
263 TFile output(outputname,"RECREATE");
268 TCanvas *cpvCanvas = new TCanvas("Cpv1","CPV analysis-I",20,20,800,600);
269 gStyle->SetOptStat(111111);
270 gStyle->SetOptFit(1);
271 gStyle->SetOptDate(1);
272 cpvCanvas->Divide(3,3);
275 gPad->SetFillColor(10);
276 hNrp->SetFillColor(16);
280 gPad->SetFillColor(10);
281 hNrpX->SetFillColor(16);
285 gPad->SetFillColor(10);
286 hNrpZ->SetFillColor(16);
290 gPad->SetFillColor(10);
291 hDx->SetFillColor(16);
296 gPad->SetFillColor(10);
297 hDz->SetFillColor(16);
302 gPad->SetFillColor(10);
303 hChi2->SetFillColor(16);
307 gPad->SetFillColor(10);
308 hEg->SetFillColor(16);
312 gPad->SetFillColor(10);
313 hEr->SetFillColor(16);
316 cpvCanvas->Write(0,kOverwrite);
321 void AliPHOSIhepAnalyze::AnalyzeEMC1(Int_t Nevents)
324 // Analyzes Emc characteristics: resolutions, cluster multiplicity,
325 // cluster lengths along Z and Phi.
326 // Author: Boris Polichtchouk, 24.08.2001
331 TH1F *hDx = new TH1F("hDx" ,"EMC x-resolution@reconstruction",100,-5. , 5.);
332 TH1F *hDz = new TH1F("hDz" ,"EMC z-resolution@reconstruction",100,-5. , 5.);
333 TH1F *hChi2 = new TH1F("hChi2" ,"Chi2/dof of one-gamma fit",30, 0. , 3.);
334 TH1S *hNrp = new TH1S("hNrp" ,"EMC rec.point multiplicity", 21,-0.5,20.5);
335 TH1S *hNrpX = new TH1S("hNrpX","EMC rec.point Phi-length" , 21,-0.5,20.5);
336 TH1S *hNrpZ = new TH1S("hNrpZ","EMC rec.point Z-length" , 21,-0.5,20.5);
338 TList * fEmcImpacts ;
339 TBranch * branchEMCimpacts;
341 AliPHOSLoader* please = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
344 AliError(Form("Could not obtain the Loader object !"));
348 const AliPHOSGeometry * fGeom = please->PHOSGeometry();
350 AliInfo(Form("Start EMC Analysis-1. Resolutions, cluster multiplicity and lengths"));
351 for ( Int_t ievent=0; ievent<Nevents; ievent++) {
355 Int_t ntracks = gAlice->GetEvent(ievent);
357 AliInfo(Form(" >>>>>>>Event %d .<<<<<<<", ievent)) ;
359 TTree* treeH = please->TreeH();
362 AliError(Form("Can not get TreeH"));
367 // Get branch of EMC impacts
368 if (! (branchEMCimpacts =treeH->GetBranch("PHOSEmcImpacts")) ) {
369 AliWarning(Form(" Couldn't find branch PHOSEmcImpacts. Exit."));
373 // Create and fill arrays of hits for each EMC module
375 Int_t nOfModules = fGeom->GetNModules();
376 TClonesArray **hitsPerModule = new TClonesArray *[nOfModules];
378 for (iModule=0; iModule < nOfModules; iModule++)
379 hitsPerModule[iModule] = new TClonesArray("AliPHOSImpact",100);
381 TClonesArray *impacts;
382 AliPHOSImpact *impact;
385 // First go through all primary tracks and fill the arrays
386 // of hits per each EMC module
388 for (Int_t itrack=0; itrack<ntracks; itrack++) {
389 branchEMCimpacts ->SetAddress(&fEmcImpacts);
390 branchEMCimpacts ->GetEntry(itrack,0);
392 for (Int_t iModule=0; iModule < nOfModules; iModule++) {
393 impacts = (TClonesArray *)fEmcImpacts->At(iModule);
394 // Do loop over impacts in the module
395 for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
396 impact=(AliPHOSImpact*)impacts->At(iImpact);
397 TClonesArray &lhits = *(TClonesArray *)hitsPerModule[iModule];
398 new(lhits[hitsPerModule[iModule]->GetEntriesFast()]) AliPHOSImpact(*impact);
401 fEmcImpacts->Clear();
403 for (iModule=0; iModule < nOfModules; iModule++) {
404 Int_t nsum = hitsPerModule[iModule]->GetEntriesFast();
405 printf("EMC module %d has %d hits\n",iModule,nsum);
409 // Then go through reconstructed points and for each find
411 // The distance from the rec.point to the closest hit
412 // gives the coordinate resolution of the EMC
414 fRunLoader->GetEvent(ievent);
415 TIter nextRP(please->EmcRecPoints()) ;
416 AliPHOSEmcRecPoint *emcRecPoint ;
417 Float_t xgen, ygen, zgen;
418 while( ( emcRecPoint = (AliPHOSEmcRecPoint *)nextRP() ) ) {
420 Float_t chi2dof = ((AliPHOSEvalRecPoint*)emcRecPoint)->Chi2Dof();
421 hChi2->Fill(chi2dof);
424 emcRecPoint->GetLocalPosition(locpos);
425 Int_t phosModule = emcRecPoint->GetPHOSMod();
426 Int_t rpMult = emcRecPoint->GetMultiplicity();
427 Int_t rpMultX, rpMultZ;
428 ((AliPHOSEvalRecPoint*)emcRecPoint)->GetClusterLengths(rpMultX,rpMultZ);
429 Float_t xrec = locpos.X();
430 Float_t zrec = locpos.Z();
431 Float_t dxmin = 1.e+10;
432 Float_t dzmin = 1.e+10;
433 Float_t r2min = 1.e+10;
436 Int_t nEMChits = (hitsPerModule[phosModule-1])->GetEntriesFast();
437 Float_t locImpX=1e10,locImpZ=1e10; // local coords of closest impact
438 Float_t gImpX=1e10, gImpZ=1e10,gImpY=1e10; // global coords of closest impact
439 for (Int_t ihit=0; ihit<nEMChits; ihit++) {
440 impact = (AliPHOSImpact*)(hitsPerModule[phosModule-1])->UncheckedAt(ihit);
446 //Transform to the local ref.frame
447 const AliPHOSGeometry* geom = please->PHOSGeometry();
448 Float_t phig = geom->GetPHOSAngle(phosModule);
449 Float_t phi = TMath::Pi()/180*phig;
450 Float_t distanceIPtoEMC = geom->GetIPtoCrystalSurface();
451 Float_t xoL,yoL,zoL ;
452 // xoL = xgen*TMath::Cos(phig)+ygen*TMath::Sin(phig) ;
453 // yoL = -xgen*TMath::Sin(phig)+ygen*TMath::Cos(phig) + distanceIPtoEMC;
454 xoL = xgen*TMath::Cos(phi)-ygen*TMath::Sin(phi) ;
455 yoL = xgen*TMath::Sin(phi)+ygen*TMath::Cos(phi) + distanceIPtoEMC;
458 r2 = TMath::Power((xoL-xrec),2) + TMath::Power((zoL-zrec),2);
470 AliInfo(Form(" Impact global (X,Z,Y) = %f %f %f", gImpX, gImpZ, gImpY));
471 AliInfo(Form(" Impact local (X,Z) = %f %f", locImpX, locImpZ));
472 AliInfo(Form(" Reconstructed (X,Z) = %f %f", xrec, zrec));
473 AliInfo(Form(" dxmin %f dzmin %f", dxmin, dzmin)) ;
476 // hDr ->Fill(TMath::Sqrt(r2min));
478 hNrpX->Fill(rpMultX);
479 hNrpZ->Fill(rpMultZ);
481 delete [] hitsPerModule;
483 AliInfo(Form("++++ Event %d : total %d impacts, %d Emc rec. points.",
484 ievent, nTotalGen, please->EmcRecPoints()->GetEntriesFast())) ;
489 Text_t outputname[80] ;
490 sprintf(outputname,"%s.analyzed",GetFileName().Data());
491 TFile output(outputname,"update");
496 TCanvas *emcCanvas = new TCanvas("Emc1","EMC analysis-I",20,20,800,400);
497 gStyle->SetOptStat(111111);
498 gStyle->SetOptFit(1);
499 gStyle->SetOptDate(1);
500 emcCanvas->Divide(3,2);
503 gPad->SetFillColor(10);
504 hNrp->SetFillColor(16);
508 gPad->SetFillColor(10);
509 hNrpX->SetFillColor(16);
513 gPad->SetFillColor(10);
514 hNrpZ->SetFillColor(16);
518 gPad->SetFillColor(10);
519 hDx->SetFillColor(16);
524 gPad->SetFillColor(10);
525 hDz->SetFillColor(16);
530 gPad->SetFillColor(10);
531 hChi2->SetFillColor(16);
534 emcCanvas->Write(0,kOverwrite);
537 //____________________________________________________________________________
538 void AliPHOSIhepAnalyze::AnalyzeCPV2(Int_t Nevents)
540 // CPV analysis - part II.
541 // Ratio of combinatoric distances between generated
542 // and reconstructed hits.
543 // Author: Boris Polichtchouk (polishchuk@mx.ihep.su)
547 TH1F* hDrijCPVr = new TH1F("Drij_cpv_r","Distance between reconstructed hits in CPV",140,0,50);
548 TH1F* hDrijCPVg = new TH1F("Drij_cpv_g","Distance between generated hits in CPV",140,0,50);
549 TH1F* hDrijCPVratio = new TH1F("Drij_cpv_ratio","R_{ij}^{rec}/R_{ij}^{gen} in CPV",140,0,50);
551 // TH1F* hT0 = new TH1F("hT0","Type of entering particle",20000,-10000,10000);
555 hDrijCPVratio->Sumw2(); //correct treatment of errors
557 TList * fCpvImpacts = new TList();
558 TBranch * branchCPVimpacts;
560 AliPHOSLoader* please = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
563 AliError(Form("Could not obtain the Loader object !"));
566 const AliPHOSGeometry * fGeom = please->PHOSGeometry();
567 fRunLoader->LoadHeader();
569 for (Int_t nev=0; nev<Nevents; nev++)
571 printf("\n=================== Event %10d ===================\n",nev);
572 fRunLoader->GetEvent(nev);
573 Int_t ntracks = fRunLoader->GetHeader()->GetNtrack();
575 Int_t nRecCPV = 0; // Reconstructed points in event
576 Int_t nGenCPV = 0; // Impacts in event
578 // Get branch of CPV impacts
579 TTree* treeH = please->TreeH();
582 AliError(Form("Can not get TreeH"));
586 if (! (branchCPVimpacts =treeH->GetBranch("PHOSCpvImpacts")) ) return;
588 // Create and fill arrays of hits for each CPV module
589 Int_t nOfModules = fGeom->GetNModules();
590 TClonesArray **hitsPerModule = new TClonesArray *[nOfModules];
592 for (iModule=0; iModule < nOfModules; iModule++)
593 hitsPerModule[iModule] = new TClonesArray("AliPHOSImpact",100);
595 TClonesArray *impacts;
596 AliPHOSImpact *impact;
598 for (Int_t itrack=0; itrack<ntracks; itrack++) {
599 branchCPVimpacts ->SetAddress(&fCpvImpacts);
600 AliInfo(Form(" branchCPVimpacts ->SetAddress(&fCpvImpacts) OK."));
601 branchCPVimpacts ->GetEntry(itrack,0);
603 for (Int_t iModule=0; iModule < nOfModules; iModule++) {
604 impacts = (TClonesArray *)fCpvImpacts->At(iModule);
605 // Do loop over impacts in the module
606 for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
607 impact=(AliPHOSImpact*)impacts->At(iImpact);
608 TClonesArray &lhits = *(TClonesArray *)hitsPerModule[iModule];
609 if(IsCharged(impact->GetPid()))
610 new(lhits[hitsPerModule[iModule]->GetEntriesFast()]) AliPHOSImpact(*impact);
613 fCpvImpacts->Clear();
616 for (iModule=0; iModule < nOfModules; iModule++) {
617 Int_t nsum = hitsPerModule[iModule]->GetEntriesFast();
618 printf("CPV module %d has %d hits\n",iModule,nsum);
620 AliPHOSImpact* genHit1;
621 AliPHOSImpact* genHit2;
623 for(irp1=0; irp1< nsum; irp1++)
625 genHit1 = (AliPHOSImpact*)((hitsPerModule[iModule])->At(irp1));
626 for(irp2 = irp1+1; irp2<nsum; irp2++)
628 genHit2 = (AliPHOSImpact*)((hitsPerModule[iModule])->At(irp2));
629 Float_t dx = genHit1->X() - genHit2->X();
630 Float_t dz = genHit1->Z() - genHit2->Z();
631 Float_t dr = TMath::Sqrt(dx*dx + dz*dz);
633 // AliInfo(Form("(dx dz dr): %f %f", dx, dz));
639 //--------- Combinatoric distance between rec. hits in CPV
641 TObjArray* cpvRecPoints = please->CpvRecPoints();
642 nRecCPV = cpvRecPoints->GetEntriesFast();
646 AliPHOSCpvRecPoint* recHit1;
647 AliPHOSCpvRecPoint* recHit2;
648 TIter nextCPVrec1(cpvRecPoints);
649 while(TObject* obj1 = nextCPVrec1() )
651 TIter nextCPVrec2(cpvRecPoints);
652 while (TObject* obj2 = nextCPVrec2())
654 if(!obj2->IsEqual(obj1))
656 recHit1 = (AliPHOSCpvRecPoint*)obj1;
657 recHit2 = (AliPHOSCpvRecPoint*)obj2;
660 recHit1->GetLocalPosition(locpos1);
661 recHit2->GetLocalPosition(locpos2);
662 Float_t dx = locpos1.X() - locpos2.X();
663 Float_t dz = locpos1.Z() - locpos2.Z();
664 Float_t dr = TMath::Sqrt(dx*dx + dz*dz);
665 if(recHit1->GetPHOSMod() == recHit2->GetPHOSMod())
672 AliInfo(Form(" Event %d . Total of %d hits, %d rec.points.",
673 nev, nGenCPV, nRecCPV)) ;
675 delete [] hitsPerModule;
677 } // End of loop over events.
680 // hDrijCPVg->Draw();
681 // hDrijCPVr->Draw();
682 hDrijCPVratio->Divide(hDrijCPVr,hDrijCPVg);
683 hDrijCPVratio->Draw();
690 void AliPHOSIhepAnalyze::CpvSingle(Int_t nevents)
692 // Distributions of coordinates and cluster lengths of rec. points
693 // in the case of single initial particle.
695 TH1F* hZr = new TH1F("Zrec","Reconstructed Z distribution",150,-5,5);
696 TH1F* hXr = new TH1F("Xrec","Reconstructed X distribution",150,-14,-2);
697 TH1F *hChi2 = new TH1F("hChi2" ,"Chi2/dof of one-gamma fit",100, 0. , 20.);
699 TH1S *hNrp = new TH1S("hNrp" ,"CPV rec.point multiplicity",21,-0.5,20.5);
700 TH1S *hNrpX = new TH1S("hNrpX","CPV rec.point Phi-length" ,21,-0.5,20.5);
701 TH1S *hNrpZ = new TH1S("hNrpZ","CPV rec.point Z-length" ,21,-0.5,20.5);
703 AliPHOSLoader* gime = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
706 AliError(Form("Could not obtain the Loader object !"));
710 for(Int_t ievent=0; ievent<nevents; ievent++)
712 fRunLoader->GetEvent(ievent);
713 if(gime->CpvRecPoints()->GetEntriesFast()>1) continue;
715 AliPHOSCpvRecPoint* pt = (AliPHOSCpvRecPoint*)(gime->CpvRecPoints())->At(0);
718 pt->GetLocalPosition(lpos);
722 Int_t rpMult = pt->GetMultiplicity();
724 Int_t rpMultX, rpMultZ;
725 pt->GetClusterLengths(rpMultX,rpMultZ);
726 hNrpX->Fill(rpMultX);
727 hNrpZ->Fill(rpMultZ);
728 hChi2->Fill(((AliPHOSEvalRecPoint*)pt)->Chi2Dof());
729 AliInfo(Form("+++++ Event %d . (Mult,MultX,MultZ) = %d %d %d +++++",
730 ievent, rpMult, rpMultX, rpMultZ)) ;
736 Text_t outputname[80] ;
737 sprintf(outputname,"%s.analyzed.single",GetFileName().Data());
738 TFile output(outputname,"RECREATE");
742 TCanvas *cpvCanvas = new TCanvas("SingleParticle","Single particle events",20,20,800,400);
743 gStyle->SetOptStat(111111);
744 gStyle->SetOptFit(1);
745 gStyle->SetOptDate(1);
746 cpvCanvas->Divide(3,2);
749 gPad->SetFillColor(10);
750 hXr->SetFillColor(16);
754 gPad->SetFillColor(10);
755 hZr->SetFillColor(16);
759 gPad->SetFillColor(10);
760 hChi2->SetFillColor(16);
764 gPad->SetFillColor(10);
765 hNrp->SetFillColor(16);
769 gPad->SetFillColor(10);
770 hNrpX->SetFillColor(16);
774 gPad->SetFillColor(10);
775 hNrpZ->SetFillColor(16);
778 cpvCanvas->Write(0,kOverwrite);
782 void AliPHOSIhepAnalyze::HitsCPV(Int_t nev)
784 // Cumulative list of charged CPV impacts in event nev.
786 TList * fCpvImpacts ;
787 TBranch * branchCPVimpacts;
789 AliPHOSLoader* please = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
792 AliError(Form("Could not obtain the Loader object !"));
795 const AliPHOSGeometry * fGeom = please->PHOSGeometry();
798 printf("\n=================== Event %10d ===================\n",nev);
799 fRunLoader->GetEvent(nev);
800 Int_t ntracks = fRunLoader->GetHeader()->GetNtrack();
802 // Int_t nRecCPV = 0; // Reconstructed points in event // 01.10.2001
803 // Int_t nGenCPV = 0; // Impacts in event
805 // Get branch of CPV impacts
806 TTree* treeH = please->TreeH();
809 AliError(Form("Can not get TreeH"));
813 if (! (branchCPVimpacts =treeH->GetBranch("PHOSCpvImpacts")) ) return;
815 // Create and fill arrays of hits for each CPV module
816 Int_t nOfModules = fGeom->GetNModules();
817 TClonesArray **hitsPerModule = new TClonesArray *[nOfModules];
819 for (iModule=0; iModule < nOfModules; iModule++)
820 hitsPerModule[iModule] = new TClonesArray("AliPHOSImpact",100);
822 TClonesArray *impacts;
823 AliPHOSImpact *impact;
825 for (Int_t itrack=0; itrack<ntracks; itrack++) {
826 branchCPVimpacts ->SetAddress(&fCpvImpacts);
827 AliInfo(Form(" branchCPVimpacts ->SetAddress(&fCpvImpacts) OK."));
828 branchCPVimpacts ->GetEntry(itrack,0);
830 for (Int_t iModule=0; iModule < nOfModules; iModule++) {
831 impacts = (TClonesArray *)fCpvImpacts->At(iModule);
832 // Do loop over impacts in the module
833 for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
834 impact=(AliPHOSImpact*)impacts->At(iImpact);
835 TClonesArray &lhits = *(TClonesArray *)hitsPerModule[iModule];
836 new(lhits[hitsPerModule[iModule]->GetEntriesFast()]) AliPHOSImpact(*impact);
839 fCpvImpacts->Clear();
842 for (iModule=0; iModule < nOfModules; iModule++) {
843 Int_t nsum = hitsPerModule[iModule]->GetEntriesFast();
844 printf("CPV module %d has %d hits\n",iModule,nsum);
847 // TList * fCpvImpacts ;
848 // TBranch * branchCPVimpacts;
849 // AliPHOSImpact* impact;
850 // TClonesArray* impacts;
852 // AliPHOSLoader * please = AliPHOSLoader::GetInstance(GetFileName().Data(),"PHOS");
853 // const AliPHOSGeometry * fGeom = please->PHOSGeometry();
855 // Int_t ntracks = gAlice->GetEvent(ievent);
856 // Int_t nOfModules = fGeom->GetNModules();
857 // AliInfo(Form(" Tracks: "<<ntracks<<" Modules: "<<nOfModules));
859 // if (! (branchCPVimpacts =gAlice->TreeH()->GetBranch("PHOSCpvImpacts")) ) return;
861 // for (Int_t itrack=0; itrack<ntracks; itrack++) {
862 // branchCPVimpacts ->SetAddress(&fCpvImpacts);
863 // Info("AnalyzeCPV1", " branchCPVimpacts ->SetAddress(&fCpvImpacts) OK.");
864 // branchCPVimpacts ->GetEntry(itrack,0);
865 // Info(Form(" branchCPVimpacts ->GetEntry(itrack,0) OK."));
867 // for (Int_t iModule=0; iModule < nOfModules; iModule++) {
868 // impacts = (TClonesArray *)fCpvImpacts->At(iModule);
869 // Info(Form(" fCpvImpacts->At(iModule) OK."));
870 // // Do loop over impacts in the module
871 // for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
872 // impact=(AliPHOSImpact*)impacts->At(iImpact);
874 // if(IsCharged(impact->GetPid()))
876 // Info(Form(" Add charged hit.."));
877 // new(hits[hits.GetEntriesFast()]) AliPHOSImpact(*impact);
878 // Info(Form("done."));
882 // fCpvImpacts->Clear();
885 // Info(Form(" PHOS event "<<ievent<<": "<<hits.GetEntries()<<" charged CPV hits."));
890 // void AliPHOSIhepAnalyze::ChargedHitsCPV(TClonesArray* hits, Int_t ievent, Int_t iModule)
892 // // Cumulative list of charged CPV hits in event ievent
893 // // in PHOS module iModule.
895 // HitsCPV(hits,ievent,iModule);
898 // while(AliPHOSCPVHit* cpvhit = (AliPHOSCPVHit*)next())
900 // if(!IsCharged(cpvhit->GetIpart()))
902 // hits->Remove(cpvhit);
908 // Info(Form(" PHOS module "<<iModule<<": "<<hits->GetEntries()<<" charged CPV hits."));
911 Bool_t AliPHOSIhepAnalyze::IsCharged(Int_t pdgCode)
914 AliInfo(Form("pdgCode %d", pdgCode));
915 if(pdgCode==211 || pdgCode==-211 || pdgCode==321 || pdgCode==-321 || pdgCode==11 || pdgCode==-11 || pdgCode==2212 || pdgCode==-2212) return kTRUE;
919 //---------------------------------------------------------------------------