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