]>
Commit | Line | Data |
---|---|---|
37c55dc0 | 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 | // This is a TTask that constructs ReconstParticles (reconstructed particles) | |
18 | // out of Digits | |
19 | // | |
20 | //-- Authors: Evgeny Karpechev(INR) and Alla Maevsksia | |
21 | ////////////////////////////////////////////////////////////////////////////// | |
22 | ||
dc8af42e | 23 | // --- ROOT system --- |
24 | #include "TTask.h" | |
25 | #include "TTree.h" | |
26 | #include "TSystem.h" | |
27 | #include "TFile.h" | |
37c55dc0 | 28 | #include "TROOT.h" |
29 | #include "TFolder.h" | |
30 | ||
dc8af42e | 31 | // --- Standard library --- |
37c55dc0 | 32 | #include <stdlib.h> |
e75fab8f | 33 | #include "Riostream.h" |
dc8af42e | 34 | |
35 | // --- AliRoot header files --- | |
36 | ||
37 | #include "AliFMDdigit.h" | |
1a42d4f2 | 38 | #include "AliFMDhit.h" |
dc8af42e | 39 | #include "AliFMDReconstParticles.h" |
40 | #include "AliFMD.h" | |
41 | #include "AliFMDv1.h" | |
42 | #include "AliFMDReconstruction.h" | |
43 | #include "AliRun.h" | |
1a42d4f2 | 44 | #include "AliHeader.h" |
45 | #include "AliGenEventHeader.h" | |
dc8af42e | 46 | ClassImp(AliFMDReconstruction) |
47 | ||
48 | ||
49 | //____________________________________________________________________________ | |
50 | ||
51 | AliFMDReconstruction::AliFMDReconstruction():TTask("AliFMDReconstruction","") | |
52 | { | |
53 | fNevents = 0 ; // Number of events to rreconnstraction, 0 means all events in current file | |
54 | // add Task to //root/Tasks folder | |
55 | TTask * roottasks = (TTask*)gROOT->GetRootFolder()->FindObject("Tasks") ; | |
56 | roottasks->Add(this) ; | |
57 | } | |
58 | //____________________________________________________________________________ | |
59 | ||
60 | AliFMDReconstruction::AliFMDReconstruction(char* HeaderFile, char *ReconstParticlesFile):TTask("AliFMDReconstruction","") | |
61 | { | |
62 | fNevents = 0 ; // Number of events to rreconnstraction, 0 means all events in current file | |
63 | fReconstParticlesFile=ReconstParticlesFile ; | |
64 | fHeadersFile=HeaderFile ; | |
65 | //add Task to //root/Tasks folder | |
66 | TTask * roottasks = (TTask*)gROOT->GetRootFolder()->FindObject("Tasks") ; | |
67 | roottasks->Add(this) ; | |
68 | } | |
69 | ||
70 | //____________________________________________________________________________ | |
71 | ||
72 | AliFMDReconstruction::~AliFMDReconstruction() | |
73 | { | |
dc8af42e | 74 | } |
cb1df35e | 75 | //---------------------------------------------------------------------------- |
dc8af42e | 76 | |
37c55dc0 | 77 | void AliFMDReconstruction::Exec(Option_t *option) |
dc8af42e | 78 | { |
a39f3926 | 79 | /* |
80 | Reconstruct nember of particles | |
81 | in given group of pads for given FMDvolume | |
82 | determine by numberOfVolume , | |
83 | numberOfMinSector,numberOfMaxSector, | |
84 | numberOfMinRing, numberOgMaxRing | |
85 | Reconstruction method choose dependence on number of empty pads | |
86 | */ | |
87 | ||
88 | ||
1a42d4f2 | 89 | cout<<"\nStart AliFMDReconstruction::Exec(...)"<<endl; |
a39f3926 | 90 | Int_t const knumVolumes=5; |
91 | Int_t const knumSectors=40; | |
92 | Int_t const knumRings=768; | |
93 | Int_t padADC[10][50][800]; | |
1a42d4f2 | 94 | Float_t eta, etain,etaout,rad,theta; |
cb1df35e | 95 | Int_t ivol, iSector, iRing; |
1a42d4f2 | 96 | Float_t rin[5]={4.2,15.4,4.2,15.4,4.2}; |
97 | Float_t rout[5]={17.4,28.4,17.4,28.4,17.4}; | |
98 | Float_t z[5]={62.8, 75.2, -83.4, -75.2, -340.}; | |
a39f3926 | 99 | Int_t numberOfRings[5]={768,384,768,384,768}; |
100 | Int_t numberOfSectors[5]= {20,40,20,40,20}; | |
101 | Int_t numberOfEtaIntervals[5]; | |
cb1df35e | 102 | // number of ring for boundary 0.1 eta |
1a42d4f2 | 103 | TBranch *brDigits=0; |
cb1df35e | 104 | |
a39f3926 | 105 | AliFMD * fFMD = (AliFMD *) gAlice->GetDetector("FMD"); |
106 | TClonesArray *fReconParticles=fFMD->ReconParticles(); | |
cb1df35e | 107 | if(fNevents == 0) fNevents=(Int_t)gAlice->TreeE()->GetEntries(); |
1a42d4f2 | 108 | cout<<" fNevents "<<fNevents<<endl; |
dc8af42e | 109 | for(Int_t ievent=0;ievent<fNevents;ievent++) |
110 | { | |
1a42d4f2 | 111 | cout<<" ievent "<<ievent<<endl; |
a39f3926 | 112 | for (Int_t i=0; i<knumVolumes; i++) |
113 | for(Int_t j=0; j<knumSectors; j++) | |
114 | for(Int_t ij=0; ij<knumRings; ij++) | |
115 | padADC[i][j][ij]=0; //zhachem ??? | |
dc8af42e | 116 | gAlice->GetEvent(ievent) ; |
117 | if(gAlice->TreeH()==0) return; | |
1a42d4f2 | 118 | if(gAlice->TreeD()==0) return; |
119 | brDigits=gAlice->TreeD()->GetBranch("FMD"); | |
120 | if (!brDigits){ | |
121 | cerr<<"EXEC Branch FMD digits not found"<<endl; | |
122 | return; | |
123 | } | |
124 | ||
dc8af42e | 125 | if(gAlice->TreeR()==0) gAlice->MakeTree("R"); |
126 | //Make branches | |
a39f3926 | 127 | fFMD->MakeBranch("R"); |
1a42d4f2 | 128 | |
37c55dc0 | 129 | |
cb1df35e | 130 | Int_t zeroADC=1; |
1a42d4f2 | 131 | |
132 | AliFMDdigit *fmdDigit; | |
a39f3926 | 133 | if (fFMD) |
dc8af42e | 134 | { |
135 | gAlice->TreeD()->GetEvent(0); | |
a39f3926 | 136 | TClonesArray *fFMDdigits=fFMD->Digits(); |
137 | Int_t nDigits=fFMDdigits->GetEntries(); | |
1a42d4f2 | 138 | cout<<" nDigits "<<nDigits<<endl; |
a39f3926 | 139 | Int_t recParticles[6]; |
dc8af42e | 140 | Int_t nRecPart=0 ; |
a39f3926 | 141 | Int_t zeroPads=0; |
142 | Int_t numberOfPads=0; //To avoid warning | |
cb1df35e | 143 | Int_t pedestal; |
144 | Float_t channelWidth=(22400*50)/1024; | |
dc8af42e | 145 | for (Int_t digit=0;digit<nDigits;digit++) |
146 | { | |
a39f3926 | 147 | fmdDigit=(AliFMDdigit*)fFMDdigits->UncheckedAt(digit); |
cb1df35e | 148 | ivol=fmdDigit->Volume(); |
149 | iSector=fmdDigit->NumberOfSector(); | |
150 | iRing=fmdDigit->NumberOfRing(); | |
151 | pedestal=Int_t(gRandom->Gaus(500,250)); | |
a39f3926 | 152 | padADC[ivol-1][iSector-1][iRing-1]= |
cb1df35e | 153 | fmdDigit->ADCsignal() |
154 | -Int_t(Float_t(pedestal)/channelWidth); | |
a39f3926 | 155 | if (padADC[ivol-1][iSector-1][iRing-1]<0) |
156 | padADC[ivol-1][iSector-1][iRing-1]=0; | |
dc8af42e | 157 | } //digit loop |
a39f3926 | 158 | Int_t rmin=0; Int_t rmax=0; //To avoid warning |
159 | Int_t smin=0; Int_t smax=0; //To avoid warning | |
1a42d4f2 | 160 | AliHeader *header = gAlice->GetHeader(); |
161 | AliGenEventHeader* genHeader = header->GenEventHeader(); | |
162 | TArrayF *o = new TArrayF(3); | |
163 | genHeader->PrimaryVertex(*o); | |
164 | Float_t zVertex=o->At(2); | |
a39f3926 | 165 | for (ivol=0; ivol<knumVolumes; ivol++) |
cb1df35e | 166 | { |
1a42d4f2 | 167 | Float_t realZ=zVertex+z[ivol]; |
168 | theta=TMath::ATan(rin[ivol]/TMath::Abs(realZ)); | |
169 | etain = - TMath::Log( TMath::Tan(theta/2.)); | |
170 | theta=TMath::ATan(rout[ivol]/TMath::Abs(realZ)); | |
171 | etaout=- TMath::Log( TMath::Tan(theta/2.)); | |
a39f3926 | 172 | numberOfEtaIntervals[ivol]=Int_t ((etain-etaout)*10); |
1a42d4f2 | 173 | eta=etain; |
a39f3926 | 174 | for (Int_t e1=0;e1<=numberOfEtaIntervals[ivol];e1++) |
cb1df35e | 175 | { |
1a42d4f2 | 176 | theta = 2.*TMath::ATan (TMath::Exp (-eta)); |
177 | Float_t radmin = TMath::Abs(realZ) * (TMath::Tan (theta)); | |
a39f3926 | 178 | rmin= Int_t ( (radmin-rin[ivol])*numberOfRings[ivol]/(rout[ivol]-rin[ivol])); |
1a42d4f2 | 179 | eta=eta-0.1; |
180 | theta = 2.*TMath::ATan (TMath::Exp (-eta)); | |
181 | rad = TMath::Abs(realZ) * (TMath::Tan (theta)); | |
a39f3926 | 182 | rmax=Int_t( (rad-rin[ivol])*numberOfRings[ivol]/(rout[ivol]-rin[ivol])); |
183 | zeroPads=0; | |
184 | smin=0; | |
185 | smax=numberOfSectors[ivol]; | |
186 | for (Int_t iring=rmin; iring<rmax; iring++) | |
cb1df35e | 187 | { |
a39f3926 | 188 | numberOfPads=(rmax-rmin)*(smax-smin); |
cb1df35e | 189 | for |
190 | (Int_t isector=1; | |
a39f3926 | 191 | isector<=numberOfSectors[ivol]; |
cb1df35e | 192 | isector++) |
a39f3926 | 193 | if(padADC[ivol][isector-1][iring-1]<zeroADC) |
194 | zeroPads++; | |
cb1df35e | 195 | } //ring |
cb1df35e | 196 | Float_t zerosRatio= |
a39f3926 | 197 | (Float_t)zeroPads/(Float_t)numberOfPads; |
198 | recParticles[0]=ivol; | |
199 | recParticles[1]=smin; | |
200 | recParticles[2]=smax; | |
201 | recParticles[3]=rmin; | |
202 | recParticles[4]=rmax; | |
203 | if (zerosRatio>0.1 ) | |
204 | recParticles[5]= | |
205 | DeterminationByPoisson | |
206 | (padADC,ivol+1,rmin,rmax,smin,smax); | |
cb1df35e | 207 | else |
a39f3926 | 208 | recParticles[5]= |
209 | DeterminationByThresholds | |
210 | (padADC,ivol+1,rmin,rmax,smin,smax); | |
cb1df35e | 211 | new((*fReconParticles)[nRecPart++]) |
a39f3926 | 212 | AliFMDReconstParticles(recParticles); |
1a42d4f2 | 213 | } // eta |
214 | } // volume | |
cb1df35e | 215 | |
cb1df35e | 216 | }//if FMD |
1a42d4f2 | 217 | gAlice->TreeR()->Reset(); |
218 | gAlice->TreeR()->Fill(); | |
219 | gAlice->TreeR()->Write(0,TObject::kOverwrite); | |
dc8af42e | 220 | } //event loop |
1a42d4f2 | 221 | cout<<"\nAliFMDReconstruction::Exec finished"<<endl; |
cb1df35e | 222 | }; |
dc8af42e | 223 | |
cb1df35e | 224 | //------------------------------------------------------------ |
a39f3926 | 225 | Int_t AliFMDReconstruction::DeterminationByThresholds |
226 | (Int_t padAdc[][50][800],Int_t volume, Int_t rmin, Int_t rmax, | |
227 | Int_t smin, Int_t smax) | |
dc8af42e | 228 | { |
a39f3926 | 229 | /* |
230 | reconstruct number of particles by | |
231 | energy deposited threshold method | |
232 | Using if number of empty pads less then 10% | |
233 | ||
234 | */ | |
1a42d4f2 | 235 | cout<<"\nStart threshold method\n"; |
cb1df35e | 236 | |
237 | Int_t thresholdInner[30]={ | |
238 | 0, 14, 28, 42, 57, | |
239 | 72, 89, 104, 124, 129, | |
240 | 164, 174, 179, 208, 228, | |
241 | 248, 268, 317, 337, 357, | |
242 | 392, 407, 416, 426, 436, | |
243 | 461, 468, 493, 506, 515}; | |
244 | ||
245 | Int_t thresholdOuter[30]={0, 18, 48, 77, 105, | |
246 | 132, 165, 198, 231, | |
247 | 264, 286, 308, 334, | |
248 | 352, 374, 418, 440, | |
249 | 462, 484, 506, 528, | |
250 | 550, 572, 594, 616}; | |
251 | Int_t threshold[30]; | |
252 | for (Int_t it=0; it<30; it++) { | |
253 | if(volume==1||volume==3||volume==5) threshold[it]=thresholdInner[it]; | |
254 | if(volume==2||volume==4) threshold[it]=thresholdOuter[it]; | |
255 | } | |
a39f3926 | 256 | Int_t thresholdArraySize = 30; |
257 | Int_t numPart=0; | |
cb1df35e | 258 | //Inner Si |
a39f3926 | 259 | for (Int_t iring=rmin; iring<rmax; iring++) |
cb1df35e | 260 | { |
a39f3926 | 261 | for (Int_t isector=smin; isector<smax; isector++) |
cb1df35e | 262 | { |
a39f3926 | 263 | for (int i=0;i<thresholdArraySize-1;i++) |
cb1df35e | 264 | { |
a39f3926 | 265 | if(padAdc[volume-1][isector][iring]>threshold[i] |
266 | &&padAdc[volume-1][isector][iring]<=threshold[i+1]) | |
cb1df35e | 267 | { |
a39f3926 | 268 | numPart+=i; |
cb1df35e | 269 | break; |
270 | }; // if in threshol interval | |
271 | }; //threshold_array_size | |
272 | }; //iring | |
273 | }; //sector | |
1a42d4f2 | 274 | cout<<"\nEnd threshol method"<<endl; |
a39f3926 | 275 | return numPart; |
dc8af42e | 276 | } |
cb1df35e | 277 | |
278 | ||
279 | //__________________________________________________________________ | |
280 | ||
a39f3926 | 281 | Int_t AliFMDReconstruction::DeterminationByPoisson |
282 | (Int_t padAdc[][50][800], | |
283 | Int_t vol, Int_t rmin, Int_t rmax, | |
284 | Int_t secmin, Int_t secmax) | |
dc8af42e | 285 | { |
a39f3926 | 286 | /* |
287 | reconstruct number of particles by Poisson statistic method | |
288 | according number of empty pad in chosen group of pads | |
289 | using if number of empty pads more then 10% | |
290 | ||
291 | */ | |
292 | ||
1a42d4f2 | 293 | // cout<<"\nStart Poisson method"; |
a39f3926 | 294 | Int_t thresholdADC=1; |
cb1df35e | 295 | Int_t zeropads=0; |
296 | for (Int_t i=rmin;i<rmax;i++) | |
297 | { | |
298 | for (Int_t j=secmin;j<secmax;j++) | |
299 | { | |
a39f3926 | 300 | if (padAdc[vol-1][j][i]<thresholdADC) zeropads++; |
cb1df35e | 301 | }; |
302 | }; | |
303 | Float_t lambda=-TMath::Log(Float_t(zeropads)/ | |
304 | ( Float_t(secmax-secmin)* | |
305 | Float_t(rmax-rmin))); //+1 zdes' ne nado | |
a39f3926 | 306 | Int_t fRecon=(Int_t)(lambda*(secmax-secmin)*(rmax-rmin)); //+1 zdes' ne nado |
cb1df35e | 307 | // cout<<"\nEnd Poisson method"<<endl; |
a39f3926 | 308 | return fRecon; |
cb1df35e | 309 | }; |
310 | ||
311 | //------------------------------------------------------------------ | |
312 | ||
313 | ||
314 | ||
315 | ||
316 | ||
dc8af42e | 317 | |
318 | ||
319 | ||
37c55dc0 | 320 | |
321 | ||
322 | ||
323 | ||
324 | ||
325 | ||
326 | ||
327 | ||
328 | ||
329 | ||
330 | ||
331 | ||
332 | ||
333 | ||
334 | ||
335 | ||
336 | ||
337 | ||
1a42d4f2 | 338 | |
dc8af42e | 339 | |
340 |