RICH/AliRICHResponseV0.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.3  2000/10/03 21:44:09  morsch
  19   Use AliSegmentation and AliHit abstract base classes.
  20
  21   Revision 1.2  2000/10/02 21:28:12  fca
  22   Removal of useless dependecies via forward declarations
  23
  24   Revision 1.1  2000/06/12 15:29:37  jbarbosa
  25   Cleaned up version.
  26
  27 */
  28
  29 #include "AliRICHResponseV0.h"
  30 #include "AliSegmentation.h"
  31 #include "AliRun.h"
  32 #include "AliMC.h"
  33
  34 #include <TMath.h>
  35 #include <TRandom.h>
  36 #include <TParticle.h>
  37 //___________________________________________
  38 ClassImp(AliRICHResponseV0)
  39
  40 Float_t AliRICHResponseV0::IntPH(Float_t eloss)
  41 {
  42     // Get number of electrons and return charge
  43
  44     Int_t nel;
  45     nel= Int_t(eloss/fEIonisation);
  46
  47     Float_t charge=0;
  48     if (nel == 0) nel=1;
  49     for (Int_t i=1;i<=nel;i++) {
  50         charge -= fChargeSlope*TMath::Log(gRandom->Rndm());
  51     }
  52     return charge;
  53 }
  54
  55 Float_t AliRICHResponseV0::IntPH()
  56 {
  57
  58 //  Get number of electrons and return charge, for a single photon
  59
  60     Float_t charge = -fChargeSlope*TMath::Log(gRandom->Rndm());
  61     return charge;
  62 }
  63
  64
  65
  66 // -------------------------------------------
  67 Float_t AliRICHResponseV0::IntXY(AliSegmentation * segmentation)
  68 {
  69
  70     const Float_t kInversePitch = 1/fPitch;
  71     Float_t response;
  72 //
  73 //  Integration limits defined by segmentation model
  74 //
  75
  76     Float_t xi1, xi2, yi1, yi2;
  77     segmentation->IntegrationLimits(xi1,xi2,yi1,yi2);
  78
  79     xi1=xi1*kInversePitch;
  80     xi2=xi2*kInversePitch;
  81     yi1=yi1*kInversePitch;
  82     yi2=yi2*kInversePitch;
  83
  84     //printf("Integration Limits: %f-%f, %f-%f\n",xi1,xi2,yi1,yi2);
  85
  86     //printf("KInversePitch:%f\n",kInversePitch);
  87
  88     //
  89 // The Mathieson function
  90     Double_t ux1=fSqrtKx3*TMath::TanH(fKx2*xi1);
  91     Double_t ux2=fSqrtKx3*TMath::TanH(fKx2*xi2);
  92
  93     Double_t uy1=fSqrtKy3*TMath::TanH(fKy2*yi1);
  94     Double_t uy2=fSqrtKy3*TMath::TanH(fKy2*yi2);
  95
  96     //printf("Integration Data: %f-%f, %f-%f\n",ux1,ux2,uy1,uy2);
  97
  98     //printf("%f %f %f %f\n",fSqrtKx3,fKx2,fKy4,fKx4);
  99
 100     response=4.*fKx4*(TMath::ATan(ux2)-TMath::ATan(ux1))*fKy4*(TMath::ATan(uy2)-TMath::ATan(uy1));
 101
 102     //printf("Response:%f\n",response);
 103
 104     return response;
 105
 106 }
 107
 108 Int_t AliRICHResponseV0::FeedBackPhotons(Float_t *source, Float_t qtot)
 109 {
 110   //
 111   // Generate FeedBack photons
 112   //
 113   Int_t j, ipart, nt;
 114
 115   Int_t sNfeed=0;
 116
 117
 118   // Local variables
 119   Float_t cthf, ranf[2], phif, enfp = 0, sthf;
 120   Int_t i, ifeed;
 121   Float_t e1[3], e2[3], e3[3];
 122   Float_t vmod, uswop;
 123   Float_t fp, random;
 124   Float_t dir[3], phi;
 125   Int_t nfp;
 126   Float_t pol[3], mom[3];
 127   TLorentzVector position;
 128   //
 129   // Determine number of feedback photons
 130
 131   //  Get weight of current particle
 132   TParticle *current = (TParticle*)
 133     (*gAlice->Particles())[gAlice->CurrentTrack()];
 134
 135   ifeed = Int_t(current->GetWeight()/100+0.5);
 136   ipart = gMC->TrackPid();
 137   fp = fAlphaFeedback * qtot;
 138   nfp = gRandom->Poisson(fp);
 139
 140   // This call to fill the time of flight
 141   gMC->TrackPosition(position);
 142   //
 143   // Generate photons
 144   for (i = 0; i <nfp; i++) {
 145
 146     // Direction
 147     gMC->Rndm(ranf, 2);
 148     cthf = ranf[0] * 2 - 1.;
 149     if (cthf < 0)  continue;
 150     sthf = TMath::Sqrt((1 - cthf) * (1 + cthf));
 151     phif = ranf[1] * 2 * TMath::Pi();
 152     //
 153     gMC->Rndm(&random, 1);
 154     if (random <= .57) {
 155       enfp = 7.5e-9;
 156     } else if (random <= .7) {
 157       enfp = 6.4e-9;
 158     } else {
 159       enfp = 7.9e-9;
 160     }
 161
 162     dir[0] = sthf * TMath::Sin(phif);
 163     dir[1] = cthf;
 164     dir[2] = sthf * TMath::Cos(phif);
 165     gMC->Gdtom(dir, mom, 2);
 166     mom[0]*=enfp;
 167     mom[1]*=enfp;
 168     mom[2]*=enfp;
 169
 170     // Polarisation
 171     e1[0] = 0;
 172     e1[1] = -dir[2];
 173     e1[2] = dir[1];
 174
 175     e2[0] = -dir[1];
 176     e2[1] = dir[0];
 177     e2[2] = 0;
 178
 179     e3[0] = dir[1];
 180     e3[1] = 0;
 181     e3[2] = -dir[0];
 182
 183     vmod=0;
 184     for(j=0;j<3;j++) vmod+=e1[j]*e1[j];
 185     if (!vmod) for(j=0;j<3;j++) {
 186       uswop=e1[j];
 187       e1[j]=e3[j];
 188       e3[j]=uswop;
 189     }
 190     vmod=0;
 191     for(j=0;j<3;j++) vmod+=e2[j]*e2[j];
 192     if (!vmod) for(j=0;j<3;j++) {
 193       uswop=e2[j];
 194       e2[j]=e3[j];
 195       e3[j]=uswop;
 196     }
 197
 198     vmod=0;
 199     for(j=0;j<3;j++) vmod+=e1[j]*e1[j];
 200     vmod=TMath::Sqrt(1/vmod);
 201     for(j=0;j<3;j++) e1[j]*=vmod;
 202
 203     vmod=0;
 204     for(j=0;j<3;j++) vmod+=e2[j]*e2[j];
 205     vmod=TMath::Sqrt(1/vmod);
 206     for(j=0;j<3;j++) e2[j]*=vmod;
 207
 208     gMC->Rndm(ranf, 1);
 209     phi = ranf[0] * 2 * TMath::Pi();
 210     for(j=0;j<3;j++) pol[j]=e1[j]*TMath::Sin(phi)+e2[j]*TMath::Cos(phi);
 211     gMC->Gdtom(pol, pol, 2);
 212
 213     // Put photon on the stack and label it as feedback (51, 52)
 214     ++sNfeed;
 215
 216     gAlice->SetTrack(Int_t(1), gAlice->CurrentTrack(), Int_t(50000051),
 217                      mom,source,pol,position[3],
 218                      kPFeedBackPhoton, nt, 1.);
 219   }
 220   return(sNfeed);
 221 }
 222
 223
 224
 225