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 }