EVGEN/AliGenHalo.cxx

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