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