]> git.uio.no Git - u/mrichter/AliRoot.git/blame - RICH/AliRICHv1.cxx
new Hits2SDigits.
[u/mrichter/AliRoot.git] / RICH / AliRICHv1.cxx
CommitLineData
c28632f0 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
88cb7938 16/* $Id$ */
c28632f0 17
88cb7938 18#include "Riostream.h"
19
c28632f0 20#include <TNode.h>
237c933d 21#include <TParticle.h>
88cb7938 22#include <TRandom.h>
23#include <TTUBE.h>
24#include <TVirtualMC.h>
d128c9d1 25#include <TPDGCode.h>
c28632f0 26
88cb7938 27#include "AliConst.h"
28#include "AliPDG.h"
29#include "AliRICHGeometry.h"
64e9b5aa 30#include "AliRICHResponse.h"
64e9b5aa 31#include "AliRICHResponseV0.h"
88cb7938 32#include "AliRICHSegmentationV1.h"
33#include "AliRICHv1.h"
c28632f0 34#include "AliRun.h"
c28632f0 35
d128c9d1 36ClassImp(AliRICHv1)
37//______________________________________________________________________________
c28632f0 38AliRICHv1::AliRICHv1(const char *name, const char *title)
d128c9d1 39 :AliRICH(name,title)
c28632f0 40{
c28632f0 41
237c933d 42// Full version of RICH with hits and diagnostics
43
64e9b5aa 44 // Version 0
45// Default Segmentation, no hits
03f3e21a 46 AliRICHSegmentationV1* segmentation = new AliRICHSegmentationV1;
64e9b5aa 47//
48// Segmentation parameters
c4384528 49 segmentation->SetPadSize(0.84,0.80);
50 segmentation->SetDAnod(0.84/2);
64e9b5aa 51//
52// Geometry parameters
53 AliRICHGeometry* geometry = new AliRICHGeometry;
54 geometry->SetGapThickness(8);
55 geometry->SetProximityGapThickness(.4);
2463d46d 56 geometry->SetQuartzLength(133);
57 geometry->SetQuartzWidth(127.9);
64e9b5aa 58 geometry->SetQuartzThickness(.5);
2463d46d 59 geometry->SetOuterFreonLength(133);
60 geometry->SetOuterFreonWidth(41.3);
61 geometry->SetInnerFreonLength(133);
62 geometry->SetInnerFreonWidth(41.3);
63 geometry->SetFreonThickness(1.5);
64e9b5aa 64//
65// Response parameters
c4384528 66 AliRICHResponseV0* response = new AliRICHResponseV0;
67 response->SetSigmaIntegration(5.);
68 response->SetChargeSlope(27.);
69 response->SetChargeSpread(0.18, 0.18);
70 response->SetMaxAdc(4096);
71 response->SetAlphaFeedback(0.036);
72 response->SetEIonisation(26.e-9);
73 response->SetSqrtKx3(0.77459667);
74 response->SetKx2(0.962);
75 response->SetKx4(0.379);
76 response->SetSqrtKy3(0.77459667);
77 response->SetKy2(0.962);
78 response->SetKy4(0.379);
79 response->SetPitch(0.25);
03f3e21a 80 response->SetWireSag(1); // 1->On, 0->Off
c4384528 81 response->SetVoltage(2150); // Should only be 2000, 2050, 2100 or 2150
82
64e9b5aa 83//
84// AliRICH *RICH = (AliRICH *) gAlice->GetDetector("RICH");
85
237c933d 86 fCkovNumber=0;
87 fFreonProd=0;
64e9b5aa 88 Int_t i=0;
c28632f0 89
64e9b5aa 90 fChambers = new TObjArray(kNCH);
91 for (i=0; i<kNCH; i++) {
92
2682e810 93 //PH (*fChambers)[i] = new AliRICHChamber();
94 fChambers->AddAt(new AliRICHChamber(), i);
64e9b5aa 95
c28632f0 96 }
64e9b5aa 97
98 for (i=0; i<kNCH; i++) {
99 SetGeometryModel(i,geometry);
c4384528 100 SetSegmentationModel(i, segmentation);
101 SetResponseModel(i, response);
64e9b5aa 102 }
103
104
c28632f0 105}
106
237c933d 107void AliRICHv1::Init()
108{
109
9e1a0ddb 110 if(fDebug) {
111 printf("%s: *********************************** RICH_INIT ***********************************\n",ClassName());
112 printf("%s: * *\n",ClassName());
113 printf("%s: * AliRICHv1 Full version started *\n",ClassName());
114 printf("%s: * *\n",ClassName());
115 }
237c933d 116
117
a2f7eaf6 118 AliSegmentation* segmentation;
237c933d 119 AliRICHGeometry* geometry;
120 AliRICHResponse* response;
121
122
123 //
124 // Initialize Tracking Chambers
125 //
2fc6393f 126 for (Int_t i=0; i<kNCH; i++) {
237c933d 127 //printf ("i:%d",i);
2682e810 128 //PH ( (AliRICHChamber*) (*fChambers)[i])->Init(i);
129 ( (AliRICHChamber*)fChambers->At(i))->Init(i);
237c933d 130 }
131
132 //
133 // Set the chamber (sensitive region) GEANT identifier
134
2682e810 135 //PH ((AliRICHChamber*)(*fChambers)[0])->SetGid(1);
136 //PH ((AliRICHChamber*)(*fChambers)[1])->SetGid(2);
137 //PH ((AliRICHChamber*)(*fChambers)[2])->SetGid(3);
138 //PH ((AliRICHChamber*)(*fChambers)[3])->SetGid(4);
139 //PH ((AliRICHChamber*)(*fChambers)[4])->SetGid(5);
140 //PH ((AliRICHChamber*)(*fChambers)[5])->SetGid(6);
141 //PH ((AliRICHChamber*)(*fChambers)[6])->SetGid(7);
142
143 ((AliRICHChamber*)fChambers->At(0))->SetGid(1);
144 ((AliRICHChamber*)fChambers->At(1))->SetGid(2);
145 ((AliRICHChamber*)fChambers->At(2))->SetGid(3);
146 ((AliRICHChamber*)fChambers->At(3))->SetGid(4);
147 ((AliRICHChamber*)fChambers->At(4))->SetGid(5);
148 ((AliRICHChamber*)fChambers->At(5))->SetGid(6);
149 ((AliRICHChamber*)fChambers->At(6))->SetGid(7);
150
237c933d 151
237c933d 152 segmentation=Chamber(0).GetSegmentationModel(0);
153 geometry=Chamber(0).GetGeometryModel();
154 response=Chamber(0).GetResponseModel();
155
449e6743 156 Float_t offset = 490 + 1.276 - geometry->GetGapThickness()/2; //distance from center of mother volume to methane
157 Float_t deltaphi = 19.5; //phi angle between center of chambers - z direction
158 Float_t deltatheta = 20; //theta angle between center of chambers - x direction
159 Float_t cosphi = TMath::Cos(deltaphi*TMath::Pi()/180);
160 Float_t sinphi = TMath::Sin(deltaphi*TMath::Pi()/180);
161 Float_t costheta = TMath::Cos(deltatheta*TMath::Pi()/180);
162 Float_t sintheta = TMath::Sin(deltatheta*TMath::Pi()/180);
163
164 Float_t pos1[3]={0. , offset*cosphi , offset*sinphi};
165 Float_t pos2[3]={offset*sintheta , offset*costheta , 0. };
166 Float_t pos3[3]={0. , offset , 0.};
167 Float_t pos4[3]={-offset*sintheta , offset*costheta , 0.};
168 Float_t pos5[3]={offset*sinphi , offset*costheta*cosphi, -offset*sinphi};
169 Float_t pos6[3]={0. , offset*cosphi , -offset*sinphi};
170 Float_t pos7[3]={ -offset*sinphi , offset*costheta*cosphi, -offset*sinphi};
171
172 Chamber(0).SetChamberTransform(pos1[0],pos1[1],pos1[2],new TRotMatrix("rot993","rot993",90., 0. , 90. - deltaphi, 90. , deltaphi, -90. ));
173 Chamber(1).SetChamberTransform(pos2[0],pos2[1],pos2[2],new TRotMatrix("rot994","rot994",90., -deltatheta , 90. , 90.- deltatheta , 0. , 0. ));
174 Chamber(2).SetChamberTransform(pos3[0],pos3[1],pos3[2],new TRotMatrix("rot995","rot995",90., 0. , 90. , 90. , 0. , 0. ));
175 Chamber(3).SetChamberTransform(pos4[0],pos4[1],pos4[2],new TRotMatrix("rot996","rot996",90., deltatheta , 90. , 90 + deltatheta , 0. , 0. ));
176 Chamber(4).SetChamberTransform(pos5[0],pos5[1],pos5[2],new TRotMatrix("rot997","rot997",90., 360. - deltatheta, 108.2 , 90.- deltatheta ,18.2 , 90 - deltatheta));
177 Chamber(5).SetChamberTransform(pos6[0],pos6[1],pos6[2],new TRotMatrix("rot998","rot998",90., 0. , 90 + deltaphi , 90. , deltaphi, 90. ));
178 Chamber(6).SetChamberTransform(pos7[0],pos7[1],pos7[2],new TRotMatrix("rot999","rot999",90., deltatheta , 108.2 , 90.+ deltatheta ,18.2 , 90 + deltatheta));
237c933d 179
d128c9d1 180 if(GetDebug()) {
9e1a0ddb 181 printf("%s: * Pads : %3dx%3d *\n",
182 ClassName(),segmentation->Npx(),segmentation->Npy());
183 printf("%s: * Pad size : %5.2f x%5.2f mm2 *\n",
184 ClassName(),segmentation->Dpx(),segmentation->Dpy());
185 printf("%s: * Gap Thickness : %5.1f cm *\n",
186 ClassName(),geometry->GetGapThickness());
187 printf("%s: * Radiator Width : %5.1f cm *\n",
188 ClassName(),geometry->GetQuartzWidth());
189 printf("%s: * Radiator Length : %5.1f cm *\n",
190 ClassName(),geometry->GetQuartzLength());
191 printf("%s: * Freon Thickness : %5.1f cm *\n",
192 ClassName(),geometry->GetFreonThickness());
193 printf("%s: * Charge Slope : %5.1f ADC *\n",
194 ClassName(),response->ChargeSlope());
195 printf("%s: * Feedback Prob. : %5.2f %% *\n",
196 ClassName(),response->AlphaFeedback()*100);
9e1a0ddb 197 printf("%s: * *\n",
198 ClassName());
199 printf("%s: *********************************************************************************\n",
200 ClassName());
201 }
237c933d 202}
c28632f0 203
d128c9d1 204//______________________________________________________________________________
205void AliRICHv1::StepManager()
206{//Full Step Manager
207
208 Int_t copy, id;
209 static Int_t idvol;
210 static Int_t vol[2];
211 Int_t ipart;
212 static Float_t hits[22];
213 static Float_t ckovData[19];
214 TLorentzVector position;
215 TLorentzVector momentum;
216 Float_t pos[3];
217 Float_t mom[4];
218 Float_t localPos[3];
219 Float_t localMom[4];
220 Float_t localTheta,localPhi;
221 Float_t theta,phi;
222 Float_t destep, step;
223 Double_t ranf[2];
224 Int_t nPads;
225 Float_t coscerenkov;
226 static Float_t eloss, xhit, yhit, tlength;
227 const Float_t kBig=1.e10;
228
229 TClonesArray &lhits = *fHits;
642f15cf 230 TParticle *current = (TParticle*)(*gAlice->Particles())[gAlice->GetCurrentTrackNumber()];
d128c9d1 231
d128c9d1 232
233
234 id=gMC->CurrentVolID(copy);
235 idvol = copy-1;
236 Float_t cherenkovLoss=0;
d128c9d1 237
238 gMC->TrackPosition(position);
239 pos[0]=position(0);
240 pos[1]=position(1);
241 pos[2]=position(2);
d128c9d1 242 ckovData[1] = pos[0]; // X-position for hit
243 ckovData[2] = pos[1]; // Y-position for hit
244 ckovData[3] = pos[2]; // Z-position for hit
245 ckovData[6] = 0; // dummy track length
d128c9d1 246
247 /********************Store production parameters for Cerenkov photons************************/
853634d3 248
249 if (gMC->TrackPid() == 50000050) { //is it a Cerenkov photon?
d128c9d1 250
251 //if (gMC->VolId("GAP ")==gMC->CurrentVolID(copy))
252 //{
253 Float_t ckovEnergy = current->Energy();
254 //energy interval for tracking
255 if (ckovEnergy > 5.6e-09 && ckovEnergy < 7.8e-09 )
256 //if (ckovEnergy > 0)
257 {
258 if (gMC->IsTrackEntering()){ //is track entering?
259 //printf("Track entered (1)\n");
260 if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
261 { //is it in freo?
262 if (gMC->IsNewTrack()){ //is it the first step?
263 //printf("I'm in!\n");
264 Int_t mother = current->GetFirstMother();
265
d128c9d1 266
267 ckovData[10] = mother;
642f15cf 268 ckovData[11] = gAlice->GetCurrentTrackNumber();
d128c9d1 269 ckovData[12] = 1; //Media where photon was produced 1->Freon, 2->Quarz
270 //printf("Produced in FREO\n");
271 fCkovNumber++;
272 fFreonProd=1;
d128c9d1 273 } //first step question
274 } //freo question
275
276 if (gMC->IsNewTrack()){ //is it first step?
277 if (gMC->VolId("QUAR")==gMC->CurrentVolID(copy)) //is it in quarz?
278 {
279 ckovData[12] = 2;
d128c9d1 280 } //quarz question
281 } //first step question
282
d128c9d1 283 } //track entering question
284
285 if (ckovData[12] == 1) //was it produced in Freon?
286 //if (fFreonProd == 1)
287 {
288 if (gMC->IsTrackEntering()){ //is track entering?
d128c9d1 289 if (gMC->VolId("META")==gMC->CurrentVolID(copy)) //is it in gap?
290 {
291 //printf("Got in META\n");
292 gMC->TrackMomentum(momentum);
293 mom[0]=momentum(0);
294 mom[1]=momentum(1);
295 mom[2]=momentum(2);
296 mom[3]=momentum(3);
297
298 gMC->Gmtod(mom,localMom,2);
299 Float_t cophi = TMath::Cos(TMath::ATan2(localMom[0], localMom[1]));
300 Float_t t = (1. - .025 / cophi) * (1. - .05 / cophi);
301 /**************** Photons lost in second grid have to be calculated by hand************/
302 gMC->GetRandom()->RndmArray(1,ranf);
303 if (ranf[0] > t) {
304 gMC->StopTrack();
305 ckovData[13] = 5;
642f15cf 306 AddCerenkov(gAlice->GetCurrentTrackNumber(),vol,ckovData);
d128c9d1 307 }
308 /**********************************************************************************/
309 } //gap
310
d128c9d1 311 if (gMC->VolId("CSI ")==gMC->CurrentVolID(copy)) //is it in csi?
312 {
313 //printf("Got in CSI\n");
314 gMC->TrackMomentum(momentum);
315 mom[0]=momentum(0);
316 mom[1]=momentum(1);
317 mom[2]=momentum(2);
318 mom[3]=momentum(3);
319
320 gMC->Gmtod(mom,localMom,2);
321 /********* Photons lost by Fresnel reflection have to be calculated by hand********/
322 /***********************Cerenkov phtons (always polarised)*************************/
323 Double_t localTc = localMom[0]*localMom[0]+localMom[2]*localMom[2];
324 Double_t localRt = TMath::Sqrt(localTc);
325 localTheta = Float_t(TMath::ATan2(localRt,Double_t(localMom[1])));
326 Double_t cotheta = TMath::Abs(cos(localTheta));
327 Float_t t = Fresnel(ckovEnergy*1e9,cotheta,1);
328 gMC->GetRandom()->RndmArray(1,ranf);
329 if (ranf[0] < t) {
330 gMC->StopTrack();
331 ckovData[13] = 6;
642f15cf 332 AddCerenkov(gAlice->GetCurrentTrackNumber(),vol,ckovData);
d128c9d1 333
334 //printf("Added One (2)!\n");
335 //printf("Lost by Fresnel\n");
336 }
337 /**********************************************************************************/
338 }
339 } //track entering?
340
341
342 /********************Evaluation of losses************************/
343 /******************still in the old fashion**********************/
344
345 TArrayI procs;
346 Int_t i1 = gMC->StepProcesses(procs); //number of physics mechanisms acting on the particle
347 for (Int_t i = 0; i < i1; ++i) {
348 // Reflection loss
349 if (procs[i] == kPLightReflection) { //was it reflected
350 ckovData[13]=10;
351 if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
352 ckovData[13]=1;
353 if (gMC->CurrentVolID(copy) == gMC->VolId("QUAR"))
354 ckovData[13]=2;
d128c9d1 355 } //reflection question
356
357 // Absorption loss
358 else if (procs[i] == kPLightAbsorption) { //was it absorbed?
359 //printf("Got in absorption\n");
360 ckovData[13]=20;
361 if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
362 ckovData[13]=11;
363 if (gMC->CurrentVolID(copy) == gMC->VolId("QUAR"))
364 ckovData[13]=12;
365 if (gMC->CurrentVolID(copy) == gMC->VolId("META"))
366 ckovData[13]=13;
367 if (gMC->CurrentVolID(copy) == gMC->VolId("GAP "))
368 ckovData[13]=13;
369
370 if (gMC->CurrentVolID(copy) == gMC->VolId("SRIC"))
371 ckovData[13]=15;
372
373 // CsI inefficiency
374 if (gMC->CurrentVolID(copy) == gMC->VolId("CSI ")) {
375 ckovData[13]=16;
376 }
377 gMC->StopTrack();
642f15cf 378 AddCerenkov(gAlice->GetCurrentTrackNumber(),vol,ckovData);
d128c9d1 379 } //absorption question
380
381
382 // Photon goes out of tracking scope
383 else if (procs[i] == kPStop) { //is it below energy treshold?
384 ckovData[13]=21;
385 gMC->StopTrack();
642f15cf 386 AddCerenkov(gAlice->GetCurrentTrackNumber(),vol,ckovData);
d128c9d1 387 } // energy treshold question
388 } //number of mechanisms cycle
389 /**********************End of evaluation************************/
390 } //freon production question
391 } //energy interval question
d128c9d1 392 } //cerenkov photon question
393
394 /**************************************End of Production Parameters Storing*********************/
395
396
397 /*******************************Treat photons that hit the CsI (Ckovs and Feedbacks)************/
398
399 if (gMC->TrackPid() == 50000050 || gMC->TrackPid() == 50000051) {
d128c9d1 400
401 if (gMC->VolId("CSI ")==gMC->CurrentVolID(copy))
402 {
d128c9d1 403 if (gMC->Edep() > 0.){
404 gMC->TrackPosition(position);
405 gMC->TrackMomentum(momentum);
406 pos[0]=position(0);
407 pos[1]=position(1);
408 pos[2]=position(2);
409 mom[0]=momentum(0);
410 mom[1]=momentum(1);
411 mom[2]=momentum(2);
412 mom[3]=momentum(3);
413 Double_t tc = mom[0]*mom[0]+mom[1]*mom[1];
414 Double_t rt = TMath::Sqrt(tc);
415 theta = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg;
416 phi = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg;
417
418 gMC->CurrentVolOffID(2,copy);
419 vol[0]=copy;
420 idvol=vol[0]-1;
421
422
423 gMC->Gmtod(pos,localPos,1);
424
d128c9d1 425
426 gMC->Gmtod(mom,localMom,2);
427
d128c9d1 428
429 gMC->CurrentVolOffID(2,copy);
430 vol[0]=copy;
431 idvol=vol[0]-1;
432
d128c9d1 433
d128c9d1 434 ((AliRICHChamber*)fChambers->At(idvol))
435 ->SigGenInit(localPos[0], localPos[2], localPos[1]);
436 if(idvol<kNCH) {
437 ckovData[0] = gMC->TrackPid(); // particle type
438 ckovData[1] = pos[0]; // X-position for hit
439 ckovData[2] = pos[1]; // Y-position for hit
440 ckovData[3] = pos[2]; // Z-position for hit
441 ckovData[4] = theta; // theta angle of incidence
442 ckovData[5] = phi; // phi angle of incidence
853634d3 443 ckovData[8] = (Float_t) fNsdigits; // first sdigit
d128c9d1 444 ckovData[9] = -1; // last pad hit
445 ckovData[13] = 4; // photon was detected
446 ckovData[14] = mom[0];
447 ckovData[15] = mom[1];
448 ckovData[16] = mom[2];
449
450 destep = gMC->Edep();
451 gMC->SetMaxStep(kBig);
452 cherenkovLoss += destep;
453 ckovData[7]=cherenkovLoss;
454
455 //nPads = Hits2SDigits(localPos[0],localPos[2],cherenkovLoss,idvol,kCerenkov);//for photons in CsI kir
456
853634d3 457 if (fNsdigits > (Int_t)ckovData[8]) {
d128c9d1 458 ckovData[8]= ckovData[8]+1;
853634d3 459 ckovData[9]= (Float_t) fNsdigits;
d128c9d1 460 }
461
d128c9d1 462
463 ckovData[17] = nPads;
853634d3 464 AliRICHhit *mipHit = (AliRICHhit*) (fHits->UncheckedAt(0));
d128c9d1 465 if (mipHit)
466 {
467 mom[0] = current->Px();
468 mom[1] = current->Py();
469 mom[2] = current->Pz();
470 Float_t mipPx = mipHit->MomX();
471 Float_t mipPy = mipHit->MomY();
472 Float_t mipPz = mipHit->MomZ();
473
474 Float_t r = mom[0]*mom[0] + mom[1]*mom[1] + mom[2]*mom[2];
475 Float_t rt = TMath::Sqrt(r);
476 Float_t mipR = mipPx*mipPx + mipPy*mipPy + mipPz*mipPz;
477 Float_t mipRt = TMath::Sqrt(mipR);
478 if ((rt*mipRt) > 0)
479 {
480 coscerenkov = (mom[0]*mipPx + mom[1]*mipPy + mom[2]*mipPz)/(rt*mipRt);
481 }
482 else
483 {
484 coscerenkov = 0;
485 }
486 Float_t cherenkov = TMath::ACos(coscerenkov);
487 ckovData[18]=cherenkov;
488 }
489 //if (sector != -1)
490 //{
642f15cf 491 AddHit(gAlice->GetCurrentTrackNumber(),vol,ckovData);
492 AddCerenkov(gAlice->GetCurrentTrackNumber(),vol,ckovData);
d128c9d1 493 //printf("Added One (5)!\n");
494 //}
495 }
496 }
497 }
498 }
499
500 /***********************************************End of photon hits*********************************************/
501
502
503 /**********************************************Charged particles treatment*************************************/
504
853634d3 505 else if (gMC->TrackCharge()){//is MIP?
d128c9d1 506 if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
507 {
508 gMC->TrackMomentum(momentum);
509 mom[0]=momentum(0);
510 mom[1]=momentum(1);
511 mom[2]=momentum(2);
512 mom[3]=momentum(3);
513 hits [19] = mom[0];
514 hits [20] = mom[1];
515 hits [21] = mom[2];
516 fFreonProd=1;
517 }
518
853634d3 519 if(gMC->VolId("GAP ")==gMC->CurrentVolID(copy)) {//is in GAP?
d128c9d1 520// Get current particle id (ipart), track position (pos) and momentum (mom)
521
522 gMC->CurrentVolOffID(3,copy);
523 vol[0]=copy;
524 idvol=vol[0]-1;
d128c9d1 525 gMC->TrackPosition(position);
526 gMC->TrackMomentum(momentum);
527 pos[0]=position(0);
528 pos[1]=position(1);
529 pos[2]=position(2);
530 mom[0]=momentum(0);
531 mom[1]=momentum(1);
532 mom[2]=momentum(2);
533 mom[3]=momentum(3);
534
535 gMC->Gmtod(pos,localPos,1);
d128c9d1 536 gMC->Gmtod(mom,localMom,2);
d128c9d1 537 ipart = gMC->TrackPid();
853634d3 538 destep = gMC->Edep();step = gMC->TrackStep();// momentum loss and steplength in last step
539 if(gMC->IsTrackEntering()){ // record hits when track enters ...
540
d128c9d1 541 Double_t tc = mom[0]*mom[0]+mom[1]*mom[1];
542 Double_t rt = TMath::Sqrt(tc);
543 theta = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg;
544 phi = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg;
545
546
547 Double_t localTc = localMom[0]*localMom[0]+localMom[2]*localMom[2];
548 Double_t localRt = TMath::Sqrt(localTc);
549 localTheta = Float_t(TMath::ATan2(localRt,Double_t(localMom[1])))*kRaddeg;
550 localPhi = Float_t(TMath::ATan2(Double_t(localMom[2]),Double_t(localMom[0])))*kRaddeg;
551
552 hits[0] = Float_t(ipart); // particle type
553 hits[1] = localPos[0]; // X-position for hit
554 hits[2] = localPos[1]; // Y-position for hit
555 hits[3] = localPos[2]; // Z-position for hit
556 hits[4] = localTheta; // theta angle of incidence
557 hits[5] = localPhi; // phi angle of incidence
853634d3 558 hits[8] = (Float_t) fNsdigits; // first sdigit
d128c9d1 559 hits[9] = -1; // last pad hit
560 hits[13] = fFreonProd; // did id hit the freon?
561 hits[14] = mom[0];
562 hits[15] = mom[1];
563 hits[16] = mom[2];
564 hits[18] = 0; // dummy cerenkov angle
565
566 tlength = 0;
567 eloss = 0;
568 fFreonProd = 0;
569
570 Chamber(idvol).LocaltoGlobal(localPos,hits+1);
571
572
573 //To make chamber coordinates x-y had to pass localPos[0], localPos[2]
574 xhit = localPos[0];
575 yhit = localPos[2];
576 // Only if not trigger chamber
577 if(idvol<kNCH) {
d128c9d1 578 // Initialize hit position (cursor) in the segmentation model
d128c9d1 579 ((AliRICHChamber*)fChambers->At(idvol))
580 ->SigGenInit(localPos[0], localPos[2], localPos[1]);
581 }
582 }
583
d128c9d1 584 // Calculate the charge induced on a pad (disintegration) in case
d128c9d1 585 // Mip left chamber ...
586 if( gMC->IsTrackExiting() || gMC->IsTrackStop() || gMC->IsTrackDisappeared()){
587 gMC->SetMaxStep(kBig);
588 eloss += destep;
589 tlength += step;
590
591
592 // Only if not trigger chamber
593 if(idvol<kNCH) {
594 if (eloss > 0)
595 {
596 if(gMC->TrackPid() == kNeutron)
597 printf("\n\n\n\n\n Neutron Making Pad Hit!!! \n\n\n\n");
598 //nPads = Hits2SDigits(xhit,yhit,eloss,idvol,kMip); //for MIP kir
599 hits[17] = nPads;
d128c9d1 600 }
601 }
602
603 hits[6]=tlength;
604 hits[7]=eloss;
853634d3 605 if (fNsdigits > (Int_t)hits[8]) {
d128c9d1 606 hits[8]= hits[8]+1;
853634d3 607 hits[9]= (Float_t) fNsdigits;
d128c9d1 608 }
609
853634d3 610 new(lhits[fNhits++]) AliRICHhit(fIshunt,gAlice->GetCurrentTrackNumber(),vol,hits);
d128c9d1 611 eloss = 0;
d128c9d1 612 // Check additional signal generation conditions
613 // defined by the segmentation
614 // model (boundary crossing conditions)
615 } else if
d128c9d1 616 (((AliRICHChamber*)fChambers->At(idvol))
617 ->SigGenCond(localPos[0], localPos[2], localPos[1]))
618 {
d128c9d1 619 ((AliRICHChamber*)fChambers->At(idvol))
620 ->SigGenInit(localPos[0], localPos[2], localPos[1]);
621 if (eloss > 0)
622 {
623 if(gMC->TrackPid() == kNeutron)
624 printf("\n\n\n\n\n Neutron Making Pad Hit!!! \n\n\n\n");
625 //nPads = Hits2SDigits(xhit,yhit,eloss,idvol,kMip);//for N kir
626 hits[17] = nPads;
627 //printf("Npads:%d",NPads);
628 }
629 xhit = localPos[0];
630 yhit = localPos[2];
631 eloss = destep;
632 tlength += step ;
633 //
634 // nothing special happened, add up energy loss
635 } else {
636 eloss += destep;
637 tlength += step ;
638 }
853634d3 639 }//is in GAP?
640 }//is MIP?
d128c9d1 641}//void AliRICHv1::StepManager()