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