]> git.uio.no Git - u/mrichter/AliRoot.git/blame - RICH/AliRICH.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Default constructor added.
[u/mrichter/AliRoot.git] / RICH / AliRICH.cxx
CommitLineData
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/*
2e5f0f7b 17 $Log$
91975aa2 18 Revision 1.39 2001/01/22 21:40:24 jbarbosa
19 Removing magic numbers
20
ee2105a4 21 Revision 1.37 2000/12/20 14:07:25 jbarbosa
22 Removed dependencies on TGeant3 (thanks to F. Carminati and I. Hrivnacova)
23
5cfcdf54 24 Revision 1.36 2000/12/18 17:45:54 jbarbosa
25 Cleaned up PadHits object.
26
176a9917 27 Revision 1.35 2000/12/15 16:49:40 jbarbosa
28 Geometry and materials updates (wire supports, pcbs, backplane supports, frame).
29
6197ac87 30 Revision 1.34 2000/11/10 18:12:12 jbarbosa
31 Bug fix for AliRICHCerenkov (thanks to P. Hristov)
32
3cc8edee 33 Revision 1.33 2000/11/02 10:09:01 jbarbosa
34 Minor bug correction (some pointers were not initialised in the default constructor)
35
9b674b76 36 Revision 1.32 2000/11/01 15:32:55 jbarbosa
37 Updated to handle both reconstruction algorithms.
38
a4622d0f 39 Revision 1.31 2000/10/26 20:18:33 jbarbosa
40 Supports for methane and freon vessels
41
683b273b 42 Revision 1.30 2000/10/24 13:19:12 jbarbosa
43 Geometry updates.
44
e293afae 45 Revision 1.29 2000/10/19 19:39:25 jbarbosa
46 Some more changes to geometry. Further correction of digitisation "per part. type"
47
53a60579 48 Revision 1.28 2000/10/17 20:50:57 jbarbosa
49 Inversed digtise by particle type (now, only the selected particle type is not digitsed).
50 Corrected several geometry minor bugs.
51 Added new parameter (opaque quartz thickness).
52
3587e146 53 Revision 1.27 2000/10/11 10:33:55 jbarbosa
54 Corrected bug introduced by earlier revisions (CerenkovData array cannot be reset to zero on wach call of StepManager)
55
8cda28a5 56 Revision 1.26 2000/10/03 21:44:08 morsch
57 Use AliSegmentation and AliHit abstract base classes.
58
a2f7eaf6 59 Revision 1.25 2000/10/02 21:28:12 fca
60 Removal of useless dependecies via forward declarations
61
94de3818 62 Revision 1.24 2000/10/02 15:43:17 jbarbosa
63 Fixed forward declarations.
64 Fixed honeycomb density.
65 Fixed cerenkov storing.
66 New electronics.
67
9dda3582 68 Revision 1.23 2000/09/13 10:42:14 hristov
69 Minor corrections for HP, DEC and Sun; strings.h included
70
e813258a 71 Revision 1.22 2000/09/12 18:11:13 fca
72 zero hits area before using
73
e51a3fa9 74 Revision 1.21 2000/07/21 10:21:07 morsch
75 fNrawch = 0; and fNrechits = 0; in the default constructor.
76
db481a6e 77 Revision 1.20 2000/07/10 15:28:39 fca
78 Correction of the inheritance scheme
79
452a64c6 80 Revision 1.19 2000/06/30 16:29:51 dibari
81 Added kDebugLevel variable to control output size on demand
82
33984590 83 Revision 1.18 2000/06/12 15:15:46 jbarbosa
84 Cleaned up version.
85
237c933d 86 Revision 1.17 2000/06/09 14:58:37 jbarbosa
87 New digitisation per particle type
88
d28b5fc2 89 Revision 1.16 2000/04/19 12:55:43 morsch
90 Newly structured and updated version (JB, AM)
91
4c039060 92*/
93
2e5f0f7b 94
ddae0931 95////////////////////////////////////////////////
96// Manager and hits classes for set:RICH //
97////////////////////////////////////////////////
fe4da5cc 98
ddae0931 99#include <TBRIK.h>
100#include <TTUBE.h>
fe4da5cc 101#include <TNode.h>
102#include <TRandom.h>
ddae0931 103#include <TObject.h>
104#include <TVector.h>
105#include <TObjArray.h>
2e5f0f7b 106#include <TArrayF.h>
237c933d 107#include <TFile.h>
108#include <TParticle.h>
94de3818 109#include <TGeometry.h>
110#include <TTree.h>
111
237c933d 112#include <iostream.h>
e813258a 113#include <strings.h>
fe4da5cc 114
fe4da5cc 115#include "AliRICH.h"
a2f7eaf6 116#include "AliSegmentation.h"
91975aa2 117#include "AliRICHSegmentationV0.h"
237c933d 118#include "AliRICHHit.h"
119#include "AliRICHCerenkov.h"
120#include "AliRICHPadHit.h"
121#include "AliRICHDigit.h"
122#include "AliRICHTransientDigit.h"
123#include "AliRICHRawCluster.h"
a4622d0f 124#include "AliRICHRecHit1D.h"
125#include "AliRICHRecHit3D.h"
237c933d 126#include "AliRICHHitMapA1.h"
2e5f0f7b 127#include "AliRICHClusterFinder.h"
fe4da5cc 128#include "AliRun.h"
ddae0931 129#include "AliMC.h"
94de3818 130#include "AliMagF.h"
452a64c6 131#include "AliConst.h"
132#include "AliPDG.h"
2e5f0f7b 133#include "AliPoints.h"
ddae0931 134#include "AliCallf77.h"
237c933d 135
ddae0931 136
137// Static variables for the pad-hit iterator routines
138static Int_t sMaxIterPad=0;
139static Int_t sCurIterPad=0;
140static TClonesArray *fClusters2;
141static TClonesArray *fHits2;
2e5f0f7b 142static TTree *TrH1;
ddae0931 143
fe4da5cc 144ClassImp(AliRICH)
ddae0931 145
146//___________________________________________
fe4da5cc 147AliRICH::AliRICH()
148{
237c933d 149// Default constructor for RICH manager class
150
ddae0931 151 fIshunt = 0;
152 fHits = 0;
2e5f0f7b 153 fPadHits = 0;
154 fNPadHits = 0;
155 fNcerenkovs = 0;
ddae0931 156 fDchambers = 0;
9b674b76 157 fRecHits1D = 0;
158 fRecHits3D = 0;
159 fRawClusters = 0;
160 fChambers = 0;
ddae0931 161 fCerenkovs = 0;
9dda3582 162 for (Int_t i=0; i<7; i++)
163 {
164 fNdch[i] = 0;
165 fNrawch[i] = 0;
a4622d0f 166 fNrechits1D[i] = 0;
167 fNrechits3D[i] = 0;
9dda3582 168 }
9b674b76 169
170 fFileName = 0;
fe4da5cc 171}
ddae0931 172
173//___________________________________________
fe4da5cc 174AliRICH::AliRICH(const char *name, const char *title)
ddae0931 175 : AliDetector(name,title)
fe4da5cc 176{
ddae0931 177//Begin_Html
178/*
179 <img src="gif/alirich.gif">
180*/
181//End_Html
182
2e5f0f7b 183 fHits = new TClonesArray("AliRICHHit",1000 );
1cedd08a 184 gAlice->AddHitList(fHits);
2e5f0f7b 185 fPadHits = new TClonesArray("AliRICHPadHit",100000);
186 fCerenkovs = new TClonesArray("AliRICHCerenkov",1000);
187 gAlice->AddHitList(fCerenkovs);
188 //gAlice->AddHitList(fHits);
189 fNPadHits = 0;
190 fNcerenkovs = 0;
191 fIshunt = 0;
ddae0931 192
9dda3582 193 //fNdch = new Int_t[kNCH];
ddae0931 194
237c933d 195 fDchambers = new TObjArray(kNCH);
2e5f0f7b 196
a4622d0f 197 fRecHits1D = new TObjArray(kNCH);
198 fRecHits3D = new TObjArray(kNCH);
ddae0931 199
200 Int_t i;
201
237c933d 202 for (i=0; i<kNCH ;i++) {
2e5f0f7b 203 (*fDchambers)[i] = new TClonesArray("AliRICHDigit",10000);
ddae0931 204 fNdch[i]=0;
205 }
2e5f0f7b 206
9dda3582 207 //fNrawch = new Int_t[kNCH];
ddae0931 208
237c933d 209 fRawClusters = new TObjArray(kNCH);
d28b5fc2 210 //printf("Created fRwClusters with adress:%p",fRawClusters);
2e5f0f7b 211
237c933d 212 for (i=0; i<kNCH ;i++) {
2e5f0f7b 213 (*fRawClusters)[i] = new TClonesArray("AliRICHRawCluster",10000);
214 fNrawch[i]=0;
215 }
216
9dda3582 217 //fNrechits = new Int_t[kNCH];
ddae0931 218
237c933d 219 for (i=0; i<kNCH ;i++) {
a4622d0f 220 (*fRecHits1D)[i] = new TClonesArray("AliRICHRecHit1D",1000);
221 }
222 for (i=0; i<kNCH ;i++) {
223 (*fRecHits3D)[i] = new TClonesArray("AliRICHRecHit3D",1000);
2e5f0f7b 224 }
d28b5fc2 225 //printf("Created fRecHits with adress:%p",fRecHits);
2e5f0f7b 226
227
ddae0931 228 SetMarkerColor(kRed);
9b674b76 229
3cc8edee 230 /*fChambers = new TObjArray(kNCH);
9b674b76 231 for (i=0; i<kNCH; i++)
3cc8edee 232 (*fChambers)[i] = new AliRICHChamber();*/
9b674b76 233
234 fFileName = 0;
235
fe4da5cc 236}
237
237c933d 238AliRICH::AliRICH(const AliRICH& RICH)
239{
240// Copy Constructor
241}
242
243
ddae0931 244//___________________________________________
fe4da5cc 245AliRICH::~AliRICH()
246{
237c933d 247
248// Destructor of RICH manager class
249
ddae0931 250 fIshunt = 0;
251 delete fHits;
2e5f0f7b 252 delete fPadHits;
ddae0931 253 delete fCerenkovs;
ee2105a4 254
255 //PH Delete TObjArrays
256 if (fChambers) {
257 fChambers->Delete();
258 delete fChambers;
259 }
260 if (fDchambers) {
261 fDchambers->Delete();
262 delete fDchambers;
263 }
264 if (fRawClusters) {
265 fRawClusters->Delete();
266 delete fRawClusters;
267 }
268 if (fRecHits1D) {
269 fRecHits1D->Delete();
270 delete fRecHits1D;
271 }
272 if (fRecHits3D) {
273 fRecHits3D->Delete();
274 delete fRecHits3D;
275 }
276
fe4da5cc 277}
278
ddae0931 279//___________________________________________
fe4da5cc 280void AliRICH::AddHit(Int_t track, Int_t *vol, Float_t *hits)
281{
237c933d 282
283//
284// Adds a hit to the Hits list
285//
286
ddae0931 287 TClonesArray &lhits = *fHits;
2e5f0f7b 288 new(lhits[fNhits++]) AliRICHHit(fIshunt,track,vol,hits);
fe4da5cc 289}
fe4da5cc 290//_____________________________________________________________________________
ddae0931 291void AliRICH::AddCerenkov(Int_t track, Int_t *vol, Float_t *cerenkovs)
292{
237c933d 293
294//
295// Adds a RICH cerenkov hit to the Cerenkov Hits list
296//
297
ddae0931 298 TClonesArray &lcerenkovs = *fCerenkovs;
299 new(lcerenkovs[fNcerenkovs++]) AliRICHCerenkov(fIshunt,track,vol,cerenkovs);
2e5f0f7b 300 //printf ("Done for Cerenkov %d\n\n\n\n",fNcerenkovs);
fe4da5cc 301}
ddae0931 302//___________________________________________
2e5f0f7b 303void AliRICH::AddPadHit(Int_t *clhits)
ddae0931 304{
237c933d 305
306//
307// Add a RICH pad hit to the list
308//
309
2e5f0f7b 310 TClonesArray &lPadHits = *fPadHits;
311 new(lPadHits[fNPadHits++]) AliRICHPadHit(clhits);
33984590 312}
fe4da5cc 313//_____________________________________________________________________________
ddae0931 314void AliRICH::AddDigits(Int_t id, Int_t *tracks, Int_t *charges, Int_t *digits)
315{
237c933d 316
317 //
318 // Add a RICH digit to the list
319 //
ddae0931 320
321 TClonesArray &ldigits = *((TClonesArray*)(*fDchambers)[id]);
2e5f0f7b 322 new(ldigits[fNdch[id]++]) AliRICHDigit(tracks,charges,digits);
fe4da5cc 323}
324
325//_____________________________________________________________________________
2e5f0f7b 326void AliRICH::AddRawCluster(Int_t id, const AliRICHRawCluster& c)
fe4da5cc 327{
ddae0931 328 //
2e5f0f7b 329 // Add a RICH digit to the list
ddae0931 330 //
2e5f0f7b 331
332 TClonesArray &lrawcl = *((TClonesArray*)(*fRawClusters)[id]);
333 new(lrawcl[fNrawch[id]++]) AliRICHRawCluster(c);
fe4da5cc 334}
335
2e5f0f7b 336//_____________________________________________________________________________
a4622d0f 337void AliRICH::AddRecHit1D(Int_t id, Float_t *rechit, Float_t *photons, Int_t *padsx, Int_t* padsy)
338{
339
340 //
341 // Add a RICH reconstructed hit to the list
342 //
343
344 TClonesArray &lrec1D = *((TClonesArray*)(*fRecHits1D)[id]);
345 new(lrec1D[fNrechits1D[id]++]) AliRICHRecHit1D(id,rechit,photons,padsx,padsy);
346}
347
348//_____________________________________________________________________________
349void AliRICH::AddRecHit3D(Int_t id, Float_t *rechit)
2e5f0f7b 350{
237c933d 351
352 //
353 // Add a RICH reconstructed hit to the list
354 //
fe4da5cc 355
a4622d0f 356 TClonesArray &lrec3D = *((TClonesArray*)(*fRecHits3D)[id]);
357 new(lrec3D[fNrechits3D[id]++]) AliRICHRecHit3D(id,rechit);
2e5f0f7b 358}
fe4da5cc 359
ddae0931 360//___________________________________________
361void AliRICH::BuildGeometry()
362
fe4da5cc 363{
237c933d 364
365 //
366 // Builds a TNode geometry for event display
367 //
53d98323 368 TNode *node, *subnode, *top;
ddae0931 369
53d98323 370 const int kColorRICH = kRed;
ddae0931 371 //
237c933d 372 top=gAlice->GetGeometry()->GetNode("alice");
53d98323 373
374 AliRICH *pRICH = (AliRICH *) gAlice->GetDetector("RICH");
91975aa2 375 AliRICHSegmentationV0* segmentation;
53d98323 376 AliRICHChamber* iChamber;
377
378 iChamber = &(pRICH->Chamber(0));
91975aa2 379 segmentation=(AliRICHSegmentationV0*) iChamber->GetSegmentationModel(0);
ddae0931 380
381 new TBRIK("S_RICH","S_RICH","void",71.09999,11.5,73.15);
53d98323 382
91975aa2 383 Float_t csi_length = segmentation->Npx()*segmentation->Dpx() + segmentation->DeadZone();
384 Float_t csi_width = segmentation->Npy()*segmentation->Dpy() + 2*segmentation->DeadZone();
385
386 Float_t padplane_width = (csi_width - 2*segmentation->DeadZone())/3;
387 Float_t padplane_length = (csi_length - segmentation->DeadZone())/2;
388
389 //if (segmentation->GetPadPlaneWidth()>0)
390 //{
391 new TBRIK("PHOTO","PHOTO","void", padplane_width/2,.1,padplane_length/2);
392 //printf("\n\n\n\n\n Padplane w: %f l: %f \n\n\n\n\n", padplane_width/2,padplane_length/2);
393 //printf("\n\n\n\n\n Padplane w: %f l: %f \n\n\n\n\n", segmentation->GetPadPlaneWidth(), segmentation->GetPadPlaneLength());
394 //}
395
237c933d 396 top->cd();
ddae0931 397 Float_t pos1[3]={0,471.8999,165.2599};
2e5f0f7b 398 //Chamber(0).SetChamberTransform(pos1[0],pos1[1],pos1[2],
399 new TRotMatrix("rot993","rot993",90,0,70.69,90,19.30999,-90);
237c933d 400 node = new TNode("RICH1","RICH1","S_RICH",pos1[0],pos1[1],pos1[2],"rot993");
237c933d 401 node->SetLineColor(kColorRICH);
91975aa2 402 node->cd();
403 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
404 subnode->SetLineColor(kGreen);
405 fNodes->Add(subnode);
406 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,padplane_length/2 + segmentation->DeadZone()/2,"");
407 subnode->SetLineColor(kGreen);
408 fNodes->Add(subnode);
409 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
410 subnode->SetLineColor(kGreen);
411 fNodes->Add(subnode);
412 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
413 subnode->SetLineColor(kGreen);
414 fNodes->Add(subnode);
415 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,-padplane_length/2 - segmentation->DeadZone()/2,"");
416 subnode->SetLineColor(kGreen);
417 fNodes->Add(subnode);
418 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
419 subnode->SetLineColor(kGreen);
420 fNodes->Add(subnode);
237c933d 421 fNodes->Add(node);
53d98323 422
423
91975aa2 424 top->cd();
ddae0931 425 Float_t pos2[3]={171,470,0};
2e5f0f7b 426 //Chamber(1).SetChamberTransform(pos2[0],pos2[1],pos2[2],
427 new TRotMatrix("rot994","rot994",90,-20,90,70,0,0);
237c933d 428 node = new TNode("RICH2","RICH2","S_RICH",pos2[0],pos2[1],pos2[2],"rot994");
237c933d 429 node->SetLineColor(kColorRICH);
91975aa2 430 node->cd();
431 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
432 subnode->SetLineColor(kGreen);
433 fNodes->Add(subnode);
434 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,padplane_length/2 + segmentation->DeadZone()/2,"");
435 subnode->SetLineColor(kGreen);
436 fNodes->Add(subnode);
437 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
438 subnode->SetLineColor(kGreen);
439 fNodes->Add(subnode);
440 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
441 subnode->SetLineColor(kGreen);
442 fNodes->Add(subnode);
443 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,-padplane_length/2 - segmentation->DeadZone()/2,"");
444 subnode->SetLineColor(kGreen);
445 fNodes->Add(subnode);
446 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
447 subnode->SetLineColor(kGreen);
448 fNodes->Add(subnode);
237c933d 449 fNodes->Add(node);
53d98323 450
451
237c933d 452 top->cd();
ddae0931 453 Float_t pos3[3]={0,500,0};
2e5f0f7b 454 //Chamber(2).SetChamberTransform(pos3[0],pos3[1],pos3[2],
455 new TRotMatrix("rot995","rot995",90,0,90,90,0,0);
237c933d 456 node = new TNode("RICH3","RICH3","S_RICH",pos3[0],pos3[1],pos3[2],"rot995");
237c933d 457 node->SetLineColor(kColorRICH);
91975aa2 458 node->cd();
459 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
460 subnode->SetLineColor(kGreen);
461 fNodes->Add(subnode);
462 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,padplane_length/2 + segmentation->DeadZone()/2,"");
463 subnode->SetLineColor(kGreen);
464 fNodes->Add(subnode);
465 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
466 subnode->SetLineColor(kGreen);
467 fNodes->Add(subnode);
468 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
469 subnode->SetLineColor(kGreen);
470 fNodes->Add(subnode);
471 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,-padplane_length/2 - segmentation->DeadZone()/2,"");
472 subnode->SetLineColor(kGreen);
473 fNodes->Add(subnode);
474 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
475 subnode->SetLineColor(kGreen);
476 fNodes->Add(subnode);
237c933d 477 fNodes->Add(node);
53d98323 478
237c933d 479 top->cd();
ddae0931 480 Float_t pos4[3]={-171,470,0};
2e5f0f7b 481 //Chamber(3).SetChamberTransform(pos4[0],pos4[1],pos4[2],
482 new TRotMatrix("rot996","rot996",90,20,90,110,0,0);
237c933d 483 node = new TNode("RICH4","RICH4","S_RICH",pos4[0],pos4[1],pos4[2],"rot996");
237c933d 484 node->SetLineColor(kColorRICH);
91975aa2 485 node->cd();
486 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
487 subnode->SetLineColor(kGreen);
488 fNodes->Add(subnode);
489 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,padplane_length/2 + segmentation->DeadZone()/2,"");
490 subnode->SetLineColor(kGreen);
491 fNodes->Add(subnode);
492 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
493 subnode->SetLineColor(kGreen);
494 fNodes->Add(subnode);
495 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
496 subnode->SetLineColor(kGreen);
497 fNodes->Add(subnode);
498 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,-padplane_length/2 - segmentation->DeadZone()/2,"");
499 subnode->SetLineColor(kGreen);
500 fNodes->Add(subnode);
501 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
502 subnode->SetLineColor(kGreen);
503 fNodes->Add(subnode);
237c933d 504 fNodes->Add(node);
53d98323 505
506
237c933d 507 top->cd();
ddae0931 508 Float_t pos5[3]={161.3999,443.3999,-165.3};
2e5f0f7b 509 //Chamber(4).SetChamberTransform(pos5[0],pos5[1],pos5[2],
510 new TRotMatrix("rot997","rot997",90,340,108.1999,70,18.2,70);
237c933d 511 node = new TNode("RICH5","RICH5","S_RICH",pos5[0],pos5[1],pos5[2],"rot997");
237c933d 512 node->SetLineColor(kColorRICH);
91975aa2 513 node->cd();
514 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
515 subnode->SetLineColor(kGreen);
516 fNodes->Add(subnode);
517 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,padplane_length/2 + segmentation->DeadZone()/2,"");
518 subnode->SetLineColor(kGreen);
519 fNodes->Add(subnode);
520 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
521 subnode->SetLineColor(kGreen);
522 fNodes->Add(subnode);
523 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
524 subnode->SetLineColor(kGreen);
525 fNodes->Add(subnode);
526 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,-padplane_length/2 - segmentation->DeadZone()/2,"");
527 subnode->SetLineColor(kGreen);
528 fNodes->Add(subnode);
529 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
530 subnode->SetLineColor(kGreen);
531 fNodes->Add(subnode);
237c933d 532 fNodes->Add(node);
53d98323 533
534
237c933d 535 top->cd();
ddae0931 536 Float_t pos6[3]={0., 471.9, -165.3,};
2e5f0f7b 537 //Chamber(5).SetChamberTransform(pos6[0],pos6[1],pos6[2],
538 new TRotMatrix("rot998","rot998",90,0,109.3099,90,19.30999,90);
237c933d 539 node = new TNode("RICH6","RICH6","S_RICH",pos6[0],pos6[1],pos6[2],"rot998");
237c933d 540 node->SetLineColor(kColorRICH);
91975aa2 541
542 fNodes->Add(node);node->cd();
543 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
544 subnode->SetLineColor(kGreen);
545 fNodes->Add(subnode);
546 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,padplane_length/2 + segmentation->DeadZone()/2,"");
547 subnode->SetLineColor(kGreen);
548 fNodes->Add(subnode);
549 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
550 subnode->SetLineColor(kGreen);
551 fNodes->Add(subnode);
552 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
553 subnode->SetLineColor(kGreen);
554 fNodes->Add(subnode);
555 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,-padplane_length/2 - segmentation->DeadZone()/2,"");
556 subnode->SetLineColor(kGreen);
557 fNodes->Add(subnode);
558 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
559 subnode->SetLineColor(kGreen);
560 fNodes->Add(subnode);
53d98323 561
562
237c933d 563 top->cd();
ddae0931 564 Float_t pos7[3]={-161.399,443.3999,-165.3};
2e5f0f7b 565 //Chamber(6).SetChamberTransform(pos7[0],pos7[1],pos7[2],
566 new TRotMatrix("rot999","rot999",90,20,108.1999,110,18.2,110);
237c933d 567 node = new TNode("RICH7","RICH7","S_RICH",pos7[0],pos7[1],pos7[2],"rot999");
568 node->SetLineColor(kColorRICH);
91975aa2 569 node->cd();
570 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
571 subnode->SetLineColor(kGreen);
572 fNodes->Add(subnode);
573 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,padplane_length/2 + segmentation->DeadZone()/2,"");
574 subnode->SetLineColor(kGreen);
575 fNodes->Add(subnode);
576 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,padplane_length/2 + segmentation->DeadZone()/2,"");
577 subnode->SetLineColor(kGreen);
578 fNodes->Add(subnode);
579 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",padplane_width + segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
580 subnode->SetLineColor(kGreen);
581 fNodes->Add(subnode);
582 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",0,0,-padplane_length/2 - segmentation->DeadZone()/2,"");
583 subnode->SetLineColor(kGreen);
584 fNodes->Add(subnode);
585 subnode = new TNode("PHOTO1","PHOTO1","PHOTO",-padplane_width - segmentation->DeadZone(),0,-padplane_length/2 - segmentation->DeadZone()/2,"");
586 subnode->SetLineColor(kGreen);
587 fNodes->Add(subnode);
237c933d 588 fNodes->Add(node);
ddae0931 589
fe4da5cc 590}
591
452a64c6 592//___________________________________________
593void AliRICH::CreateGeometry()
594{
595 //
596 // Create the geometry for RICH version 1
597 //
598 // Modified by: N. Colonna (INFN - BARI, Nicola.Colonna@ba.infn.it)
599 // R.A. Fini (INFN - BARI, Rosanna.Fini@ba.infn.it)
600 // R.A. Loconsole (Bari University, loco@riscom.ba.infn.it)
601 //
602 //Begin_Html
603 /*
604 <img src="picts/AliRICHv1.gif">
605 */
606 //End_Html
607 //Begin_Html
608 /*
609 <img src="picts/AliRICHv1Tree.gif">
610 */
611 //End_Html
612
613 AliRICH *pRICH = (AliRICH *) gAlice->GetDetector("RICH");
a2f7eaf6 614 AliSegmentation* segmentation;
452a64c6 615 AliRICHGeometry* geometry;
616 AliRICHChamber* iChamber;
617
618 iChamber = &(pRICH->Chamber(0));
619 segmentation=iChamber->GetSegmentationModel(0);
620 geometry=iChamber->GetGeometryModel();
621
622 Float_t distance;
623 distance = geometry->GetFreonThickness()/2 + geometry->GetQuartzThickness() + geometry->GetGapThickness();
624 geometry->SetRadiatorToPads(distance);
625
3587e146 626 //Opaque quartz thickness
e293afae 627 Float_t oqua_thickness = .5;
683b273b 628 //CsI dimensions
629
53d98323 630 //Float_t csi_length = 160*.8 + 2.6;
631 //Float_t csi_width = 144*.84 + 2*2.6;
632
633 Float_t deadzone=2.6;
634
635 Float_t csi_length = segmentation->Npx()*segmentation->Dpx() + deadzone;
636 Float_t csi_width = segmentation->Npy()*segmentation->Dpy() + 2*deadzone;
637
638
639 Int_t *idtmed = fIdtmed->GetArray()-999;
452a64c6 640
641 Int_t i;
642 Float_t zs;
643 Int_t idrotm[1099];
644 Float_t par[3];
645
646 // --- Define the RICH detector
647 // External aluminium box
e293afae 648 par[0] = 68.8;
6197ac87 649 par[1] = 13; //Original Settings
e293afae 650 par[2] = 70.86;
452a64c6 651 /*par[0] = 73.15;
652 par[1] = 11.5;
653 par[2] = 71.1;*/
654 gMC->Gsvolu("RICH", "BOX ", idtmed[1009], par, 3);
655
683b273b 656 // Air
657 par[0] = 66.3;
6197ac87 658 par[1] = 13; //Original Settings
683b273b 659 par[2] = 68.35;
452a64c6 660 /*par[0] = 66.55;
661 par[1] = 11.5;
662 par[2] = 64.8;*/
663 gMC->Gsvolu("SRIC", "BOX ", idtmed[1000], par, 3);
664
6197ac87 665 // Air 2 (cutting the lower part of the box)
666
667 par[0] = 1.25;
668 par[1] = 3; //Original Settings
669 par[2] = 70.86;
670 gMC->Gsvolu("AIR2", "BOX ", idtmed[1000], par, 3);
671
672 // Air 3 (cutting the lower part of the box)
673
674 par[0] = 66.3;
675 par[1] = 3; //Original Settings
676 par[2] = 1.2505;
677 gMC->Gsvolu("AIR3", "BOX ", idtmed[1000], par, 3);
678
452a64c6 679 // Honeycomb
e293afae 680 par[0] = 66.3;
452a64c6 681 par[1] = .188; //Original Settings
e293afae 682 par[2] = 68.35;
452a64c6 683 /*par[0] = 66.55;
684 par[1] = .188;
685 par[2] = 63.1;*/
686 gMC->Gsvolu("HONE", "BOX ", idtmed[1001], par, 3);
687
688 // Aluminium sheet
e293afae 689 par[0] = 66.3;
452a64c6 690 par[1] = .025; //Original Settings
e293afae 691 par[2] = 68.35;
452a64c6 692 /*par[0] = 66.5;
693 par[1] = .025;
694 par[2] = 63.1;*/
695 gMC->Gsvolu("ALUM", "BOX ", idtmed[1009], par, 3);
696
697 // Quartz
698 par[0] = geometry->GetQuartzWidth()/2;
699 par[1] = geometry->GetQuartzThickness()/2;
700 par[2] = geometry->GetQuartzLength()/2;
701 /*par[0] = 63.1;
702 par[1] = .25; //Original Settings
703 par[2] = 65.5;*/
704 /*par[0] = geometry->GetQuartzWidth()/2;
705 par[1] = geometry->GetQuartzThickness()/2;
706 par[2] = geometry->GetQuartzLength()/2;*/
707 //printf("\n\n\n\n\n\n\n\\n\n\n\n Gap Thickness: %f %f %f\n\n\n\n\n\n\n\n\n\n\n\n\n\n",par[0],par[1],par[2]);
708 gMC->Gsvolu("QUAR", "BOX ", idtmed[1002], par, 3);
709
710 // Spacers (cylinders)
711 par[0] = 0.;
712 par[1] = .5;
713 par[2] = geometry->GetFreonThickness()/2;
714 gMC->Gsvolu("SPAC", "TUBE", idtmed[1002], par, 3);
715
683b273b 716 // Feet (freon slabs supports)
717
718 par[0] = .7;
719 par[1] = .3;
720 par[2] = 1.9;
721 gMC->Gsvolu("FOOT", "BOX", idtmed[1009], par, 3);
722
452a64c6 723 // Opaque quartz
e293afae 724 par[0] = geometry->GetQuartzWidth()/2;
725 par[1] = .2;
726 par[2] = geometry->GetQuartzLength()/2;
727 /*par[0] = 61.95;
452a64c6 728 par[1] = .2; //Original Settings
e293afae 729 par[2] = 66.5;*/
452a64c6 730 /*par[0] = 66.5;
731 par[1] = .2;
732 par[2] = 61.95;*/
733 gMC->Gsvolu("OQUA", "BOX ", idtmed[1007], par, 3);
734
735 // Frame of opaque quartz
53a60579 736 par[0] = geometry->GetOuterFreonWidth()/2;
737 //+ oqua_thickness;
452a64c6 738 par[1] = geometry->GetFreonThickness()/2;
53a60579 739 par[2] = geometry->GetOuterFreonLength()/2;
740 //+ oqua_thickness;
452a64c6 741 /*par[0] = 20.65;
742 par[1] = .5; //Original Settings
743 par[2] = 66.5;*/
744 /*par[0] = 66.5;
745 par[1] = .5;
746 par[2] = 20.65;*/
747 gMC->Gsvolu("OQF1", "BOX ", idtmed[1007], par, 3);
748
e293afae 749 par[0] = geometry->GetInnerFreonWidth()/2;
452a64c6 750 par[1] = geometry->GetFreonThickness()/2;
e293afae 751 par[2] = geometry->GetInnerFreonLength()/2;
452a64c6 752 gMC->Gsvolu("OQF2", "BOX ", idtmed[1007], par, 3);
753
754 // Little bar of opaque quartz
e293afae 755 //par[0] = .275;
756 //par[1] = geometry->GetQuartzThickness()/2;
53a60579 757 //par[2] = geometry->GetInnerFreonLength()/2 - 2.4;
e293afae 758 //par[2] = geometry->GetInnerFreonLength()/2;
53a60579 759 //+ oqua_thickness;
452a64c6 760 /*par[0] = .275;
761 par[1] = .25; //Original Settings
762 par[2] = 63.1;*/
763 /*par[0] = 63.1;
764 par[1] = .25;
765 par[2] = .275;*/
e293afae 766 //gMC->Gsvolu("BARR", "BOX ", idtmed[1007], par, 3);
452a64c6 767
768 // Freon
e293afae 769 par[0] = geometry->GetOuterFreonWidth()/2 - oqua_thickness;
452a64c6 770 par[1] = geometry->GetFreonThickness()/2;
e293afae 771 par[2] = geometry->GetOuterFreonLength()/2 - 2*oqua_thickness;
452a64c6 772 /*par[0] = 20.15;
773 par[1] = .5; //Original Settings
774 par[2] = 65.5;*/
775 /*par[0] = 65.5;
776 par[1] = .5;
777 par[2] = 20.15;*/
778 gMC->Gsvolu("FRE1", "BOX ", idtmed[1003], par, 3);
779
e293afae 780 par[0] = geometry->GetInnerFreonWidth()/2 - oqua_thickness;
452a64c6 781 par[1] = geometry->GetFreonThickness()/2;
e293afae 782 par[2] = geometry->GetInnerFreonLength()/2 - 2*oqua_thickness;
452a64c6 783 gMC->Gsvolu("FRE2", "BOX ", idtmed[1003], par, 3);
784
785 // Methane
e293afae 786 //par[0] = 64.8;
683b273b 787 par[0] = csi_width/2;
452a64c6 788 par[1] = geometry->GetGapThickness()/2;
789 //printf("\n\n\n\n\n\n\n\\n\n\n\n Gap Thickness: %f\n\n\n\n\n\n\n\n\n\n\n\n\n\n",par[1]);
e293afae 790 //par[2] = 64.8;
683b273b 791 par[2] = csi_length/2;
452a64c6 792 gMC->Gsvolu("META", "BOX ", idtmed[1004], par, 3);
793
794 // Methane gap
e293afae 795 //par[0] = 64.8;
683b273b 796 par[0] = csi_width/2;
452a64c6 797 par[1] = geometry->GetProximityGapThickness()/2;
798 //printf("\n\n\n\n\n\n\n\\n\n\n\n Gap Thickness: %f\n\n\n\n\n\n\n\n\n\n\n\n\n\n",par[1]);
e293afae 799 //par[2] = 64.8;
683b273b 800 par[2] = csi_length/2;
452a64c6 801 gMC->Gsvolu("GAP ", "BOX ", idtmed[1008], par, 3);
802
803 // CsI photocathode
e293afae 804 //par[0] = 64.8;
683b273b 805 par[0] = csi_width/2;
452a64c6 806 par[1] = .25;
e293afae 807 //par[2] = 64.8;
683b273b 808 par[2] = csi_length/2;
452a64c6 809 gMC->Gsvolu("CSI ", "BOX ", idtmed[1005], par, 3);
810
811 // Anode grid
812 par[0] = 0.;
813 par[1] = .001;
814 par[2] = 20.;
815 gMC->Gsvolu("GRID", "TUBE", idtmed[1006], par, 3);
683b273b 816
6197ac87 817 // Wire supports
818 // Bar of metal
819
820 par[0] = csi_width/2;
821 par[1] = 1.05;
822 par[2] = 1.05;
823 gMC->Gsvolu("WSMe", "BOX ", idtmed[1009], par, 3);
824
825 // Ceramic pick up (base)
826
827 par[0] = csi_width/2;
828 par[1] = .25;
829 par[2] = 1.05;
830 gMC->Gsvolu("WSG1", "BOX ", idtmed[1010], par, 3);
831
832 // Ceramic pick up (head)
833
834 par[0] = csi_width/2;
835 par[1] = .1;
836 par[2] = .1;
837 gMC->Gsvolu("WSG2", "BOX ", idtmed[1010], par, 3);
838
683b273b 839 // Aluminium supports for methane and CsI
840 // Short bar
841
842 par[0] = csi_width/2;
843 par[1] = geometry->GetGapThickness()/2 + .25;
844 par[2] = (68.35 - csi_length/2)/2;
845 gMC->Gsvolu("SMSH", "BOX", idtmed[1009], par, 3);
452a64c6 846
683b273b 847 // Long bar
848
849 par[0] = (66.3 - csi_width/2)/2;
850 par[1] = geometry->GetGapThickness()/2 + .25;
851 par[2] = csi_length/2 + 68.35 - csi_length/2;
852 gMC->Gsvolu("SMLG", "BOX", idtmed[1009], par, 3);
853
854 // Aluminium supports for freon
855 // Short bar
856
857 par[0] = geometry->GetQuartzWidth()/2;
858 par[1] = .3;
859 par[2] = (68.35 - geometry->GetQuartzLength()/2)/2;
860 gMC->Gsvolu("SFSH", "BOX", idtmed[1009], par, 3);
861
862 // Long bar
863
864 par[0] = (66.3 - geometry->GetQuartzWidth()/2)/2;
865 par[1] = .3;
866 par[2] = geometry->GetQuartzLength()/2 + 68.35 - geometry->GetQuartzLength()/2;
867 gMC->Gsvolu("SFLG", "BOX", idtmed[1009], par, 3);
868
6197ac87 869 // PCB backplane
683b273b 870
6197ac87 871 par[0] = csi_width/2;
872 par[1] = .25;
873 par[2] = csi_length/4 -.5025;
874 gMC->Gsvolu("PCB ", "BOX", idtmed[1011], par, 3);
683b273b 875
6197ac87 876
877 // Backplane supports
878
879 // Aluminium slab
880
881 par[0] = 33.15;
882 par[1] = 2;
883 par[2] = 21.65;
884 gMC->Gsvolu("BACK", "BOX", idtmed[1009], par, 3);
885
886 // Big hole
887
888 par[0] = 9.05;
889 par[1] = 2;
890 par[2] = 4.4625;
891 gMC->Gsvolu("BKHL", "BOX", idtmed[1000], par, 3);
892
893 // Small hole
894
895 par[0] = 5.7;
896 par[1] = 2;
897 par[2] = 4.4625;
898 gMC->Gsvolu("BKHS", "BOX", idtmed[1000], par, 3);
899
900 // Place holes inside backplane support
901
902 gMC->Gspos("BKHS", 1, "BACK", .8 + 5.7,0., .6 + 4.4625, 0, "ONLY");
903 gMC->Gspos("BKHS", 2, "BACK", -.8 - 5.7,0., .6 + 4.4625, 0, "ONLY");
904 gMC->Gspos("BKHS", 3, "BACK", .8 + 5.7,0., -.6 - 4.4625, 0, "ONLY");
905 gMC->Gspos("BKHS", 4, "BACK", -.8 - 5.7,0., -.6 - 4.4625, 0, "ONLY");
906 gMC->Gspos("BKHS", 5, "BACK", .8 + 5.7,0., .6 + 8.925 + 1.2 + 4.4625, 0, "ONLY");
907 gMC->Gspos("BKHS", 6, "BACK", -.8 - 5.7,0., .6 + 8.925 + 1.2 + 4.4625, 0, "ONLY");
908 gMC->Gspos("BKHS", 7, "BACK", .8 + 5.7,0., -.6 - 8.925 - 1.2 - 4.4625, 0, "ONLY");
909 gMC->Gspos("BKHS", 8, "BACK", -.8 - 5.7,0., -.6 - 8.925 - 1.2 - 4.4625, 0, "ONLY");
910 gMC->Gspos("BKHL", 1, "BACK", .8 + 11.4 + 1.6 + 9.05, 0., .6 + 4.4625, 0, "ONLY");
911 gMC->Gspos("BKHL", 2, "BACK", -.8 - 11.4 - 1.6 - 9.05, 0., .6 + 4.4625, 0, "ONLY");
912 gMC->Gspos("BKHL", 3, "BACK", .8 + 11.4 + 1.6 + 9.05, 0., -.6 - 4.4625, 0, "ONLY");
913 gMC->Gspos("BKHL", 4, "BACK", -.8 - 11.4 - 1.6 - 9.05, 0., -.6 - 4.4625, 0, "ONLY");
914 gMC->Gspos("BKHL", 5, "BACK", .8 + 11.4+ 1.6 + 9.05, 0., .6 + 8.925 + 1.2 + 4.4625, 0, "ONLY");
915 gMC->Gspos("BKHL", 6, "BACK", -.8 - 11.4 - 1.6 - 9.05, 0., .6 + 8.925 + 1.2 + 4.4625, 0, "ONLY");
916 gMC->Gspos("BKHL", 7, "BACK", .8 + 11.4 + 1.6 + 9.05, 0., -.6 - 8.925 - 1.2 - 4.4625, 0, "ONLY");
917 gMC->Gspos("BKHL", 8, "BACK", -.8 - 11.4 - 1.6 - 9.05, 0., -.6 - 8.925 - 1.2 - 4.4625, 0, "ONLY");
918
919
920
452a64c6 921 // --- Places the detectors defined with GSVOLU
922 // Place material inside RICH
6197ac87 923 gMC->Gspos("SRIC", 1, "RICH", 0.,0., 0., 0, "ONLY");
924 gMC->Gspos("AIR2", 1, "RICH", 66.3 + 1.2505, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .05 - .376 -.5 - 3.5, 0., 0, "ONLY");
925 gMC->Gspos("AIR2", 2, "RICH", -66.3 - 1.2505, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .05 - .376 -.5 - 3.5, 0., 0, "ONLY");
926 gMC->Gspos("AIR3", 1, "RICH", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .05 - .376 -.5 - 3.5, -68.35 - 1.25, 0, "ONLY");
927 gMC->Gspos("AIR3", 2, "RICH", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .05 - .376 -.5 - 3.5, 68.35 + 1.25, 0, "ONLY");
928
683b273b 929
930 gMC->Gspos("ALUM", 1, "SRIC", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .05 - .376 -.025, 0., 0, "ONLY");
931 gMC->Gspos("HONE", 1, "SRIC", 0., 1.276- geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .05 - .188, 0., 0, "ONLY");
932 gMC->Gspos("ALUM", 2, "SRIC", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .025, 0., 0, "ONLY");
933 gMC->Gspos("FOOT", 1, "SRIC", 64.95, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, 36.9, 0, "ONLY");
934 gMC->Gspos("FOOT", 2, "SRIC", 21.65, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3 , 36.9, 0, "ONLY");
935 gMC->Gspos("FOOT", 3, "SRIC", -21.65, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, 36.9, 0, "ONLY");
936 gMC->Gspos("FOOT", 4, "SRIC", -64.95, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, 36.9, 0, "ONLY");
937 gMC->Gspos("FOOT", 5, "SRIC", 64.95, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, -36.9, 0, "ONLY");
938 gMC->Gspos("FOOT", 6, "SRIC", 21.65, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, -36.9, 0, "ONLY");
939 gMC->Gspos("FOOT", 7, "SRIC", -21.65, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, -36.9, 0, "ONLY");
940 gMC->Gspos("FOOT", 8, "SRIC", -64.95, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, -36.9, 0, "ONLY");
452a64c6 941 gMC->Gspos("OQUA", 1, "SRIC", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .2, 0., 0, "ONLY");
942
683b273b 943 // Supports placing
944
945 // Methane supports
946 gMC->Gspos("SMLG", 1, "SRIC", csi_width/2 + (66.3 - csi_width/2)/2, 1.276 + .25, 0., 0, "ONLY");
947 gMC->Gspos("SMLG", 2, "SRIC", - csi_width/2 - (66.3 - csi_width/2)/2, 1.276 + .25, 0., 0, "ONLY");
948 gMC->Gspos("SMSH", 1, "SRIC", 0., 1.276 + .25, csi_length/2 + (68.35 - csi_length/2)/2, 0, "ONLY");
949 gMC->Gspos("SMSH", 2, "SRIC", 0., 1.276 + .25, - csi_length/2 - (68.35 - csi_length/2)/2, 0, "ONLY");
950
951 //Freon supports
952
953 Float_t supp_y = 1.276 - geometry->GetGapThickness()/2- geometry->GetQuartzThickness() -geometry->GetFreonThickness() - .2 + .3; //y position of freon supports
954
955 gMC->Gspos("SFLG", 1, "SRIC", geometry->GetQuartzWidth()/2 + (66.3 - geometry->GetQuartzWidth()/2)/2, supp_y, 0., 0, "ONLY");
956 gMC->Gspos("SFLG", 2, "SRIC", - geometry->GetQuartzWidth()/2 - (66.3 - geometry->GetQuartzWidth()/2)/2, supp_y, 0., 0, "ONLY");
957 gMC->Gspos("SFSH", 1, "SRIC", 0., supp_y, geometry->GetQuartzLength()/2 + (68.35 - geometry->GetQuartzLength()/2)/2, 0, "ONLY");
958 gMC->Gspos("SFSH", 2, "SRIC", 0., supp_y, - geometry->GetQuartzLength()/2 - (68.35 - geometry->GetQuartzLength()/2)/2, 0, "ONLY");
959
452a64c6 960 AliMatrix(idrotm[1019], 0., 0., 90., 0., 90., 90.);
961
e293afae 962 //Int_t nspacers = (Int_t)(TMath::Abs(geometry->GetInnerFreonLength()/14.4));
963 Int_t nspacers = 9;
452a64c6 964 //printf("\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n Spacers:%d\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n",nspacers);
965
966 //printf("Nspacers: %d", nspacers);
967
968 //for (i = 1; i <= 9; ++i) {
969 //zs = (5 - i) * 14.4; //Original settings
970 for (i = 0; i < nspacers; i++) {
971 zs = (TMath::Abs(nspacers/2) - i) * 14.4;
972 gMC->Gspos("SPAC", i, "FRE1", 6.7, 0., zs, idrotm[1019], "ONLY"); //Original settings
973 //gMC->Gspos("SPAC", i, "FRE1", zs, 0., 6.7, idrotm[1019], "ONLY");
974 }
975 //for (i = 10; i <= 18; ++i) {
976 //zs = (14 - i) * 14.4; //Original settings
977 for (i = nspacers; i < nspacers*2; ++i) {
978 zs = (nspacers + TMath::Abs(nspacers/2) - i) * 14.4;
979 gMC->Gspos("SPAC", i, "FRE1", -6.7, 0., zs, idrotm[1019], "ONLY"); //Original settings
980 //gMC->Gspos("SPAC", i, "FRE1", zs, 0., -6.7, idrotm[1019], "ONLY");
981 }
982
983 //for (i = 1; i <= 9; ++i) {
984 //zs = (5 - i) * 14.4; //Original settings
985 for (i = 0; i < nspacers; i++) {
986 zs = (TMath::Abs(nspacers/2) - i) * 14.4;
987 gMC->Gspos("SPAC", i, "FRE2", 6.7, 0., zs, idrotm[1019], "ONLY"); //Original settings
988 //gMC->Gspos("SPAC", i, "FRE2", zs, 0., 6.7, idrotm[1019], "ONLY");
989 }
990 //for (i = 10; i <= 18; ++i) {
991 //zs = (5 - i) * 14.4; //Original settings
992 for (i = nspacers; i < nspacers*2; ++i) {
993 zs = (nspacers + TMath::Abs(nspacers/2) - i) * 14.4;
994 gMC->Gspos("SPAC", i, "FRE2", -6.7, 0., zs, idrotm[1019], "ONLY"); //Original settings
995 //gMC->Gspos("SPAC", i, "FRE2", zs, 0., -6.7, idrotm[1019], "ONLY");
996 }
997
998 /*gMC->Gspos("FRE1", 1, "OQF1", 0., 0., 0., 0, "ONLY");
999 gMC->Gspos("FRE2", 1, "OQF2", 0., 0., 0., 0, "ONLY");
1000 gMC->Gspos("OQF1", 1, "SRIC", 31.3, -4.724, 41.3, 0, "ONLY");
1001 gMC->Gspos("OQF2", 2, "SRIC", 0., -4.724, 0., 0, "ONLY");
1002 gMC->Gspos("OQF1", 3, "SRIC", -31.3, -4.724, -41.3, 0, "ONLY");
1003 gMC->Gspos("BARR", 1, "QUAR", -21.65, 0., 0., 0, "ONLY"); //Original settings
1004 gMC->Gspos("BARR", 2, "QUAR", 21.65, 0., 0., 0, "ONLY"); //Original settings
1005 gMC->Gspos("QUAR", 1, "SRIC", 0., -3.974, 0., 0, "ONLY");
1006 gMC->Gspos("GAP ", 1, "META", 0., 4.8, 0., 0, "ONLY");
1007 gMC->Gspos("META", 1, "SRIC", 0., 1.276, 0., 0, "ONLY");
1008 gMC->Gspos("CSI ", 1, "SRIC", 0., 6.526, 0., 0, "ONLY");*/
1009
1010
1011 gMC->Gspos("FRE1", 1, "OQF1", 0., 0., 0., 0, "ONLY");
1012 gMC->Gspos("FRE2", 1, "OQF2", 0., 0., 0., 0, "ONLY");
683b273b 1013 gMC->Gspos("OQF1", 1, "SRIC", geometry->GetOuterFreonWidth()/2 + geometry->GetInnerFreonWidth()/2 + 2, 1.276 - geometry->GetGapThickness()/2- geometry->GetQuartzThickness() -geometry->GetFreonThickness()/2, 0., 0, "ONLY"); //Original settings (31.3)
452a64c6 1014 gMC->Gspos("OQF2", 2, "SRIC", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()/2, 0., 0, "ONLY"); //Original settings
683b273b 1015 gMC->Gspos("OQF1", 3, "SRIC", - (geometry->GetOuterFreonWidth()/2 + geometry->GetInnerFreonWidth()/2) - 2, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()/2, 0., 0, "ONLY"); //Original settings (-31.3)
e293afae 1016 //gMC->Gspos("BARR", 1, "QUAR", - geometry->GetInnerFreonWidth()/2 - oqua_thickness, 0., 0., 0, "ONLY"); //Original settings (-21.65)
1017 //gMC->Gspos("BARR", 2, "QUAR", geometry->GetInnerFreonWidth()/2 + oqua_thickness, 0., 0., 0, "ONLY"); //Original settings (21.65)
452a64c6 1018 gMC->Gspos("QUAR", 1, "SRIC", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness()/2, 0., 0, "ONLY");
1019 gMC->Gspos("GAP ", 1, "META", 0., geometry->GetGapThickness()/2 - geometry->GetProximityGapThickness()/2 - 0.0001, 0., 0, "ONLY");
1020 gMC->Gspos("META", 1, "SRIC", 0., 1.276, 0., 0, "ONLY");
1021 gMC->Gspos("CSI ", 1, "SRIC", 0., 1.276 + geometry->GetGapThickness()/2 + .25, 0., 0, "ONLY");
6197ac87 1022
1023 // Wire support placing
1024
1025 gMC->Gspos("WSG2", 1, "GAP ", 0., geometry->GetProximityGapThickness()/2 - .1, 0., 0, "ONLY");
1026 gMC->Gspos("WSG1", 1, "CSI ", 0., 0., 0., 0, "ONLY");
1027 gMC->Gspos("WSMe", 1, "SRIC ", 0., 1.276 + geometry->GetGapThickness()/2 + .5 + 1.05, 0., 0, "ONLY");
1028
1029 // Backplane placing
1030
1031 gMC->Gspos("BACK", 1, "SRIC ", -33.15, 1.276 + geometry->GetGapThickness()/2 + .5 + 2.1 + 2, 43.3, 0, "ONLY");
1032 gMC->Gspos("BACK", 2, "SRIC ", 33.15, 1.276 + geometry->GetGapThickness()/2 + .5 + 2.1 + 2 , 43.3, 0, "ONLY");
1033 gMC->Gspos("BACK", 3, "SRIC ", -33.15, 1.276 + geometry->GetGapThickness()/2 + .5 + 2.1 + 2, 0., 0, "ONLY");
1034 gMC->Gspos("BACK", 4, "SRIC ", 33.15, 1.276 + geometry->GetGapThickness()/2 + .5 + 2.1 + 2, 0., 0, "ONLY");
1035 gMC->Gspos("BACK", 5, "SRIC ", 33.15, 1.276 + geometry->GetGapThickness()/2 + .5 + 2.1 + 2, -43.3, 0, "ONLY");
1036 gMC->Gspos("BACK", 6, "SRIC ", -33.15, 1.276 + geometry->GetGapThickness()/2 + .5 + 2.1 + 2, -43.3, 0, "ONLY");
1037
1038 // PCB placing
1039
1040 gMC->Gspos("PCB ", 1, "SRIC ", 0., 1.276 + geometry->GetGapThickness()/2 + .5 + 1.05, csi_width/4 + .5025 + 1.2, 0, "ONLY");
1041 gMC->Gspos("PCB ", 2, "SRIC ", 0., 1.276 + geometry->GetGapThickness()/2 + .5 + 1.05, -csi_width/4 - .5025 - 1.2, 0, "ONLY");
1042
1043
452a64c6 1044
1045 //printf("Position of the gap: %f to %f\n", 1.276 + geometry->GetGapThickness()/2 - geometry->GetProximityGapThickness()/2 - .2, 1.276 + geometry->GetGapThickness()/2 - geometry->GetProximityGapThickness()/2 + .2);
1046
1047 // Place RICH inside ALICE apparatus
1048
1049 AliMatrix(idrotm[1000], 90., 0., 70.69, 90., 19.31, -90.);
1050 AliMatrix(idrotm[1001], 90., -20., 90., 70., 0., 0.);
1051 AliMatrix(idrotm[1002], 90., 0., 90., 90., 0., 0.);
1052 AliMatrix(idrotm[1003], 90., 20., 90., 110., 0., 0.);
1053 AliMatrix(idrotm[1004], 90., 340., 108.2, 70., 18.2, 70.);
1054 AliMatrix(idrotm[1005], 90., 0., 109.31, 90., 19.31, 90.);
1055 AliMatrix(idrotm[1006], 90., 20., 108.2, 110., 18.2, 110.);
1056
1057 gMC->Gspos("RICH", 1, "ALIC", 0., 471.9, 165.26, idrotm[1000], "ONLY");
1058 gMC->Gspos("RICH", 2, "ALIC", 171., 470., 0., idrotm[1001], "ONLY");
1059 gMC->Gspos("RICH", 3, "ALIC", 0., 500., 0., idrotm[1002], "ONLY");
1060 gMC->Gspos("RICH", 4, "ALIC", -171., 470., 0., idrotm[1003], "ONLY");
1061 gMC->Gspos("RICH", 5, "ALIC", 161.4, 443.4, -165.3, idrotm[1004], "ONLY");
1062 gMC->Gspos("RICH", 6, "ALIC", 0., 471.9, -165.3, idrotm[1005], "ONLY");
1063 gMC->Gspos("RICH", 7, "ALIC", -161.4, 443.4, -165.3, idrotm[1006], "ONLY");
1064
1065}
1066
1067
1068//___________________________________________
1069void AliRICH::CreateMaterials()
1070{
1071 //
1072 // *** DEFINITION OF AVAILABLE RICH MATERIALS ***
1073 // ORIGIN : NICK VAN EIJNDHOVEN
1074 // Modified by: N. Colonna (INFN - BARI, Nicola.Colonna@ba.infn.it)
1075 // R.A. Fini (INFN - BARI, Rosanna.Fini@ba.infn.it)
1076 // R.A. Loconsole (Bari University, loco@riscom.ba.infn.it)
1077 //
1078 Int_t isxfld = gAlice->Field()->Integ();
1079 Float_t sxmgmx = gAlice->Field()->Max();
1080 Int_t i;
1081
1082 /************************************Antonnelo's Values (14-vectors)*****************************************/
1083 /*
1084 Float_t ppckov[14] = { 5.63e-9,5.77e-9,5.9e-9,6.05e-9,6.2e-9,6.36e-9,6.52e-9,
1085 6.7e-9,6.88e-9,7.08e-9,7.3e-9,7.51e-9,7.74e-9,8e-9 };
1086 Float_t rIndexQuarz[14] = { 1.528309,1.533333,
1087 1.538243,1.544223,1.550568,1.55777,
1088 1.565463,1.574765,1.584831,1.597027,
1089 1.611858,1.6277,1.6472,1.6724 };
1090 Float_t rIndexOpaqueQuarz[14] = { 1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1. };
1091 Float_t rIndexMethane[14] = { 1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1. };
1092 Float_t rIndexGrid[14] = { 1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1. };
1093 Float_t abscoFreon[14] = { 179.0987,179.0987,
1094 179.0987,179.0987,179.0987,142.92,56.65,13.95,10.43,7.07,2.03,.5773,.33496,0. };
1095 //Float_t abscoFreon[14] = { 1e-5,1e-5,1e-5,1e-5,1e-5,1e-5,1e-5,1e-5,1e-5,
1096 // 1e-5,1e-5,1e-5,1e-5,1e-5 };
1097 Float_t abscoQuarz[14] = { 64.035,39.98,35.665,31.262,27.527,22.815,21.04,17.52,
1098 14.177,9.282,4.0925,1.149,.3627,.10857 };
1099 Float_t abscoOpaqueQuarz[14] = { 1e-5,1e-5,1e-5,1e-5,1e-5,1e-5,1e-5,1e-5,1e-5,
1100 1e-5,1e-5,1e-5,1e-5,1e-5 };
1101 Float_t abscoCsI[14] = { 1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,
1102 1e-4,1e-4,1e-4,1e-4 };
1103 Float_t abscoMethane[14] = { 1e6,1e6,1e6,1e6,1e6,1e6,1e6,1e6,1e6,1e6,1e6,
1104 1e6,1e6,1e6 };
1105 Float_t abscoGrid[14] = { 1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,
1106 1e-4,1e-4,1e-4,1e-4 };
1107 Float_t efficAll[14] = { 1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1. };
1108 Float_t efficCsI[14] = { 6e-4,.005,.0075,.01125,.045,.117,.135,.16575,
1109 .17425,.1785,.1836,.1904,.1938,.221 };
1110 Float_t efficGrid[14] = { 1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1. };
1111 */
1112
1113
1114 /**********************************End of Antonnelo's Values**********************************/
1115
1116 /**********************************Values from rich_media.f (31-vectors)**********************************/
1117
1118
1119 //Photons energy intervals
1120 Float_t ppckov[26];
1121 for (i=0;i<26;i++)
1122 {
1123 ppckov[i] = (Float_t(i)*0.1+5.5)*1e-9;
1124 //printf ("Energy intervals: %e\n",ppckov[i]);
1125 }
1126
1127
1128 //Refraction index for quarz
1129 Float_t rIndexQuarz[26];
1130 Float_t e1= 10.666;
1131 Float_t e2= 18.125;
1132 Float_t f1= 46.411;
1133 Float_t f2= 228.71;
1134 for (i=0;i<26;i++)
1135 {
1136 Float_t ene=ppckov[i]*1e9;
1137 Float_t a=f1/(e1*e1 - ene*ene);
1138 Float_t b=f2/(e2*e2 - ene*ene);
1139 rIndexQuarz[i] = TMath::Sqrt(1. + a + b );
1140 //printf ("rIndexQuarz: %e\n",rIndexQuarz[i]);
1141 }
1142
1143 //Refraction index for opaque quarz, methane and grid
1144 Float_t rIndexOpaqueQuarz[26];
1145 Float_t rIndexMethane[26];
1146 Float_t rIndexGrid[26];
1147 for (i=0;i<26;i++)
1148 {
1149 rIndexOpaqueQuarz[i]=1;
1150 rIndexMethane[i]=1.000444;
1151 rIndexGrid[i]=1;
1152 //printf ("rIndexOpaqueQuarz , etc: %e, %e, %e\n",rIndexOpaqueQuarz[i], rIndexMethane[i], rIndexGrid[i]=1);
1153 }
1154
1155 //Absorption index for freon
1156 Float_t abscoFreon[26] = {179.0987, 179.0987, 179.0987, 179.0987, 179.0987, 179.0987, 179.0987, 179.0987,
1157 179.0987, 142.9206, 56.64957, 25.58622, 13.95293, 12.03905, 10.42953, 8.804196,
1158 7.069031, 4.461292, 2.028366, 1.293013, .577267, .40746, .334964, 0., 0., 0.};
1159
1160 //Absorption index for quarz
1161 /*Float_t Qzt [21] = {.0,.0,.005,.04,.35,.647,.769,.808,.829,.844,.853,.858,.869,.887,.903,.902,.902,
1162 .906,.907,.907,.907};
1163 Float_t Wavl2[] = {150.,155.,160.0,165.0,170.0,175.0,180.0,185.0,190.0,195.0,200.0,205.0,210.0,
1164 215.0,220.0,225.0,230.0,235.0,240.0,245.0,250.0};
1165 Float_t abscoQuarz[31];
1166 for (Int_t i=0;i<31;i++)
1167 {
1168 Float_t Xlam = 1237.79 / (ppckov[i]*1e9);
1169 if (Xlam <= 160) abscoQuarz[i] = 0;
1170 if (Xlam > 250) abscoQuarz[i] = 1;
1171 else
1172 {
1173 for (Int_t j=0;j<21;j++)
1174 {
1175 //printf ("Passed\n");
1176 if (Xlam > Wavl2[j] && Xlam < Wavl2[j+1])
1177 {
1178 Float_t Dabs = (Qzt[j+1] - Qzt[j])/(Wavl2[j+1] - Wavl2[j]);
1179 Float_t Abso = Qzt[j] + Dabs*(Xlam - Wavl2[j]);
1180 abscoQuarz[i] = -5.0/(TMath::Log(Abso));
1181 }
1182 }
1183 }
1184 printf ("abscoQuarz: %e abscoFreon: %e for energy: %e\n",abscoQuarz[i],abscoFreon[i],ppckov[i]);
1185 }*/
1186
1187 /*Float_t abscoQuarz[31] = {49.64211, 48.41296, 47.46989, 46.50492, 45.13682, 44.47883, 43.1929 , 41.30922, 40.5943 ,
1188 39.82956, 38.98623, 38.6247 , 38.43448, 37.41084, 36.22575, 33.74852, 30.73901, 24.25086,
1189 17.94531, 11.88753, 5.99128, 3.83503, 2.36661, 1.53155, 1.30582, 1.08574, .8779708,
1190 .675275, 0., 0., 0.};
1191
1192 for (Int_t i=0;i<31;i++)
1193 {
1194 abscoQuarz[i] = abscoQuarz[i]/10;
1195 }*/
1196
1197 Float_t abscoQuarz [26] = {105.8, 65.52, 48.58, 42.85, 35.79, 31.262, 28.598, 27.527, 25.007, 22.815, 21.004,
1198 19.266, 17.525, 15.878, 14.177, 11.719, 9.282, 6.62, 4.0925, 2.601, 1.149, .667, .3627,
1199 .192, .1497, .10857};
1200
1201 //Absorption index for methane
1202 Float_t abscoMethane[26];
1203 for (i=0;i<26;i++)
1204 {
1205 abscoMethane[i]=AbsoCH4(ppckov[i]*1e9);
1206 //printf("abscoMethane: %e for energy: %e\n", abscoMethane[i],ppckov[i]*1e9);
1207 }
1208
1209 //Absorption index for opaque quarz, csi and grid, efficiency for all and grid
1210 Float_t abscoOpaqueQuarz[26];
1211 Float_t abscoCsI[26];
1212 Float_t abscoGrid[26];
1213 Float_t efficAll[26];
1214 Float_t efficGrid[26];
1215 for (i=0;i<26;i++)
1216 {
1217 abscoOpaqueQuarz[i]=1e-5;
1218 abscoCsI[i]=1e-4;
1219 abscoGrid[i]=1e-4;
1220 efficAll[i]=1;
1221 efficGrid[i]=1;
1222 //printf ("All must be 1: %e, %e, %e, %e, %e\n",abscoOpaqueQuarz[i],abscoCsI[i],abscoGrid[i],efficAll[i],efficGrid[i]);
1223 }
1224
1225 //Efficiency for csi
1226
1227 Float_t efficCsI[26] = {0.000199999995, 0.000600000028, 0.000699999975, 0.00499999989, 0.00749999983, 0.010125,
1228 0.0242999997, 0.0405000001, 0.0688500032, 0.105299994, 0.121500008, 0.141749993, 0.157949999,
1229 0.162, 0.166050002, 0.167669997, 0.174299985, 0.176789999, 0.179279998, 0.182599992, 0.18592,
1230 0.187579989, 0.189239994, 0.190899998, 0.207499996, 0.215799987};
1231
1232
1233
1234 //FRESNEL LOSS CORRECTION FOR PERPENDICULAR INCIDENCE AND
1235 //UNPOLARIZED PHOTONS
1236
1237 for (i=0;i<26;i++)
1238 {
1239 efficCsI[i] = efficCsI[i]/(1.-Fresnel(ppckov[i]*1e9,1.,0));
1240 //printf ("Fresnel result: %e for energy: %e\n",Fresnel(ppckov[i]*1e9,1.,0),ppckov[i]*1e9);
1241 }
1242
1243 /*******************************************End of rich_media.f***************************************/
1244
1245
1246
1247
1248
1249
1250 Float_t afre[2], agri, amet[2], aqua[2], ahon, zfre[2], zgri, zhon,
1251 zmet[2], zqua[2];
1252 Int_t nlmatfre;
1253 Float_t densquao;
1254 Int_t nlmatmet, nlmatqua;
1255 Float_t wmatquao[2], rIndexFreon[26];
1256 Float_t aquao[2], epsil, stmin, zquao[2];
1257 Int_t nlmatquao;
1258 Float_t radlal, densal, tmaxfd, deemax, stemax;
1259 Float_t aal, zal, radlgri, densfre, radlhon, densgri, denshon,densqua, densmet, wmatfre[2], wmatmet[2], wmatqua[2];
1260
1261 Int_t *idtmed = fIdtmed->GetArray()-999;
1262
452a64c6 1263 // --- Photon energy (GeV)
1264 // --- Refraction indexes
1265 for (i = 0; i < 26; ++i) {
1266 rIndexFreon[i] = ppckov[i] * .0172 * 1e9 + 1.177;
1267 //rIndexFreon[i] = 1;
1268 //printf ("rIndexFreon: %e \n efficCsI: %e for energy: %e\n",rIndexFreon[i], efficCsI[i], ppckov[i]);
1269 }
1270
1271 // --- Detection efficiencies (quantum efficiency for CsI)
1272 // --- Define parameters for honeycomb.
1273 // Used carbon of equivalent rad. lenght
1274
1275 ahon = 12.01;
1276 zhon = 6.;
9dda3582 1277 denshon = 0.1;
452a64c6 1278 radlhon = 18.8;
1279
1280 // --- Parameters to include in GSMIXT, relative to Quarz (SiO2)
1281
1282 aqua[0] = 28.09;
1283 aqua[1] = 16.;
1284 zqua[0] = 14.;
1285 zqua[1] = 8.;
1286 densqua = 2.64;
1287 nlmatqua = -2;
1288 wmatqua[0] = 1.;
1289 wmatqua[1] = 2.;
1290
1291 // --- Parameters to include in GSMIXT, relative to opaque Quarz (SiO2)
1292
1293 aquao[0] = 28.09;
1294 aquao[1] = 16.;
1295 zquao[0] = 14.;
1296 zquao[1] = 8.;
1297 densquao = 2.64;
1298 nlmatquao = -2;
1299 wmatquao[0] = 1.;
1300 wmatquao[1] = 2.;
1301
1302 // --- Parameters to include in GSMIXT, relative to Freon (C6F14)
1303
1304 afre[0] = 12.;
1305 afre[1] = 19.;
1306 zfre[0] = 6.;
1307 zfre[1] = 9.;
1308 densfre = 1.7;
1309 nlmatfre = -2;
1310 wmatfre[0] = 6.;
1311 wmatfre[1] = 14.;
1312
1313 // --- Parameters to include in GSMIXT, relative to methane (CH4)
1314
1315 amet[0] = 12.01;
1316 amet[1] = 1.;
1317 zmet[0] = 6.;
1318 zmet[1] = 1.;
1319 densmet = 7.17e-4;
1320 nlmatmet = -2;
1321 wmatmet[0] = 1.;
1322 wmatmet[1] = 4.;
1323
1324 // --- Parameters to include in GSMIXT, relative to anode grid (Cu)
1325
1326 agri = 63.54;
1327 zgri = 29.;
1328 densgri = 8.96;
1329 radlgri = 1.43;
1330
1331 // --- Parameters to include in GSMATE related to aluminium sheet
1332
1333 aal = 26.98;
1334 zal = 13.;
1335 densal = 2.7;
1336 radlal = 8.9;
6197ac87 1337
1338 // --- Glass parameters
1339
1340 Float_t aglass[5]={12.01, 28.09, 16., 10.8, 23.};
1341 Float_t zglass[5]={ 6., 14., 8., 5., 11.};
1342 Float_t wglass[5]={ 0.5, 0.105, 0.355, 0.03, 0.01};
1343 Float_t dglass=1.74;
1344
452a64c6 1345
1346 AliMaterial(1, "Air $", 14.61, 7.3, .001205, 30420., 67500);
1347 AliMaterial(6, "HON", ahon, zhon, denshon, radlhon, 0);
1348 AliMaterial(16, "CSI", ahon, zhon, denshon, radlhon, 0);
1349 AliMixture(20, "QUA", aqua, zqua, densqua, nlmatqua, wmatqua);
1350 AliMixture(21, "QUAO", aquao, zquao, densquao, nlmatquao, wmatquao);
1351 AliMixture(30, "FRE", afre, zfre, densfre, nlmatfre, wmatfre);
1352 AliMixture(40, "MET", amet, zmet, densmet, nlmatmet, wmatmet);
1353 AliMixture(41, "METG", amet, zmet, densmet, nlmatmet, wmatmet);
1354 AliMaterial(11, "GRI", agri, zgri, densgri, radlgri, 0);
1355 AliMaterial(50, "ALUM", aal, zal, densal, radlal, 0);
6197ac87 1356 AliMixture(32, "GLASS",aglass, zglass, dglass, 5, wglass);
1357 AliMaterial(31, "COPPER$", 63.54, 29., 8.96, 1.4, 0.);
452a64c6 1358
1359 tmaxfd = -10.;
1360 stemax = -.1;
1361 deemax = -.2;
1362 epsil = .001;
1363 stmin = -.001;
1364
1365 AliMedium(1, "DEFAULT MEDIUM AIR$", 1, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
1366 AliMedium(2, "HONEYCOMB$", 6, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
1367 AliMedium(3, "QUARZO$", 20, 1, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
1368 AliMedium(4, "FREON$", 30, 1, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
1369 AliMedium(5, "METANO$", 40, 1, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
1370 AliMedium(6, "CSI$", 16, 1, isxfld, sxmgmx,tmaxfd, stemax, deemax, epsil, stmin);
1371 AliMedium(7, "GRIGLIA$", 11, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
1372 AliMedium(8, "QUARZOO$", 21, 1, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
1373 AliMedium(9, "GAP$", 41, 1, isxfld, sxmgmx,tmaxfd, .1, -deemax, epsil, -stmin);
1374 AliMedium(10, "ALUMINUM$", 50, 1, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
6197ac87 1375 AliMedium(11, "GLASS", 32, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
1376 AliMedium(12, "PCB_COPPER", 31, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
452a64c6 1377
1378
5cfcdf54 1379 gMC->SetCerenkov(idtmed[1000], 26, ppckov, abscoMethane, efficAll, rIndexMethane);
1380 gMC->SetCerenkov(idtmed[1001], 26, ppckov, abscoMethane, efficAll, rIndexMethane);
1381 gMC->SetCerenkov(idtmed[1002], 26, ppckov, abscoQuarz, efficAll,rIndexQuarz);
1382 gMC->SetCerenkov(idtmed[1003], 26, ppckov, abscoFreon, efficAll,rIndexFreon);
1383 gMC->SetCerenkov(idtmed[1004], 26, ppckov, abscoMethane, efficAll, rIndexMethane);
1384 gMC->SetCerenkov(idtmed[1005], 26, ppckov, abscoCsI, efficCsI, rIndexMethane);
1385 gMC->SetCerenkov(idtmed[1006], 26, ppckov, abscoGrid, efficGrid, rIndexGrid);
1386 gMC->SetCerenkov(idtmed[1007], 26, ppckov, abscoOpaqueQuarz, efficAll, rIndexOpaqueQuarz);
1387 gMC->SetCerenkov(idtmed[1008], 26, ppckov, abscoMethane, efficAll, rIndexMethane);
1388 gMC->SetCerenkov(idtmed[1009], 26, ppckov, abscoGrid, efficGrid, rIndexGrid);
1389 gMC->SetCerenkov(idtmed[1010], 26, ppckov, abscoOpaqueQuarz, efficAll, rIndexOpaqueQuarz);
452a64c6 1390}
1391
1392//___________________________________________
1393
1394Float_t AliRICH::Fresnel(Float_t ene,Float_t pdoti, Bool_t pola)
1395{
1396
1397 //ENE(EV), PDOTI=COS(INC.ANG.), PDOTR=COS(POL.PLANE ROT.ANG.)
1398
1399 Float_t en[36] = {5.0,5.1,5.2,5.3,5.4,5.5,5.6,5.7,5.8,5.9,6.0,6.1,6.2,
1400 6.3,6.4,6.5,6.6,6.7,6.8,6.9,7.0,7.1,7.2,7.3,7.4,7.5,7.6,7.7,
1401 7.8,7.9,8.0,8.1,8.2,8.3,8.4,8.5};
1402
1403
1404 Float_t csin[36] = {2.14,2.21,2.33,2.48,2.76,2.97,2.99,2.59,2.81,3.05,
1405 2.86,2.53,2.55,2.66,2.79,2.96,3.18,3.05,2.84,2.81,2.38,2.11,
1406 2.01,2.13,2.39,2.73,3.08,3.15,2.95,2.73,2.56,2.41,2.12,1.95,
1407 1.72,1.53};
1408
1409 Float_t csik[36] = {0.,0.,0.,0.,0.,0.196,0.408,0.208,0.118,0.49,0.784,0.543,
1410 0.424,0.404,0.371,0.514,0.922,1.102,1.139,1.376,1.461,1.253,0.878,
1411 0.69,0.612,0.649,0.824,1.347,1.571,1.678,1.763,1.857,1.824,1.824,
1412 1.714,1.498};
1413 Float_t xe=ene;
1414 Int_t j=Int_t(xe*10)-49;
1415 Float_t cn=csin[j]+((csin[j+1]-csin[j])/0.1)*(xe-en[j]);
1416 Float_t ck=csik[j]+((csik[j+1]-csik[j])/0.1)*(xe-en[j]);
1417
1418 //FORMULAE FROM HANDBOOK OF OPTICS, 33.23 OR
1419 //W.R. HUNTER, J.O.S.A. 54 (1964),15 , J.O.S.A. 55(1965),1197
1420
1421 Float_t sinin=TMath::Sqrt(1-pdoti*pdoti);
1422 Float_t tanin=sinin/pdoti;
1423
1424 Float_t c1=cn*cn-ck*ck-sinin*sinin;
1425 Float_t c2=4*cn*cn*ck*ck;
1426 Float_t aO=TMath::Sqrt(0.5*(TMath::Sqrt(c1*c1+c2)+c1));
1427 Float_t b2=0.5*(TMath::Sqrt(c1*c1+c2)-c1);
1428
1429 Float_t rs=((aO-pdoti)*(aO-pdoti)+b2)/((aO+pdoti)*(aO+pdoti)+b2);
1430 Float_t rp=rs*((aO-sinin*tanin)*(aO-sinin*tanin)+b2)/((aO+sinin*tanin)*(aO+sinin*tanin)+b2);
1431
1432
1433 //CORRECTION FACTOR FOR SURFACE ROUGHNESS
1434 //B.J. STAGG APPLIED OPTICS, 30(1991),4113
1435
1436 Float_t sigraf=18.;
1437 Float_t lamb=1240/ene;
1438 Float_t fresn;
1439
1440 Float_t rO=TMath::Exp(-(4*TMath::Pi()*pdoti*sigraf/lamb)*(4*TMath::Pi()*pdoti*sigraf/lamb));
1441
1442 if(pola)
1443 {
1444 Float_t pdotr=0.8; //DEGREE OF POLARIZATION : 1->P , -1->S
1445 fresn=0.5*(rp*(1+pdotr)+rs*(1-pdotr));
1446 }
1447 else
1448 fresn=0.5*(rp+rs);
1449
1450 fresn = fresn*rO;
1451 return(fresn);
1452}
1453
1454//__________________________________________
1455Float_t AliRICH::AbsoCH4(Float_t x)
1456{
1457
1458 //KLOSCH,SCH4(9),WL(9),EM(9),ALENGTH(31)
1459 Float_t sch4[9] = {.12,.16,.23,.38,.86,2.8,7.9,28.,80.}; //MB X 10^22
1460 //Float_t wl[9] = {153.,152.,151.,150.,149.,148.,147.,146.,145};
1461 Float_t em[9] = {8.1,8.158,8.212,8.267,8.322,8.378,8.435,8.493,8.55};
1462 const Float_t kLosch=2.686763E19; // LOSCHMIDT NUMBER IN CM-3
1463 const Float_t kIgas1=100, kIgas2=0, kOxy=10., kWater=5., kPressure=750.,kTemperature=283.;
1464 Float_t pn=kPressure/760.;
1465 Float_t tn=kTemperature/273.16;
1466
1467
1468// ------- METHANE CROSS SECTION -----------------
1469// ASTROPH. J. 214, L47 (1978)
1470
1471 Float_t sm=0;
1472 if (x<7.75)
1473 sm=.06e-22;
1474
1475 if(x>=7.75 && x<=8.1)
1476 {
1477 Float_t c0=-1.655279e-1;
1478 Float_t c1=6.307392e-2;
1479 Float_t c2=-8.011441e-3;
1480 Float_t c3=3.392126e-4;
1481 sm=(c0+c1*x+c2*x*x+c3*x*x*x)*1.e-18;
1482 }
1483
1484 if (x> 8.1)
1485 {
1486 Int_t j=0;
1487 while (x<=em[j] && x>=em[j+1])
1488 {
1489 j++;
1490 Float_t a=(sch4[j+1]-sch4[j])/(em[j+1]-em[j]);
1491 sm=(sch4[j]+a*(x-em[j]))*1e-22;
1492 }
1493 }
1494
1495 Float_t dm=(kIgas1/100.)*(1.-((kOxy+kWater)/1.e6))*kLosch*pn/tn;
1496 Float_t abslm=1./sm/dm;
1497
1498// ------- ISOBUTHANE CROSS SECTION --------------
1499// i-C4H10 (ai) abs. length from curves in
1500// Lu-McDonald paper for BARI RICH workshop .
1501// -----------------------------------------------------------
1502
1503 Float_t ai;
1504 Float_t absli;
1505 if (kIgas2 != 0)
1506 {
1507 if (x<7.25)
1508 ai=100000000.;
1509
1510 if(x>=7.25 && x<7.375)
1511 ai=24.3;
1512
1513 if(x>=7.375)
1514 ai=.0000000001;
1515
1516 Float_t si = 1./(ai*kLosch*273.16/293.); // ISOB. CRO.SEC.IN CM2
1517 Float_t di=(kIgas2/100.)*(1.-((kOxy+kWater)/1.e6))*kLosch*pn/tn;
1518 absli =1./si/di;
1519 }
1520 else
1521 absli=1.e18;
1522// ---------------------------------------------------------
1523//
1524// transmission of O2
1525//
1526// y= path in cm, x=energy in eV
1527// so= cross section for UV absorption in cm2
1528// do= O2 molecular density in cm-3
1529// ---------------------------------------------------------
1530
1531 Float_t abslo;
1532 Float_t so=0;
1533 if(x>=6.0)
1534 {
1535 if(x>=6.0 && x<6.5)
1536 {
1537 so=3.392709e-13 * TMath::Exp(2.864104 *x);
1538 so=so*1e-18;
1539 }
1540
1541 if(x>=6.5 && x<7.0)
1542 {
1543 so=2.910039e-34 * TMath::Exp(10.3337*x);
1544 so=so*1e-18;
1545 }
1546
1547
1548 if (x>=7.0)
1549 {
1550 Float_t a0=-73770.76;
1551 Float_t a1=46190.69;
1552 Float_t a2=-11475.44;
1553 Float_t a3=1412.611;
1554 Float_t a4=-86.07027;
1555 Float_t a5=2.074234;
1556 so= a0+(a1*x)+(a2*x*x)+(a3*x*x*x)+(a4*x*x*x*x)+(a5*x*x*x*x*x);
1557 so=so*1e-18;
1558 }
1559
1560 Float_t dox=(kOxy/1e6)*kLosch*pn/tn;
1561 abslo=1./so/dox;
1562 }
1563 else
1564 abslo=1.e18;
1565// ---------------------------------------------------------
1566//
1567// transmission of H2O
1568//
1569// y= path in cm, x=energy in eV
1570// sw= cross section for UV absorption in cm2
1571// dw= H2O molecular density in cm-3
1572// ---------------------------------------------------------
1573
1574 Float_t abslw;
1575
1576 Float_t b0=29231.65;
1577 Float_t b1=-15807.74;
1578 Float_t b2=3192.926;
1579 Float_t b3=-285.4809;
1580 Float_t b4=9.533944;
1581
1582 if(x>6.75)
1583 {
1584 Float_t sw= b0+(b1*x)+(b2*x*x)+(b3*x*x*x)+(b4*x*x*x*x);
1585 sw=sw*1e-18;
1586 Float_t dw=(kWater/1e6)*kLosch*pn/tn;
1587 abslw=1./sw/dw;
1588 }
1589 else
1590 abslw=1.e18;
1591
1592// ---------------------------------------------------------
1593
1594 Float_t alength=1./(1./abslm+1./absli+1./abslo+1./abslw);
1595 return (alength);
1596}
1597
1598
1599
ddae0931 1600//___________________________________________
1601Int_t AliRICH::DistancetoPrimitive(Int_t , Int_t )
1602{
237c933d 1603
1604// Default value
1605
ddae0931 1606 return 9999;
1607}
fe4da5cc 1608
ddae0931 1609//___________________________________________
1610void AliRICH::MakeBranch(Option_t* option)
fe4da5cc 1611{
237c933d 1612 // Create Tree branches for the RICH.
fe4da5cc 1613
237c933d 1614 const Int_t kBufferSize = 4000;
ddae0931 1615 char branchname[20];
fe4da5cc 1616
ddae0931 1617
1618 AliDetector::MakeBranch(option);
1619 sprintf(branchname,"%sCerenkov",GetName());
1620 if (fCerenkovs && gAlice->TreeH()) {
237c933d 1621 gAlice->TreeH()->Branch(branchname,&fCerenkovs, kBufferSize);
ddae0931 1622 printf("Making Branch %s for Cerenkov Hits\n",branchname);
fe4da5cc 1623 }
ddae0931 1624
2e5f0f7b 1625 sprintf(branchname,"%sPadHits",GetName());
1626 if (fPadHits && gAlice->TreeH()) {
237c933d 1627 gAlice->TreeH()->Branch(branchname,&fPadHits, kBufferSize);
2e5f0f7b 1628 printf("Making Branch %s for PadHits\n",branchname);
fe4da5cc 1629 }
1630
ddae0931 1631// one branch for digits per chamber
1632 Int_t i;
1633
237c933d 1634 for (i=0; i<kNCH ;i++) {
ddae0931 1635 sprintf(branchname,"%sDigits%d",GetName(),i+1);
1636
1637 if (fDchambers && gAlice->TreeD()) {
237c933d 1638 gAlice->TreeD()->Branch(branchname,&((*fDchambers)[i]), kBufferSize);
ddae0931 1639 printf("Making Branch %s for digits in chamber %d\n",branchname,i+1);
1640 }
fe4da5cc 1641 }
2e5f0f7b 1642
1643// one branch for raw clusters per chamber
237c933d 1644 for (i=0; i<kNCH ;i++) {
2e5f0f7b 1645 sprintf(branchname,"%sRawClusters%d",GetName(),i+1);
1646
1647 if (fRawClusters && gAlice->TreeR()) {
237c933d 1648 gAlice->TreeR()->Branch(branchname,&((*fRawClusters)[i]), kBufferSize);
2e5f0f7b 1649 printf("Making Branch %s for raw clusters in chamber %d\n",branchname,i+1);
1650 }
1651 }
1652
1653 // one branch for rec hits per chamber
237c933d 1654 for (i=0; i<kNCH ;i++) {
a4622d0f 1655 sprintf(branchname,"%sRecHits1D%d",GetName(),i+1);
1656
1657 if (fRecHits1D && gAlice->TreeR()) {
1658 gAlice->TreeR()->Branch(branchname,&((*fRecHits1D)[i]), kBufferSize);
1659 printf("Making Branch %s for 1D rec. hits in chamber %d\n",branchname,i+1);
1660 }
1661 }
1662 for (i=0; i<kNCH ;i++) {
1663 sprintf(branchname,"%sRecHits3D%d",GetName(),i+1);
2e5f0f7b 1664
a4622d0f 1665 if (fRecHits3D && gAlice->TreeR()) {
1666 gAlice->TreeR()->Branch(branchname,&((*fRecHits3D)[i]), kBufferSize);
1667 printf("Making Branch %s for 3D rec. hits in chamber %d\n",branchname,i+1);
2e5f0f7b 1668 }
1669 }
a4622d0f 1670
fe4da5cc 1671}
1672
ddae0931 1673//___________________________________________
1674void AliRICH::SetTreeAddress()
fe4da5cc 1675{
237c933d 1676 // Set branch address for the Hits and Digits Tree.
2e5f0f7b 1677 char branchname[20];
1678 Int_t i;
1679
ddae0931 1680 AliDetector::SetTreeAddress();
1681
1682 TBranch *branch;
1683 TTree *treeH = gAlice->TreeH();
1684 TTree *treeD = gAlice->TreeD();
2e5f0f7b 1685 TTree *treeR = gAlice->TreeR();
ddae0931 1686
1687 if (treeH) {
2e5f0f7b 1688 if (fPadHits) {
1689 branch = treeH->GetBranch("RICHPadHits");
1690 if (branch) branch->SetAddress(&fPadHits);
fe4da5cc 1691 }
ddae0931 1692 if (fCerenkovs) {
1693 branch = treeH->GetBranch("RICHCerenkov");
1694 if (branch) branch->SetAddress(&fCerenkovs);
fe4da5cc 1695 }
ddae0931 1696 }
1697
1698 if (treeD) {
237c933d 1699 for (int i=0; i<kNCH; i++) {
ddae0931 1700 sprintf(branchname,"%sDigits%d",GetName(),i+1);
1701 if (fDchambers) {
1702 branch = treeD->GetBranch(branchname);
1703 if (branch) branch->SetAddress(&((*fDchambers)[i]));
1704 }
fe4da5cc 1705 }
fe4da5cc 1706 }
2e5f0f7b 1707 if (treeR) {
237c933d 1708 for (i=0; i<kNCH; i++) {
2e5f0f7b 1709 sprintf(branchname,"%sRawClusters%d",GetName(),i+1);
1710 if (fRawClusters) {
1711 branch = treeR->GetBranch(branchname);
1712 if (branch) branch->SetAddress(&((*fRawClusters)[i]));
1713 }
1714 }
1715
237c933d 1716 for (i=0; i<kNCH; i++) {
a4622d0f 1717 sprintf(branchname,"%sRecHits1D%d",GetName(),i+1);
1718 if (fRecHits1D) {
2e5f0f7b 1719 branch = treeR->GetBranch(branchname);
a4622d0f 1720 if (branch) branch->SetAddress(&((*fRecHits1D)[i]));
2e5f0f7b 1721 }
1722 }
1723
a4622d0f 1724 for (i=0; i<kNCH; i++) {
1725 sprintf(branchname,"%sRecHits3D%d",GetName(),i+1);
1726 if (fRecHits3D) {
1727 branch = treeR->GetBranch(branchname);
1728 if (branch) branch->SetAddress(&((*fRecHits3D)[i]));
1729 }
1730 }
1731
2e5f0f7b 1732 }
ddae0931 1733}
1734//___________________________________________
1735void AliRICH::ResetHits()
1736{
237c933d 1737 // Reset number of clusters and the cluster array for this detector
ddae0931 1738 AliDetector::ResetHits();
2e5f0f7b 1739 fNPadHits = 0;
1740 fNcerenkovs = 0;
1741 if (fPadHits) fPadHits->Clear();
ddae0931 1742 if (fCerenkovs) fCerenkovs->Clear();
fe4da5cc 1743}
1744
2e5f0f7b 1745
ddae0931 1746//____________________________________________
1747void AliRICH::ResetDigits()
1748{
237c933d 1749 //
1750 // Reset number of digits and the digits array for this detector
1751 //
1752 for ( int i=0;i<kNCH;i++ ) {
ddae0931 1753 if ((*fDchambers)[i]) (*fDchambers)[i]->Clear();
1754 if (fNdch) fNdch[i]=0;
1755 }
1756}
2e5f0f7b 1757
ddae0931 1758//____________________________________________
2e5f0f7b 1759void AliRICH::ResetRawClusters()
fe4da5cc 1760{
237c933d 1761 //
1762 // Reset number of raw clusters and the raw clust array for this detector
1763 //
1764 for ( int i=0;i<kNCH;i++ ) {
2e5f0f7b 1765 if ((*fRawClusters)[i]) ((TClonesArray*)(*fRawClusters)[i])->Clear();
1766 if (fNrawch) fNrawch[i]=0;
fe4da5cc 1767 }
ddae0931 1768}
fe4da5cc 1769
2e5f0f7b 1770//____________________________________________
a4622d0f 1771void AliRICH::ResetRecHits1D()
fe4da5cc 1772{
237c933d 1773 //
1774 // Reset number of raw clusters and the raw clust array for this detector
1775 //
2e5f0f7b 1776
237c933d 1777 for ( int i=0;i<kNCH;i++ ) {
a4622d0f 1778 if ((*fRecHits1D)[i]) ((TClonesArray*)(*fRecHits1D)[i])->Clear();
1779 if (fNrechits1D) fNrechits1D[i]=0;
1780 }
1781}
1782
1783//____________________________________________
1784void AliRICH::ResetRecHits3D()
1785{
1786 //
1787 // Reset number of raw clusters and the raw clust array for this detector
1788 //
1789
1790 for ( int i=0;i<kNCH;i++ ) {
1791 if ((*fRecHits3D)[i]) ((TClonesArray*)(*fRecHits3D)[i])->Clear();
1792 if (fNrechits3D) fNrechits3D[i]=0;
2e5f0f7b 1793 }
fe4da5cc 1794}
1795
ddae0931 1796//___________________________________________
2e5f0f7b 1797void AliRICH::SetGeometryModel(Int_t id, AliRICHGeometry *geometry)
fe4da5cc 1798{
237c933d 1799
1800//
1801// Setter for the RICH geometry model
1802//
1803
1804
2e5f0f7b 1805 ((AliRICHChamber*) (*fChambers)[id])->GeometryModel(geometry);
fe4da5cc 1806}
1807
ddae0931 1808//___________________________________________
a2f7eaf6 1809void AliRICH::SetSegmentationModel(Int_t id, AliSegmentation *segmentation)
ddae0931 1810{
237c933d 1811
1812//
1813// Setter for the RICH segmentation model
1814//
1815
a2f7eaf6 1816 ((AliRICHChamber*) (*fChambers)[id])->SetSegmentationModel(segmentation);
ddae0931 1817}
fe4da5cc 1818
ddae0931 1819//___________________________________________
2e5f0f7b 1820void AliRICH::SetResponseModel(Int_t id, AliRICHResponse *response)
fe4da5cc 1821{
237c933d 1822
1823//
1824// Setter for the RICH response model
1825//
1826
2e5f0f7b 1827 ((AliRICHChamber*) (*fChambers)[id])->ResponseModel(response);
fe4da5cc 1828}
1829
2e5f0f7b 1830void AliRICH::SetReconstructionModel(Int_t id, AliRICHClusterFinder *reconst)
ddae0931 1831{
237c933d 1832
1833//
1834// Setter for the RICH reconstruction model (clusters)
1835//
1836
a2f7eaf6 1837 ((AliRICHChamber*) (*fChambers)[id])->SetReconstructionModel(reconst);
fe4da5cc 1838}
1839
ddae0931 1840void AliRICH::SetNsec(Int_t id, Int_t nsec)
1841{
237c933d 1842
1843//
1844// Sets the number of padplanes
1845//
1846
2e5f0f7b 1847 ((AliRICHChamber*) (*fChambers)[id])->SetNsec(nsec);
ddae0931 1848}
1849
1850
1851//___________________________________________
ddae0931 1852void AliRICH::StepManager()
fe4da5cc 1853{
237c933d 1854
452a64c6 1855// Full Step Manager
1856
1857 Int_t copy, id;
1858 static Int_t idvol;
1859 static Int_t vol[2];
1860 Int_t ipart;
1861 static Float_t hits[18];
1862 static Float_t ckovData[19];
1863 TLorentzVector position;
1864 TLorentzVector momentum;
1865 Float_t pos[3];
1866 Float_t mom[4];
1867 Float_t localPos[3];
1868 Float_t localMom[4];
1869 Float_t localTheta,localPhi;
1870 Float_t theta,phi;
1871 Float_t destep, step;
1872 Float_t ranf[2];
1873 Int_t nPads;
1874 Float_t coscerenkov;
1875 static Float_t eloss, xhit, yhit, tlength;
1876 const Float_t kBig=1.e10;
1877
1878 TClonesArray &lhits = *fHits;
452a64c6 1879 TParticle *current = (TParticle*)(*gAlice->Particles())[gAlice->CurrentTrack()];
237c933d 1880
452a64c6 1881 //if (current->Energy()>1)
1882 //{
1883
1884 // Only gas gap inside chamber
1885 // Tag chambers and record hits when track enters
1886
1887 idvol=-1;
1888 id=gMC->CurrentVolID(copy);
1889 Float_t cherenkovLoss=0;
1890 //gAlice->KeepTrack(gAlice->CurrentTrack());
1891
1892 gMC->TrackPosition(position);
1893 pos[0]=position(0);
1894 pos[1]=position(1);
1895 pos[2]=position(2);
8cda28a5 1896 //bzero((char *)ckovData,sizeof(ckovData)*19);
452a64c6 1897 ckovData[1] = pos[0]; // X-position for hit
1898 ckovData[2] = pos[1]; // Y-position for hit
1899 ckovData[3] = pos[2]; // Z-position for hit
3cc8edee 1900 ckovData[6] = 0; // dummy track length
452a64c6 1901 //ckovData[11] = gAlice->CurrentTrack();
8cda28a5 1902
1903 //printf("\n+++++++++++\nTrack: %d\n++++++++++++\n",gAlice->CurrentTrack());
452a64c6 1904
1905 //AliRICH *RICH = (AliRICH *) gAlice->GetDetector("RICH");
1906
1907 /********************Store production parameters for Cerenkov photons************************/
1908//is it a Cerenkov photon?
8cda28a5 1909 if (gMC->TrackPid() == 50000050) {
452a64c6 1910
1911 //if (gMC->VolId("GAP ")==gMC->CurrentVolID(copy))
1912 //{
1913 Float_t ckovEnergy = current->Energy();
1914 //energy interval for tracking
1915 if (ckovEnergy > 5.6e-09 && ckovEnergy < 7.8e-09 )
1916 //if (ckovEnergy > 0)
1917 {
8cda28a5 1918 if (gMC->IsTrackEntering()){ //is track entering?
1919 //printf("Track entered (1)\n");
452a64c6 1920 if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
1921 { //is it in freo?
5cfcdf54 1922 if (gMC->IsNewTrack()){ //is it the first step?
8cda28a5 1923 //printf("I'm in!\n");
452a64c6 1924 Int_t mother = current->GetFirstMother();
1925
1926 //printf("Second Mother:%d\n",current->GetSecondMother());
1927
1928 ckovData[10] = mother;
1929 ckovData[11] = gAlice->CurrentTrack();
1930 ckovData[12] = 1; //Media where photon was produced 1->Freon, 2->Quarz
8cda28a5 1931 //printf("Produced in FREO\n");
452a64c6 1932 fCkovNumber++;
1933 fFreonProd=1;
1934 //printf("Index: %d\n",fCkovNumber);
1935 } //first step question
1936 } //freo question
1937
5cfcdf54 1938 if (gMC->IsNewTrack()){ //is it first step?
452a64c6 1939 if (gMC->VolId("QUAR")==gMC->CurrentVolID(copy)) //is it in quarz?
1940 {
1941 ckovData[12] = 2;
8cda28a5 1942 //printf("Produced in QUAR\n");
452a64c6 1943 } //quarz question
1944 } //first step question
1945
1946 //printf("Before %d\n",fFreonProd);
1947 } //track entering question
1948
1949 if (ckovData[12] == 1) //was it produced in Freon?
1950 //if (fFreonProd == 1)
1951 {
1952 if (gMC->IsTrackEntering()){ //is track entering?
8cda28a5 1953 //printf("Track entered (2)\n");
1954 //printf("Current volume (should be META): %s\n",gMC->CurrentVolName());
1955 //printf("VolId: %d, CurrentVolID: %d\n",gMC->VolId("META"),gMC->CurrentVolID(copy));
452a64c6 1956 if (gMC->VolId("META")==gMC->CurrentVolID(copy)) //is it in gap?
1957 {
8cda28a5 1958 //printf("Got in META\n");
452a64c6 1959 gMC->TrackMomentum(momentum);
1960 mom[0]=momentum(0);
1961 mom[1]=momentum(1);
1962 mom[2]=momentum(2);
1963 mom[3]=momentum(3);
1964 // Z-position for hit
1965
1966
1967 /**************** Photons lost in second grid have to be calculated by hand************/
1968
1969 Float_t cophi = TMath::Cos(TMath::ATan2(mom[0], mom[1]));
1970 Float_t t = (1. - .025 / cophi) * (1. - .05 / cophi);
1971 gMC->Rndm(ranf, 1);
1972 //printf("grid calculation:%f\n",t);
1973 if (ranf[0] > t) {
5cfcdf54 1974 gMC->StopTrack();
452a64c6 1975 ckovData[13] = 5;
1976 AddCerenkov(gAlice->CurrentTrack(),vol,ckovData);
8cda28a5 1977 //printf("Added One (1)!\n");
452a64c6 1978 //printf("Lost one in grid\n");
1979 }
1980 /**********************************************************************************/
1981 } //gap
1982
8cda28a5 1983 //printf("Current volume (should be CSI) (1): %s\n",gMC->CurrentVolName());
1984 //printf("VolId: %d, CurrentVolID: %d\n",gMC->VolId("CSI "),gMC->CurrentVolID(copy));
452a64c6 1985 if (gMC->VolId("CSI ")==gMC->CurrentVolID(copy)) //is it in csi?
1986 {
8cda28a5 1987 //printf("Got in CSI\n");
452a64c6 1988 gMC->TrackMomentum(momentum);
1989 mom[0]=momentum(0);
1990 mom[1]=momentum(1);
1991 mom[2]=momentum(2);
1992 mom[3]=momentum(3);
1993
1994 /********* Photons lost by Fresnel reflection have to be calculated by hand********/
1995 /***********************Cerenkov phtons (always polarised)*************************/
1996
1997 Float_t cophi = TMath::Cos(TMath::ATan2(mom[0], mom[1]));
1998 Float_t t = Fresnel(ckovEnergy*1e9,cophi,1);
1999 gMC->Rndm(ranf, 1);
2000 if (ranf[0] < t) {
5cfcdf54 2001 gMC->StopTrack();
452a64c6 2002 ckovData[13] = 6;
2003 AddCerenkov(gAlice->CurrentTrack(),vol,ckovData);
8cda28a5 2004 //printf("Added One (2)!\n");
452a64c6 2005 //printf("Lost by Fresnel\n");
2006 }
2007 /**********************************************************************************/
2008 }
2009 } //track entering?
2010
2011
2012 /********************Evaluation of losses************************/
2013 /******************still in the old fashion**********************/
2014
5cfcdf54 2015 TArrayI procs;
2016 Int_t i1 = gMC->StepProcesses(procs); //number of physics mechanisms acting on the particle
452a64c6 2017 for (Int_t i = 0; i < i1; ++i) {
2018 // Reflection loss
5cfcdf54 2019 if (procs[i] == kPLightReflection) { //was it reflected
452a64c6 2020 ckovData[13]=10;
2021 if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
2022 ckovData[13]=1;
2023 if (gMC->CurrentVolID(copy) == gMC->VolId("QUAR"))
2024 ckovData[13]=2;
5cfcdf54 2025 //gMC->StopTrack();
9dda3582 2026 //AddCerenkov(gAlice->CurrentTrack(),vol,ckovData);
452a64c6 2027 } //reflection question
8cda28a5 2028
452a64c6 2029 // Absorption loss
5cfcdf54 2030 else if (procs[i] == kPLightAbsorption) { //was it absorbed?
8cda28a5 2031 //printf("Got in absorption\n");
452a64c6 2032 ckovData[13]=20;
2033 if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
2034 ckovData[13]=11;
2035 if (gMC->CurrentVolID(copy) == gMC->VolId("QUAR"))
2036 ckovData[13]=12;
2037 if (gMC->CurrentVolID(copy) == gMC->VolId("META"))
2038 ckovData[13]=13;
2039 if (gMC->CurrentVolID(copy) == gMC->VolId("GAP "))
2040 ckovData[13]=13;
2041
2042 if (gMC->CurrentVolID(copy) == gMC->VolId("SRIC"))
2043 ckovData[13]=15;
2044
2045 // CsI inefficiency
2046 if (gMC->CurrentVolID(copy) == gMC->VolId("CSI ")) {
2047 ckovData[13]=16;
2048 }
5cfcdf54 2049 gMC->StopTrack();
452a64c6 2050 AddCerenkov(gAlice->CurrentTrack(),vol,ckovData);
8cda28a5 2051 //printf("Added One (3)!\n");
452a64c6 2052 //printf("Added cerenkov %d\n",fCkovNumber);
2053 } //absorption question
2054
2055
2056 // Photon goes out of tracking scope
5cfcdf54 2057 else if (procs[i] == kPStop) { //is it below energy treshold?
452a64c6 2058 ckovData[13]=21;
5cfcdf54 2059 gMC->StopTrack();
452a64c6 2060 AddCerenkov(gAlice->CurrentTrack(),vol,ckovData);
8cda28a5 2061 //printf("Added One (4)!\n");
452a64c6 2062 } // energy treshold question
2063 } //number of mechanisms cycle
2064 /**********************End of evaluation************************/
2065 } //freon production question
2066 } //energy interval question
2067 //}//inside the proximity gap question
2068 } //cerenkov photon question
2069
2070 /**************************************End of Production Parameters Storing*********************/
2071
2072
2073 /*******************************Treat photons that hit the CsI (Ckovs and Feedbacks)************/
2074
2075 if (gMC->TrackPid() == 50000050 || gMC->TrackPid() == 50000051) {
2076 //printf("Cerenkov\n");
2077 if (gMC->VolId("CSI ")==gMC->CurrentVolID(copy))
2078 {
8cda28a5 2079 //printf("Current volume (should be CSI) (2): %s\n",gMC->CurrentVolName());
2080 //printf("VolId: %d, CurrentVolID: %d\n",gMC->VolId("CSI "),gMC->CurrentVolID(copy));
2081 //printf("Got in CSI\n");
452a64c6 2082 if (gMC->Edep() > 0.){
2083 gMC->TrackPosition(position);
2084 gMC->TrackMomentum(momentum);
2085 pos[0]=position(0);
2086 pos[1]=position(1);
2087 pos[2]=position(2);
2088 mom[0]=momentum(0);
2089 mom[1]=momentum(1);
2090 mom[2]=momentum(2);
2091 mom[3]=momentum(3);
2092 Double_t tc = mom[0]*mom[0]+mom[1]*mom[1];
2093 Double_t rt = TMath::Sqrt(tc);
2094 theta = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg;
2095 phi = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg;
2096 gMC->Gmtod(pos,localPos,1);
2097 gMC->Gmtod(mom,localMom,2);
2098
2099 gMC->CurrentVolOffID(2,copy);
2100 vol[0]=copy;
2101 idvol=vol[0]-1;
2102
2103 //Int_t sector=((AliRICHChamber*) (*fChambers)[idvol])
2104 //->Sector(localPos[0], localPos[2]);
2105 //printf("Sector:%d\n",sector);
2106
2107 /*if (gMC->TrackPid() == 50000051){
2108 fFeedbacks++;
2109 printf("Feedbacks:%d\n",fFeedbacks);
2110 }*/
2111
2112 ((AliRICHChamber*) (*fChambers)[idvol])
2113 ->SigGenInit(localPos[0], localPos[2], localPos[1]);
2114 if(idvol<kNCH) {
2115 ckovData[0] = gMC->TrackPid(); // particle type
2116 ckovData[1] = pos[0]; // X-position for hit
2117 ckovData[2] = pos[1]; // Y-position for hit
2118 ckovData[3] = pos[2]; // Z-position for hit
2119 ckovData[4] = theta; // theta angle of incidence
2120 ckovData[5] = phi; // phi angle of incidence
2121 ckovData[8] = (Float_t) fNPadHits; // first padhit
2122 ckovData[9] = -1; // last pad hit
2123 ckovData[13] = 4; // photon was detected
2124 ckovData[14] = mom[0];
2125 ckovData[15] = mom[1];
2126 ckovData[16] = mom[2];
2127
2128 destep = gMC->Edep();
2129 gMC->SetMaxStep(kBig);
2130 cherenkovLoss += destep;
2131 ckovData[7]=cherenkovLoss;
2132
2133 nPads = MakePadHits(localPos[0],localPos[2],cherenkovLoss,idvol,kCerenkov);
2134 if (fNPadHits > (Int_t)ckovData[8]) {
2135 ckovData[8]= ckovData[8]+1;
2136 ckovData[9]= (Float_t) fNPadHits;
2137 }
2138
2139 ckovData[17] = nPads;
2140 //printf("nPads:%d",nPads);
2141
2142 //TClonesArray *Hits = RICH->Hits();
2143 AliRICHHit *mipHit = (AliRICHHit*) (fHits->UncheckedAt(0));
2144 if (mipHit)
2145 {
2146 mom[0] = current->Px();
2147 mom[1] = current->Py();
2148 mom[2] = current->Pz();
2149 Float_t mipPx = mipHit->fMomX;
2150 Float_t mipPy = mipHit->fMomY;
2151 Float_t mipPz = mipHit->fMomZ;
2152
2153 Float_t r = mom[0]*mom[0] + mom[1]*mom[1] + mom[2]*mom[2];
2154 Float_t rt = TMath::Sqrt(r);
2155 Float_t mipR = mipPx*mipPx + mipPy*mipPy + mipPz*mipPz;
2156 Float_t mipRt = TMath::Sqrt(mipR);
2157 if ((rt*mipRt) > 0)
2158 {
2159 coscerenkov = (mom[0]*mipPx + mom[1]*mipPy + mom[2]*mipPz)/(rt*mipRt);
2160 }
2161 else
2162 {
2163 coscerenkov = 0;
2164 }
2165 Float_t cherenkov = TMath::ACos(coscerenkov);
2166 ckovData[18]=cherenkov;
2167 }
2168 //if (sector != -1)
2169 //{
2170 AddHit(gAlice->CurrentTrack(),vol,ckovData);
2171 AddCerenkov(gAlice->CurrentTrack(),vol,ckovData);
8cda28a5 2172 //printf("Added One (5)!\n");
452a64c6 2173 //}
2174 }
2175 }
2176 }
2177 }
2178
2179 /***********************************************End of photon hits*********************************************/
2180
2181
2182 /**********************************************Charged particles treatment*************************************/
2183
2184 else if (gMC->TrackCharge())
2185 //else if (1 == 1)
2186 {
2187//If MIP
2188 /*if (gMC->IsTrackEntering())
2189 {
2190 hits[13]=20;//is track entering?
2191 }*/
2192 if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
2193 {
2194 fFreonProd=1;
2195 }
2196
2197 if (gMC->VolId("GAP ")== gMC->CurrentVolID(copy)) {
2198// Get current particle id (ipart), track position (pos) and momentum (mom)
2199
2200 gMC->CurrentVolOffID(3,copy);
2201 vol[0]=copy;
2202 idvol=vol[0]-1;
2203
2204 //Int_t sector=((AliRICHChamber*) (*fChambers)[idvol])
2205 //->Sector(localPos[0], localPos[2]);
2206 //printf("Sector:%d\n",sector);
2207
2208 gMC->TrackPosition(position);
2209 gMC->TrackMomentum(momentum);
2210 pos[0]=position(0);
2211 pos[1]=position(1);
2212 pos[2]=position(2);
2213 mom[0]=momentum(0);
2214 mom[1]=momentum(1);
2215 mom[2]=momentum(2);
2216 mom[3]=momentum(3);
2217 gMC->Gmtod(pos,localPos,1);
2218 gMC->Gmtod(mom,localMom,2);
2219
2220 ipart = gMC->TrackPid();
2221 //
2222 // momentum loss and steplength in last step
2223 destep = gMC->Edep();
2224 step = gMC->TrackStep();
2225
2226 //
2227 // record hits when track enters ...
2228 if( gMC->IsTrackEntering()) {
2229// gMC->SetMaxStep(fMaxStepGas);
2230 Double_t tc = mom[0]*mom[0]+mom[1]*mom[1];
2231 Double_t rt = TMath::Sqrt(tc);
2232 theta = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg;
2233 phi = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg;
2234
2235
2236 Double_t localTc = localMom[0]*localMom[0]+localMom[2]*localMom[2];
2237 Double_t localRt = TMath::Sqrt(localTc);
2238 localTheta = Float_t(TMath::ATan2(localRt,Double_t(localMom[1])))*kRaddeg;
2239 localPhi = Float_t(TMath::ATan2(Double_t(localMom[2]),Double_t(localMom[0])))*kRaddeg;
2240
2241 hits[0] = Float_t(ipart); // particle type
2242 hits[1] = localPos[0]; // X-position for hit
2243 hits[2] = localPos[1]; // Y-position for hit
2244 hits[3] = localPos[2]; // Z-position for hit
2245 hits[4] = localTheta; // theta angle of incidence
2246 hits[5] = localPhi; // phi angle of incidence
2247 hits[8] = (Float_t) fNPadHits; // first padhit
2248 hits[9] = -1; // last pad hit
2249 hits[13] = fFreonProd; // did id hit the freon?
2250 hits[14] = mom[0];
2251 hits[15] = mom[1];
2252 hits[16] = mom[2];
2253
2254 tlength = 0;
2255 eloss = 0;
2256 fFreonProd = 0;
2257
2258 Chamber(idvol).LocaltoGlobal(localPos,hits+1);
2259
2260
2261 //To make chamber coordinates x-y had to pass localPos[0], localPos[2]
2262 xhit = localPos[0];
2263 yhit = localPos[2];
2264 // Only if not trigger chamber
2265 if(idvol<kNCH) {
2266 //
2267 // Initialize hit position (cursor) in the segmentation model
2268 ((AliRICHChamber*) (*fChambers)[idvol])
2269 ->SigGenInit(localPos[0], localPos[2], localPos[1]);
2270 }
2271 }
2272
2273 //
2274 // Calculate the charge induced on a pad (disintegration) in case
2275 //
2276 // Mip left chamber ...
2277 if( gMC->IsTrackExiting() || gMC->IsTrackStop() || gMC->IsTrackDisappeared()){
2278 gMC->SetMaxStep(kBig);
2279 eloss += destep;
2280 tlength += step;
2281
2282
2283 // Only if not trigger chamber
2284 if(idvol<kNCH) {
2285 if (eloss > 0)
2286 {
2287 if(gMC->TrackPid() == kNeutron)
2288 printf("\n\n\n\n\n Neutron Making Pad Hit!!! \n\n\n\n");
2289 nPads = MakePadHits(xhit,yhit,eloss,idvol,kMip);
2290 hits[17] = nPads;
2291 //printf("nPads:%d",nPads);
2292 }
2293 }
2294
2295 hits[6]=tlength;
2296 hits[7]=eloss;
2297 if (fNPadHits > (Int_t)hits[8]) {
2298 hits[8]= hits[8]+1;
2299 hits[9]= (Float_t) fNPadHits;
2300 }
2301
2302 //if(sector !=-1)
2303 new(lhits[fNhits++]) AliRICHHit(fIshunt,gAlice->CurrentTrack(),vol,hits);
2304 eloss = 0;
2305 //
2306 // Check additional signal generation conditions
2307 // defined by the segmentation
2308 // model (boundary crossing conditions)
2309 } else if
2310 (((AliRICHChamber*) (*fChambers)[idvol])
2311 ->SigGenCond(localPos[0], localPos[2], localPos[1]))
2312 {
2313 ((AliRICHChamber*) (*fChambers)[idvol])
2314 ->SigGenInit(localPos[0], localPos[2], localPos[1]);
2315 if (eloss > 0)
2316 {
2317 if(gMC->TrackPid() == kNeutron)
2318 printf("\n\n\n\n\n Neutron Making Pad Hit!!! \n\n\n\n");
2319 nPads = MakePadHits(xhit,yhit,eloss,idvol,kMip);
2320 hits[17] = nPads;
2321 //printf("Npads:%d",NPads);
2322 }
2323 xhit = localPos[0];
2324 yhit = localPos[2];
2325 eloss = destep;
2326 tlength += step ;
2327 //
2328 // nothing special happened, add up energy loss
2329 } else {
2330 eloss += destep;
2331 tlength += step ;
2332 }
2333 }
2334 }
2335 /*************************************************End of MIP treatment**************************************/
2336 //}
ddae0931 2337}
2e5f0f7b 2338
237c933d 2339void AliRICH::FindClusters(Int_t nev,Int_t lastEntry)
ddae0931 2340{
2e5f0f7b 2341
ddae0931 2342//
2343// Loop on chambers and on cathode planes
2344//
2e5f0f7b 2345 for (Int_t icat=1;icat<2;icat++) {
2346 gAlice->ResetDigits();
2347 gAlice->TreeD()->GetEvent(1); // spurious +1 ...
237c933d 2348 for (Int_t ich=0;ich<kNCH;ich++) {
2e5f0f7b 2349 AliRICHChamber* iChamber=(AliRICHChamber*) (*fChambers)[ich];
237c933d 2350 TClonesArray *pRICHdigits = this->DigitsAddress(ich);
2351 if (pRICHdigits == 0)
2e5f0f7b 2352 continue;
2353 //
2354 // Get ready the current chamber stuff
2355 //
2356 AliRICHResponse* response = iChamber->GetResponseModel();
a2f7eaf6 2357 AliSegmentation* seg = iChamber->GetSegmentationModel();
2e5f0f7b 2358 AliRICHClusterFinder* rec = iChamber->GetReconstructionModel();
2359 if (seg) {
2360 rec->SetSegmentation(seg);
2361 rec->SetResponse(response);
237c933d 2362 rec->SetDigits(pRICHdigits);
2e5f0f7b 2363 rec->SetChamber(ich);
2364 if (nev==0) rec->CalibrateCOG();
2365 rec->FindRawClusters();
2366 }
2367 TClonesArray *fRch;
2368 fRch=RawClustAddress(ich);
2369 fRch->Sort();
2370 } // for ich
2371
2372 gAlice->TreeR()->Fill();
2373 TClonesArray *fRch;
237c933d 2374 for (int i=0;i<kNCH;i++) {
2e5f0f7b 2375 fRch=RawClustAddress(i);
2376 int nraw=fRch->GetEntriesFast();
2377 printf ("Chamber %d, raw clusters %d\n",i,nraw);
2378 }
2379
2380 ResetRawClusters();
2381
2382 } // for icat
2383
2384 char hname[30];
2385 sprintf(hname,"TreeR%d",nev);
33984590 2386 gAlice->TreeR()->Write(hname,kOverwrite,0);
2e5f0f7b 2387 gAlice->TreeR()->Reset();
2388
2389 //gObjectTable->Print();
ddae0931 2390}
2391
2392
2393//______________________________________________________________________________
2394void AliRICH::Streamer(TBuffer &R__b)
2395{
2396 // Stream an object of class AliRICH.
2e5f0f7b 2397 AliRICHChamber *iChamber;
a2f7eaf6 2398 AliSegmentation *segmentation;
2e5f0f7b 2399 AliRICHResponse *response;
ddae0931 2400 TClonesArray *digitsaddress;
2e5f0f7b 2401 TClonesArray *rawcladdress;
a4622d0f 2402 TClonesArray *rechitaddress1D;
2403 TClonesArray *rechitaddress3D;
ddae0931 2404
2405 if (R__b.IsReading()) {
2406 Version_t R__v = R__b.ReadVersion(); if (R__v) { }
2407 AliDetector::Streamer(R__b);
2e5f0f7b 2408 R__b >> fNPadHits;
2409 R__b >> fPadHits; // diff
2410 R__b >> fNcerenkovs;
2411 R__b >> fCerenkovs; // diff
ddae0931 2412 R__b >> fDchambers;
2e5f0f7b 2413 R__b >> fRawClusters;
a4622d0f 2414 R__b >> fRecHits1D; //diff
2415 R__b >> fRecHits3D; //diff
33984590 2416 R__b >> fDebugLevel; //diff
9dda3582 2417 R__b.ReadStaticArray(fNdch);
2418 R__b.ReadStaticArray(fNrawch);
a4622d0f 2419 R__b.ReadStaticArray(fNrechits1D);
2420 R__b.ReadStaticArray(fNrechits3D);
2e5f0f7b 2421//
ddae0931 2422 R__b >> fChambers;
2423// Stream chamber related information
237c933d 2424 for (Int_t i =0; i<kNCH; i++) {
2e5f0f7b 2425 iChamber=(AliRICHChamber*) (*fChambers)[i];
ddae0931 2426 iChamber->Streamer(R__b);
2e5f0f7b 2427 segmentation=iChamber->GetSegmentationModel();
2428 segmentation->Streamer(R__b);
2429 response=iChamber->GetResponseModel();
ddae0931 2430 response->Streamer(R__b);
2e5f0f7b 2431 rawcladdress=(TClonesArray*) (*fRawClusters)[i];
2432 rawcladdress->Streamer(R__b);
a4622d0f 2433 rechitaddress1D=(TClonesArray*) (*fRecHits1D)[i];
2434 rechitaddress1D->Streamer(R__b);
2435 rechitaddress3D=(TClonesArray*) (*fRecHits3D)[i];
2436 rechitaddress3D->Streamer(R__b);
ddae0931 2437 digitsaddress=(TClonesArray*) (*fDchambers)[i];
2438 digitsaddress->Streamer(R__b);
2439 }
452a64c6 2440 R__b >> fDebugLevel;
2441 R__b >> fCkovNumber;
2442 R__b >> fCkovQuarz;
2443 R__b >> fCkovGap;
2444 R__b >> fCkovCsi;
2445 R__b >> fLostRfreo;
2446 R__b >> fLostRquar;
2447 R__b >> fLostAfreo;
2448 R__b >> fLostAquarz;
2449 R__b >> fLostAmeta;
2450 R__b >> fLostCsi;
2451 R__b >> fLostWires;
2452 R__b >> fFreonProd;
2453 R__b >> fMipx;
2454 R__b >> fMipy;
2455 R__b >> fFeedbacks;
2456 R__b >> fLostFresnel;
ddae0931 2457
2458 } else {
2459 R__b.WriteVersion(AliRICH::IsA());
2460 AliDetector::Streamer(R__b);
2e5f0f7b 2461 R__b << fNPadHits;
2462 R__b << fPadHits; // diff
2463 R__b << fNcerenkovs;
2464 R__b << fCerenkovs; // diff
ddae0931 2465 R__b << fDchambers;
2e5f0f7b 2466 R__b << fRawClusters;
a4622d0f 2467 R__b << fRecHits1D; //diff
2468 R__b << fRecHits3D; //diff
33984590 2469 R__b << fDebugLevel; //diff
237c933d 2470 R__b.WriteArray(fNdch, kNCH);
2471 R__b.WriteArray(fNrawch, kNCH);
a4622d0f 2472 R__b.WriteArray(fNrechits1D, kNCH);
2473 R__b.WriteArray(fNrechits3D, kNCH);
2474//
ddae0931 2475 R__b << fChambers;
2476// Stream chamber related information
237c933d 2477 for (Int_t i =0; i<kNCH; i++) {
2e5f0f7b 2478 iChamber=(AliRICHChamber*) (*fChambers)[i];
ddae0931 2479 iChamber->Streamer(R__b);
2e5f0f7b 2480 segmentation=iChamber->GetSegmentationModel();
2481 segmentation->Streamer(R__b);
2482 response=iChamber->GetResponseModel();
ddae0931 2483 response->Streamer(R__b);
2e5f0f7b 2484 rawcladdress=(TClonesArray*) (*fRawClusters)[i];
2485 rawcladdress->Streamer(R__b);
a4622d0f 2486 rechitaddress1D=(TClonesArray*) (*fRecHits1D)[i];
2487 rechitaddress1D->Streamer(R__b);
2488 rechitaddress3D=(TClonesArray*) (*fRecHits3D)[i];
2489 rechitaddress3D->Streamer(R__b);
ddae0931 2490 digitsaddress=(TClonesArray*) (*fDchambers)[i];
2491 digitsaddress->Streamer(R__b);
2492 }
452a64c6 2493 R__b << fDebugLevel;
2494 R__b << fCkovNumber;
2495 R__b << fCkovQuarz;
2496 R__b << fCkovGap;
2497 R__b << fCkovCsi;
2498 R__b << fLostRfreo;
2499 R__b << fLostRquar;
2500 R__b << fLostAfreo;
2501 R__b << fLostAquarz;
2502 R__b << fLostAmeta;
2503 R__b << fLostCsi;
2504 R__b << fLostWires;
2505 R__b << fFreonProd;
2506 R__b << fMipx;
2507 R__b << fMipy;
2508 R__b << fFeedbacks;
2509 R__b << fLostFresnel;
fe4da5cc 2510 }
ddae0931 2511}
2e5f0f7b 2512AliRICHPadHit* AliRICH::FirstPad(AliRICHHit* hit,TClonesArray *clusters )
ddae0931 2513{
2514//
2515 // Initialise the pad iterator
2516 // Return the address of the first padhit for hit
2517 TClonesArray *theClusters = clusters;
2518 Int_t nclust = theClusters->GetEntriesFast();
2519 if (nclust && hit->fPHlast > 0) {
2e5f0f7b 2520 sMaxIterPad=Int_t(hit->fPHlast);
2521 sCurIterPad=Int_t(hit->fPHfirst);
2522 return (AliRICHPadHit*) clusters->UncheckedAt(sCurIterPad-1);
ddae0931 2523 } else {
2524 return 0;
fe4da5cc 2525 }
ddae0931 2526
fe4da5cc 2527}
2528
2e5f0f7b 2529AliRICHPadHit* AliRICH::NextPad(TClonesArray *clusters)
fe4da5cc 2530{
237c933d 2531
2532 // Iterates over pads
2533
ddae0931 2534 sCurIterPad++;
2535 if (sCurIterPad <= sMaxIterPad) {
2e5f0f7b 2536 return (AliRICHPadHit*) clusters->UncheckedAt(sCurIterPad-1);
ddae0931 2537 } else {
2538 return 0;
2539 }
fe4da5cc 2540}
2541
fe4da5cc 2542
d28b5fc2 2543void AliRICH::Digitise(Int_t nev, Int_t flag, Option_t *option,Text_t *filename)
fe4da5cc 2544{
ddae0931 2545 // keep galice.root for signal and name differently the file for
2546 // background when add! otherwise the track info for signal will be lost !
2e5f0f7b 2547
c90dd3e2 2548 static Bool_t first=kTRUE;
237c933d 2549 static TFile *pFile;
2550 char *addBackground = strstr(option,"Add");
53a60579 2551 Int_t particle;
ddae0931 2552
2e5f0f7b 2553 FILE* points; //these will be the digits...
2554
2555 points=fopen("points.dat","w");
2556
2557 AliRICHChamber* iChamber;
a2f7eaf6 2558 AliSegmentation* segmentation;
ddae0931 2559
53a60579 2560 Int_t digitise=0;
ddae0931 2561 Int_t trk[50];
2562 Int_t chtrk[50];
2563 TObjArray *list=new TObjArray;
237c933d 2564 static TClonesArray *pAddress=0;
2565 if(!pAddress) pAddress=new TClonesArray("TVector",1000);
2e5f0f7b 2566 Int_t digits[5];
ddae0931 2567
237c933d 2568 AliRICH *pRICH = (AliRICH *) gAlice->GetDetector("RICH");
a2f7eaf6 2569 AliHitMap* pHitMap[10];
2e5f0f7b 2570 Int_t i;
237c933d 2571 for (i=0; i<10; i++) {pHitMap[i]=0;}
2572 if (addBackground ) {
ddae0931 2573 if(first) {
2574 fFileName=filename;
2575 cout<<"filename"<<fFileName<<endl;
237c933d 2576 pFile=new TFile(fFileName);
ddae0931 2577 cout<<"I have opened "<<fFileName<<" file "<<endl;
2e5f0f7b 2578 fHits2 = new TClonesArray("AliRICHHit",1000 );
2579 fClusters2 = new TClonesArray("AliRICHPadHit",10000);
c90dd3e2 2580 first=kFALSE;
ddae0931 2581 }
237c933d 2582 pFile->cd();
ddae0931 2583 // Get Hits Tree header from file
2584 if(fHits2) fHits2->Clear();
2585 if(fClusters2) fClusters2->Clear();
2e5f0f7b 2586 if(TrH1) delete TrH1;
2587 TrH1=0;
2588
ddae0931 2589 char treeName[20];
2590 sprintf(treeName,"TreeH%d",nev);
2e5f0f7b 2591 TrH1 = (TTree*)gDirectory->Get(treeName);
2592 if (!TrH1) {
ddae0931 2593 printf("ERROR: cannot find Hits Tree for event:%d\n",nev);
2594 }
2595 // Set branch addresses
2596 TBranch *branch;
2597 char branchname[20];
2598 sprintf(branchname,"%s",GetName());
2e5f0f7b 2599 if (TrH1 && fHits2) {
2600 branch = TrH1->GetBranch(branchname);
ddae0931 2601 if (branch) branch->SetAddress(&fHits2);
2602 }
2e5f0f7b 2603 if (TrH1 && fClusters2) {
2604 branch = TrH1->GetBranch("RICHCluster");
ddae0931 2605 if (branch) branch->SetAddress(&fClusters2);
2606 }
2607 }
33984590 2608
a2f7eaf6 2609 AliHitMap* hm;
2e5f0f7b 2610 Int_t countadr=0;
33984590 2611 Int_t counter=0;
2612 for (i =0; i<kNCH; i++) {
2613 iChamber=(AliRICHChamber*) (*fChambers)[i];
2614 segmentation=iChamber->GetSegmentationModel(1);
2615 pHitMap[i] = new AliRICHHitMapA1(segmentation, list);
2616 }
2617 //
2618 // Loop over tracks
2619 //
2620
2621 TTree *treeH = gAlice->TreeH();
2622 Int_t ntracks =(Int_t) treeH->GetEntries();
2623 for (Int_t track=0; track<ntracks; track++) {
2624 gAlice->ResetHits();
2625 treeH->GetEvent(track);
2626 //
2627 // Loop over hits
2628 for(AliRICHHit* mHit=(AliRICHHit*)pRICH->FirstHit(-1);
2629 mHit;
2630 mHit=(AliRICHHit*)pRICH->NextHit())
2631 {
2632
33984590 2633 Int_t nch = mHit->fChamber-1; // chamber number
53a60579 2634 Int_t index = mHit->Track();
33984590 2635 if (nch >kNCH) continue;
2636 iChamber = &(pRICH->Chamber(nch));
2637
53a60579 2638 TParticle *current = (TParticle*)(*gAlice->Particles())[index];
33984590 2639
53a60579 2640 if (current->GetPdgCode() >= 50000050)
2641 {
2642 TParticle *motherofcurrent = (TParticle*)(*gAlice->Particles())[current->GetFirstMother()];
2643 particle = motherofcurrent->GetPdgCode();
2644 }
2645 else
2646 {
2647 particle = current->GetPdgCode();
2648 }
2649
33984590 2650
2651 //printf("Flag:%d\n",flag);
2652 //printf("Track:%d\n",track);
2653 //printf("Particle:%d\n",particle);
2654
53a60579 2655 digitise=1;
33984590 2656
2657 if (flag == 1)
53a60579 2658 if(TMath::Abs(particle)==211 || TMath::Abs(particle)==111)
2659 digitise=0;
33984590 2660
2661 if (flag == 2)
2662 if(TMath::Abs(particle)==321 || TMath::Abs(particle)==130 || TMath::Abs(particle)==310
2663 || TMath::Abs(particle)==311)
53a60579 2664 digitise=0;
33984590 2665
2666 if (flag == 3 && TMath::Abs(particle)==2212)
53a60579 2667 digitise=0;
33984590 2668
2669 if (flag == 4 && TMath::Abs(particle)==13)
53a60579 2670 digitise=0;
33984590 2671
2672 if (flag == 5 && TMath::Abs(particle)==11)
53a60579 2673 digitise=0;
33984590 2674
2675 if (flag == 6 && TMath::Abs(particle)==2112)
53a60579 2676 digitise=0;
33984590 2677
2678
53a60579 2679 //printf ("Particle: %d, Flag: %d, Digitise: %d\n",particle,flag,digitise);
33984590 2680
2681
53a60579 2682 if (digitise)
ddae0931 2683 {
d28b5fc2 2684
33984590 2685 //
2686 // Loop over pad hits
2687 for (AliRICHPadHit* mPad=
2688 (AliRICHPadHit*)pRICH->FirstPad(mHit,fPadHits);
2689 mPad;
2690 mPad=(AliRICHPadHit*)pRICH->NextPad(fPadHits))
d28b5fc2 2691 {
33984590 2692 Int_t ipx = mPad->fPadX; // pad number on X
2693 Int_t ipy = mPad->fPadY; // pad number on Y
2694 Int_t iqpad = mPad->fQpad; // charge per pad
2695 //
2696 //
2697 //printf("X:%d, Y:%d, Q:%d\n",ipx,ipy,iqpad);
2698
a2f7eaf6 2699 Float_t thex, they, thez;
176a9917 2700 segmentation=iChamber->GetSegmentationModel(0);
a2f7eaf6 2701 segmentation->GetPadC(ipx,ipy,thex,they,thez);
33984590 2702 new((*pAddress)[countadr++]) TVector(2);
2703 TVector &trinfo=*((TVector*) (*pAddress)[countadr-1]);
2704 trinfo(0)=(Float_t)track;
2705 trinfo(1)=(Float_t)iqpad;
2706
2707 digits[0]=ipx;
2708 digits[1]=ipy;
2709 digits[2]=iqpad;
2710
2711 AliRICHTransientDigit* pdigit;
2712 // build the list of fired pads and update the info
2713 if (!pHitMap[nch]->TestHit(ipx, ipy)) {
2714 list->AddAtAndExpand(new AliRICHTransientDigit(nch,digits),counter);
2715 pHitMap[nch]->SetHit(ipx, ipy, counter);
2716 counter++;
2717 pdigit=(AliRICHTransientDigit*)list->At(list->GetLast());
2718 // list of tracks
2719 TObjArray *trlist=(TObjArray*)pdigit->TrackList();
2720 trlist->Add(&trinfo);
2721 } else {
2722 pdigit=(AliRICHTransientDigit*) pHitMap[nch]->GetHit(ipx, ipy);
2723 // update charge
2724 (*pdigit).fSignal+=iqpad;
2725 // update list of tracks
2726 TObjArray* trlist=(TObjArray*)pdigit->TrackList();
2727 Int_t lastEntry=trlist->GetLast();
2728 TVector *ptrkP=(TVector*)trlist->At(lastEntry);
2729 TVector &ptrk=*ptrkP;
2730 Int_t lastTrack=Int_t(ptrk(0));
2731 Int_t lastCharge=Int_t(ptrk(1));
2732 if (lastTrack==track) {
2733 lastCharge+=iqpad;
2734 trlist->RemoveAt(lastEntry);
2735 trinfo(0)=lastTrack;
2736 trinfo(1)=lastCharge;
2737 trlist->AddAt(&trinfo,lastEntry);
2738 } else {
2739 trlist->Add(&trinfo);
2740 }
2741 // check the track list
2742 Int_t nptracks=trlist->GetEntriesFast();
2743 if (nptracks > 2) {
2744 printf("Attention - tracks: %d (>2)\n",nptracks);
2745 //printf("cat,nch,ix,iy %d %d %d %d \n",icat+1,nch,ipx,ipy);
2746 for (Int_t tr=0;tr<nptracks;tr++) {
2747 TVector *pptrkP=(TVector*)trlist->At(tr);
2748 TVector &pptrk=*pptrkP;
2749 trk[tr]=Int_t(pptrk(0));
2750 chtrk[tr]=Int_t(pptrk(1));
2751 }
2752 } // end if nptracks
2753 } // end if pdigit
2754 } //end loop over clusters
2755 }// track type condition
2756 } // hit loop
2757 } // track loop
2758
2759 // open the file with background
2760
2761 if (addBackground ) {
2762 ntracks =(Int_t)TrH1->GetEntries();
2763 //printf("background - icat,ntracks1 %d %d\n",icat,ntracks);
2764 //printf("background - Start loop over tracks \n");
2765 //
2766 // Loop over tracks
2767 //
2768 for (Int_t trak=0; trak<ntracks; trak++) {
2769 if (fHits2) fHits2->Clear();
2770 if (fClusters2) fClusters2->Clear();
2771 TrH1->GetEvent(trak);
2772 //
2773 // Loop over hits
2774 AliRICHHit* mHit;
2775 for(int j=0;j<fHits2->GetEntriesFast();++j)
2776 {
2777 mHit=(AliRICHHit*) (*fHits2)[j];
2778 Int_t nch = mHit->fChamber-1; // chamber number
2779 if (nch >6) continue;
2780 iChamber = &(pRICH->Chamber(nch));
2781 Int_t rmin = (Int_t)iChamber->RInner();
2782 Int_t rmax = (Int_t)iChamber->ROuter();
2783 //
2784 // Loop over pad hits
2785 for (AliRICHPadHit* mPad=
2786 (AliRICHPadHit*)pRICH->FirstPad(mHit,fClusters2);
2787 mPad;
2788 mPad=(AliRICHPadHit*)pRICH->NextPad(fClusters2))
2789 {
33984590 2790 Int_t ipx = mPad->fPadX; // pad number on X
2791 Int_t ipy = mPad->fPadY; // pad number on Y
2792 Int_t iqpad = mPad->fQpad; // charge per pad
2e5f0f7b 2793
a2f7eaf6 2794 Float_t thex, they, thez;
176a9917 2795 segmentation=iChamber->GetSegmentationModel(0);
a2f7eaf6 2796 segmentation->GetPadC(ipx,ipy,thex,they,thez);
33984590 2797 Float_t rpad=TMath::Sqrt(thex*thex+they*they);
2798 if (rpad < rmin || iqpad ==0 || rpad > rmax) continue;
2799 new((*pAddress)[countadr++]) TVector(2);
2800 TVector &trinfo=*((TVector*) (*pAddress)[countadr-1]);
2801 trinfo(0)=-1; // tag background
2802 trinfo(1)=-1;
2803 digits[0]=ipx;
2804 digits[1]=ipy;
2805 digits[2]=iqpad;
2806 if (trak <4 && nch==0)
2807 printf("bgr - pHitMap[nch]->TestHit(ipx, ipy),trak %d %d\n",
2808 pHitMap[nch]->TestHit(ipx, ipy),trak);
2809 AliRICHTransientDigit* pdigit;
2810 // build the list of fired pads and update the info
2811 if (!pHitMap[nch]->TestHit(ipx, ipy)) {
2812 list->AddAtAndExpand(new AliRICHTransientDigit(nch,digits),counter);
2813
2814 pHitMap[nch]->SetHit(ipx, ipy, counter);
2815 counter++;
2816 printf("bgr new elem in list - counter %d\n",counter);
2817
2818 pdigit=(AliRICHTransientDigit*)list->At(list->GetLast());
2819 // list of tracks
2820 TObjArray *trlist=(TObjArray*)pdigit->TrackList();
2821 trlist->Add(&trinfo);
2822 } else {
2823 pdigit=(AliRICHTransientDigit*) pHitMap[nch]->GetHit(ipx, ipy);
2824 // update charge
2825 (*pdigit).fSignal+=iqpad;
2826 // update list of tracks
2827 TObjArray* trlist=(TObjArray*)pdigit->TrackList();
2828 Int_t lastEntry=trlist->GetLast();
2829 TVector *ptrkP=(TVector*)trlist->At(lastEntry);
2830 TVector &ptrk=*ptrkP;
2831 Int_t lastTrack=Int_t(ptrk(0));
2832 if (lastTrack==-1) {
2833 continue;
2834 } else {
2835 trlist->Add(&trinfo);
2836 }
2837 // check the track list
2838 Int_t nptracks=trlist->GetEntriesFast();
2839 if (nptracks > 0) {
2840 for (Int_t tr=0;tr<nptracks;tr++) {
2841 TVector *pptrkP=(TVector*)trlist->At(tr);
2842 TVector &pptrk=*pptrkP;
2843 trk[tr]=Int_t(pptrk(0));
2844 chtrk[tr]=Int_t(pptrk(1));
2845 }
2846 } // end if nptracks
2847 } // end if pdigit
2848 } //end loop over clusters
2849 } // hit loop
2850 } // track loop
ddae0931 2851 TTree *fAli=gAlice->TreeK();
237c933d 2852 if (fAli) pFile =fAli->GetCurrentFile();
2853 pFile->cd();
33984590 2854 } // if Add
2855
2856 Int_t tracks[10];
2857 Int_t charges[10];
2858 //cout<<"Start filling digits \n "<<endl;
2859 Int_t nentries=list->GetEntriesFast();
2860 //printf(" \n \n nentries %d \n",nentries);
2861
2862 // start filling the digits
2863
2864 for (Int_t nent=0;nent<nentries;nent++) {
2865 AliRICHTransientDigit *address=(AliRICHTransientDigit*)list->At(nent);
2866 if (address==0) continue;
2867
2868 Int_t ich=address->fChamber;
2869 Int_t q=address->fSignal;
2870 iChamber=(AliRICHChamber*) (*fChambers)[ich];
2871 AliRICHResponse * response=iChamber->GetResponseModel();
2872 Int_t adcmax= (Int_t) response->MaxAdc();
2873
2874
2875 // add white noise and do zero-suppression and signal truncation (new electronics,old electronics gaus 1.2,0.2)
2876 //printf("Treshold: %d\n",iChamber->fTresh->GetHitIndex(address->fPadX,address->fPadY));
2877 Int_t pedestal = iChamber->fTresh->GetHitIndex(address->fPadX,address->fPadY);
ddae0931 2878
33984590 2879 //printf("Pedestal:%d\n",pedestal);
2880 //Int_t pedestal=0;
9dda3582 2881 Float_t treshold = (pedestal + 4*2.2);
33984590 2882
9dda3582 2883 Float_t meanNoise = gRandom->Gaus(2.2, 0.3);
33984590 2884 Float_t noise = gRandom->Gaus(0, meanNoise);
2885 q+=(Int_t)(noise + pedestal);
2886 //q+=(Int_t)(noise);
2887 // magic number to be parametrised !!!
2888 if ( q <= treshold)
2889 {
2890 q = q - pedestal;
2891 continue;
2e5f0f7b 2892 }
33984590 2893 q = q - pedestal;
2894 if ( q >= adcmax) q=adcmax;
2895 digits[0]=address->fPadX;
2896 digits[1]=address->fPadY;
2897 digits[2]=q;
2898
2899 TObjArray* trlist=(TObjArray*)address->TrackList();
2900 Int_t nptracks=trlist->GetEntriesFast();
2901
2902 // this was changed to accomodate the real number of tracks
2903 if (nptracks > 10) {
2904 cout<<"Attention - tracks > 10 "<<nptracks<<endl;
2905 nptracks=10;
2906 }
2907 if (nptracks > 2) {
2908 printf("Attention - tracks > 2 %d \n",nptracks);
2909 //printf("cat,ich,ix,iy,q %d %d %d %d %d \n",
2910 //icat,ich,digits[0],digits[1],q);
2911 }
2912 for (Int_t tr=0;tr<nptracks;tr++) {
2913 TVector *ppP=(TVector*)trlist->At(tr);
2914 TVector &pp =*ppP;
2915 tracks[tr]=Int_t(pp(0));
2916 charges[tr]=Int_t(pp(1));
2917 } //end loop over list of tracks for one pad
2918 if (nptracks < 10 ) {
2919 for (Int_t t=nptracks; t<10; t++) {
2920 tracks[t]=0;
2921 charges[t]=0;
ddae0931 2922 }
33984590 2923 }
2924 //write file
2925 if (ich==2)
2926 fprintf(points,"%4d, %4d, %4d\n",digits[0],digits[1],digits[2]);
2927
2928 // fill digits
2929 pRICH->AddDigits(ich,tracks,charges,digits);
2930 }
2931 gAlice->TreeD()->Fill();
2932
2933 list->Delete();
2934 for(Int_t ii=0;ii<kNCH;++ii) {
2935 if (pHitMap[ii]) {
2936 hm=pHitMap[ii];
2937 delete hm;
2938 pHitMap[ii]=0;
2939 }
2940 }
2941
2942 //TTree *TD=gAlice->TreeD();
2943 //Stat_t ndig=TD->GetEntries();
2944 //cout<<"number of digits "<<ndig<<endl;
2945 TClonesArray *fDch;
2946 for (int k=0;k<kNCH;k++) {
2947 fDch= pRICH->DigitsAddress(k);
2948 int ndigit=fDch->GetEntriesFast();
2949 printf ("Chamber %d digits %d \n",k,ndigit);
2950 }
2951 pRICH->ResetDigits();
ddae0931 2952 char hname[30];
2953 sprintf(hname,"TreeD%d",nev);
33984590 2954 gAlice->TreeD()->Write(hname,kOverwrite,0);
2955
2956 // reset tree
2957 // gAlice->TreeD()->Reset();
2e5f0f7b 2958 delete list;
237c933d 2959 pAddress->Clear();
33984590 2960 // gObjectTable->Print();
2e5f0f7b 2961}
ddae0931 2962
237c933d 2963AliRICH& AliRICH::operator=(const AliRICH& rhs)
2964{
2965// Assignment operator
2966 return *this;
2967
2968}
ddae0931 2969
237c933d 2970Int_t AliRICH::MakePadHits(Float_t xhit,Float_t yhit,Float_t eloss, Int_t idvol, ResponseType res)
2971{
2972//
2973// Calls the charge disintegration method of the current chamber and adds
2974// the simulated cluster to the root treee
2975//
176a9917 2976 Int_t clhits[5];
2977 Float_t newclust[4][500];
237c933d 2978 Int_t nnew;
2979
2980//
2981// Integrated pulse height on chamber
2982
2983 clhits[0]=fNhits+1;
2984
2985 ((AliRICHChamber*) (*fChambers)[idvol])->DisIntegration(eloss, xhit, yhit, nnew, newclust, res);
2986 Int_t ic=0;
2987
2988//
2989// Add new clusters
2990 for (Int_t i=0; i<nnew; i++) {
176a9917 2991 if (Int_t(newclust[0][i]) > 0) {
237c933d 2992 ic++;
237c933d 2993// Cluster Charge
176a9917 2994 clhits[1] = Int_t(newclust[0][i]);
237c933d 2995// Pad: ix
176a9917 2996 clhits[2] = Int_t(newclust[1][i]);
237c933d 2997// Pad: iy
176a9917 2998 clhits[3] = Int_t(newclust[2][i]);
237c933d 2999// Pad: chamber sector
176a9917 3000 clhits[4] = Int_t(newclust[3][i]);
237c933d 3001
3002 AddPadHit(clhits);
3003 }
3004 }
3005return nnew;
3006}
ddae0931 3007