1 // **************************************************************************
2 // * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 // * Author: The ALICE Off-line Project. *
5 // * Contributors are mentioned in the code where appropriate. *
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 // **************************************************************************
17 #include "AliRICHv1.h"
18 #include "AliRICHParam.h"
19 #include "AliRICHChamber.h"
20 #include <TParticle.h>
22 #include <TVirtualMC.h>
31 //__________________________________________________________________________________________________
32 void AliRICHv1::StepManager()
35 // 3- Ckovs absorbed in Collection electrods
36 // 5- Ckovs absorbed in Cathode wires
37 // 6- Ckovs absorbed in Anod wires
40 static Int_t iCurrentChamber;
41 //history of Cerenkovs
42 if(gMC->TrackPid()==kCerenkov){
43 if(gMC->IsNewTrack() && (gMC->CurrentVolID(copy)==gMC->VolId("FRE1")||gMC->CurrentVolID(copy)==gMC->VolId("FRE2"))) fCounters(0)++;// 0- Ckovs produced in Freon
44 if(!gMC->IsTrackAlive() && (gMC->CurrentVolID(copy)==gMC->VolId("FRE1")||gMC->CurrentVolID(copy)==gMC->VolId("FRE2"))) fCounters(1)++;// 1- Ckovs absorbed in Freon
45 if(!gMC->IsTrackAlive() && gMC->CurrentVolID(copy)==gMC->VolId("QUAR")) fCounters(2)++;//2- Ckovs absorbed in Quartz
46 if(!gMC->IsTrackAlive() && gMC->CurrentVolID(copy)==gMC->VolId("META")) fCounters(4)++;//4- Ckovs absorbed in CH4
50 static TLorentzVector cerX4;
51 if((gMC->TrackPid()==kCerenkov||gMC->TrackPid()==kFeedback)&&gMC->CurrentVolID(copy)==gMC->VolId("CSI ")){//photon in CSI
52 if(gMC->Edep()>0){//CF+CSI+DE
53 if(IsLostByFresnel()){
54 if(gMC->TrackPid()==kCerenkov) fCounters(7)++;// 7- Ckovs reflected on CsI
58 gMC->TrackPosition(cerX4); gMC->CurrentVolOffID(2,iCurrentChamber);
60 AddHit(iCurrentChamber,gAlice->GetMCApp()->GetCurrentTrackNumber(),cerX4.Vect(),cerX4.Vect());//HIT for PHOTON in conditions CF+CSI+DE
61 fCounters(8)++;//4- Ckovs converted in CsI
62 GenerateFeedbacks(iCurrentChamber);
66 //Treat charged particles
68 static TLorentzVector mipInX4,mipOutX4;
69 if(gMC->TrackCharge() && gMC->CurrentVolID(copy)==gMC->VolId("GAP ")){//MIP in GAP
70 gMC->CurrentVolOffID(3,iCurrentChamber);
71 if(gMC->IsTrackEntering()||gMC->IsNewTrack()) {//MIP in GAP entering or newly created
73 gMC->TrackPosition(mipInX4);
74 }else if(gMC->IsTrackExiting()||gMC->IsTrackStop()||gMC->IsTrackDisappeared()){//MIP in GAP exiting or disappeared
75 eloss+=gMC->Edep();//take into account last step dEdX
76 gMC->TrackPosition(mipOutX4);
77 AddHit(iCurrentChamber,gAlice->GetMCApp()->GetCurrentTrackNumber(),mipInX4.Vect(),mipOutX4.Vect(),eloss);//HIT for MIP: MIP in GAP Exiting
78 GenerateFeedbacks(iCurrentChamber,eloss);//MIP+GAP+Exit
79 }else//MIP in GAP going inside
83 //__________________________________________________________________________________________________
84 Bool_t AliRICHv1::IsLostByFresnel()
86 // Calculate probability for the photon to be lost by Fresnel reflection.
88 Double_t mom[3],localMom[3];
89 gMC->TrackMomentum(p4); mom[0]=p4(1); mom[1]=p4(2); mom[2]=p4(3);
90 localMom[0]=0; localMom[1]=0; localMom[2]=0;
91 gMC->Gmtod(mom,localMom,2);
92 Double_t localTc = localMom[0]*localMom[0]+localMom[2]*localMom[2];
93 Double_t localTheta = TMath::ATan2(TMath::Sqrt(localTc),localMom[1]);
94 Double_t cotheta = TMath::Abs(TMath::Cos(localTheta));
95 if(gMC->GetRandom()->Rndm() < Fresnel(p4.E()*1e9,cotheta,1)){
96 if(GetDebug()) Info("IsLostByFresnel","Photon lost");
101 //__________________________________________________________________________________________________
102 void AliRICHv1::GenerateFeedbacks(Int_t iChamber,Float_t eloss)
104 // Generate FeedBack photons for the current particle. To be invoked from StepManager().
105 // eloss=0 means photon so only pulse height distribution is to be analysed. This one is done in AliRICHParam::TotQdc()
108 gMC->TrackPosition(x4);
109 TVector2 x2=C(iChamber)->Glob2Loc(x4);
110 Int_t sector=P()->Loc2Sec(x2); if(sector==kBad) return; //hit in dead zone, nothing to produce
111 Int_t iTotQdc=P()->TotQdc(x2,eloss);
112 Int_t iNphotons=gMC->GetRandom()->Poisson(P()->AlphaFeedback(sector)*iTotQdc);
113 if(GetDebug())Info("GenerateFeedbacks","N photons=%i",iNphotons);
115 Float_t cthf, phif, enfp = 0, sthf, e1[3], e2[3], e3[3], vmod, uswop,dir[3], phi,pol[3], mom[4];
117 for(Int_t i=0;i<iNphotons;i++){//feedbacks loop
119 gMC->GetRandom()->RndmArray(2,ranf); //Sample direction
122 sthf = TMath::Sqrt((1 - cthf) * (1 + cthf));
123 phif = ranf[1] * 2 * TMath::Pi();
125 if(Double_t randomNumber=gMC->GetRandom()->Rndm()<=0.57)
127 else if(randomNumber<=0.7)
133 dir[0] = sthf * TMath::Sin(phif); dir[1] = cthf; dir[2] = sthf * TMath::Cos(phif);
134 gMC->Gdtom(dir, mom, 2);
135 mom[0]*=enfp; mom[1]*=enfp; mom[2]*=enfp;
136 mom[3] = TMath::Sqrt(mom[0]*mom[0]+mom[1]*mom[1]+mom[2]*mom[2]);
139 e1[0]= 0; e1[1]=-dir[2]; e1[2]= dir[1];
140 e2[0]=-dir[1]; e2[1]= dir[0]; e2[2]= 0;
141 e3[0]= dir[1]; e3[1]= 0; e3[2]=-dir[0];
144 for(j=0;j<3;j++) vmod+=e1[j]*e1[j];
145 if (!vmod) for(j=0;j<3;j++) {
151 for(j=0;j<3;j++) vmod+=e2[j]*e2[j];
152 if (!vmod) for(j=0;j<3;j++) {
158 vmod=0; for(j=0;j<3;j++) vmod+=e1[j]*e1[j]; vmod=TMath::Sqrt(1/vmod); for(j=0;j<3;j++) e1[j]*=vmod;
159 vmod=0; for(j=0;j<3;j++) vmod+=e2[j]*e2[j]; vmod=TMath::Sqrt(1/vmod); for(j=0;j<3;j++) e2[j]*=vmod;
161 phi = gMC->GetRandom()->Rndm()* 2 * TMath::Pi();
162 for(j=0;j<3;j++) pol[j]=e1[j]*TMath::Sin(phi)+e2[j]*TMath::Cos(phi);
163 gMC->Gdtom(pol, pol, 2);
164 Int_t outputNtracksStored;
165 gAlice->GetMCApp()->PushTrack(1, //do not transport
166 gAlice->GetMCApp()->GetCurrentTrackNumber(),//parent track
168 mom[0],mom[1],mom[2],mom[3], //track momentum
169 x4.X(),x4.Y(),x4.Z(),x4.T(), //track origin
170 pol[0],pol[1],pol[2], //polarization
175 if(GetDebug()) Info("GenerateFeedbacks","Stop.");
176 }//GenerateFeedbacks()