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 //____________________________________________________________________________
62 AliPHOSIhepAnalyze::AliPHOSIhepAnalyze():
68 //____________________________________________________________________________
69 AliPHOSIhepAnalyze::AliPHOSIhepAnalyze(Text_t * name) :
73 // Constructor: open a header file
74 fRunLoader = AliRunLoader::Open(fFileName);
75 if (fRunLoader == 0x0)
77 AliFatal(Form("Can not load event from file %s",name));
81 //____________________________________________________________________________
82 void AliPHOSIhepAnalyze::AnalyzeCPV1(Int_t Nevents)
85 // Analyzes CPV characteristics: resolutions, cluster multiplicity,
86 // cluster lengths along Z and Phi.
87 // Author: Yuri Kharlov
89 // Modified by Boris Polichtchouk, 3.07.2001
94 TH1F *hDx = new TH1F("hDx" ,"CPV x-resolution@reconstruction",100,-5. , 5.);
95 TH1F *hDz = new TH1F("hDz" ,"CPV z-resolution@reconstruction",100,-5. , 5.);
96 TH1F *hChi2 = new TH1F("hChi2" ,"Chi2/dof of one-gamma fit",30, 0. , 10.);
97 TH1S *hNrp = new TH1S("hNrp" ,"CPV rec.point multiplicity", 21,-0.5,20.5);
98 TH1S *hNrpX = new TH1S("hNrpX","CPV rec.point Phi-length" , 21,-0.5,20.5);
99 TH1S *hNrpZ = new TH1S("hNrpZ","CPV rec.point Z-length" , 21,-0.5,20.5);
101 TH1F *hEg = new TH1F("hEg" ,"Energies of impacts",30,0.,6.);
102 TH1F *hEr = new TH1F("hEr" ,"Amplitudes of rec. points",50,0.,20.);
104 TList * fCpvImpacts ;
105 TBranch * branchCPVimpacts;
109 AliPHOSLoader* please = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
112 AliError(Form("Could not obtain the Loader object !"));
116 const AliPHOSGeometry * fGeom = please->PHOSGeometry();
118 AliInfo(Form("Start CPV Analysis-1. Resolutions, cluster multiplicity and lengths")) ;
119 for ( Int_t ievent=0; ievent<Nevents; ievent++) {
122 Int_t nChargedGen = 0;
124 Int_t ntracks = gAlice->GetEvent(ievent);
125 AliInfo(Form(">>>>>>>Event %d .<<<<<<<", ievent)) ;
127 /******************************************************************/
128 TTree* treeH = please->TreeH();
131 AliError(Form("Can not get TreeH"));
134 /******************************************************************/
136 // Get branch of CPV impacts
137 if (! (branchCPVimpacts =treeH->GetBranch("PHOSCpvImpacts")) ) {
138 AliWarning(Form("Couldn't find branch PHOSCpvImpacts. Exit.")) ;
142 // Create and fill arrays of hits for each CPV module
144 Int_t nOfModules = fGeom->GetNModules();
145 TClonesArray **hitsPerModule = new TClonesArray *[nOfModules];
147 for (iModule=0; iModule < nOfModules; iModule++)
148 hitsPerModule[iModule] = new TClonesArray("AliPHOSImpact",100);
150 TClonesArray *impacts;
151 AliPHOSImpact *impact;
154 // First go through all primary tracks and fill the arrays
155 // of hits per each CPV module
157 for (Int_t itrack=0; itrack<ntracks; itrack++) {
158 branchCPVimpacts ->SetAddress(&fCpvImpacts);
159 branchCPVimpacts ->GetEntry(itrack,0);
161 for (Int_t iModule=0; iModule < nOfModules; iModule++) {
162 impacts = (TClonesArray *)fCpvImpacts->At(iModule);
163 // Do loop over impacts in the module
164 for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
165 impact=(AliPHOSImpact*)impacts->At(iImpact);
166 hEg->Fill(impact->GetMomentum().E());
167 TClonesArray &lhits = *(TClonesArray *)hitsPerModule[iModule];
168 if(IsCharged(impact->GetPid())) nChargedGen++;
169 new(lhits[hitsPerModule[iModule]->GetEntriesFast()]) AliPHOSImpact(*impact);
172 fCpvImpacts->Clear();
174 for (iModule=0; iModule < nOfModules; iModule++) {
175 Int_t nsum = hitsPerModule[iModule]->GetEntriesFast();
176 printf("CPV module %d has %d impacts\n",iModule,nsum);
180 // Then go through reconstructed points and for each find
182 // The distance from the rec.point to the closest hit
183 // gives the coordinate resolution of the CPV
185 fRunLoader->GetEvent(ievent);
186 TIter nextRP(please->CpvRecPoints()) ;
187 AliPHOSCpvRecPoint *cpvRecPoint ;
188 Float_t xgen, ygen, zgen;
189 while( ( cpvRecPoint = (AliPHOSCpvRecPoint *)nextRP() ) ) {
191 Float_t chi2dof = ((AliPHOSEvalRecPoint*)cpvRecPoint)->Chi2Dof();
192 hChi2->Fill(chi2dof);
193 hEr->Fill(cpvRecPoint->GetEnergy());
196 cpvRecPoint->GetLocalPosition(locpos);
197 Int_t phosModule = cpvRecPoint->GetPHOSMod();
198 Int_t rpMult = cpvRecPoint->GetMultiplicity();
199 Int_t rpMultX, rpMultZ;
200 cpvRecPoint->GetClusterLengths(rpMultX,rpMultZ);
201 Float_t xrec = locpos.X();
202 Float_t zrec = locpos.Z();
203 Float_t dxmin = 1.e+10;
204 Float_t dzmin = 1.e+10;
205 Float_t r2min = 1.e+10;
208 Int_t nCPVhits = (hitsPerModule[phosModule-1])->GetEntriesFast();
209 Float_t locImpX=1e10,locImpZ=1e10; // local coords of closest impact
210 Float_t gImpX=1e10, gImpZ=1e10,gImpY=1e10; // global coords of closest impact
211 for (Int_t ihit=0; ihit<nCPVhits; ihit++) {
212 impact = (AliPHOSImpact*)(hitsPerModule[phosModule-1])->UncheckedAt(ihit);
217 //Transform to the local ref.frame
218 const AliPHOSGeometry* geom = please->PHOSGeometry();
219 Float_t phig = geom->GetPHOSAngle(phosModule);
220 Float_t phi = TMath::Pi()/180*phig;
221 Float_t distanceIPtoCPV = geom->GetIPtoOuterCoverDistance() -
222 (geom->GetFTPosition(1)+
223 geom->GetFTPosition(2)+
224 geom->GetCPVTextoliteThickness()
226 Float_t xoL,yoL,zoL ;
227 // xoL = xgen*TMath::Cos(phig)+ygen*TMath::Sin(phig) ;
228 // yoL = -xgen*TMath::Sin(phig)+ygen*TMath::Cos(phig) + distanceIPtoCPV;
229 xoL = xgen*TMath::Cos(phi)-ygen*TMath::Sin(phi) ;
230 yoL = xgen*TMath::Sin(phi)+ygen*TMath::Cos(phi) + distanceIPtoCPV;
233 r2 = TMath::Power((xoL-xrec),2) + TMath::Power((zoL-zrec),2);
245 AliInfo(Form("Impact global (X,Z,Y) = %f %f %f", gImpX, gImpZ, gImpY));
246 AliInfo(Form("Impact local (X,Z) = %f %f", locImpX, locImpZ));
247 AliInfo(Form("Reconstructed (X,Z) = %f %f", xrec, zrec));
248 AliInfo(Form("dxmin %f dzmin %f", dxmin, dzmin));
251 // hDr ->Fill(TMath::Sqrt(r2min));
253 hNrpX->Fill(rpMultX);
254 hNrpZ->Fill(rpMultZ);
256 delete [] hitsPerModule;
258 AliInfo(Form("++++ Event %d : total %d impacts, %d charged impacts and %d rec. points.",
259 ievent, nTotalGen, nChargedGen, please->CpvRecPoints()->GetEntries())) ;
263 Text_t outputname[80] ;
264 sprintf(outputname,"%s.analyzed",GetFileName().Data());
265 TFile output(outputname,"RECREATE");
270 TCanvas *cpvCanvas = new TCanvas("Cpv1","CPV analysis-I",20,20,800,600);
271 gStyle->SetOptStat(111111);
272 gStyle->SetOptFit(1);
273 gStyle->SetOptDate(1);
274 cpvCanvas->Divide(3,3);
277 gPad->SetFillColor(10);
278 hNrp->SetFillColor(16);
282 gPad->SetFillColor(10);
283 hNrpX->SetFillColor(16);
287 gPad->SetFillColor(10);
288 hNrpZ->SetFillColor(16);
292 gPad->SetFillColor(10);
293 hDx->SetFillColor(16);
298 gPad->SetFillColor(10);
299 hDz->SetFillColor(16);
304 gPad->SetFillColor(10);
305 hChi2->SetFillColor(16);
309 gPad->SetFillColor(10);
310 hEg->SetFillColor(16);
314 gPad->SetFillColor(10);
315 hEr->SetFillColor(16);
318 cpvCanvas->Write(0,kOverwrite);
323 void AliPHOSIhepAnalyze::AnalyzeEMC1(Int_t Nevents)
326 // Analyzes Emc characteristics: resolutions, cluster multiplicity,
327 // cluster lengths along Z and Phi.
328 // Author: Boris Polichtchouk, 24.08.2001
333 TH1F *hDx = new TH1F("hDx" ,"EMC x-resolution@reconstruction",100,-5. , 5.);
334 TH1F *hDz = new TH1F("hDz" ,"EMC z-resolution@reconstruction",100,-5. , 5.);
335 TH1F *hChi2 = new TH1F("hChi2" ,"Chi2/dof of one-gamma fit",30, 0. , 3.);
336 TH1S *hNrp = new TH1S("hNrp" ,"EMC rec.point multiplicity", 21,-0.5,20.5);
337 TH1S *hNrpX = new TH1S("hNrpX","EMC rec.point Phi-length" , 21,-0.5,20.5);
338 TH1S *hNrpZ = new TH1S("hNrpZ","EMC rec.point Z-length" , 21,-0.5,20.5);
340 TList * fEmcImpacts ;
341 TBranch * branchEMCimpacts;
343 AliPHOSLoader* please = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
346 AliError(Form("Could not obtain the Loader object !"));
350 const AliPHOSGeometry * fGeom = please->PHOSGeometry();
352 AliInfo(Form("Start EMC Analysis-1. Resolutions, cluster multiplicity and lengths"));
353 for ( Int_t ievent=0; ievent<Nevents; ievent++) {
357 Int_t ntracks = gAlice->GetEvent(ievent);
359 AliInfo(Form(" >>>>>>>Event %d .<<<<<<<", ievent)) ;
361 TTree* treeH = please->TreeH();
364 AliError(Form("Can not get TreeH"));
369 // Get branch of EMC impacts
370 if (! (branchEMCimpacts =treeH->GetBranch("PHOSEmcImpacts")) ) {
371 AliWarning(Form(" Couldn't find branch PHOSEmcImpacts. Exit."));
375 // Create and fill arrays of hits for each EMC module
377 Int_t nOfModules = fGeom->GetNModules();
378 TClonesArray **hitsPerModule = new TClonesArray *[nOfModules];
380 for (iModule=0; iModule < nOfModules; iModule++)
381 hitsPerModule[iModule] = new TClonesArray("AliPHOSImpact",100);
383 TClonesArray *impacts;
384 AliPHOSImpact *impact;
387 // First go through all primary tracks and fill the arrays
388 // of hits per each EMC module
390 for (Int_t itrack=0; itrack<ntracks; itrack++) {
391 branchEMCimpacts ->SetAddress(&fEmcImpacts);
392 branchEMCimpacts ->GetEntry(itrack,0);
394 for (Int_t iModule=0; iModule < nOfModules; iModule++) {
395 impacts = (TClonesArray *)fEmcImpacts->At(iModule);
396 // Do loop over impacts in the module
397 for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
398 impact=(AliPHOSImpact*)impacts->At(iImpact);
399 TClonesArray &lhits = *(TClonesArray *)hitsPerModule[iModule];
400 new(lhits[hitsPerModule[iModule]->GetEntriesFast()]) AliPHOSImpact(*impact);
403 fEmcImpacts->Clear();
405 for (iModule=0; iModule < nOfModules; iModule++) {
406 Int_t nsum = hitsPerModule[iModule]->GetEntriesFast();
407 printf("EMC module %d has %d hits\n",iModule,nsum);
411 // Then go through reconstructed points and for each find
413 // The distance from the rec.point to the closest hit
414 // gives the coordinate resolution of the EMC
416 fRunLoader->GetEvent(ievent);
417 TIter nextRP(please->EmcRecPoints()) ;
418 AliPHOSEmcRecPoint *emcRecPoint ;
419 Float_t xgen, ygen, zgen;
420 while( ( emcRecPoint = (AliPHOSEmcRecPoint *)nextRP() ) ) {
422 Float_t chi2dof = ((AliPHOSEvalRecPoint*)emcRecPoint)->Chi2Dof();
423 hChi2->Fill(chi2dof);
426 emcRecPoint->GetLocalPosition(locpos);
427 Int_t phosModule = emcRecPoint->GetPHOSMod();
428 Int_t rpMult = emcRecPoint->GetMultiplicity();
429 Int_t rpMultX, rpMultZ;
430 ((AliPHOSEvalRecPoint*)emcRecPoint)->GetClusterLengths(rpMultX,rpMultZ);
431 Float_t xrec = locpos.X();
432 Float_t zrec = locpos.Z();
433 Float_t dxmin = 1.e+10;
434 Float_t dzmin = 1.e+10;
435 Float_t r2min = 1.e+10;
438 Int_t nEMChits = (hitsPerModule[phosModule-1])->GetEntriesFast();
439 Float_t locImpX=1e10,locImpZ=1e10; // local coords of closest impact
440 Float_t gImpX=1e10, gImpZ=1e10,gImpY=1e10; // global coords of closest impact
441 for (Int_t ihit=0; ihit<nEMChits; ihit++) {
442 impact = (AliPHOSImpact*)(hitsPerModule[phosModule-1])->UncheckedAt(ihit);
448 //Transform to the local ref.frame
449 const AliPHOSGeometry* geom = please->PHOSGeometry();
450 Float_t phig = geom->GetPHOSAngle(phosModule);
451 Float_t phi = TMath::Pi()/180*phig;
452 Float_t distanceIPtoEMC = geom->GetIPtoCrystalSurface();
453 Float_t xoL,yoL,zoL ;
454 // xoL = xgen*TMath::Cos(phig)+ygen*TMath::Sin(phig) ;
455 // yoL = -xgen*TMath::Sin(phig)+ygen*TMath::Cos(phig) + distanceIPtoEMC;
456 xoL = xgen*TMath::Cos(phi)-ygen*TMath::Sin(phi) ;
457 yoL = xgen*TMath::Sin(phi)+ygen*TMath::Cos(phi) + distanceIPtoEMC;
460 r2 = TMath::Power((xoL-xrec),2) + TMath::Power((zoL-zrec),2);
472 AliInfo(Form(" Impact global (X,Z,Y) = %f %f %f", gImpX, gImpZ, gImpY));
473 AliInfo(Form(" Impact local (X,Z) = %f %f", locImpX, locImpZ));
474 AliInfo(Form(" Reconstructed (X,Z) = %f %f", xrec, zrec));
475 AliInfo(Form(" dxmin %f dzmin %f", dxmin, dzmin)) ;
478 // hDr ->Fill(TMath::Sqrt(r2min));
480 hNrpX->Fill(rpMultX);
481 hNrpZ->Fill(rpMultZ);
483 delete [] hitsPerModule;
485 AliInfo(Form("++++ Event %d : total %d impacts, %d Emc rec. points.",
486 ievent, nTotalGen, please->EmcRecPoints()->GetEntriesFast())) ;
491 Text_t outputname[80] ;
492 sprintf(outputname,"%s.analyzed",GetFileName().Data());
493 TFile output(outputname,"update");
498 TCanvas *emcCanvas = new TCanvas("Emc1","EMC analysis-I",20,20,800,400);
499 gStyle->SetOptStat(111111);
500 gStyle->SetOptFit(1);
501 gStyle->SetOptDate(1);
502 emcCanvas->Divide(3,2);
505 gPad->SetFillColor(10);
506 hNrp->SetFillColor(16);
510 gPad->SetFillColor(10);
511 hNrpX->SetFillColor(16);
515 gPad->SetFillColor(10);
516 hNrpZ->SetFillColor(16);
520 gPad->SetFillColor(10);
521 hDx->SetFillColor(16);
526 gPad->SetFillColor(10);
527 hDz->SetFillColor(16);
532 gPad->SetFillColor(10);
533 hChi2->SetFillColor(16);
536 emcCanvas->Write(0,kOverwrite);
539 //____________________________________________________________________________
540 void AliPHOSIhepAnalyze::AnalyzeCPV2(Int_t Nevents)
542 // CPV analysis - part II.
543 // Ratio of combinatoric distances between generated
544 // and reconstructed hits.
545 // Author: Boris Polichtchouk (polishchuk@mx.ihep.su)
549 TH1F* hDrijCPVr = new TH1F("Drij_cpv_r","Distance between reconstructed hits in CPV",140,0,50);
550 TH1F* hDrijCPVg = new TH1F("Drij_cpv_g","Distance between generated hits in CPV",140,0,50);
551 TH1F* hDrijCPVratio = new TH1F("Drij_cpv_ratio","R_{ij}^{rec}/R_{ij}^{gen} in CPV",140,0,50);
553 // TH1F* hT0 = new TH1F("hT0","Type of entering particle",20000,-10000,10000);
557 hDrijCPVratio->Sumw2(); //correct treatment of errors
559 TList * fCpvImpacts = new TList();
560 TBranch * branchCPVimpacts;
562 AliPHOSLoader* please = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
565 AliError(Form("Could not obtain the Loader object !"));
568 const AliPHOSGeometry * fGeom = please->PHOSGeometry();
569 fRunLoader->LoadHeader();
571 for (Int_t nev=0; nev<Nevents; nev++)
573 printf("\n=================== Event %10d ===================\n",nev);
574 fRunLoader->GetEvent(nev);
575 Int_t ntracks = fRunLoader->GetHeader()->GetNtrack();
577 Int_t nRecCPV = 0; // Reconstructed points in event
578 Int_t nGenCPV = 0; // Impacts in event
580 // Get branch of CPV impacts
581 TTree* treeH = please->TreeH();
584 AliError(Form("Can not get TreeH"));
588 if (! (branchCPVimpacts =treeH->GetBranch("PHOSCpvImpacts")) ) return;
590 // Create and fill arrays of hits for each CPV module
591 Int_t nOfModules = fGeom->GetNModules();
592 TClonesArray **hitsPerModule = new TClonesArray *[nOfModules];
594 for (iModule=0; iModule < nOfModules; iModule++)
595 hitsPerModule[iModule] = new TClonesArray("AliPHOSImpact",100);
597 TClonesArray *impacts;
598 AliPHOSImpact *impact;
600 for (Int_t itrack=0; itrack<ntracks; itrack++) {
601 branchCPVimpacts ->SetAddress(&fCpvImpacts);
602 AliInfo(Form(" branchCPVimpacts ->SetAddress(&fCpvImpacts) OK."));
603 branchCPVimpacts ->GetEntry(itrack,0);
605 for (Int_t iModule=0; iModule < nOfModules; iModule++) {
606 impacts = (TClonesArray *)fCpvImpacts->At(iModule);
607 // Do loop over impacts in the module
608 for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
609 impact=(AliPHOSImpact*)impacts->At(iImpact);
610 TClonesArray &lhits = *(TClonesArray *)hitsPerModule[iModule];
611 if(IsCharged(impact->GetPid()))
612 new(lhits[hitsPerModule[iModule]->GetEntriesFast()]) AliPHOSImpact(*impact);
615 fCpvImpacts->Clear();
618 for (iModule=0; iModule < nOfModules; iModule++) {
619 Int_t nsum = hitsPerModule[iModule]->GetEntriesFast();
620 printf("CPV module %d has %d hits\n",iModule,nsum);
622 AliPHOSImpact* genHit1;
623 AliPHOSImpact* genHit2;
625 for(irp1=0; irp1< nsum; irp1++)
627 genHit1 = (AliPHOSImpact*)((hitsPerModule[iModule])->At(irp1));
628 for(irp2 = irp1+1; irp2<nsum; irp2++)
630 genHit2 = (AliPHOSImpact*)((hitsPerModule[iModule])->At(irp2));
631 Float_t dx = genHit1->X() - genHit2->X();
632 Float_t dz = genHit1->Z() - genHit2->Z();
633 Float_t dr = TMath::Sqrt(dx*dx + dz*dz);
635 // AliInfo(Form("(dx dz dr): %f %f", dx, dz));
641 //--------- Combinatoric distance between rec. hits in CPV
643 TObjArray* cpvRecPoints = please->CpvRecPoints();
644 nRecCPV = cpvRecPoints->GetEntriesFast();
648 AliPHOSCpvRecPoint* recHit1;
649 AliPHOSCpvRecPoint* recHit2;
650 TIter nextCPVrec1(cpvRecPoints);
651 while(TObject* obj1 = nextCPVrec1() )
653 TIter nextCPVrec2(cpvRecPoints);
654 while (TObject* obj2 = nextCPVrec2())
656 if(!obj2->IsEqual(obj1))
658 recHit1 = (AliPHOSCpvRecPoint*)obj1;
659 recHit2 = (AliPHOSCpvRecPoint*)obj2;
662 recHit1->GetLocalPosition(locpos1);
663 recHit2->GetLocalPosition(locpos2);
664 Float_t dx = locpos1.X() - locpos2.X();
665 Float_t dz = locpos1.Z() - locpos2.Z();
666 Float_t dr = TMath::Sqrt(dx*dx + dz*dz);
667 if(recHit1->GetPHOSMod() == recHit2->GetPHOSMod())
674 AliInfo(Form(" Event %d . Total of %d hits, %d rec.points.",
675 nev, nGenCPV, nRecCPV)) ;
677 delete [] hitsPerModule;
679 } // End of loop over events.
682 // hDrijCPVg->Draw();
683 // hDrijCPVr->Draw();
684 hDrijCPVratio->Divide(hDrijCPVr,hDrijCPVg);
685 hDrijCPVratio->Draw();
692 void AliPHOSIhepAnalyze::CpvSingle(Int_t nevents)
694 // Distributions of coordinates and cluster lengths of rec. points
695 // in the case of single initial particle.
697 TH1F* hZr = new TH1F("Zrec","Reconstructed Z distribution",150,-5,5);
698 TH1F* hXr = new TH1F("Xrec","Reconstructed X distribution",150,-14,-2);
699 TH1F *hChi2 = new TH1F("hChi2" ,"Chi2/dof of one-gamma fit",100, 0. , 20.);
701 TH1S *hNrp = new TH1S("hNrp" ,"CPV rec.point multiplicity",21,-0.5,20.5);
702 TH1S *hNrpX = new TH1S("hNrpX","CPV rec.point Phi-length" ,21,-0.5,20.5);
703 TH1S *hNrpZ = new TH1S("hNrpZ","CPV rec.point Z-length" ,21,-0.5,20.5);
705 AliPHOSLoader* gime = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
708 AliError(Form("Could not obtain the Loader object !"));
712 for(Int_t ievent=0; ievent<nevents; ievent++)
714 fRunLoader->GetEvent(ievent);
715 if(gime->CpvRecPoints()->GetEntriesFast()>1) continue;
717 AliPHOSCpvRecPoint* pt = (AliPHOSCpvRecPoint*)(gime->CpvRecPoints())->At(0);
720 pt->GetLocalPosition(lpos);
724 Int_t rpMult = pt->GetMultiplicity();
726 Int_t rpMultX, rpMultZ;
727 pt->GetClusterLengths(rpMultX,rpMultZ);
728 hNrpX->Fill(rpMultX);
729 hNrpZ->Fill(rpMultZ);
730 hChi2->Fill(((AliPHOSEvalRecPoint*)pt)->Chi2Dof());
731 AliInfo(Form("+++++ Event %d . (Mult,MultX,MultZ) = %d %d %d +++++",
732 ievent, rpMult, rpMultX, rpMultZ)) ;
738 Text_t outputname[80] ;
739 sprintf(outputname,"%s.analyzed.single",GetFileName().Data());
740 TFile output(outputname,"RECREATE");
744 TCanvas *cpvCanvas = new TCanvas("SingleParticle","Single particle events",20,20,800,400);
745 gStyle->SetOptStat(111111);
746 gStyle->SetOptFit(1);
747 gStyle->SetOptDate(1);
748 cpvCanvas->Divide(3,2);
751 gPad->SetFillColor(10);
752 hXr->SetFillColor(16);
756 gPad->SetFillColor(10);
757 hZr->SetFillColor(16);
761 gPad->SetFillColor(10);
762 hChi2->SetFillColor(16);
766 gPad->SetFillColor(10);
767 hNrp->SetFillColor(16);
771 gPad->SetFillColor(10);
772 hNrpX->SetFillColor(16);
776 gPad->SetFillColor(10);
777 hNrpZ->SetFillColor(16);
780 cpvCanvas->Write(0,kOverwrite);
784 void AliPHOSIhepAnalyze::HitsCPV(Int_t nev)
786 // Cumulative list of charged CPV impacts in event nev.
788 TList * fCpvImpacts ;
789 TBranch * branchCPVimpacts;
791 AliPHOSLoader* please = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
794 AliError(Form("Could not obtain the Loader object !"));
797 const AliPHOSGeometry * fGeom = please->PHOSGeometry();
800 printf("\n=================== Event %10d ===================\n",nev);
801 fRunLoader->GetEvent(nev);
802 Int_t ntracks = fRunLoader->GetHeader()->GetNtrack();
804 // Int_t nRecCPV = 0; // Reconstructed points in event // 01.10.2001
805 // Int_t nGenCPV = 0; // Impacts in event
807 // Get branch of CPV impacts
808 TTree* treeH = please->TreeH();
811 AliError(Form("Can not get TreeH"));
815 if (! (branchCPVimpacts =treeH->GetBranch("PHOSCpvImpacts")) ) return;
817 // Create and fill arrays of hits for each CPV module
818 Int_t nOfModules = fGeom->GetNModules();
819 TClonesArray **hitsPerModule = new TClonesArray *[nOfModules];
821 for (iModule=0; iModule < nOfModules; iModule++)
822 hitsPerModule[iModule] = new TClonesArray("AliPHOSImpact",100);
824 TClonesArray *impacts;
825 AliPHOSImpact *impact;
827 for (Int_t itrack=0; itrack<ntracks; itrack++) {
828 branchCPVimpacts ->SetAddress(&fCpvImpacts);
829 AliInfo(Form(" branchCPVimpacts ->SetAddress(&fCpvImpacts) OK."));
830 branchCPVimpacts ->GetEntry(itrack,0);
832 for (Int_t iModule=0; iModule < nOfModules; iModule++) {
833 impacts = (TClonesArray *)fCpvImpacts->At(iModule);
834 // Do loop over impacts in the module
835 for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
836 impact=(AliPHOSImpact*)impacts->At(iImpact);
837 TClonesArray &lhits = *(TClonesArray *)hitsPerModule[iModule];
838 new(lhits[hitsPerModule[iModule]->GetEntriesFast()]) AliPHOSImpact(*impact);
841 fCpvImpacts->Clear();
844 for (iModule=0; iModule < nOfModules; iModule++) {
845 Int_t nsum = hitsPerModule[iModule]->GetEntriesFast();
846 printf("CPV module %d has %d hits\n",iModule,nsum);
849 // TList * fCpvImpacts ;
850 // TBranch * branchCPVimpacts;
851 // AliPHOSImpact* impact;
852 // TClonesArray* impacts;
854 // AliPHOSLoader * please = AliPHOSLoader::GetInstance(GetFileName().Data(),"PHOS");
855 // const AliPHOSGeometry * fGeom = please->PHOSGeometry();
857 // Int_t ntracks = gAlice->GetEvent(ievent);
858 // Int_t nOfModules = fGeom->GetNModules();
859 // AliInfo(Form(" Tracks: "<<ntracks<<" Modules: "<<nOfModules));
861 // if (! (branchCPVimpacts =gAlice->TreeH()->GetBranch("PHOSCpvImpacts")) ) return;
863 // for (Int_t itrack=0; itrack<ntracks; itrack++) {
864 // branchCPVimpacts ->SetAddress(&fCpvImpacts);
865 // Info("AnalyzeCPV1", " branchCPVimpacts ->SetAddress(&fCpvImpacts) OK.");
866 // branchCPVimpacts ->GetEntry(itrack,0);
867 // Info(Form(" branchCPVimpacts ->GetEntry(itrack,0) OK."));
869 // for (Int_t iModule=0; iModule < nOfModules; iModule++) {
870 // impacts = (TClonesArray *)fCpvImpacts->At(iModule);
871 // Info(Form(" fCpvImpacts->At(iModule) OK."));
872 // // Do loop over impacts in the module
873 // for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
874 // impact=(AliPHOSImpact*)impacts->At(iImpact);
876 // if(IsCharged(impact->GetPid()))
878 // Info(Form(" Add charged hit.."));
879 // new(hits[hits.GetEntriesFast()]) AliPHOSImpact(*impact);
880 // Info(Form("done."));
884 // fCpvImpacts->Clear();
887 // Info(Form(" PHOS event "<<ievent<<": "<<hits.GetEntries()<<" charged CPV hits."));
892 // void AliPHOSIhepAnalyze::ChargedHitsCPV(TClonesArray* hits, Int_t ievent, Int_t iModule)
894 // // Cumulative list of charged CPV hits in event ievent
895 // // in PHOS module iModule.
897 // HitsCPV(hits,ievent,iModule);
900 // while(AliPHOSCPVHit* cpvhit = (AliPHOSCPVHit*)next())
902 // if(!IsCharged(cpvhit->GetIpart()))
904 // hits->Remove(cpvhit);
910 // Info(Form(" PHOS module "<<iModule<<": "<<hits->GetEntries()<<" charged CPV hits."));
913 Bool_t AliPHOSIhepAnalyze::IsCharged(Int_t pdgCode)
916 AliInfo(Form("pdgCode %d", pdgCode));
917 if(pdgCode==211 || pdgCode==-211 || pdgCode==321 || pdgCode==-321 || pdgCode==11 || pdgCode==-11 || pdgCode==2212 || pdgCode==-2212) return kTRUE;
921 //---------------------------------------------------------------------------