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