]>
Commit | Line | Data |
---|---|---|
ddae0931 | 1 | //////////////////////////////////////////////// |
2 | // Manager and hits classes for set:RICH // | |
3 | //////////////////////////////////////////////// | |
fe4da5cc | 4 | |
ddae0931 | 5 | #include <TBRIK.h> |
6 | #include <TTUBE.h> | |
fe4da5cc | 7 | #include <TNode.h> |
8 | #include <TRandom.h> | |
ddae0931 | 9 | #include <TObject.h> |
10 | #include <TVector.h> | |
11 | #include <TObjArray.h> | |
fe4da5cc | 12 | |
fe4da5cc | 13 | #include "AliRICH.h" |
ddae0931 | 14 | #include "AliRICHHitMap.h" |
fe4da5cc | 15 | #include "AliRun.h" |
ddae0931 | 16 | #include "AliMC.h" |
17 | #include "iostream.h" | |
18 | #include "AliCallf77.h" | |
19 | ||
20 | // Static variables for the pad-hit iterator routines | |
21 | static Int_t sMaxIterPad=0; | |
22 | static Int_t sCurIterPad=0; | |
23 | static TClonesArray *fClusters2; | |
24 | static TClonesArray *fHits2; | |
25 | ||
fe4da5cc | 26 | ClassImp(AliRICH) |
ddae0931 | 27 | |
28 | //___________________________________________ | |
fe4da5cc | 29 | AliRICH::AliRICH() |
30 | { | |
ddae0931 | 31 | fIshunt = 0; |
32 | fHits = 0; | |
33 | fClusters = 0; | |
34 | fNclusters = 0; | |
35 | fDchambers = 0; | |
36 | fRecClusters= 0; | |
37 | fCerenkovs = 0; | |
38 | fNdch = 0; | |
fe4da5cc | 39 | } |
ddae0931 | 40 | |
41 | //___________________________________________ | |
fe4da5cc | 42 | AliRICH::AliRICH(const char *name, const char *title) |
ddae0931 | 43 | : AliDetector(name,title) |
fe4da5cc | 44 | { |
ddae0931 | 45 | //Begin_Html |
46 | /* | |
47 | <img src="gif/alirich.gif"> | |
48 | */ | |
49 | //End_Html | |
50 | ||
51 | fHits = new TClonesArray("AliRICHhit",1000 ); | |
52 | fClusters = new TClonesArray("AliRICHcluster",10000); | |
53 | fCerenkovs = new TClonesArray("AliRICHCerenkov",20000); | |
54 | fNclusters = 0; | |
55 | fIshunt = 0; | |
56 | ||
57 | fNdch = new Int_t[7]; | |
58 | ||
59 | fDchambers = new TObjArray(7); | |
60 | ||
61 | Int_t i; | |
62 | ||
63 | for (i=0; i<7 ;i++) { | |
64 | (*fDchambers)[i] = new TClonesArray("AliRICHdigit",10000); | |
65 | fNdch[i]=0; | |
66 | } | |
67 | ||
68 | fRecClusters=new TObjArray(7); | |
69 | for (i=0; i<7;i++) | |
70 | (*fRecClusters)[i] = new TObjArray(1000); | |
71 | ||
72 | // | |
73 | // Transport angular cut | |
74 | fAccCut=0; | |
75 | fAccMin=2; | |
76 | fAccMax=9; | |
77 | ||
78 | SetMarkerColor(kRed); | |
fe4da5cc | 79 | } |
80 | ||
ddae0931 | 81 | //___________________________________________ |
fe4da5cc | 82 | AliRICH::~AliRICH() |
83 | { | |
ddae0931 | 84 | fIshunt = 0; |
85 | delete fHits; | |
86 | delete fClusters; | |
87 | delete fCerenkovs; | |
88 | for (Int_t i=0;i<7;i++) | |
89 | delete (*fRecClusters)[i]; | |
90 | delete fRecClusters; | |
91 | ||
fe4da5cc | 92 | } |
93 | ||
ddae0931 | 94 | //___________________________________________ |
fe4da5cc | 95 | void AliRICH::AddHit(Int_t track, Int_t *vol, Float_t *hits) |
96 | { | |
ddae0931 | 97 | TClonesArray &lhits = *fHits; |
98 | new(lhits[fNhits++]) AliRICHhit(fIshunt,track,vol,hits); | |
fe4da5cc | 99 | } |
fe4da5cc | 100 | //_____________________________________________________________________________ |
ddae0931 | 101 | void AliRICH::AddCerenkov(Int_t track, Int_t *vol, Float_t *cerenkovs) |
102 | { | |
103 | TClonesArray &lcerenkovs = *fCerenkovs; | |
104 | new(lcerenkovs[fNcerenkovs++]) AliRICHCerenkov(fIshunt,track,vol,cerenkovs); | |
fe4da5cc | 105 | } |
ddae0931 | 106 | //___________________________________________ |
107 | void AliRICH::AddCluster(Int_t *clhits) | |
108 | { | |
109 | TClonesArray &lclusters = *fClusters; | |
110 | new(lclusters[fNclusters++]) AliRICHcluster(clhits); | |
fe4da5cc | 111 | } |
fe4da5cc | 112 | //_____________________________________________________________________________ |
ddae0931 | 113 | void AliRICH::AddDigits(Int_t id, Int_t *tracks, Int_t *charges, Int_t *digits) |
114 | { | |
115 | // | |
116 | // Add a RICH digit to the list | |
117 | // | |
118 | ||
119 | TClonesArray &ldigits = *((TClonesArray*)(*fDchambers)[id]); | |
120 | new(ldigits[fNdch[id]++]) AliRICHdigit(tracks,charges,digits); | |
fe4da5cc | 121 | } |
122 | ||
ddae0931 | 123 | |
fe4da5cc | 124 | //_____________________________________________________________________________ |
ddae0931 | 125 | void AliRICH::AddRecCluster(Int_t iCh, Int_t iCat, AliRICHRecCluster* Cluster) |
fe4da5cc | 126 | { |
ddae0931 | 127 | // |
128 | // Add a RICH reconstructed cluster to the list | |
129 | // | |
130 | TObjArray* ClusterList = RecClusters(iCh,iCat); | |
131 | ClusterList->Add(Cluster); | |
fe4da5cc | 132 | } |
133 | ||
fe4da5cc | 134 | |
135 | ||
ddae0931 | 136 | //___________________________________________ |
137 | void AliRICH::BuildGeometry() | |
138 | ||
fe4da5cc | 139 | { |
ddae0931 | 140 | // |
141 | // Builds a TNode geometry for event display | |
142 | // | |
143 | TNode *Node, *Top; | |
144 | ||
145 | const int kColorRICH = kGreen; | |
146 | // | |
147 | Top=gAlice->GetGeometry()->GetNode("alice"); | |
148 | ||
149 | ||
150 | new TBRIK("S_RICH","S_RICH","void",71.09999,11.5,73.15); | |
151 | ||
152 | Top->cd(); | |
153 | Float_t pos1[3]={0,471.8999,165.2599}; | |
154 | Chamber(0).SetChamberTransform(pos1[0],pos1[1],pos1[2],new TRotMatrix("rot993","rot993",90,0,70.69,90,19.30999,-90)); | |
155 | Node = new TNode("RICH1","RICH1","S_RICH",pos1[0],pos1[1],pos1[2],"rot993"); | |
156 | ||
157 | ||
158 | Node->SetLineColor(kColorRICH); | |
159 | fNodes->Add(Node); | |
160 | Top->cd(); | |
161 | ||
162 | Float_t pos2[3]={171,470,0}; | |
163 | Chamber(1).SetChamberTransform(pos2[0],pos2[1],pos2[2],new TRotMatrix("rot994","rot994",90,-20,90,70,0,0)); | |
164 | Node = new TNode("RICH2","RICH2","S_RICH",pos2[0],pos2[1],pos2[2],"rot994"); | |
165 | ||
166 | ||
167 | Node->SetLineColor(kColorRICH); | |
168 | fNodes->Add(Node); | |
169 | Top->cd(); | |
170 | Float_t pos3[3]={0,500,0}; | |
171 | Chamber(2).SetChamberTransform(pos3[0],pos3[1],pos3[2],new TRotMatrix("rot995","rot995",90,0,90,90,0,0)); | |
172 | Node = new TNode("RICH3","RICH3","S_RICH",pos3[0],pos3[1],pos3[2],"rot995"); | |
173 | ||
174 | ||
175 | Node->SetLineColor(kColorRICH); | |
176 | fNodes->Add(Node); | |
177 | Top->cd(); | |
178 | Float_t pos4[3]={-171,470,0}; | |
179 | Chamber(3).SetChamberTransform(pos4[0],pos4[1],pos4[2], new TRotMatrix("rot996","rot996",90,20,90,110,0,0)); | |
180 | Node = new TNode("RICH4","RICH4","S_RICH",pos4[0],pos4[1],pos4[2],"rot996"); | |
181 | ||
182 | ||
183 | Node->SetLineColor(kColorRICH); | |
184 | fNodes->Add(Node); | |
185 | Top->cd(); | |
186 | Float_t pos5[3]={161.3999,443.3999,-165.3}; | |
187 | Chamber(4).SetChamberTransform(pos5[0],pos5[1],pos5[2],new TRotMatrix("rot997","rot997",90,340,108.1999,70,18.2,70)); | |
188 | Node = new TNode("RICH5","RICH5","S_RICH",pos5[0],pos5[1],pos5[2],"rot997"); | |
189 | ||
190 | Node->SetLineColor(kColorRICH); | |
191 | fNodes->Add(Node); | |
192 | Top->cd(); | |
193 | Float_t pos6[3]={0., 471.9, -165.3,}; | |
194 | Chamber(5).SetChamberTransform(pos6[0],pos6[1],pos6[2],new TRotMatrix("rot998","rot998",90,0,109.3099,90,19.30999,90)); | |
195 | Node = new TNode("RICH6","RICH6","S_RICH",pos6[0],pos6[1],pos6[2],"rot998"); | |
196 | ||
197 | ||
198 | Node->SetLineColor(kColorRICH); | |
199 | fNodes->Add(Node); | |
200 | Top->cd(); | |
201 | Float_t pos7[3]={-161.399,443.3999,-165.3}; | |
202 | Chamber(6).SetChamberTransform(pos7[0],pos7[1],pos7[2],new TRotMatrix("rot999","rot999",90,20,108.1999,110,18.2,110)); | |
203 | Node = new TNode("RICH7","RICH7","S_RICH",pos7[0],pos7[1],pos7[2],"rot999"); | |
204 | Node->SetLineColor(kColorRICH); | |
205 | fNodes->Add(Node); | |
206 | ||
fe4da5cc | 207 | } |
208 | ||
ddae0931 | 209 | //___________________________________________ |
210 | Int_t AliRICH::DistancetoPrimitive(Int_t , Int_t ) | |
211 | { | |
212 | return 9999; | |
213 | } | |
fe4da5cc | 214 | |
ddae0931 | 215 | //___________________________________________ |
216 | void AliRICH::MakeBranch(Option_t* option) | |
fe4da5cc | 217 | { |
ddae0931 | 218 | // Create Tree branches for the RICH. |
fe4da5cc | 219 | |
ddae0931 | 220 | const Int_t buffersize = 4000; |
221 | char branchname[20]; | |
fe4da5cc | 222 | |
ddae0931 | 223 | |
224 | AliDetector::MakeBranch(option); | |
225 | sprintf(branchname,"%sCerenkov",GetName()); | |
226 | if (fCerenkovs && gAlice->TreeH()) { | |
227 | gAlice->TreeH()->Branch(branchname,&fCerenkovs, buffersize); | |
228 | printf("Making Branch %s for Cerenkov Hits\n",branchname); | |
fe4da5cc | 229 | } |
ddae0931 | 230 | |
231 | sprintf(branchname,"%sCluster",GetName()); | |
232 | if (fClusters && gAlice->TreeH()) { | |
233 | gAlice->TreeH()->Branch(branchname,&fClusters, buffersize); | |
234 | printf("Making Branch %s for clusters\n",branchname); | |
fe4da5cc | 235 | } |
236 | ||
ddae0931 | 237 | // one branch for digits per chamber |
238 | Int_t i; | |
239 | ||
240 | for (i=0; i<7 ;i++) { | |
241 | sprintf(branchname,"%sDigits%d",GetName(),i+1); | |
242 | ||
243 | if (fDchambers && gAlice->TreeD()) { | |
244 | gAlice->TreeD()->Branch(branchname,&((*fDchambers)[i]), buffersize); | |
245 | printf("Making Branch %s for digits in chamber %d\n",branchname,i+1); | |
246 | } | |
fe4da5cc | 247 | } |
ddae0931 | 248 | // one branch for rec clusters |
249 | for (i=0; i<7; i++) { | |
250 | sprintf(branchname,"%sRecClus%d",GetName(),i+1); | |
251 | if (fRecClusters && gAlice->TreeD()) { | |
252 | gAlice->TreeR() | |
253 | ->Branch(branchname,"TObjArray", | |
254 | &((*fRecClusters)[i]), buffersize,0); | |
255 | printf("Making Branch %s for clusters in chamber %d\n", | |
256 | branchname,i+1); | |
fe4da5cc | 257 | } |
ddae0931 | 258 | } |
fe4da5cc | 259 | } |
260 | ||
ddae0931 | 261 | //___________________________________________ |
262 | void AliRICH::SetTreeAddress() | |
fe4da5cc | 263 | { |
ddae0931 | 264 | // Set branch address for the Hits and Digits Tree. |
265 | char branchname[20]; | |
266 | AliDetector::SetTreeAddress(); | |
267 | ||
268 | TBranch *branch; | |
269 | TTree *treeH = gAlice->TreeH(); | |
270 | TTree *treeD = gAlice->TreeD(); | |
271 | ||
272 | if (treeH) { | |
273 | if (fClusters) { | |
274 | branch = treeH->GetBranch("RICHCluster"); | |
275 | if (branch) branch->SetAddress(&fClusters); | |
fe4da5cc | 276 | } |
ddae0931 | 277 | if (fCerenkovs) { |
278 | branch = treeH->GetBranch("RICHCerenkov"); | |
279 | if (branch) branch->SetAddress(&fCerenkovs); | |
fe4da5cc | 280 | } |
ddae0931 | 281 | |
282 | } | |
283 | ||
284 | if (treeD) { | |
285 | for (int i=0; i<7; i++) { | |
286 | sprintf(branchname,"%sDigits%d",GetName(),i+1); | |
287 | if (fDchambers) { | |
288 | branch = treeD->GetBranch(branchname); | |
289 | if (branch) branch->SetAddress(&((*fDchambers)[i])); | |
290 | } | |
fe4da5cc | 291 | } |
fe4da5cc | 292 | } |
ddae0931 | 293 | } |
294 | //___________________________________________ | |
295 | void AliRICH::ResetHits() | |
296 | { | |
297 | // Reset number of clusters and the cluster array for this detector | |
298 | AliDetector::ResetHits(); | |
299 | fNclusters = 0; | |
300 | if (fClusters) fClusters->Clear(); | |
301 | if (fCerenkovs) fCerenkovs->Clear(); | |
fe4da5cc | 302 | } |
303 | ||
ddae0931 | 304 | //____________________________________________ |
305 | void AliRICH::ResetDigits() | |
306 | { | |
307 | // | |
308 | // Reset number of digits and the digits array for this detector | |
309 | // | |
310 | for ( int i=0;i<7;i++ ) { | |
311 | if ((*fDchambers)[i]) (*fDchambers)[i]->Clear(); | |
312 | if (fNdch) fNdch[i]=0; | |
313 | } | |
314 | } | |
315 | //____________________________________________ | |
316 | void AliRICH::ResetRecClusters() | |
fe4da5cc | 317 | { |
ddae0931 | 318 | // |
319 | // Reset the rec clusters | |
320 | // | |
321 | for ( int i=0;i<7;i++ ) { | |
322 | if ((*fRecClusters)[i]) (*fRecClusters)[i]->Clear(); | |
fe4da5cc | 323 | } |
ddae0931 | 324 | } |
325 | //___________________________________________ | |
fe4da5cc | 326 | |
ddae0931 | 327 | void AliRICH::SetPADSIZ(Int_t id, Int_t isec, Float_t p1, Float_t p2) |
fe4da5cc | 328 | { |
ddae0931 | 329 | Int_t i=2*(id-1); |
330 | ((AliRICHchamber*) (*fChambers)[i]) ->SetPADSIZ(isec,p1,p2); | |
331 | ((AliRICHchamber*) (*fChambers)[i+1])->SetPADSIZ(isec,p1,p2); | |
332 | } | |
fe4da5cc | 333 | |
ddae0931 | 334 | //___________________________________________ |
335 | void AliRICH::SetMUCHSP(Int_t id, Float_t p1) | |
fe4da5cc | 336 | { |
ddae0931 | 337 | Int_t i=2*(id-1); |
338 | ((AliRICHchamber*) (*fChambers)[i])->SetMUCHSP(p1); | |
339 | ((AliRICHchamber*) (*fChambers)[i+1])->SetMUCHSP(p1); | |
340 | } | |
fe4da5cc | 341 | |
ddae0931 | 342 | //___________________________________________ |
343 | void AliRICH::SetMUSIGM(Int_t id, Float_t p1, Float_t p2) | |
344 | { | |
345 | Int_t i=2*(id-1); | |
346 | ((AliRICHchamber*) (*fChambers)[i])->SetMUSIGM(p1,p2); | |
347 | ((AliRICHchamber*) (*fChambers)[i+1])->SetMUSIGM(p1,p2); | |
348 | } | |
fe4da5cc | 349 | |
ddae0931 | 350 | //___________________________________________ |
351 | void AliRICH::SetRSIGM(Int_t id, Float_t p1) | |
352 | { | |
353 | Int_t i=2*(id-1); | |
354 | ((AliRICHchamber*) (*fChambers)[i])->SetRSIGM(p1); | |
355 | ((AliRICHchamber*) (*fChambers)[i+1])->SetRSIGM(p1); | |
356 | } | |
fe4da5cc | 357 | |
ddae0931 | 358 | //___________________________________________ |
359 | void AliRICH::SetMAXADC(Int_t id, Float_t p1) | |
fe4da5cc | 360 | { |
ddae0931 | 361 | Int_t i=2*(id-1); |
362 | ((AliRICHchamber*) (*fChambers)[i])->SetMAXADC(p1); | |
363 | ((AliRICHchamber*) (*fChambers)[i+1])->SetMAXADC(p1); | |
fe4da5cc | 364 | } |
365 | ||
ddae0931 | 366 | //___________________________________________ |
367 | void AliRICH::SetSMAXAR(Float_t p1) | |
fe4da5cc | 368 | { |
ddae0931 | 369 | fMaxStepGas=p1; |
fe4da5cc | 370 | } |
371 | ||
ddae0931 | 372 | //___________________________________________ |
373 | void AliRICH::SetSMAXAL(Float_t p1) | |
374 | { | |
375 | fMaxStepAlu=p1; | |
376 | } | |
fe4da5cc | 377 | |
ddae0931 | 378 | //___________________________________________ |
379 | void AliRICH::SetDMAXAR(Float_t p1) | |
fe4da5cc | 380 | { |
ddae0931 | 381 | fMaxDestepGas=p1; |
fe4da5cc | 382 | } |
383 | ||
ddae0931 | 384 | //___________________________________________ |
385 | void AliRICH::SetDMAXAL(Float_t p1) | |
386 | { | |
387 | fMaxDestepAlu=p1; | |
388 | } | |
389 | //___________________________________________ | |
390 | void AliRICH::SetRICHACC(Bool_t acc, Float_t angmin, Float_t angmax) | |
391 | { | |
392 | fAccCut=acc; | |
393 | fAccMin=angmin; | |
394 | fAccMax=angmax; | |
395 | } | |
396 | //___________________________________________ | |
397 | void AliRICH::SetSegmentationModel(Int_t id, Int_t isec, AliRICHsegmentation *segmentation) | |
398 | { | |
399 | ((AliRICHchamber*) (*fChambers)[id])->SegmentationModel(isec, segmentation); | |
fe4da5cc | 400 | |
ddae0931 | 401 | } |
402 | //___________________________________________ | |
403 | void AliRICH::SetResponseModel(Int_t id, Response_t res, AliRICHresponse *response) | |
fe4da5cc | 404 | { |
ddae0931 | 405 | ((AliRICHchamber*) (*fChambers)[id])->ResponseModel(res, response); |
fe4da5cc | 406 | } |
407 | ||
ddae0931 | 408 | void AliRICH::SetNsec(Int_t id, Int_t nsec) |
409 | { | |
410 | ((AliRICHchamber*) (*fChambers)[id])->SetNsec(nsec); | |
411 | } | |
412 | ||
413 | ||
414 | //___________________________________________ | |
415 | ||
416 | void AliRICH::StepManager() | |
fe4da5cc | 417 | { |
ddae0931 | 418 | printf("Dummy version of RICH step -- it should never happen!!\n"); |
419 | const Float_t kRaddeg = 180/TMath::Pi(); | |
420 | AliMC* pMC = AliMC::GetMC(); | |
421 | Int_t nsec, ipart; | |
422 | Float_t x[4], p[4]; | |
423 | Float_t pt, th0, th1; | |
424 | char proc[5]; | |
425 | if(fAccCut) { | |
426 | if((nsec=pMC->NSecondaries())>0) { | |
427 | pMC->ProdProcess(proc); | |
428 | if((pMC->TrackPid()==113 || pMC->TrackPid()==114) && !strcmp(proc,"DCAY")) { | |
429 | ||
430 | // Check angular acceptance | |
431 | //* --- and have muons from resonance decays in the wanted window --- | |
432 | if(nsec != 2) { | |
433 | printf(" AliRICH::StepManager: Strange resonance Decay into %d particles\n",nsec); | |
434 | pMC->StopEvent(); | |
435 | } else { | |
436 | pMC->GetSecondary(0,ipart,x,p); | |
437 | pt = TMath::Sqrt(p[0]*p[0]+p[1]*p[1]); | |
438 | th0 = TMath::ATan2(pt,p[2])*kRaddeg; | |
439 | pMC->GetSecondary(1,ipart,x,p); | |
440 | pt = TMath::Sqrt(p[0]*p[0]+p[1]*p[1]); | |
441 | th1 = TMath::ATan2(pt,p[2])*kRaddeg; | |
442 | if(!(fAccMin < th0 && th0 < fAccMax) || | |
443 | !(fAccMin < th1 && th1 < fAccMax)) | |
444 | pMC->StopEvent(); | |
445 | } | |
446 | } | |
447 | } | |
fe4da5cc | 448 | } |
ddae0931 | 449 | } |
450 | void AliRICH::ReconstructClusters() | |
451 | { | |
452 | // | |
453 | // Initialize the necessary correspondance table | |
454 | // | |
f91473f6 | 455 | static const Int_t kMaxNPadx = 600; |
456 | static const Int_t kMaxNPady = 600; | |
457 | Int_t elem[kMaxNPadx*2][kMaxNPady*2]; | |
ddae0931 | 458 | // |
459 | // Loop on chambers and on cathode planes | |
460 | // | |
461 | for (Int_t ich=0;ich<7;ich++) | |
462 | for (Int_t icat=0;icat<1;icat++) { | |
463 | // | |
464 | // Get ready the current chamber stuff | |
465 | // | |
466 | ||
467 | printf ("Olarilole"); | |
468 | AliRICHchamber* iChamber= &(this->Chamber(ich)); | |
469 | AliRICHsegmentation* segmentation = | |
470 | iChamber->GetSegmentationModel(icat+1); | |
471 | if (!segmentation) | |
472 | continue; | |
473 | TClonesArray *RICHdigits = this->DigitsAddress(ich); | |
474 | if (RICHdigits == 0) | |
475 | continue; | |
476 | cout << "Npx " << segmentation->Npx() | |
477 | << " Npy " << segmentation->Npy() << endl; | |
478 | // | |
479 | // Ready the digits | |
480 | // | |
481 | gAlice->ResetDigits(); | |
482 | gAlice->TreeD()->GetEvent(icat+1); // spurious +1 ... | |
483 | Int_t ndigits = RICHdigits->GetEntriesFast(); | |
484 | if (ndigits == 0) | |
485 | continue; | |
486 | printf("Found %d digits for cathode %d in chamber %d \n", | |
487 | ndigits,icat,ich+1); | |
488 | AliRICHdigit *mdig; | |
489 | AliRICHRecCluster *Cluster; | |
490 | // | |
491 | // Build the correspondance table | |
492 | // | |
f91473f6 | 493 | memset(elem,0,sizeof(Int_t)*kMaxNPadx*kMaxNPady*4); |
ddae0931 | 494 | Int_t digit; |
495 | for (digit=0; digit<ndigits; digit++) | |
496 | { | |
497 | mdig = (AliRICHdigit*)RICHdigits->UncheckedAt(digit); | |
f91473f6 | 498 | elem[kMaxNPadx+mdig->fPadX][kMaxNPady+mdig->fPadY] = digit+1; |
ddae0931 | 499 | // because default is 0 |
500 | } | |
501 | // | |
502 | // Declare some useful variables | |
503 | // | |
504 | Int_t Xlist[10]; | |
505 | Int_t Ylist[10]; | |
506 | Int_t Nlist; | |
507 | Int_t nclust=0; | |
508 | // | |
509 | // loop over digits | |
510 | // | |
511 | for (digit=0;digit<ndigits;digit++) { | |
512 | mdig = (AliRICHdigit*)RICHdigits->UncheckedAt(digit); | |
513 | // | |
514 | // if digit still available, start clustering | |
515 | // | |
f91473f6 | 516 | if (elem[kMaxNPadx+mdig->fPadX][kMaxNPady+mdig->fPadY]) { |
ddae0931 | 517 | Cluster = new AliRICHRecCluster(digit, ich, icat); |
f91473f6 | 518 | elem[kMaxNPadx+mdig->fPadX][kMaxNPady+mdig->fPadY]=0; |
ddae0931 | 519 | // |
520 | // loop over the current list of digits | |
521 | // and look for neighbours | |
522 | // | |
523 | for(Int_t clusDigit=Cluster->FirstDigitIndex(); | |
524 | clusDigit!=Cluster->InvalidDigitIndex(); | |
525 | clusDigit=Cluster->NextDigitIndex()) { | |
526 | AliRICHdigit* pDigit=(AliRICHdigit*)RICHdigits | |
527 | ->UncheckedAt(clusDigit); | |
528 | segmentation->Neighbours(pDigit->fPadX,pDigit->fPadY, | |
529 | &Nlist, Xlist, Ylist); | |
530 | for (Int_t Ilist=0;Ilist<Nlist;Ilist++) { | |
f91473f6 | 531 | if (elem[kMaxNPadx+Xlist[Ilist]][kMaxNPady |
ddae0931 | 532 | +Ylist[Ilist]]) { |
533 | // | |
534 | // Add the digit at the end and reset the table | |
535 | // | |
f91473f6 | 536 | Cluster->AddDigit(elem[kMaxNPadx+Xlist[Ilist]][kMaxNPady+Ylist[Ilist]]-1); |
537 | elem[kMaxNPadx+Xlist[Ilist]][kMaxNPady | |
ddae0931 | 538 | +Ylist[Ilist]] =0; |
539 | } // if elem | |
540 | } // for Ilist | |
541 | } // for pDigit | |
542 | // | |
543 | // Store the cluster (good time to do Cluster polishing) | |
544 | // | |
545 | segmentation->FitXY(Cluster,RICHdigits); | |
546 | nclust ++; | |
547 | AddRecCluster(ich,icat,Cluster); | |
548 | } | |
549 | } | |
550 | printf("===> %d Clusters\n",nclust); | |
551 | } // for icat | |
552 | } | |
553 | ||
554 | ||
555 | //______________________________________________________________________________ | |
556 | void AliRICH::Streamer(TBuffer &R__b) | |
557 | { | |
558 | // Stream an object of class AliRICH. | |
559 | AliRICHchamber *iChamber; | |
560 | AliRICHsegmentation *segmentation; | |
561 | AliRICHresponse *response; | |
562 | TClonesArray *digitsaddress; | |
563 | ||
564 | if (R__b.IsReading()) { | |
565 | Version_t R__v = R__b.ReadVersion(); if (R__v) { } | |
566 | AliDetector::Streamer(R__b); | |
567 | R__b >> fNclusters; | |
568 | R__b >> fClusters; // diff | |
569 | R__b >> fDchambers; | |
570 | R__b.ReadArray(fNdch); | |
571 | // | |
572 | R__b >> fAccCut; | |
573 | R__b >> fAccMin; | |
574 | R__b >> fAccMax; | |
575 | // | |
576 | R__b >> fChambers; | |
577 | // Stream chamber related information | |
578 | for (Int_t i =0; i<7; i++) { | |
579 | iChamber=(AliRICHchamber*) (*fChambers)[i]; | |
580 | iChamber->Streamer(R__b); | |
581 | if (iChamber->Nsec()==1) { | |
582 | segmentation=iChamber->GetSegmentationModel(1); | |
583 | segmentation->Streamer(R__b); | |
584 | } else { | |
585 | segmentation=iChamber->GetSegmentationModel(1); | |
586 | segmentation->Streamer(R__b); | |
587 | segmentation=iChamber->GetSegmentationModel(2); | |
588 | segmentation->Streamer(R__b); | |
589 | } | |
590 | response=iChamber->GetResponseModel(mip); | |
591 | response->Streamer(R__b); | |
592 | response=iChamber->GetResponseModel(cerenkov); | |
593 | response->Streamer(R__b); | |
594 | ||
595 | digitsaddress=(TClonesArray*) (*fDchambers)[i]; | |
596 | digitsaddress->Streamer(R__b); | |
597 | } | |
598 | ||
599 | } else { | |
600 | R__b.WriteVersion(AliRICH::IsA()); | |
601 | AliDetector::Streamer(R__b); | |
602 | R__b << fNclusters; | |
603 | R__b << fClusters; // diff | |
604 | R__b << fDchambers; | |
605 | R__b.WriteArray(fNdch, 7); | |
606 | // | |
607 | R__b << fAccCut; | |
608 | R__b << fAccMin; | |
609 | R__b << fAccMax; | |
610 | // | |
611 | R__b << fChambers; | |
612 | // Stream chamber related information | |
613 | for (Int_t i =0; i<7; i++) { | |
614 | iChamber=(AliRICHchamber*) (*fChambers)[i]; | |
615 | iChamber->Streamer(R__b); | |
616 | if (iChamber->Nsec()==1) { | |
617 | segmentation=iChamber->GetSegmentationModel(1); | |
618 | segmentation->Streamer(R__b); | |
619 | } else { | |
620 | segmentation=iChamber->GetSegmentationModel(1); | |
621 | segmentation->Streamer(R__b); | |
622 | segmentation=iChamber->GetSegmentationModel(2); | |
623 | segmentation->Streamer(R__b); | |
624 | } | |
625 | response=iChamber->GetResponseModel(mip); | |
626 | response->Streamer(R__b); | |
627 | response=iChamber->GetResponseModel(cerenkov); | |
628 | response->Streamer(R__b); | |
629 | ||
630 | digitsaddress=(TClonesArray*) (*fDchambers)[i]; | |
631 | digitsaddress->Streamer(R__b); | |
632 | } | |
fe4da5cc | 633 | } |
ddae0931 | 634 | } |
635 | AliRICHcluster* AliRICH::FirstPad(AliRICHhit* hit,TClonesArray *clusters ) | |
636 | { | |
637 | // | |
638 | // Initialise the pad iterator | |
639 | // Return the address of the first padhit for hit | |
640 | TClonesArray *theClusters = clusters; | |
641 | Int_t nclust = theClusters->GetEntriesFast(); | |
642 | if (nclust && hit->fPHlast > 0) { | |
643 | sMaxIterPad=hit->fPHlast; | |
644 | sCurIterPad=hit->fPHfirst; | |
645 | return (AliRICHcluster*) clusters->UncheckedAt(sCurIterPad-1); | |
646 | } else { | |
647 | return 0; | |
fe4da5cc | 648 | } |
ddae0931 | 649 | |
fe4da5cc | 650 | } |
651 | ||
ddae0931 | 652 | AliRICHcluster* AliRICH::NextPad(TClonesArray *clusters) |
fe4da5cc | 653 | { |
ddae0931 | 654 | sCurIterPad++; |
655 | if (sCurIterPad <= sMaxIterPad) { | |
656 | return (AliRICHcluster*) clusters->UncheckedAt(sCurIterPad-1); | |
657 | } else { | |
658 | return 0; | |
659 | } | |
fe4da5cc | 660 | } |
661 | ||
ddae0931 | 662 | ClassImp(AliRICHcluster) |
663 | ||
664 | //___________________________________________ | |
665 | AliRICHcluster::AliRICHcluster(Int_t *clhits) | |
666 | { | |
667 | fHitNumber=clhits[0]; | |
668 | fCathode=clhits[1]; | |
669 | fQ=clhits[2]; | |
670 | fPadX=clhits[3]; | |
671 | fPadY=clhits[4]; | |
672 | fQpad=clhits[5]; | |
673 | fRSec=clhits[6]; | |
674 | } | |
675 | ClassImp(AliRICHdigit) | |
fe4da5cc | 676 | //_____________________________________________________________________________ |
ddae0931 | 677 | AliRICHdigit::AliRICHdigit(Int_t *digits) |
fe4da5cc | 678 | { |
ddae0931 | 679 | // |
680 | // Creates a RICH digit object to be updated | |
681 | // | |
682 | fPadX = digits[0]; | |
683 | fPadY = digits[1]; | |
684 | fSignal = digits[2]; | |
685 | ||
fe4da5cc | 686 | } |
fe4da5cc | 687 | //_____________________________________________________________________________ |
ddae0931 | 688 | AliRICHdigit::AliRICHdigit(Int_t *tracks, Int_t *charges, Int_t *digits) |
fe4da5cc | 689 | { |
ddae0931 | 690 | // |
691 | // Creates a RICH digit object | |
692 | // | |
693 | fPadX = digits[0]; | |
694 | fPadY = digits[1]; | |
695 | fSignal = digits[2]; | |
696 | for(Int_t i=0; i<10; i++) { | |
697 | fTcharges[i] = charges[i]; | |
698 | fTracks[i] = tracks[i]; | |
699 | } | |
fe4da5cc | 700 | } |
701 | ||
ddae0931 | 702 | ClassImp(AliRICHlist) |
703 | ||
704 | //____________________________________________________________________________ | |
705 | AliRICHlist::AliRICHlist(Int_t ich, Int_t *digits): | |
706 | AliRICHdigit(digits) | |
fe4da5cc | 707 | { |
ddae0931 | 708 | // |
709 | // Creates a RICH digit list object | |
710 | // | |
711 | ||
712 | fChamber = ich; | |
713 | fTrackList = new TObjArray; | |
714 | ||
fe4da5cc | 715 | } |
ddae0931 | 716 | //_____________________________________________________________________________ |
fe4da5cc | 717 | |
fe4da5cc | 718 | |
ddae0931 | 719 | ClassImp(AliRICHhit) |
720 | ||
721 | //___________________________________________ | |
722 | AliRICHhit::AliRICHhit(Int_t shunt, Int_t track, Int_t *vol, Float_t *hits): | |
723 | AliHit(shunt, track) | |
fe4da5cc | 724 | { |
ddae0931 | 725 | fChamber=vol[0]; |
726 | fParticle=hits[0]; | |
727 | fX=hits[1]; | |
728 | fY=hits[2]; | |
729 | fZ=hits[3]; | |
730 | fTheta=hits[4]; | |
731 | fPhi=hits[5]; | |
732 | fTlength=hits[6]; | |
733 | fEloss=hits[7]; | |
734 | fPHfirst=(Int_t) hits[8]; | |
735 | fPHlast=(Int_t) hits[9]; | |
fe4da5cc | 736 | } |
ddae0931 | 737 | ClassImp(AliRICHreccluster) |
738 | ||
739 | ClassImp(AliRICHCerenkov) | |
740 | //___________________________________________ | |
741 | AliRICHCerenkov::AliRICHCerenkov(Int_t shunt, Int_t track, Int_t *vol, Float_t *hits): | |
742 | AliHit(shunt, track) | |
fe4da5cc | 743 | { |
ddae0931 | 744 | fChamber=vol[0]; |
745 | fX=hits[1]; | |
746 | fY=hits[2]; | |
747 | fZ=hits[3]; | |
748 | fTheta=hits[4]; | |
749 | fPhi=hits[5]; | |
750 | fTlength=hits[6]; | |
751 | fPHfirst=(Int_t) hits[8]; | |
752 | fPHlast=(Int_t) hits[9]; | |
fe4da5cc | 753 | } |
754 | ||
ddae0931 | 755 | ClassImp(AliRICHRecCluster) |
756 | ||
757 | //_____________________________________________________________________ | |
758 | AliRICHRecCluster::AliRICHRecCluster() | |
fe4da5cc | 759 | { |
ddae0931 | 760 | fDigits=0; |
761 | fNdigit=-1; | |
fe4da5cc | 762 | } |
ddae0931 | 763 | |
764 | AliRICHRecCluster::AliRICHRecCluster(Int_t FirstDigit,Int_t Ichamber, Int_t Icathod) | |
fe4da5cc | 765 | { |
ddae0931 | 766 | fX = 0.; |
767 | fY = 0.; | |
768 | fDigits = new TArrayI(10); | |
769 | fNdigit=0; | |
770 | AddDigit(FirstDigit); | |
771 | fChamber=Ichamber; | |
772 | fCathod=Icathod; | |
fe4da5cc | 773 | } |
774 | ||
ddae0931 | 775 | void AliRICHRecCluster::AddDigit(Int_t Digit) |
fe4da5cc | 776 | { |
ddae0931 | 777 | if (fNdigit==fDigits->GetSize()) { |
778 | //enlarge the list by hand! | |
779 | Int_t *array= new Int_t[fNdigit*2]; | |
780 | for(Int_t i=0;i<fNdigit;i++) | |
781 | array[i] = fDigits->At(i); | |
782 | fDigits->Adopt(fNdigit*2,array); | |
783 | } | |
784 | fDigits->AddAt(Digit,fNdigit); | |
785 | fNdigit++; | |
fe4da5cc | 786 | } |
787 | ||
ddae0931 | 788 | |
789 | AliRICHRecCluster::~AliRICHRecCluster() | |
fe4da5cc | 790 | { |
ddae0931 | 791 | if (fDigits) |
792 | delete fDigits; | |
fe4da5cc | 793 | } |
794 | ||
ddae0931 | 795 | Int_t AliRICHRecCluster::FirstDigitIndex() |
fe4da5cc | 796 | { |
ddae0931 | 797 | fCurrentDigit=0; |
798 | return fDigits->At(fCurrentDigit); | |
fe4da5cc | 799 | } |
800 | ||
ddae0931 | 801 | Int_t AliRICHRecCluster::NextDigitIndex() |
802 | { | |
803 | fCurrentDigit++; | |
804 | if (fCurrentDigit<fNdigit) | |
805 | return fDigits->At(fCurrentDigit); | |
806 | else | |
807 | return InvalidDigitIndex(); | |
808 | } | |
fe4da5cc | 809 | |
ddae0931 | 810 | Int_t AliRICHRecCluster::NDigits() |
fe4da5cc | 811 | { |
ddae0931 | 812 | return fNdigit; |
fe4da5cc | 813 | } |
814 | ||
ddae0931 | 815 | void AliRICHRecCluster::Finish() |
fe4da5cc | 816 | { |
ddae0931 | 817 | // In order to reconstruct coordinates, one has to |
818 | // get back to the digits which is not trivial here, | |
819 | // because we don't know where digits are stored! | |
820 | // Center Of Gravity, or other method should be | |
821 | // a property of AliRICH class! | |
fe4da5cc | 822 | } |
823 | ||
fe4da5cc | 824 | |
825 | ||
ddae0931 | 826 | void AliRICH::Digitise(Int_t nev,Option_t *option,Text_t *filename) |
fe4da5cc | 827 | { |
ddae0931 | 828 | // keep galice.root for signal and name differently the file for |
829 | // background when add! otherwise the track info for signal will be lost ! | |
830 | ||
c90dd3e2 | 831 | static Bool_t first=kTRUE; |
ddae0931 | 832 | static TTree *TH1; |
833 | static TFile *File; | |
c90dd3e2 | 834 | Int_t i; |
ddae0931 | 835 | char *Add = strstr(option,"Add"); |
836 | ||
837 | AliRICHchamber* iChamber; | |
838 | AliRICHsegmentation* segmentation; | |
839 | ||
840 | ||
841 | Int_t trk[50]; | |
842 | Int_t chtrk[50]; | |
843 | TObjArray *list=new TObjArray; | |
844 | Int_t digits[3]; | |
845 | ||
846 | AliRICH *RICH = (AliRICH *) gAlice->GetDetector("RICH"); | |
847 | AliRICHHitMap* HitMap[10]; | |
c90dd3e2 | 848 | for (i=0; i<10; i++) {HitMap[i]=0;} |
ddae0931 | 849 | if (Add ) { |
850 | if(first) { | |
851 | fFileName=filename; | |
852 | cout<<"filename"<<fFileName<<endl; | |
853 | File=new TFile(fFileName); | |
854 | cout<<"I have opened "<<fFileName<<" file "<<endl; | |
855 | fHits2 = new TClonesArray("AliRICHhit",1000 ); | |
856 | fClusters2 = new TClonesArray("AliRICHcluster",10000); | |
c90dd3e2 | 857 | first=kFALSE; |
ddae0931 | 858 | } |
859 | File->cd(); | |
860 | File->ls(); | |
861 | // Get Hits Tree header from file | |
862 | if(fHits2) fHits2->Clear(); | |
863 | if(fClusters2) fClusters2->Clear(); | |
864 | if(TH1) delete TH1; | |
865 | TH1=0; | |
866 | // | |
867 | char treeName[20]; | |
868 | sprintf(treeName,"TreeH%d",nev); | |
869 | TH1 = (TTree*)gDirectory->Get(treeName); | |
870 | if (!TH1) { | |
871 | printf("ERROR: cannot find Hits Tree for event:%d\n",nev); | |
872 | } | |
873 | // Set branch addresses | |
874 | TBranch *branch; | |
875 | char branchname[20]; | |
876 | sprintf(branchname,"%s",GetName()); | |
877 | if (TH1 && fHits2) { | |
878 | branch = TH1->GetBranch(branchname); | |
879 | if (branch) branch->SetAddress(&fHits2); | |
880 | } | |
881 | if (TH1 && fClusters2) { | |
882 | branch = TH1->GetBranch("RICHCluster"); | |
883 | if (branch) branch->SetAddress(&fClusters2); | |
884 | } | |
885 | } | |
886 | // | |
887 | // loop over cathodes | |
888 | // | |
889 | AliRICHHitMap* hm; | |
890 | ||
f91473f6 | 891 | for (int icat=0; icat<1; icat++) { |
c90dd3e2 | 892 | for (i=0; i<7; i++) { |
ddae0931 | 893 | if (HitMap[i]) { |
894 | hm=HitMap[i]; | |
895 | delete hm; | |
896 | HitMap[i]=0; | |
897 | } | |
898 | } | |
899 | Int_t counter=0; | |
c90dd3e2 | 900 | for (i =0; i<7; i++) { |
ddae0931 | 901 | iChamber=(AliRICHchamber*) (*fChambers)[i]; |
902 | if (iChamber->Nsec()==1 && icat==1) { | |
903 | continue; | |
904 | } else { | |
905 | segmentation=iChamber->GetSegmentationModel(icat+1); | |
906 | } | |
907 | HitMap[i] = new AliRICHHitMapA1(segmentation, list); | |
908 | } | |
909 | printf("Start loop over tracks \n"); | |
910 | // | |
911 | // Loop over tracks | |
912 | // | |
913 | ||
914 | TTree *TH = gAlice->TreeH(); | |
915 | Int_t ntracks =(Int_t) TH->GetEntries(); | |
916 | for (Int_t track=0; track<ntracks; track++) { | |
917 | gAlice->ResetHits(); | |
918 | TH->GetEvent(track); | |
919 | // | |
920 | // Loop over hits | |
921 | for(AliRICHhit* mHit=(AliRICHhit*)RICH->FirstHit(-1); | |
922 | mHit; | |
923 | mHit=(AliRICHhit*)RICH->NextHit()) | |
924 | { | |
925 | Int_t nch = mHit->fChamber-1; // chamber number | |
926 | if (nch >7) continue; | |
927 | iChamber = &(RICH->Chamber(nch)); | |
928 | ||
929 | // | |
930 | // Loop over pad hits | |
931 | for (AliRICHcluster* mPad= | |
932 | (AliRICHcluster*)RICH->FirstPad(mHit,fClusters); | |
933 | mPad; | |
934 | mPad=(AliRICHcluster*)RICH->NextPad(fClusters)) | |
935 | { | |
936 | Int_t cathode = mPad->fCathode; // cathode number | |
937 | Int_t ipx = mPad->fPadX; // pad number on X | |
938 | Int_t ipy = mPad->fPadY; // pad number on Y | |
939 | Int_t iqpad = mPad->fQpad; // charge per pad | |
940 | // | |
941 | // | |
942 | ||
943 | if (cathode != (icat+1)) continue; | |
944 | // fill the info array | |
945 | Float_t thex, they; | |
946 | segmentation=iChamber->GetSegmentationModel(cathode); | |
947 | segmentation->GetPadCxy(ipx,ipy,thex,they); | |
948 | TVector *trinfo_p= new TVector(2); | |
949 | TVector &trinfo = *trinfo_p; | |
950 | trinfo(0)=(Float_t)track; | |
951 | trinfo(1)=(Float_t)iqpad; | |
952 | ||
953 | digits[0]=ipx; | |
954 | digits[1]=ipy; | |
955 | digits[2]=iqpad; | |
956 | ||
957 | AliRICHlist* pdigit; | |
958 | // build the list of fired pads and update the info | |
959 | if (!HitMap[nch]->TestHit(ipx, ipy)) { | |
960 | list->AddAtAndExpand( | |
961 | new AliRICHlist(nch,digits),counter); | |
962 | HitMap[nch]->SetHit(ipx, ipy, counter); | |
963 | counter++; | |
964 | pdigit=(AliRICHlist*)list->At(list->GetLast()); | |
965 | // list of tracks | |
966 | TObjArray *trlist=(TObjArray*)pdigit->TrackList(); | |
967 | trlist->Add(&trinfo); | |
968 | } else { | |
969 | pdigit=(AliRICHlist*) HitMap[nch]->GetHit(ipx, ipy); | |
970 | // update charge | |
971 | (*pdigit).fSignal+=iqpad; | |
972 | // update list of tracks | |
973 | TObjArray* trlist=(TObjArray*)pdigit->TrackList(); | |
974 | Int_t last_entry=trlist->GetLast(); | |
975 | TVector *ptrk_p=(TVector*)trlist->At(last_entry); | |
976 | TVector &ptrk=*ptrk_p; | |
977 | Int_t last_track=Int_t(ptrk(0)); | |
978 | Int_t last_charge=Int_t(ptrk(1)); | |
979 | if (last_track==track) { | |
980 | last_charge+=iqpad; | |
981 | trlist->RemoveAt(last_entry); | |
982 | trinfo(0)=last_track; | |
983 | trinfo(1)=last_charge; | |
984 | trlist->AddAt(&trinfo,last_entry); | |
985 | } else { | |
986 | trlist->Add(&trinfo); | |
987 | } | |
988 | // check the track list | |
989 | Int_t nptracks=trlist->GetEntriesFast(); | |
990 | if (nptracks > 2) { | |
991 | printf("Attention - nptracks > 2 %d \n",nptracks); | |
992 | printf("cat,nch,ix,iy %d %d %d %d \n",icat+1,nch,ipx,ipy); | |
993 | for (Int_t tr=0;tr<nptracks;tr++) { | |
994 | TVector *pptrk_p=(TVector*)trlist->At(tr); | |
995 | TVector &pptrk=*pptrk_p; | |
996 | trk[tr]=Int_t(pptrk(0)); | |
997 | chtrk[tr]=Int_t(pptrk(1)); | |
998 | } | |
999 | } // end if nptracks | |
1000 | } // end if pdigit | |
1001 | } //end loop over clusters | |
1002 | } // hit loop | |
1003 | } // track loop | |
1004 | ||
1005 | Int_t nentr1=list->GetEntriesFast(); | |
1006 | printf(" \n counter, nentr1 %d %d\n",counter,nentr1); | |
1007 | ||
1008 | // open the file with background | |
1009 | ||
1010 | if (Add ) { | |
1011 | ntracks =(Int_t)TH1->GetEntries(); | |
1012 | printf("background - icat,ntracks1 %d %d\n",icat,ntracks); | |
1013 | printf("background - Start loop over tracks \n"); | |
1014 | // | |
1015 | // Loop over tracks | |
1016 | // | |
f91473f6 | 1017 | for (Int_t trak=0; trak<ntracks; trak++) { |
ddae0931 | 1018 | |
1019 | if (fHits2) fHits2->Clear(); | |
1020 | if (fClusters2) fClusters2->Clear(); | |
1021 | ||
f91473f6 | 1022 | TH1->GetEvent(trak); |
ddae0931 | 1023 | // |
1024 | // Loop over hits | |
1025 | AliRICHhit* mHit; | |
f91473f6 | 1026 | for(int j=0;j<fHits2->GetEntriesFast();++j) |
ddae0931 | 1027 | { |
f91473f6 | 1028 | mHit=(AliRICHhit*) (*fHits2)[j]; |
ddae0931 | 1029 | Int_t nch = mHit->fChamber-1; // chamber number |
1030 | if (nch >9) continue; | |
1031 | iChamber = &(RICH->Chamber(nch)); | |
1032 | Int_t rmin = (Int_t)iChamber->RInner(); | |
1033 | Int_t rmax = (Int_t)iChamber->ROuter(); | |
1034 | // | |
1035 | // Loop over pad hits | |
1036 | for (AliRICHcluster* mPad= | |
1037 | (AliRICHcluster*)RICH->FirstPad(mHit,fClusters2); | |
1038 | mPad; | |
1039 | mPad=(AliRICHcluster*)RICH->NextPad(fClusters2)) | |
1040 | { | |
1041 | Int_t cathode = mPad->fCathode; // cathode number | |
1042 | Int_t ipx = mPad->fPadX; // pad number on X | |
1043 | Int_t ipy = mPad->fPadY; // pad number on Y | |
1044 | Int_t iqpad = mPad->fQpad; // charge per pad | |
f91473f6 | 1045 | if (trak==3 && nch==0 && icat==0) printf("bgr - trak,iqpad,ipx,ipy %d %d %d %d\n",trak,iqpad,ipx,ipy); |
ddae0931 | 1046 | // |
1047 | // | |
1048 | if (cathode != (icat+1)) continue; | |
1049 | // fill the info array | |
1050 | Float_t thex, they; | |
1051 | segmentation=iChamber->GetSegmentationModel(cathode); | |
1052 | segmentation->GetPadCxy(ipx,ipy,thex,they); | |
1053 | Float_t rpad=TMath::Sqrt(thex*thex+they*they); | |
1054 | if (rpad < rmin || iqpad ==0 || rpad > rmax) continue; | |
1055 | ||
1056 | TVector *trinfo_p; | |
1057 | trinfo_p = new TVector(2); | |
1058 | TVector &trinfo = *trinfo_p; | |
1059 | trinfo(0)=-1; // tag background | |
1060 | trinfo(1)=-1; | |
1061 | ||
1062 | digits[0]=ipx; | |
1063 | digits[1]=ipy; | |
1064 | digits[2]=iqpad; | |
1065 | ||
1066 | ||
f91473f6 | 1067 | if (trak <4 && icat==0 && nch==0) |
1068 | printf("bgr - HitMap[nch]->TestHit(ipx, ipy),trak %d %d\n", | |
1069 | HitMap[nch]->TestHit(ipx, ipy),trak); | |
ddae0931 | 1070 | AliRICHlist* pdigit; |
1071 | // build the list of fired pads and update the info | |
1072 | if (!HitMap[nch]->TestHit(ipx, ipy)) { | |
1073 | list->AddAtAndExpand(new AliRICHlist(nch,digits),counter); | |
1074 | ||
1075 | HitMap[nch]->SetHit(ipx, ipy, counter); | |
1076 | counter++; | |
1077 | printf("bgr new elem in list - counter %d\n",counter); | |
1078 | ||
1079 | pdigit=(AliRICHlist*)list->At(list->GetLast()); | |
1080 | // list of tracks | |
1081 | TObjArray *trlist=(TObjArray*)pdigit->TrackList(); | |
1082 | trlist->Add(&trinfo); | |
1083 | } else { | |
1084 | pdigit=(AliRICHlist*) HitMap[nch]->GetHit(ipx, ipy); | |
1085 | // update charge | |
1086 | (*pdigit).fSignal+=iqpad; | |
1087 | // update list of tracks | |
1088 | TObjArray* trlist=(TObjArray*)pdigit->TrackList(); | |
1089 | Int_t last_entry=trlist->GetLast(); | |
1090 | TVector *ptrk_p=(TVector*)trlist->At(last_entry); | |
1091 | TVector &ptrk=*ptrk_p; | |
1092 | Int_t last_track=Int_t(ptrk(0)); | |
1093 | if (last_track==-1) { | |
1094 | continue; | |
1095 | } else { | |
1096 | trlist->Add(&trinfo); | |
1097 | } | |
1098 | // check the track list | |
1099 | Int_t nptracks=trlist->GetEntriesFast(); | |
1100 | if (nptracks > 0) { | |
1101 | for (Int_t tr=0;tr<nptracks;tr++) { | |
1102 | TVector *pptrk_p=(TVector*)trlist->At(tr); | |
1103 | TVector &pptrk=*pptrk_p; | |
1104 | trk[tr]=Int_t(pptrk(0)); | |
1105 | chtrk[tr]=Int_t(pptrk(1)); | |
1106 | } | |
1107 | } // end if nptracks | |
1108 | } // end if pdigit | |
1109 | } //end loop over clusters | |
1110 | } // hit loop | |
1111 | } // track loop | |
1112 | Int_t nentr2=list->GetEntriesFast(); | |
1113 | printf(" \n counter2, nentr2 %d %d \n",counter,nentr2); | |
1114 | TTree *fAli=gAlice->TreeK(); | |
1115 | if (fAli) File =fAli->GetCurrentFile(); | |
1116 | File->cd(); | |
1117 | } // if Add | |
1118 | ||
1119 | Int_t tracks[10]; | |
1120 | Int_t charges[10]; | |
1121 | cout<<"start filling digits \n "<<endl; | |
1122 | Int_t nentries=list->GetEntriesFast(); | |
1123 | printf(" \n \n nentries %d \n",nentries); | |
1124 | ||
1125 | // start filling the digits | |
1126 | ||
1127 | for (Int_t nent=0;nent<nentries;nent++) { | |
1128 | AliRICHlist *address=(AliRICHlist*)list->At(nent); | |
1129 | if (address==0) continue; | |
1130 | Int_t ich=address->fChamber; | |
1131 | Int_t q=address->fSignal; | |
1132 | iChamber=(AliRICHchamber*) (*fChambers)[ich]; | |
1133 | // add white noise and do zero-suppression and signal truncation | |
1134 | Float_t MeanNoise = gRandom->Gaus(1, 0.2); | |
1135 | Float_t ZeroSupp=5*MeanNoise; | |
1136 | Float_t Noise = gRandom->Gaus(0, MeanNoise); | |
1137 | q+=(Int_t)Noise; | |
1138 | if ( q <= ZeroSupp) continue; | |
1139 | digits[0]=address->fPadX; | |
1140 | digits[1]=address->fPadY; | |
1141 | digits[2]=q; | |
1142 | ||
1143 | TObjArray* trlist=(TObjArray*)address->TrackList(); | |
1144 | Int_t nptracks=trlist->GetEntriesFast(); | |
1145 | ||
1146 | // this was changed to accomodate the real number of tracks | |
1147 | if (nptracks > 10) { | |
1148 | cout<<"Attention - nptracks > 10 "<<nptracks<<endl; | |
1149 | nptracks=10; | |
1150 | } | |
1151 | if (nptracks > 2) { | |
1152 | printf("Attention - nptracks > 2 %d \n",nptracks); | |
1153 | printf("cat,ich,ix,iy,q %d %d %d %d %d \n",icat,ich,digits[0],digits[1],q); | |
1154 | } | |
1155 | for (Int_t tr=0;tr<nptracks;tr++) { | |
1156 | TVector *pp_p=(TVector*)trlist->At(tr); | |
1157 | TVector &pp =*pp_p; | |
1158 | tracks[tr]=Int_t(pp(0)); | |
1159 | charges[tr]=Int_t(pp(1)); | |
1160 | } //end loop over list of tracks for one pad | |
1161 | if (nptracks < 10 ) { | |
f91473f6 | 1162 | for (Int_t t=nptracks; t<10; t++) { |
1163 | tracks[t]=0; | |
1164 | charges[t]=0; | |
ddae0931 | 1165 | } |
1166 | } | |
1167 | // fill digits | |
1168 | RICH->AddDigits(ich,tracks,charges,digits); | |
1169 | ||
1170 | delete address; | |
1171 | } | |
1172 | cout<<"I'm out of the loops for digitisation"<<endl; | |
1173 | gAlice->TreeD()->Fill(); | |
1174 | TTree *TD=gAlice->TreeD(); | |
1175 | Stat_t ndig=TD->GetEntries(); | |
1176 | cout<<"number of digits "<<ndig<<endl; | |
1177 | TClonesArray *fDch; | |
f91473f6 | 1178 | for (int k=0;k<7;k++) { |
1179 | fDch= RICH->DigitsAddress(k); | |
c90dd3e2 | 1180 | int ndigit=fDch->GetEntriesFast(); |
1181 | printf (" k, ndigits %d %d \n",k,ndigit); | |
ddae0931 | 1182 | } |
1183 | RICH->ResetDigits(); | |
1184 | ||
1185 | list->Clear(); | |
1186 | ||
1187 | } //end loop over cathodes | |
1188 | char hname[30]; | |
1189 | sprintf(hname,"TreeD%d",nev); | |
1190 | gAlice->TreeD()->Write(hname); | |
fe4da5cc | 1191 | } |
1192 | ||
1193 | ||
ddae0931 | 1194 | |
1195 | ||
1196 | ||
1197 | ||
1198 | ||
1199 | ||
1200 |