Removed setters.
[u/mrichter/AliRoot.git] / RICH / AliRICHv1.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /*
17   $Log$
18   Revision 1.4  2000/06/13 13:13:40  jbarbosa
19   Correcting previous correction...
20
21   Revision 1.3  2000/06/13 13:06:38  jbarbosa
22   Fixed compiling error for HP (multiple declaration)
23
24   Revision 1.2  2000/06/12 15:36:16  jbarbosa
25   Cleaned up version.
26
27   Revision 1.1  2000/06/09 15:00:31  jbarbosa
28   New full version. All parameters configurable.
29
30   Revision 1.9  2000/05/31 08:19:38  jbarbosa
31   Fixed bug in StepManager
32
33   Revision 1.8  2000/05/26 17:30:08  jbarbosa
34   Cerenkov angle now stored within cerenkov data structure.
35
36   Revision 1.7  2000/05/18 10:31:36  jbarbosa
37   Fixed positioning of spacers inside freon.
38   Fixed positioning of proximity gap
39   inside methane.
40   Fixed cut on neutral particles in the StepManager.
41
42   Revision 1.6  2000/04/28 11:51:58  morsch
43    Dimensions of arrays hits and Ckov_data corrected.
44
45   Revision 1.5  2000/04/19 13:28:46  morsch
46   Major changes in geometry (parametrised), materials (updated) and
47   step manager (diagnostics) (JB, AM)
48
49 */
50
51
52
53 //////////////////////////////////////////////////////////
54 //  Manager and hits classes for set: RICH full version //
55 //////////////////////////////////////////////////////////
56
57 #include <TTUBE.h>
58 #include <TNode.h> 
59 #include <TRandom.h> 
60 #include <TParticle.h> 
61
62 #include "AliRICHv1.h"
63 #include "AliRICHHit.h"
64 #include "AliRICHSegmentation.h"
65 #include "AliRICHResponse.h"
66 #include "AliRICHSegmentationV0.h"
67 #include "AliRICHResponseV0.h"
68 #include "AliRICHGeometry.h"
69 #include "AliRun.h"
70 #include "AliMC.h"
71 #include "iostream.h"
72 #include "AliCallf77.h"
73 #include "AliConst.h" 
74 #include "AliPDG.h" 
75 #include "TGeant3.h"
76
77 ClassImp(AliRICHv1)
78     
79 //___________________________________________
80 AliRICHv1::AliRICHv1() : AliRICHv0()
81 {
82
83 // Default constructor fo AliRICHvv1 (full version)
84
85     //fChambers = 0;
86 }
87
88 //___________________________________________
89 AliRICHv1::AliRICHv1(const char *name, const char *title)
90     : AliRICHv0(name,title)
91 {
92
93 // Full version of RICH with hits and diagnostics
94
95   // Version 0
96 // Default Segmentation, no hits
97     AliRICHSegmentationV0* segmentationV0 = new AliRICHSegmentationV0;
98 //
99 //  Segmentation parameters
100     segmentationV0->SetPadSize(0.84,0.80);
101     segmentationV0->SetDAnod(0.84/2);
102 //
103 //  Geometry parameters
104     AliRICHGeometry* geometry = new AliRICHGeometry;
105     geometry->SetGapThickness(8);
106     geometry->SetProximityGapThickness(.4);
107     geometry->SetQuartzLength(131);
108     geometry->SetQuartzWidth(126.2);
109     geometry->SetQuartzThickness(.5);
110     geometry->SetOuterFreonLength(131);
111     geometry->SetOuterFreonWidth(40.3);
112     geometry->SetInnerFreonLength(131);
113     geometry->SetInnerFreonWidth(40.3);
114     geometry->SetFreonThickness(1);
115 //
116 //  Response parameters
117     AliRICHResponseV0*  responseV0   = new AliRICHResponseV0;
118     responseV0->SetSigmaIntegration(5.);
119     responseV0->SetChargeSlope(40.);
120     responseV0->SetChargeSpread(0.18, 0.18);
121     responseV0->SetMaxAdc(1024);
122     responseV0->SetAlphaFeedback(0.05);
123     responseV0->SetEIonisation(26.e-9);
124     responseV0->SetSqrtKx3(0.77459667);
125     responseV0->SetKx2(0.962);
126     responseV0->SetKx4(0.379);
127     responseV0->SetSqrtKy3(0.77459667);
128     responseV0->SetKy2(0.962);
129     responseV0->SetKy4(0.379);
130     responseV0->SetPitch(0.25);
131 //
132 //
133 //    AliRICH *RICH = (AliRICH *) gAlice->GetDetector("RICH"); 
134     
135     fCkovNumber=0;
136     fFreonProd=0;
137     Int_t i=0;
138     
139     fChambers = new TObjArray(kNCH);
140     for (i=0; i<kNCH; i++) {
141       
142       (*fChambers)[i] = new AliRICHChamber();  
143       
144     }
145   
146     for (i=0; i<kNCH; i++) {
147       SetGeometryModel(i,geometry);
148       SetSegmentationModel(i, segmentationV0);
149       SetResponseModel(i, responseV0);
150       SetNsec(i,1);
151       SetDebugLevel(0);
152     }
153
154
155 }
156
157 void AliRICHv1::Init()
158 {
159
160   printf("*********************************** RICH_INIT ***********************************\n");
161   printf("*                                                                               *\n");
162   printf("*                        AliRICHv1 Full version started                         *\n");
163   printf("*                                                                               *\n");
164
165   
166   AliRICHSegmentation*  segmentation;
167   AliRICHGeometry*  geometry;
168   AliRICHResponse*  response;
169
170
171     // 
172     // Initialize Tracking Chambers
173     //
174     for (Int_t i=1; i<kNCH; i++) {
175         //printf ("i:%d",i);
176         ( (AliRICHChamber*) (*fChambers)[i])->Init();  
177     }  
178     
179     //
180     // Set the chamber (sensitive region) GEANT identifier
181     
182     ((AliRICHChamber*)(*fChambers)[0])->SetGid(1);  
183     ((AliRICHChamber*)(*fChambers)[1])->SetGid(2);  
184     ((AliRICHChamber*)(*fChambers)[2])->SetGid(3);  
185     ((AliRICHChamber*)(*fChambers)[3])->SetGid(4);  
186     ((AliRICHChamber*)(*fChambers)[4])->SetGid(5);  
187     ((AliRICHChamber*)(*fChambers)[5])->SetGid(6);  
188     ((AliRICHChamber*)(*fChambers)[6])->SetGid(7); 
189
190     Float_t pos1[3]={0,471.8999,165.2599};
191     Chamber(0).SetChamberTransform(pos1[0],pos1[1],pos1[2],new TRotMatrix("rot993","rot993",90,0,70.69,90,19.30999,-90));
192
193     Float_t pos2[3]={171,470,0};
194     Chamber(1).SetChamberTransform(pos2[0],pos2[1],pos2[2],new TRotMatrix("rot994","rot994",90,-20,90,70,0,0));
195
196     Float_t pos3[3]={0,500,0};
197     Chamber(2).SetChamberTransform(pos3[0],pos3[1],pos3[2],new TRotMatrix("rot995","rot995",90,0,90,90,0,0));
198     
199     Float_t pos4[3]={-171,470,0};
200     Chamber(3).SetChamberTransform(pos4[0],pos4[1],pos4[2], new TRotMatrix("rot996","rot996",90,20,90,110,0,0));  
201
202     Float_t pos5[3]={161.3999,443.3999,-165.3};
203     Chamber(4).SetChamberTransform(pos5[0],pos5[1],pos5[2],new TRotMatrix("rot997","rot997",90,340,108.1999,70,18.2,70));
204
205     Float_t pos6[3]={0., 471.9, -165.3,};
206     Chamber(5).SetChamberTransform(pos6[0],pos6[1],pos6[2],new TRotMatrix("rot998","rot998",90,0,109.3099,90,19.30999,90));
207
208     Float_t pos7[3]={-161.399,443.3999,-165.3};
209     Chamber(6).SetChamberTransform(pos7[0],pos7[1],pos7[2],new TRotMatrix("rot999","rot999",90,20,108.1999,110,18.2,110));
210     
211     segmentation=Chamber(0).GetSegmentationModel(0);
212     geometry=Chamber(0).GetGeometryModel();
213     response=Chamber(0).GetResponseModel();
214     
215      
216     printf("*                            Pads            : %3dx%3d                          *\n",segmentation->Npx(),segmentation->Npy());
217     printf("*                            Pad size        : %5.2f x%5.2f mm2                 *\n",segmentation->Dpx(),segmentation->Dpy()); 
218     printf("*                            Gap Thickness   : %5.1f cm                         *\n",geometry->GetGapThickness());
219     printf("*                            Radiator Width  : %5.1f cm                         *\n",geometry->GetQuartzWidth());
220     printf("*                            Radiator Length : %5.1f cm                         *\n",geometry->GetQuartzLength());
221     printf("*                            Freon Thickness : %5.1f cm                         *\n",geometry->GetFreonThickness());
222     printf("*                            Charge Slope    : %5.1f ADC                        *\n",response->ChargeSlope());
223     printf("*                            Feedback Prob.  : %5.2f %%                         *\n",response->AlphaFeedback()*100);
224     printf("*                            Debug Level     : %3d                              *\n",GetDebugLevel());
225     printf("*                                                                               *\n");
226     printf("*                                   Success!                                    *\n");
227     printf("*                                                                               *\n");
228     printf("*********************************************************************************\n");
229
230 }
231
232 //___________________________________________
233 void AliRICHv1::StepManager()
234 {
235
236 // Full Step Manager
237
238     Int_t          copy, id;
239     static Int_t   idvol;
240     static Int_t   vol[2];
241     Int_t          ipart;
242     static Float_t hits[18];
243     static Float_t ckovData[19];
244     TLorentzVector position;
245     TLorentzVector momentum;
246     Float_t        pos[3];
247     Float_t        mom[4];
248     Float_t        localPos[3];
249     Float_t        localMom[4];
250     Float_t        localTheta,localPhi;
251     Float_t        theta,phi;
252     Float_t        destep, step;
253     Float_t        ranf[2];
254     Int_t          nPads;
255     Float_t        coscerenkov;
256     static Float_t eloss, xhit, yhit, tlength;
257     const  Float_t kBig=1.e10;
258        
259     TClonesArray &lhits = *fHits;
260     TGeant3 *geant3 = (TGeant3*) gMC;
261     TParticle *current = (TParticle*)(*gAlice->Particles())[gAlice->CurrentTrack()];
262
263  //if (current->Energy()>1)
264    //{
265         
266     // Only gas gap inside chamber
267     // Tag chambers and record hits when track enters 
268     
269     idvol=-1;
270     id=gMC->CurrentVolID(copy);
271     Float_t cherenkovLoss=0;
272     //gAlice->KeepTrack(gAlice->CurrentTrack());
273     
274     gMC->TrackPosition(position);
275     pos[0]=position(0);
276     pos[1]=position(1);
277     pos[2]=position(2);
278     ckovData[1] = pos[0];                 // X-position for hit
279     ckovData[2] = pos[1];                 // Y-position for hit
280     ckovData[3] = pos[2];                 // Z-position for hit
281     //ckovData[11] = gAlice->CurrentTrack();
282
283     //AliRICH *RICH = (AliRICH *) gAlice->GetDetector("RICH"); 
284     
285     /********************Store production parameters for Cerenkov photons************************/ 
286 //is it a Cerenkov photon? 
287     if (gMC->TrackPid() == 50000050) {          
288
289       //if (gMC->VolId("GAP ")==gMC->CurrentVolID(copy))
290         //{                    
291           Float_t ckovEnergy = current->Energy();
292           //energy interval for tracking
293           if  (ckovEnergy > 5.6e-09 && ckovEnergy < 7.8e-09 )       
294             //if (ckovEnergy > 0)
295             {
296               if (gMC->IsTrackEntering()){                                     //is track entering?
297                 if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
298                   {                                                          //is it in freo?
299                     if (geant3->Gctrak()->nstep<1){                          //is it the first step?
300                       Int_t mother = current->GetFirstMother(); 
301                       
302                       //printf("Second Mother:%d\n",current->GetSecondMother());
303                       
304                       ckovData[10] = mother;
305                       ckovData[11] = gAlice->CurrentTrack();
306                       ckovData[12] = 1;             //Media where photon was produced 1->Freon, 2->Quarz
307                       fCkovNumber++;
308                       fFreonProd=1;
309                       //printf("Index: %d\n",fCkovNumber);
310                     }    //first step question
311                   }        //freo question
312                 
313                 if (geant3->Gctrak()->nstep<1){                                  //is it first step?
314                   if (gMC->VolId("QUAR")==gMC->CurrentVolID(copy))             //is it in quarz?
315                     {
316                       ckovData[12] = 2;
317                     }    //quarz question
318                 }        //first step question
319                 
320                 //printf("Before %d\n",fFreonProd);
321               }   //track entering question
322               
323               if (ckovData[12] == 1)                                        //was it produced in Freon?
324                 //if (fFreonProd == 1)
325                 {
326                   if (gMC->IsTrackEntering()){                                     //is track entering?
327                     //printf("Got in");
328                     if (gMC->VolId("META")==gMC->CurrentVolID(copy))                //is it in gap?      
329                       {
330                         //printf("Got in\n");
331                         gMC->TrackMomentum(momentum);
332                         mom[0]=momentum(0);
333                         mom[1]=momentum(1);
334                         mom[2]=momentum(2);
335                         mom[3]=momentum(3);
336                         // Z-position for hit
337                         
338                         
339                         /**************** Photons lost in second grid have to be calculated by hand************/ 
340                         
341                         Float_t cophi = TMath::Cos(TMath::ATan2(mom[0], mom[1]));
342                         Float_t t = (1. - .025 / cophi) * (1. - .05 /  cophi);
343                         gMC->Rndm(ranf, 1);
344                         //printf("grid calculation:%f\n",t);
345                         if (ranf[0] > t) {
346                           //geant3->StopTrack();
347                           ckovData[13] = 5;
348                           AddCerenkov(gAlice->CurrentTrack(),vol,ckovData);
349                           //printf("Lost one in grid\n");
350                         }
351                         /**********************************************************************************/
352                       }    //gap
353                     
354                     if (gMC->VolId("CSI ")==gMC->CurrentVolID(copy))             //is it in csi?      
355                       {
356                         gMC->TrackMomentum(momentum);
357                         mom[0]=momentum(0);
358                         mom[1]=momentum(1);
359                         mom[2]=momentum(2);
360                         mom[3]=momentum(3);
361                         
362                         /********* Photons lost by Fresnel reflection have to be calculated by hand********/ 
363                         /***********************Cerenkov phtons (always polarised)*************************/
364                         
365                         Float_t cophi = TMath::Cos(TMath::ATan2(mom[0], mom[1]));
366                         Float_t t = Fresnel(ckovEnergy*1e9,cophi,1);
367                         gMC->Rndm(ranf, 1);
368                         if (ranf[0] < t) {
369                           //geant3->StopTrack();
370                           ckovData[13] = 6;
371                           AddCerenkov(gAlice->CurrentTrack(),vol,ckovData);
372                           //printf("Lost by Fresnel\n");
373                         }
374                         /**********************************************************************************/
375                       }
376                   } //track entering?
377                   
378                   
379                   /********************Evaluation of losses************************/
380                   /******************still in the old fashion**********************/
381                   
382                   Int_t i1 = geant3->Gctrak()->nmec;            //number of physics mechanisms acting on the particle
383                   for (Int_t i = 0; i < i1; ++i) {
384                     //        Reflection loss 
385                     if (geant3->Gctrak()->lmec[i] == 106) {        //was it reflected
386                       ckovData[13]=10;
387                       if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy)) 
388                         ckovData[13]=1;
389                       if (gMC->CurrentVolID(copy) == gMC->VolId("QUAR")) 
390                         ckovData[13]=2;
391                       //geant3->StopTrack();
392                       AddCerenkov(gAlice->CurrentTrack(),vol,ckovData);
393                     } //reflection question
394                     
395                     
396                     //        Absorption loss 
397                     else if (geant3->Gctrak()->lmec[i] == 101) {              //was it absorbed?
398                       ckovData[13]=20;
399                       if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy)) 
400                         ckovData[13]=11;
401                       if (gMC->CurrentVolID(copy) == gMC->VolId("QUAR")) 
402                         ckovData[13]=12;
403                       if (gMC->CurrentVolID(copy) == gMC->VolId("META")) 
404                         ckovData[13]=13;
405                       if (gMC->CurrentVolID(copy) == gMC->VolId("GAP ")) 
406                         ckovData[13]=13;
407                       
408                       if (gMC->CurrentVolID(copy) == gMC->VolId("SRIC")) 
409                         ckovData[13]=15;
410                       
411                       //        CsI inefficiency 
412                       if (gMC->CurrentVolID(copy) == gMC->VolId("CSI ")) {
413                         ckovData[13]=16;
414                       }
415                       //geant3->StopTrack();
416                       AddCerenkov(gAlice->CurrentTrack(),vol,ckovData);
417                       //printf("Added cerenkov %d\n",fCkovNumber);
418                     } //absorption question 
419                     
420                     
421                     //        Photon goes out of tracking scope 
422                     else if (geant3->Gctrak()->lmec[i] == 30) {                 //is it below energy treshold?
423                       ckovData[13]=21;
424                       //geant3->StopTrack();
425                       AddCerenkov(gAlice->CurrentTrack(),vol,ckovData);
426                     }   // energy treshold question         
427                   }  //number of mechanisms cycle
428                   /**********************End of evaluation************************/
429                 } //freon production question
430             } //energy interval question
431         //}//inside the proximity gap question
432     } //cerenkov photon question
433       
434     /**************************************End of Production Parameters Storing*********************/ 
435     
436     
437     /*******************************Treat photons that hit the CsI (Ckovs and Feedbacks)************/ 
438     
439     if (gMC->TrackPid() == 50000050 || gMC->TrackPid() == 50000051) {
440       //printf("Cerenkov\n");
441         if (gMC->VolId("CSI ")==gMC->CurrentVolID(copy))
442         {
443             
444           if (gMC->Edep() > 0.){
445                 gMC->TrackPosition(position);
446                 gMC->TrackMomentum(momentum);
447                 pos[0]=position(0);
448                 pos[1]=position(1);
449                 pos[2]=position(2);
450                 mom[0]=momentum(0);
451                 mom[1]=momentum(1);
452                 mom[2]=momentum(2);
453                 mom[3]=momentum(3);
454                 Double_t tc = mom[0]*mom[0]+mom[1]*mom[1];
455                 Double_t rt = TMath::Sqrt(tc);
456                 theta   = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg;
457                 phi     = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg;
458                 gMC->Gmtod(pos,localPos,1);                                                                    
459                 gMC->Gmtod(mom,localMom,2);
460                 
461                 gMC->CurrentVolOffID(2,copy);
462                 vol[0]=copy;
463                 idvol=vol[0]-1;
464
465                 //Int_t sector=((AliRICHChamber*) (*fChambers)[idvol])
466                         //->Sector(localPos[0], localPos[2]);
467                 //printf("Sector:%d\n",sector);
468
469                 /*if (gMC->TrackPid() == 50000051){
470                   fFeedbacks++;
471                   printf("Feedbacks:%d\n",fFeedbacks);
472                 }*/     
473                 
474                 ((AliRICHChamber*) (*fChambers)[idvol])
475                     ->SigGenInit(localPos[0], localPos[2], localPos[1]);
476                 if(idvol<kNCH) {        
477                     ckovData[0] = gMC->TrackPid();        // particle type
478                     ckovData[1] = pos[0];                 // X-position for hit
479                     ckovData[2] = pos[1];                 // Y-position for hit
480                     ckovData[3] = pos[2];                 // Z-position for hit
481                     ckovData[4] = theta;                      // theta angle of incidence
482                     ckovData[5] = phi;                      // phi angle of incidence 
483                     ckovData[8] = (Float_t) fNPadHits;      // first padhit
484                     ckovData[9] = -1;                       // last pad hit
485                     ckovData[13] = 4;                       // photon was detected
486                     ckovData[14] = mom[0];
487                     ckovData[15] = mom[1];
488                     ckovData[16] = mom[2];
489                     
490                     destep = gMC->Edep();
491                     gMC->SetMaxStep(kBig);
492                     cherenkovLoss  += destep;
493                     ckovData[7]=cherenkovLoss;
494                     
495                     nPads = MakePadHits(localPos[0],localPos[2],cherenkovLoss,idvol,kCerenkov);
496                     if (fNPadHits > (Int_t)ckovData[8]) {
497                         ckovData[8]= ckovData[8]+1;
498                         ckovData[9]= (Float_t) fNPadHits;
499                     }
500
501                     ckovData[17] = nPads;
502                     //printf("nPads:%d",nPads);
503                     
504                     //TClonesArray *Hits = RICH->Hits();
505                     AliRICHHit *mipHit =  (AliRICHHit*) (fHits->UncheckedAt(0));
506                     if (mipHit)
507                       {
508                         mom[0] = current->Px();
509                         mom[1] = current->Py();
510                         mom[2] = current->Pz();
511                         Float_t mipPx = mipHit->fMomX;
512                         Float_t mipPy = mipHit->fMomY;
513                         Float_t mipPz = mipHit->fMomZ;
514                         
515                         Float_t r = mom[0]*mom[0] + mom[1]*mom[1] + mom[2]*mom[2];
516                         Float_t rt = TMath::Sqrt(r);
517                         Float_t mipR = mipPx*mipPx + mipPy*mipPy + mipPz*mipPz; 
518                         Float_t mipRt = TMath::Sqrt(mipR);
519                         if ((rt*mipRt) > 0)
520                           {
521                             coscerenkov = (mom[0]*mipPx + mom[1]*mipPy + mom[2]*mipPz)/(rt*mipRt);
522                           }
523                         else
524                           {
525                             coscerenkov = 0;
526                           }
527                         Float_t cherenkov = TMath::ACos(coscerenkov);
528                         ckovData[18]=cherenkov;
529                       }
530                     //if (sector != -1)
531                     //{
532                     AddHit(gAlice->CurrentTrack(),vol,ckovData);
533                     AddCerenkov(gAlice->CurrentTrack(),vol,ckovData);
534                     //}
535                 }
536             }
537         }
538     }
539     
540     /***********************************************End of photon hits*********************************************/
541     
542
543     /**********************************************Charged particles treatment*************************************/
544
545     else if (gMC->TrackCharge())
546     //else if (1 == 1)
547       {
548 //If MIP
549         /*if (gMC->IsTrackEntering())
550           {                
551             hits[13]=20;//is track entering?
552           }*/
553         if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
554           {
555             fFreonProd=1;
556           }
557
558         if (gMC->VolId("GAP ")== gMC->CurrentVolID(copy)) {
559 // Get current particle id (ipart), track position (pos)  and momentum (mom)
560             
561             gMC->CurrentVolOffID(3,copy);
562             vol[0]=copy;
563             idvol=vol[0]-1;
564
565             //Int_t sector=((AliRICHChamber*) (*fChambers)[idvol])
566                         //->Sector(localPos[0], localPos[2]);
567             //printf("Sector:%d\n",sector);
568             
569             gMC->TrackPosition(position);
570             gMC->TrackMomentum(momentum);
571             pos[0]=position(0);
572             pos[1]=position(1);
573             pos[2]=position(2);
574             mom[0]=momentum(0);
575             mom[1]=momentum(1);
576             mom[2]=momentum(2);
577             mom[3]=momentum(3);
578             gMC->Gmtod(pos,localPos,1);                                                                    
579             gMC->Gmtod(mom,localMom,2);
580             
581             ipart  = gMC->TrackPid();
582             //
583             // momentum loss and steplength in last step
584             destep = gMC->Edep();
585             step   = gMC->TrackStep();
586   
587             //
588             // record hits when track enters ...
589             if( gMC->IsTrackEntering()) {
590 //              gMC->SetMaxStep(fMaxStepGas);
591                 Double_t tc = mom[0]*mom[0]+mom[1]*mom[1];
592                 Double_t rt = TMath::Sqrt(tc);
593                 theta   = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg;
594                 phi     = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg;
595                 
596
597                 Double_t localTc = localMom[0]*localMom[0]+localMom[2]*localMom[2];
598                 Double_t localRt = TMath::Sqrt(localTc);
599                 localTheta   = Float_t(TMath::ATan2(localRt,Double_t(localMom[1])))*kRaddeg;                       
600                 localPhi     = Float_t(TMath::ATan2(Double_t(localMom[2]),Double_t(localMom[0])))*kRaddeg;    
601                 
602                 hits[0] = Float_t(ipart);         // particle type
603                 hits[1] = localPos[0];                 // X-position for hit
604                 hits[2] = localPos[1];                 // Y-position for hit
605                 hits[3] = localPos[2];                 // Z-position for hit
606                 hits[4] = localTheta;                  // theta angle of incidence
607                 hits[5] = localPhi;                    // phi angle of incidence 
608                 hits[8] = (Float_t) fNPadHits;    // first padhit
609                 hits[9] = -1;                     // last pad hit
610                 hits[13] = fFreonProd;           // did id hit the freon?
611                 hits[14] = mom[0];
612                 hits[15] = mom[1];
613                 hits[16] = mom[2];
614
615                 tlength = 0;
616                 eloss   = 0;
617                 fFreonProd = 0;
618         
619                 Chamber(idvol).LocaltoGlobal(localPos,hits+1);
620            
621                 
622                 //To make chamber coordinates x-y had to pass localPos[0], localPos[2]
623                 xhit    = localPos[0];
624                 yhit    = localPos[2];
625                 // Only if not trigger chamber
626                 if(idvol<kNCH) {
627                     //
628                     //  Initialize hit position (cursor) in the segmentation model 
629                     ((AliRICHChamber*) (*fChambers)[idvol])
630                         ->SigGenInit(localPos[0], localPos[2], localPos[1]);
631                 }
632             }
633             
634             // 
635             // Calculate the charge induced on a pad (disintegration) in case 
636             //
637             // Mip left chamber ...
638             if( gMC->IsTrackExiting() || gMC->IsTrackStop() || gMC->IsTrackDisappeared()){
639                 gMC->SetMaxStep(kBig);
640                 eloss   += destep;
641                 tlength += step;
642                 
643                                 
644                 // Only if not trigger chamber
645                 if(idvol<kNCH) {
646                   if (eloss > 0) 
647                     {
648                       if(gMC->TrackPid() == kNeutron)
649                         printf("\n\n\n\n\n Neutron Making Pad Hit!!! \n\n\n\n");
650                       nPads = MakePadHits(xhit,yhit,eloss,idvol,kMip);
651                       hits[17] = nPads;
652                       //printf("nPads:%d",nPads);
653                     }
654                 }
655                 
656                 hits[6]=tlength;
657                 hits[7]=eloss;
658                 if (fNPadHits > (Int_t)hits[8]) {
659                     hits[8]= hits[8]+1;
660                     hits[9]= (Float_t) fNPadHits;
661                 }
662                 
663                 //if(sector !=-1)
664                 new(lhits[fNhits++]) AliRICHHit(fIshunt,gAlice->CurrentTrack(),vol,hits);
665                 eloss = 0; 
666                 //
667                 // Check additional signal generation conditions 
668                 // defined by the segmentation
669                 // model (boundary crossing conditions) 
670             } else if 
671                 (((AliRICHChamber*) (*fChambers)[idvol])
672                  ->SigGenCond(localPos[0], localPos[2], localPos[1]))
673             {
674                 ((AliRICHChamber*) (*fChambers)[idvol])
675                     ->SigGenInit(localPos[0], localPos[2], localPos[1]);
676                 if (eloss > 0) 
677                   {
678                     if(gMC->TrackPid() == kNeutron)
679                       printf("\n\n\n\n\n Neutron Making Pad Hit!!! \n\n\n\n");
680                     nPads = MakePadHits(xhit,yhit,eloss,idvol,kMip);
681                     hits[17] = nPads;
682                     //printf("Npads:%d",NPads);
683                   }
684                 xhit     = localPos[0];
685                 yhit     = localPos[2]; 
686                 eloss    = destep;
687                 tlength += step ;
688                 //
689                 // nothing special  happened, add up energy loss
690             } else {        
691                 eloss   += destep;
692                 tlength += step ;
693             }
694         }
695       }
696     /*************************************************End of MIP treatment**************************************/
697    //}
698 }