EVGEN/AliGenHalo.cxx

   1 #include "AliGenHalo.h"
   2 #include "AliGenMUONlib.h"
   3 #include "AliMC.h"
   4 #include "AliRun.h"
   5 #include "AliPDG.h"
   6
   7 #include <TDirectory.h>
   8 #include <TDatabasePDG.h>
   9 #include <TFile.h>
  10 #include <TTree.h>
  11 #include <stdlib.h>
  12  ClassImp(AliGenHalo)
  13      AliGenHalo::AliGenHalo()
  14          :AliGenerator(-1)
  15 {
  16     fName="Halo";
  17     fTitle="Halo from LHC Tunnel";
  18     // Set the default file
  19     fFileName="~/marsip/marsip5.mu";
  20 //
  21 //  Read all particles
  22     fNpart=-1;
  23     fp=0;
  24 }
  25
  26 AliGenHalo::AliGenHalo(Int_t npart)
  27     :AliGenerator(npart)
  28 {
  29     fName="Halo";
  30     fTitle="Halo from LHC Tunnel";
  31     // Set the default file
  32     fFileName="~/marsip/marsip5.mu";
  33 //
  34 //  Read all particles
  35     fNpart=-1;
  36     fp=0;
  37 }
  38
  39 //____________________________________________________________
  40 AliGenHalo::~AliGenHalo()
  41 {
  42 }
  43
  44 //____________________________________________________________
  45 void AliGenHalo::Init()
  46 {}
  47
  48 //____________________________________________________________
  49 void AliGenHalo::Generate()
  50 {
  51     FILE *fp = fopen(fFileName,"r");
  52     if (fp) {
  53         printf("\n File %s opened for reading ! \n ", fFileName);
  54     } else {
  55         printf("\n Opening of file %s failed ! \n ", fFileName);
  56     }
  57 //
  58 // MARS particle codes
  59     // const Int_t imars[12]={0,14, 13, 8, 9, 11, 12, 5, 6, 1, 3, 2};
  60   const Int_t imars[12]={0,kProton,kNeutron,kPiPlus,kPiMinus,kKPlus,kKMinus,
  61                          kMuonPlus,kMuonMinus,kGamma,kElectron,kPositron};
  62
  63   Float_t polar[3]= {0,0,0};
  64   Float_t origin[3];
  65   Float_t p[3], p0;
  66   Float_t ekin, wgt, tx, ty, tz, txy;
  67   Float_t amass;
  68   //
  69   Int_t ipart, ncols, nt;
  70
  71   Int_t nread=0;
  72   origin[2]=2650;
  73
  74   while(1) {
  75       ncols = fscanf(fp,"%i %f %f %f %f %f %f",
  76                      &ipart, &ekin, &wgt,
  77                      &origin[0], &origin[1],
  78                      &tx, &ty);
  79       if (ncols < 0) break;
  80       nread++;
  81       if (fNpart !=-1 && nread > fNpart) break;
  82       ipart = imars[ipart];
  83       amass = TDatabasePDG::Instance()->GetParticle(ipart)->Mass();
  84       p0=sqrt(ekin*ekin + 2.*amass);
  85
  86       txy=TMath::Sqrt(tx*tx+ty*ty);
  87       if (txy == 1.) {
  88           tz=0;
  89       } else {
  90           tz=-TMath::Sqrt(1.-txy);
  91       }
  92       p[0]=p0*tx;
  93       p[1]=p0*ty;
  94       p[2]=p0*tz;
  95       fParentWeight=wgt;
  96       gAlice->SetTrack(1,-1,ipart,p,origin,polar,0,"Halo+",nt,fParentWeight);
  97       origin[2]=-origin[2];
  98       p[2]=-p[2];
  99       gAlice->SetTrack(1,-1,ipart,p,origin,polar,0,"Halo-",nt,fParentWeight);
 100       origin[2]=-origin[2];
 101       p[2]=-p[2];
 102   }
 103 }
 104
 105
 106
 107
 108
 109
 110
 111
 112