Line data Source code
1 : /**************************************************************************
2 : * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 : * *
4 : * Author: The ALICE Off-line Project. *
5 : * Contributors are mentioned in the code where appropriate. *
6 : * *
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 : **************************************************************************/
15 :
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
21 : //*--
22 : //*-- Author: B. Polichtchouk (IHEP)
23 : //////////////////////////////////////////////////////////////////////////////
24 :
25 : // --- ROOT system ---
26 :
27 : #include "TFile.h"
28 : #include "TH1.h"
29 : #include "TPad.h"
30 : #include "TH2.h"
31 : #include "TParticle.h"
32 : #include "TClonesArray.h"
33 : #include "TTree.h"
34 : #include "TMath.h"
35 : #include "TCanvas.h"
36 : #include "TStyle.h"
37 :
38 : // --- Standard library ---
39 :
40 : // --- AliRoot header files ---
41 :
42 : #include "AliRunLoader.h"
43 : #include "AliHeader.h"
44 :
45 : // --- PHOS header files ---
46 : #include "AliLog.h"
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"
55 : #include "AliRun.h"
56 : #include "AliPHOSGeometry.h"
57 : #include "AliPHOSEvalRecPoint.h"
58 :
59 20 : ClassImp(AliPHOSIhepAnalyze)
60 :
61 : //____________________________________________________________________________
62 0 : AliPHOSIhepAnalyze::AliPHOSIhepAnalyze():
63 0 : fRunLoader(0),
64 0 : fFileName()
65 0 : {
66 0 : }
67 :
68 : //____________________________________________________________________________
69 0 : AliPHOSIhepAnalyze::AliPHOSIhepAnalyze(Text_t * name) :
70 0 : fRunLoader(0),
71 0 : fFileName(name)
72 0 : {
73 : // Constructor: open a header file
74 0 : fRunLoader = AliRunLoader::Open(fFileName);
75 0 : if (fRunLoader == 0x0)
76 : {
77 0 : AliFatal(Form("Can not load event from file %s",name));
78 : }
79 0 : }
80 :
81 : //____________________________________________________________________________
82 : void AliPHOSIhepAnalyze::AnalyzeCPV1(Int_t Nevents)
83 : {
84 : //
85 : // Analyzes CPV characteristics: resolutions, cluster multiplicity,
86 : // cluster lengths along Z and Phi.
87 : // Author: Yuri Kharlov
88 : // 9 October 2000
89 : // Modified by Boris Polichtchouk, 3.07.2001
90 : //
91 :
92 : // Book histograms
93 :
94 0 : TH1F *hDx = new TH1F("hDx" ,"CPV x-resolution@reconstruction",100,-5. , 5.);
95 0 : TH1F *hDz = new TH1F("hDz" ,"CPV z-resolution@reconstruction",100,-5. , 5.);
96 0 : TH1F *hChi2 = new TH1F("hChi2" ,"Chi2/dof of one-gamma fit",30, 0. , 10.);
97 0 : TH1S *hNrp = new TH1S("hNrp" ,"CPV rec.point multiplicity", 21,-0.5,20.5);
98 0 : TH1S *hNrpX = new TH1S("hNrpX","CPV rec.point Phi-length" , 21,-0.5,20.5);
99 0 : TH1S *hNrpZ = new TH1S("hNrpZ","CPV rec.point Z-length" , 21,-0.5,20.5);
100 :
101 0 : TH1F *hEg = new TH1F("hEg" ,"Energies of impacts",30,0.,6.);
102 0 : TH1F *hEr = new TH1F("hEr" ,"Amplitudes of rec. points",50,0.,20.);
103 :
104 0 : TList * fCpvImpacts ;
105 : TBranch * branchCPVimpacts;
106 :
107 :
108 :
109 0 : AliPHOSLoader* please = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
110 0 : if ( please == 0 )
111 : {
112 0 : AliError(Form("Could not obtain the Loader object !"));
113 0 : return ;
114 : }
115 :
116 0 : AliPHOSGeometry * phosgeom = AliPHOSGeometry::GetInstance() ;
117 :
118 0 : AliInfo(Form("Start CPV Analysis-1. Resolutions, cluster multiplicity and lengths")) ;
119 0 : for ( Int_t ievent=0; ievent<Nevents; ievent++) {
120 :
121 : Int_t nTotalGen = 0;
122 : Int_t nChargedGen = 0;
123 :
124 0 : Int_t ntracks = gAlice->GetEvent(ievent);
125 0 : AliInfo(Form(">>>>>>>Event %d .<<<<<<<", ievent)) ;
126 :
127 : /******************************************************************/
128 0 : TTree* treeH = please->TreeH();
129 0 : if (treeH == 0x0)
130 : {
131 0 : AliError(Form("Can not get TreeH"));
132 0 : return;
133 : }
134 : /******************************************************************/
135 :
136 : // Get branch of CPV impacts
137 0 : if (! (branchCPVimpacts =treeH->GetBranch("PHOSCpvImpacts")) ) {
138 0 : AliWarning(Form("Couldn't find branch PHOSCpvImpacts. Exit.")) ;
139 0 : return;
140 : }
141 :
142 : // Create and fill arrays of hits for each CPV module
143 :
144 0 : Int_t nOfModules = phosgeom->GetNModules();
145 0 : TClonesArray **hitsPerModule = new TClonesArray *[nOfModules];
146 : Int_t iModule = 0;
147 0 : for (iModule=0; iModule < nOfModules; iModule++)
148 0 : hitsPerModule[iModule] = new TClonesArray("AliPHOSImpact",100);
149 :
150 : TClonesArray *impacts;
151 : AliPHOSImpact *impact;
152 0 : TLorentzVector p;
153 :
154 : // First go through all primary tracks and fill the arrays
155 : // of hits per each CPV module
156 :
157 0 : for (Int_t itrack=0; itrack<ntracks; itrack++) {
158 0 : branchCPVimpacts ->SetAddress(&fCpvImpacts);
159 0 : branchCPVimpacts ->GetEntry(itrack,0);
160 :
161 0 : for (iModule=0; iModule < nOfModules; iModule++) {
162 0 : impacts = (TClonesArray *)fCpvImpacts->At(iModule);
163 : // Do loop over impacts in the module
164 0 : for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
165 0 : impact=(AliPHOSImpact*)impacts->At(iImpact);
166 0 : hEg->Fill(impact->GetMomentum().E());
167 0 : TClonesArray &lhits = *(TClonesArray *)hitsPerModule[iModule];
168 0 : if(IsCharged(impact->GetPid())) nChargedGen++;
169 0 : new(lhits[hitsPerModule[iModule]->GetEntriesFast()]) AliPHOSImpact(*impact);
170 : }
171 : }
172 0 : fCpvImpacts->Clear();
173 : }
174 0 : for (iModule=0; iModule < nOfModules; iModule++) {
175 0 : Int_t nsum = hitsPerModule[iModule]->GetEntriesFast();
176 0 : printf("CPV module %d has %d impacts\n",iModule,nsum);
177 0 : nTotalGen += nsum;
178 : }
179 :
180 : // Then go through reconstructed points and for each find
181 : // the closeset hit
182 : // The distance from the rec.point to the closest hit
183 : // gives the coordinate resolution of the CPV
184 :
185 0 : fRunLoader->GetEvent(ievent);
186 0 : TIter nextRP(please->CpvRecPoints()) ;
187 : AliPHOSCpvRecPoint *cpvRecPoint ;
188 : Float_t xgen, ygen, zgen;
189 0 : while( ( cpvRecPoint = (AliPHOSCpvRecPoint *)nextRP() ) ) {
190 :
191 0 : Float_t chi2dof = ((AliPHOSEvalRecPoint*)cpvRecPoint)->Chi2Dof();
192 0 : hChi2->Fill(chi2dof);
193 0 : hEr->Fill(cpvRecPoint->GetEnergy());
194 :
195 0 : TVector3 locpos;
196 0 : cpvRecPoint->GetLocalPosition(locpos);
197 0 : Int_t phosModule = cpvRecPoint->GetPHOSMod();
198 0 : Int_t rpMult = cpvRecPoint->GetMultiplicity();
199 0 : Int_t rpMultX, rpMultZ;
200 0 : cpvRecPoint->GetClusterLengths(rpMultX,rpMultZ);
201 0 : Float_t xrec = locpos.X();
202 0 : 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;
206 : Float_t r2;
207 :
208 0 : 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 0 : for (Int_t ihit=0; ihit<nCPVhits; ihit++) {
212 0 : impact = (AliPHOSImpact*)(hitsPerModule[phosModule-1])->UncheckedAt(ihit);
213 0 : xgen = impact->X();
214 0 : zgen = impact->Z();
215 0 : ygen = impact->Y();
216 :
217 : //Transform to the local ref.frame
218 0 : Float_t phig = phosgeom->GetPHOSAngle(phosModule);
219 0 : Float_t phi = TMath::Pi()/180*phig;
220 0 : Float_t distanceIPtoCPV = phosgeom->GetIPtoOuterCoverDistance() -
221 0 : (phosgeom->GetFTPosition(1)+
222 0 : phosgeom->GetFTPosition(2)+
223 0 : phosgeom->GetCPVTextoliteThickness())/2;
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 0 : xoL = xgen*TMath::Cos(phi)-ygen*TMath::Sin(phi) ;
228 0 : yoL = xgen*TMath::Sin(phi)+ygen*TMath::Cos(phi) + distanceIPtoCPV;
229 : zoL = zgen;
230 :
231 0 : r2 = TMath::Power((xoL-xrec),2) + TMath::Power((zoL-zrec),2);
232 0 : if ( r2 < r2min ) {
233 : r2min = r2;
234 : dxmin = xoL - xrec;
235 : dzmin = zoL - zrec;
236 : locImpX = xoL;
237 : locImpZ = zoL;
238 : gImpX = xgen;
239 : gImpZ = zgen;
240 : gImpY = ygen;
241 0 : }
242 : }
243 0 : AliInfo(Form("Impact global (X,Z,Y) = %f %f %f", gImpX, gImpZ, gImpY));
244 0 : AliInfo(Form("Impact local (X,Z) = %f %f", locImpX, locImpZ));
245 0 : AliInfo(Form("Reconstructed (X,Z) = %f %f", xrec, zrec));
246 0 : AliInfo(Form("dxmin %f dzmin %f", dxmin, dzmin));
247 0 : hDx ->Fill(dxmin);
248 0 : hDz ->Fill(dzmin);
249 : // hDr ->Fill(TMath::Sqrt(r2min));
250 0 : hNrp ->Fill(rpMult);
251 0 : hNrpX->Fill(rpMultX);
252 0 : hNrpZ->Fill(rpMultZ);
253 0 : }
254 0 : delete [] hitsPerModule;
255 :
256 0 : AliInfo(Form("++++ Event %d : total %d impacts, %d charged impacts and %d rec. points.",
257 : ievent, nTotalGen, nChargedGen, please->CpvRecPoints()->GetEntries())) ;
258 0 : }
259 : // Save histograms
260 :
261 0 : Text_t outputname[80] ;
262 0 : snprintf(outputname,80,"%s.analyzed",GetFileName().Data());
263 0 : TFile output(outputname,"RECREATE");
264 0 : output.cd();
265 :
266 : // Plot histograms
267 :
268 0 : TCanvas *cpvCanvas = new TCanvas("Cpv1","CPV analysis-I",20,20,800,600);
269 0 : gStyle->SetOptStat(111111);
270 0 : gStyle->SetOptFit(1);
271 0 : gStyle->SetOptDate(1);
272 0 : cpvCanvas->Divide(3,3);
273 :
274 0 : cpvCanvas->cd(1);
275 0 : gPad->SetFillColor(10);
276 0 : hNrp->SetFillColor(16);
277 0 : hNrp->Draw();
278 :
279 0 : cpvCanvas->cd(2);
280 0 : gPad->SetFillColor(10);
281 0 : hNrpX->SetFillColor(16);
282 0 : hNrpX->Draw();
283 :
284 0 : cpvCanvas->cd(3);
285 0 : gPad->SetFillColor(10);
286 0 : hNrpZ->SetFillColor(16);
287 0 : hNrpZ->Draw();
288 :
289 0 : cpvCanvas->cd(4);
290 0 : gPad->SetFillColor(10);
291 0 : hDx->SetFillColor(16);
292 0 : hDx->Fit("gaus");
293 0 : hDx->Draw();
294 :
295 0 : cpvCanvas->cd(5);
296 0 : gPad->SetFillColor(10);
297 0 : hDz->SetFillColor(16);
298 0 : hDz->Fit("gaus");
299 0 : hDz->Draw();
300 :
301 0 : cpvCanvas->cd(6);
302 0 : gPad->SetFillColor(10);
303 0 : hChi2->SetFillColor(16);
304 0 : hChi2->Draw();
305 :
306 0 : cpvCanvas->cd(7);
307 0 : gPad->SetFillColor(10);
308 0 : hEg->SetFillColor(16);
309 0 : hEg->Draw();
310 :
311 0 : cpvCanvas->cd(8);
312 0 : gPad->SetFillColor(10);
313 0 : hEr->SetFillColor(16);
314 0 : hEr->Draw();
315 :
316 0 : cpvCanvas->Write(0,kOverwrite);
317 :
318 0 : }
319 :
320 :
321 : void AliPHOSIhepAnalyze::AnalyzeEMC1(Int_t Nevents)
322 : {
323 : //
324 : // Analyzes Emc characteristics: resolutions, cluster multiplicity,
325 : // cluster lengths along Z and Phi.
326 : // Author: Boris Polichtchouk, 24.08.2001
327 : //
328 :
329 : // Book histograms
330 :
331 0 : TH1F *hDx = new TH1F("hDx" ,"EMC x-resolution@reconstruction",100,-5. , 5.);
332 0 : TH1F *hDz = new TH1F("hDz" ,"EMC z-resolution@reconstruction",100,-5. , 5.);
333 0 : TH1F *hChi2 = new TH1F("hChi2" ,"Chi2/dof of one-gamma fit",30, 0. , 3.);
334 0 : TH1S *hNrp = new TH1S("hNrp" ,"EMC rec.point multiplicity", 21,-0.5,20.5);
335 0 : TH1S *hNrpX = new TH1S("hNrpX","EMC rec.point Phi-length" , 21,-0.5,20.5);
336 0 : TH1S *hNrpZ = new TH1S("hNrpZ","EMC rec.point Z-length" , 21,-0.5,20.5);
337 :
338 0 : TList * fEmcImpacts ;
339 : TBranch * branchEMCimpacts;
340 :
341 0 : AliPHOSLoader* please = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
342 0 : if ( please == 0 )
343 : {
344 0 : AliError(Form("Could not obtain the Loader object !"));
345 0 : return ;
346 : }
347 :
348 0 : AliPHOSGeometry * phosgeom = AliPHOSGeometry::GetInstance() ;
349 :
350 0 : AliInfo(Form("Start EMC Analysis-1. Resolutions, cluster multiplicity and lengths"));
351 0 : for ( Int_t ievent=0; ievent<Nevents; ievent++) {
352 :
353 : Int_t nTotalGen = 0;
354 :
355 0 : Int_t ntracks = gAlice->GetEvent(ievent);
356 :
357 0 : AliInfo(Form(" >>>>>>>Event %d .<<<<<<<", ievent)) ;
358 :
359 0 : TTree* treeH = please->TreeH();
360 0 : if (treeH == 0x0)
361 : {
362 0 : AliError(Form("Can not get TreeH"));
363 0 : return;
364 : }
365 :
366 :
367 : // Get branch of EMC impacts
368 0 : if (! (branchEMCimpacts =treeH->GetBranch("PHOSEmcImpacts")) ) {
369 0 : AliWarning(Form(" Couldn't find branch PHOSEmcImpacts. Exit."));
370 0 : return;
371 : }
372 :
373 : // Create and fill arrays of hits for each EMC module
374 :
375 0 : Int_t nOfModules = phosgeom->GetNModules();
376 0 : TClonesArray **hitsPerModule = new TClonesArray *[nOfModules];
377 : Int_t iModule = 0;
378 0 : for (iModule=0; iModule < nOfModules; iModule++)
379 0 : hitsPerModule[iModule] = new TClonesArray("AliPHOSImpact",100);
380 :
381 : TClonesArray *impacts;
382 : AliPHOSImpact *impact;
383 0 : TLorentzVector p;
384 :
385 : // First go through all primary tracks and fill the arrays
386 : // of hits per each EMC module
387 :
388 0 : for (Int_t itrack=0; itrack<ntracks; itrack++) {
389 0 : branchEMCimpacts ->SetAddress(&fEmcImpacts);
390 0 : branchEMCimpacts ->GetEntry(itrack,0);
391 :
392 0 : for (iModule=0; iModule < nOfModules; iModule++) {
393 0 : impacts = (TClonesArray *)fEmcImpacts->At(iModule);
394 : // Do loop over impacts in the module
395 0 : for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
396 0 : impact=(AliPHOSImpact*)impacts->At(iImpact);
397 0 : TClonesArray &lhits = *(TClonesArray *)hitsPerModule[iModule];
398 0 : new(lhits[hitsPerModule[iModule]->GetEntriesFast()]) AliPHOSImpact(*impact);
399 : }
400 : }
401 0 : fEmcImpacts->Clear();
402 : }
403 0 : for (iModule=0; iModule < nOfModules; iModule++) {
404 0 : Int_t nsum = hitsPerModule[iModule]->GetEntriesFast();
405 0 : printf("EMC module %d has %d hits\n",iModule,nsum);
406 0 : nTotalGen += nsum;
407 : }
408 :
409 : // Then go through reconstructed points and for each find
410 : // the closeset hit
411 : // The distance from the rec.point to the closest hit
412 : // gives the coordinate resolution of the EMC
413 :
414 0 : fRunLoader->GetEvent(ievent);
415 0 : TIter nextRP(please->EmcRecPoints()) ;
416 : AliPHOSEmcRecPoint *emcRecPoint ;
417 : Float_t xgen, ygen, zgen;
418 0 : while( ( emcRecPoint = (AliPHOSEmcRecPoint *)nextRP() ) ) {
419 :
420 0 : Float_t chi2dof = ((AliPHOSEvalRecPoint*)emcRecPoint)->Chi2Dof();
421 0 : hChi2->Fill(chi2dof);
422 :
423 0 : TVector3 locpos;
424 0 : emcRecPoint->GetLocalPosition(locpos);
425 0 : Int_t phosModule = emcRecPoint->GetPHOSMod();
426 0 : Int_t rpMult = emcRecPoint->GetMultiplicity();
427 0 : Int_t rpMultX, rpMultZ;
428 0 : ((AliPHOSEvalRecPoint*)emcRecPoint)->GetClusterLengths(rpMultX,rpMultZ);
429 0 : Float_t xrec = locpos.X();
430 0 : 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;
434 : Float_t r2;
435 :
436 0 : 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 0 : for (Int_t ihit=0; ihit<nEMChits; ihit++) {
440 0 : impact = (AliPHOSImpact*)(hitsPerModule[phosModule-1])->UncheckedAt(ihit);
441 0 : xgen = impact->X();
442 0 : zgen = impact->Z();
443 0 : ygen = impact->Y();
444 :
445 :
446 : //Transform to the local ref.frame
447 0 : Float_t phig = phosgeom->GetPHOSAngle(phosModule);
448 0 : Float_t phi = TMath::Pi()/180*phig;
449 0 : Float_t distanceIPtoEMC = phosgeom->GetIPtoCrystalSurface();
450 : Float_t xoL,yoL,zoL ;
451 : // xoL = xgen*TMath::Cos(phig)+ygen*TMath::Sin(phig) ;
452 : // yoL = -xgen*TMath::Sin(phig)+ygen*TMath::Cos(phig) + distanceIPtoEMC;
453 0 : xoL = xgen*TMath::Cos(phi)-ygen*TMath::Sin(phi) ;
454 0 : yoL = xgen*TMath::Sin(phi)+ygen*TMath::Cos(phi) + distanceIPtoEMC;
455 : zoL = zgen;
456 :
457 0 : r2 = TMath::Power((xoL-xrec),2) + TMath::Power((zoL-zrec),2);
458 0 : if ( r2 < r2min ) {
459 : r2min = r2;
460 : dxmin = xoL - xrec;
461 : dzmin = zoL - zrec;
462 : locImpX = xoL;
463 : locImpZ = zoL;
464 : gImpX = xgen;
465 : gImpZ = zgen;
466 : gImpY = ygen;
467 0 : }
468 : }
469 0 : AliInfo(Form(" Impact global (X,Z,Y) = %f %f %f", gImpX, gImpZ, gImpY));
470 0 : AliInfo(Form(" Impact local (X,Z) = %f %f", locImpX, locImpZ));
471 0 : AliInfo(Form(" Reconstructed (X,Z) = %f %f", xrec, zrec));
472 0 : AliInfo(Form(" dxmin %f dzmin %f", dxmin, dzmin)) ;
473 0 : hDx ->Fill(dxmin);
474 0 : hDz ->Fill(dzmin);
475 : // hDr ->Fill(TMath::Sqrt(r2min));
476 0 : hNrp ->Fill(rpMult);
477 0 : hNrpX->Fill(rpMultX);
478 0 : hNrpZ->Fill(rpMultZ);
479 0 : }
480 0 : delete [] hitsPerModule;
481 :
482 0 : AliInfo(Form("++++ Event %d : total %d impacts, %d Emc rec. points.",
483 : ievent, nTotalGen, please->EmcRecPoints()->GetEntriesFast())) ;
484 :
485 0 : }
486 : // Save histograms
487 :
488 0 : Text_t outputname[80] ;
489 0 : snprintf(outputname,80,"%s.analyzed",GetFileName().Data());
490 0 : TFile output(outputname,"update");
491 0 : output.cd();
492 :
493 : // Plot histograms
494 :
495 0 : TCanvas *emcCanvas = new TCanvas("Emc1","EMC analysis-I",20,20,800,400);
496 0 : gStyle->SetOptStat(111111);
497 0 : gStyle->SetOptFit(1);
498 0 : gStyle->SetOptDate(1);
499 0 : emcCanvas->Divide(3,2);
500 :
501 0 : emcCanvas->cd(1);
502 0 : gPad->SetFillColor(10);
503 0 : hNrp->SetFillColor(16);
504 0 : hNrp->Draw();
505 :
506 0 : emcCanvas->cd(2);
507 0 : gPad->SetFillColor(10);
508 0 : hNrpX->SetFillColor(16);
509 0 : hNrpX->Draw();
510 :
511 0 : emcCanvas->cd(3);
512 0 : gPad->SetFillColor(10);
513 0 : hNrpZ->SetFillColor(16);
514 0 : hNrpZ->Draw();
515 :
516 0 : emcCanvas->cd(4);
517 0 : gPad->SetFillColor(10);
518 0 : hDx->SetFillColor(16);
519 0 : hDx->Fit("gaus");
520 0 : hDx->Draw();
521 :
522 0 : emcCanvas->cd(5);
523 0 : gPad->SetFillColor(10);
524 0 : hDz->SetFillColor(16);
525 0 : hDz->Fit("gaus");
526 0 : hDz->Draw();
527 :
528 0 : emcCanvas->cd(6);
529 0 : gPad->SetFillColor(10);
530 0 : hChi2->SetFillColor(16);
531 0 : hChi2->Draw();
532 :
533 0 : emcCanvas->Write(0,kOverwrite);
534 0 : }
535 :
536 : //____________________________________________________________________________
537 : void AliPHOSIhepAnalyze::AnalyzeCPV2(Int_t Nevents)
538 : {
539 : // CPV analysis - part II.
540 : // Ratio of combinatoric distances between generated
541 : // and reconstructed hits.
542 : // Author: Boris Polichtchouk (polishchuk@mx.ihep.su)
543 : // 24 March 2001
544 :
545 :
546 0 : TH1F* hDrijCPVr = new TH1F("Drij_cpv_r","Distance between reconstructed hits in CPV",140,0,50);
547 0 : TH1F* hDrijCPVg = new TH1F("Drij_cpv_g","Distance between generated hits in CPV",140,0,50);
548 0 : TH1F* hDrijCPVratio = new TH1F("Drij_cpv_ratio","R_{ij}^{rec}/R_{ij}^{gen} in CPV",140,0,50);
549 :
550 : // TH1F* hT0 = new TH1F("hT0","Type of entering particle",20000,-10000,10000);
551 :
552 0 : hDrijCPVr->Sumw2();
553 0 : hDrijCPVg->Sumw2();
554 0 : hDrijCPVratio->Sumw2(); //correct treatment of errors
555 :
556 0 : TList * fCpvImpacts = new TList();
557 : TBranch * branchCPVimpacts;
558 :
559 0 : AliPHOSLoader* please = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
560 0 : if ( please == 0 )
561 : {
562 0 : AliError(Form("Could not obtain the Loader object !"));
563 0 : return ;
564 : }
565 0 : AliPHOSGeometry * phosgeom = AliPHOSGeometry::GetInstance() ;
566 0 : fRunLoader->LoadHeader();
567 :
568 0 : for (Int_t nev=0; nev<Nevents; nev++)
569 : {
570 0 : printf("\n=================== Event %10d ===================\n",nev);
571 0 : fRunLoader->GetEvent(nev);
572 0 : Int_t ntracks = fRunLoader->GetHeader()->GetNtrack();
573 :
574 : Int_t nRecCPV = 0; // Reconstructed points in event
575 : Int_t nGenCPV = 0; // Impacts in event
576 :
577 : // Get branch of CPV impacts
578 0 : TTree* treeH = please->TreeH();
579 0 : if (treeH == 0x0)
580 : {
581 0 : AliError(Form("Can not get TreeH"));
582 0 : return;
583 : }
584 :
585 0 : if (! (branchCPVimpacts =treeH->GetBranch("PHOSCpvImpacts")) ) return;
586 :
587 : // Create and fill arrays of hits for each CPV module
588 0 : Int_t nOfModules = phosgeom->GetNModules();
589 0 : TClonesArray **hitsPerModule = new TClonesArray *[nOfModules];
590 : Int_t iModule = 0;
591 0 : for (iModule=0; iModule < nOfModules; iModule++)
592 0 : hitsPerModule[iModule] = new TClonesArray("AliPHOSImpact",100);
593 :
594 : TClonesArray *impacts;
595 : AliPHOSImpact *impact;
596 :
597 0 : for (Int_t itrack=0; itrack<ntracks; itrack++) {
598 0 : branchCPVimpacts ->SetAddress(&fCpvImpacts);
599 0 : AliInfo(Form(" branchCPVimpacts ->SetAddress(&fCpvImpacts) OK."));
600 0 : branchCPVimpacts ->GetEntry(itrack,0);
601 :
602 0 : for (iModule=0; iModule < nOfModules; iModule++) {
603 0 : impacts = (TClonesArray *)fCpvImpacts->At(iModule);
604 : // Do loop over impacts in the module
605 0 : for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
606 0 : impact=(AliPHOSImpact*)impacts->At(iImpact);
607 0 : TClonesArray &lhits = *(TClonesArray *)hitsPerModule[iModule];
608 0 : if(IsCharged(impact->GetPid()))
609 0 : new(lhits[hitsPerModule[iModule]->GetEntriesFast()]) AliPHOSImpact(*impact);
610 : }
611 : }
612 0 : fCpvImpacts->Clear();
613 : }
614 :
615 0 : for (iModule=0; iModule < nOfModules; iModule++) {
616 0 : Int_t nsum = hitsPerModule[iModule]->GetEntriesFast();
617 0 : printf("CPV module %d has %d hits\n",iModule,nsum);
618 :
619 : AliPHOSImpact* genHit1;
620 : AliPHOSImpact* genHit2;
621 : Int_t irp1,irp2;
622 0 : for(irp1=0; irp1< nsum; irp1++)
623 : {
624 0 : genHit1 = (AliPHOSImpact*)((hitsPerModule[iModule])->At(irp1));
625 0 : for(irp2 = irp1+1; irp2<nsum; irp2++)
626 : {
627 0 : genHit2 = (AliPHOSImpact*)((hitsPerModule[iModule])->At(irp2));
628 0 : Float_t dx = genHit1->X() - genHit2->X();
629 0 : Float_t dz = genHit1->Z() - genHit2->Z();
630 0 : Float_t dr = TMath::Sqrt(dx*dx + dz*dz);
631 0 : hDrijCPVg->Fill(dr);
632 : // AliInfo(Form("(dx dz dr): %f %f", dx, dz));
633 : }
634 : }
635 : }
636 :
637 :
638 : //--------- Combinatoric distance between rec. hits in CPV
639 :
640 0 : TObjArray* cpvRecPoints = please->CpvRecPoints();
641 0 : nRecCPV = cpvRecPoints->GetEntriesFast();
642 :
643 0 : if(nRecCPV)
644 : {
645 : AliPHOSCpvRecPoint* recHit1;
646 : AliPHOSCpvRecPoint* recHit2;
647 0 : TIter nextCPVrec1(cpvRecPoints);
648 0 : while(TObject* obj1 = nextCPVrec1() )
649 : {
650 0 : TIter nextCPVrec2(cpvRecPoints);
651 0 : while (TObject* obj2 = nextCPVrec2())
652 : {
653 0 : if(!obj2->IsEqual(obj1))
654 : {
655 0 : recHit1 = (AliPHOSCpvRecPoint*)obj1;
656 0 : recHit2 = (AliPHOSCpvRecPoint*)obj2;
657 0 : TVector3 locpos1;
658 0 : TVector3 locpos2;
659 0 : recHit1->GetLocalPosition(locpos1);
660 0 : recHit2->GetLocalPosition(locpos2);
661 0 : Float_t dx = locpos1.X() - locpos2.X();
662 0 : Float_t dz = locpos1.Z() - locpos2.Z();
663 0 : Float_t dr = TMath::Sqrt(dx*dx + dz*dz);
664 0 : if(recHit1->GetPHOSMod() == recHit2->GetPHOSMod())
665 0 : hDrijCPVr->Fill(dr);
666 0 : }
667 0 : }
668 0 : }
669 0 : }
670 :
671 0 : AliInfo(Form(" Event %d . Total of %d hits, %d rec.points.",
672 : nev, nGenCPV, nRecCPV)) ;
673 :
674 0 : delete [] hitsPerModule;
675 :
676 0 : } // End of loop over events.
677 :
678 :
679 : // hDrijCPVg->Draw();
680 : // hDrijCPVr->Draw();
681 0 : hDrijCPVratio->Divide(hDrijCPVr,hDrijCPVg);
682 0 : hDrijCPVratio->Draw();
683 :
684 : // hT0->Draw();
685 :
686 0 : }
687 :
688 :
689 : void AliPHOSIhepAnalyze::CpvSingle(Int_t nevents)
690 : {
691 : // Distributions of coordinates and cluster lengths of rec. points
692 : // in the case of single initial particle.
693 :
694 0 : TH1F* hZr = new TH1F("Zrec","Reconstructed Z distribution",150,-5,5);
695 0 : TH1F* hXr = new TH1F("Xrec","Reconstructed X distribution",150,-14,-2);
696 0 : TH1F *hChi2 = new TH1F("hChi2" ,"Chi2/dof of one-gamma fit",100, 0. , 20.);
697 :
698 0 : TH1S *hNrp = new TH1S("hNrp" ,"CPV rec.point multiplicity",21,-0.5,20.5);
699 0 : TH1S *hNrpX = new TH1S("hNrpX","CPV rec.point Phi-length" ,21,-0.5,20.5);
700 0 : TH1S *hNrpZ = new TH1S("hNrpZ","CPV rec.point Z-length" ,21,-0.5,20.5);
701 :
702 0 : AliPHOSLoader* gime = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
703 0 : if ( gime == 0 )
704 : {
705 0 : AliError(Form("Could not obtain the Loader object !"));
706 0 : return ;
707 : }
708 :
709 0 : for(Int_t ievent=0; ievent<nevents; ievent++)
710 : {
711 0 : fRunLoader->GetEvent(ievent);
712 0 : if(gime->CpvRecPoints()->GetEntriesFast()>1) continue;
713 :
714 0 : AliPHOSCpvRecPoint* pt = (AliPHOSCpvRecPoint*)(gime->CpvRecPoints())->At(0);
715 0 : if(pt) {
716 0 : TVector3 lpos;
717 0 : pt->GetLocalPosition(lpos);
718 0 : hXr->Fill(lpos.X());
719 0 : hZr->Fill(lpos.Z());
720 :
721 0 : Int_t rpMult = pt->GetMultiplicity();
722 0 : hNrp->Fill(rpMult);
723 0 : Int_t rpMultX, rpMultZ;
724 0 : pt->GetClusterLengths(rpMultX,rpMultZ);
725 0 : hNrpX->Fill(rpMultX);
726 0 : hNrpZ->Fill(rpMultZ);
727 0 : hChi2->Fill(((AliPHOSEvalRecPoint*)pt)->Chi2Dof());
728 0 : AliInfo(Form("+++++ Event %d . (Mult,MultX,MultZ) = %d %d %d +++++",
729 : ievent, rpMult, rpMultX, rpMultZ)) ;
730 :
731 0 : }
732 :
733 0 : }
734 :
735 0 : Text_t outputname[80] ;
736 0 : snprintf(outputname,80,"%s.analyzed.single",GetFileName().Data());
737 0 : TFile output(outputname,"RECREATE");
738 0 : output.cd();
739 :
740 : // Plot histograms
741 0 : TCanvas *cpvCanvas = new TCanvas("SingleParticle","Single particle events",20,20,800,400);
742 0 : gStyle->SetOptStat(111111);
743 0 : gStyle->SetOptFit(1);
744 0 : gStyle->SetOptDate(1);
745 0 : cpvCanvas->Divide(3,2);
746 :
747 0 : cpvCanvas->cd(1);
748 0 : gPad->SetFillColor(10);
749 0 : hXr->SetFillColor(16);
750 0 : hXr->Draw();
751 :
752 0 : cpvCanvas->cd(2);
753 0 : gPad->SetFillColor(10);
754 0 : hZr->SetFillColor(16);
755 0 : hZr->Draw();
756 :
757 0 : cpvCanvas->cd(3);
758 0 : gPad->SetFillColor(10);
759 0 : hChi2->SetFillColor(16);
760 0 : hChi2->Draw();
761 :
762 0 : cpvCanvas->cd(4);
763 0 : gPad->SetFillColor(10);
764 0 : hNrp->SetFillColor(16);
765 0 : hNrp->Draw();
766 :
767 0 : cpvCanvas->cd(5);
768 0 : gPad->SetFillColor(10);
769 0 : hNrpX->SetFillColor(16);
770 0 : hNrpX->Draw();
771 :
772 0 : cpvCanvas->cd(6);
773 0 : gPad->SetFillColor(10);
774 0 : hNrpZ->SetFillColor(16);
775 0 : hNrpZ->Draw();
776 :
777 0 : cpvCanvas->Write(0,kOverwrite);
778 :
779 0 : }
780 :
781 : void AliPHOSIhepAnalyze::HitsCPV(Int_t nev)
782 : {
783 : // Cumulative list of charged CPV impacts in event nev.
784 :
785 0 : AliPHOSLoader* please = dynamic_cast<AliPHOSLoader*>(fRunLoader->GetLoader("PHOSLoader"));
786 0 : if ( please == 0 )
787 : {
788 0 : AliError(Form("Could not obtain the Loader object !"));
789 0 : return ;
790 : }
791 :
792 :
793 0 : printf("\n=================== Event %10d ===================\n",nev);
794 : //16.03.2011: DP. Code below seems to be obsollete
795 : //Comment it to sutisfy Coverity
796 : /*
797 : TList * fCpvImpacts ;
798 : TBranch * branchCPVimpacts;
799 :
800 : AliPHOSGeometry * phosgeom = AliPHOSGeometry::GetInstance() ;
801 :
802 : fRunLoader->GetEvent(nev);
803 : Int_t ntracks = fRunLoader->GetHeader()->GetNtrack();
804 :
805 : // Int_t nRecCPV = 0; // Reconstructed points in event // 01.10.2001
806 : // Int_t nGenCPV = 0; // Impacts in event
807 :
808 : // Get branch of CPV impacts
809 : TTree* treeH = please->TreeH();
810 : if (treeH == 0x0)
811 : {
812 : AliError(Form("Can not get TreeH"));
813 : return;
814 : }
815 :
816 : if (! (branchCPVimpacts =treeH->GetBranch("PHOSCpvImpacts")) ) return;
817 :
818 : // Create and fill arrays of hits for each CPV module
819 : Int_t nOfModules = phosgeom->GetNModules();
820 : TClonesArray **hitsPerModule = new TClonesArray *[nOfModules];
821 : Int_t iModule = 0;
822 : for (iModule=0; iModule < nOfModules; iModule++)
823 : hitsPerModule[iModule] = new TClonesArray("AliPHOSImpact",100);
824 :
825 : TClonesArray *impacts;
826 : AliPHOSImpact *impact;
827 :
828 : for (Int_t itrack=0; itrack<ntracks; itrack++) {
829 : branchCPVimpacts ->SetAddress(&fCpvImpacts);
830 : AliInfo(Form(" branchCPVimpacts ->SetAddress(&fCpvImpacts) OK."));
831 : branchCPVimpacts ->GetEntry(itrack,0);
832 :
833 : for (iModule=0; iModule < nOfModules; iModule++) {
834 : impacts = (TClonesArray *)fCpvImpacts->At(iModule);
835 : // Do loop over impacts in the module
836 : for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
837 : impact=(AliPHOSImpact*)impacts->At(iImpact);
838 : TClonesArray &lhits = *(TClonesArray *)hitsPerModule[iModule];
839 : new(lhits[hitsPerModule[iModule]->GetEntriesFast()]) AliPHOSImpact(*impact);
840 : }
841 : }
842 : fCpvImpacts->Clear();
843 : }
844 :
845 : for (iModule=0; iModule < nOfModules; iModule++) {
846 : Int_t nsum = hitsPerModule[iModule]->GetEntriesFast();
847 : printf("CPV module %d has %d hits\n",iModule,nsum);
848 : }
849 :
850 : */
851 :
852 :
853 : // TList * fCpvImpacts ;
854 : // TBranch * branchCPVimpacts;
855 : // AliPHOSImpact* impact;
856 : // TClonesArray* impacts;
857 :
858 : // AliPHOSLoader * please = AliPHOSLoader::GetInstance(GetFileName().Data(),"PHOS");
859 : // const AliPHOSGeometry * fGeom = please->PHOSGeometry();
860 :
861 : // Int_t ntracks = gAlice->GetEvent(ievent);
862 : // Int_t nOfModules = fGeom->GetNModules();
863 : // AliInfo(Form(" Tracks: "<<ntracks<<" Modules: "<<nOfModules));
864 :
865 : // if (! (branchCPVimpacts =gAlice->TreeH()->GetBranch("PHOSCpvImpacts")) ) return;
866 :
867 : // for (Int_t itrack=0; itrack<ntracks; itrack++) {
868 : // branchCPVimpacts ->SetAddress(&fCpvImpacts);
869 : // Info("AnalyzeCPV1", " branchCPVimpacts ->SetAddress(&fCpvImpacts) OK.");
870 : // branchCPVimpacts ->GetEntry(itrack,0);
871 : // Info(Form(" branchCPVimpacts ->GetEntry(itrack,0) OK."));
872 :
873 : // for (Int_t iModule=0; iModule < nOfModules; iModule++) {
874 : // impacts = (TClonesArray *)fCpvImpacts->At(iModule);
875 : // Info(Form(" fCpvImpacts->At(iModule) OK."));
876 : // // Do loop over impacts in the module
877 : // for (Int_t iImpact=0; iImpact<impacts->GetEntries(); iImpact++) {
878 : // impact=(AliPHOSImpact*)impacts->At(iImpact);
879 : // impact->Print();
880 : // if(IsCharged(impact->GetPid()))
881 : // {
882 : // Info(Form(" Add charged hit.."));
883 : // new(hits[hits.GetEntriesFast()]) AliPHOSImpact(*impact);
884 : // Info(Form("done."));
885 : // }
886 : // }
887 : // }
888 : // fCpvImpacts->Clear();
889 : // }
890 :
891 : // Info(Form(" PHOS event "<<ievent<<": "<<hits.GetEntries()<<" charged CPV hits."));
892 :
893 0 : }
894 :
895 :
896 : // void AliPHOSIhepAnalyze::ChargedHitsCPV(TClonesArray* hits, Int_t ievent, Int_t iModule)
897 : // {
898 : // // Cumulative list of charged CPV hits in event ievent
899 : // // in PHOS module iModule.
900 :
901 : // HitsCPV(hits,ievent,iModule);
902 : // TIter next(hits);
903 :
904 : // while(AliPHOSCPVHit* cpvhit = (AliPHOSCPVHit*)next())
905 : // {
906 : // if(!IsCharged(cpvhit->GetIpart()))
907 : // {
908 : // hits->Remove(cpvhit);
909 : // delete cpvhit;
910 : // hits->Compress();
911 : // }
912 : // }
913 :
914 : // Info(Form(" PHOS module "<<iModule<<": "<<hits->GetEntries()<<" charged CPV hits."));
915 : // }
916 :
917 : Bool_t AliPHOSIhepAnalyze::IsCharged(Int_t pdgCode)
918 : {
919 : // For HIJING
920 0 : AliInfo(Form("pdgCode %d", pdgCode));
921 0 : if(pdgCode==211 || pdgCode==-211 || pdgCode==321 || pdgCode==-321 || pdgCode==11 || pdgCode==-11 || pdgCode==2212 || pdgCode==-2212) return kTRUE;
922 : else
923 0 : return kFALSE;
924 0 : }
925 : //---------------------------------------------------------------------------
926 :
927 :
928 :
929 :
930 :
931 :
|
