]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - EVGEN/AliGenFLUKAsource.cxx
Wrong comments removed from file
[u/mrichter/AliRoot.git] / EVGEN / AliGenFLUKAsource.cxx
index 5a5ade01edca76f3e6f478d37e3e0e3bc8fb9008..cd12e166a6f0f5a27f71a61f971a357b28b6546c 100644 (file)
@@ -1,7 +1,33 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/*
+$Log$
+Revision 1.7  1999/11/03 17:43:20  fca
+New version from G.Martinez & A.Morsch
+
+Revision 1.6  1999/09/29 09:24:12  fca
+Introduction of the Copyright and cvs Log
+
+*/
+
 #include "AliGenFLUKAsource.h"
 #include "AliGenMUONlib.h"
 #include "AliMC.h"
 #include "AliRun.h"
+#include "AliPDG.h"
 #include <TDirectory.h>
 #include <TFile.h>
 #include <TTree.h>
@@ -29,6 +55,7 @@
     fFileName="flukasource.root";
 
     fTreeFluka=0;
+    fTreeChain = new TChain("h1");
 //
 //  Read all particles
     fNpart=-1;
@@ -57,37 +84,28 @@ AliGenFLUKAsource::AliGenFLUKAsource(Int_t npart)
     fFileName="flukasource.root";
 
     fTreeFluka=0;
-
+    fTreeChain = new TChain("h1"); 
 }
 
 //____________________________________________________________
 AliGenFLUKAsource::~AliGenFLUKAsource()
 {
-    delete fTreeFluka;
+ if (fTreeFluka) delete fTreeFluka;
+ if (fTreeChain) delete fTreeChain;
+ // if (fFileName)  delete fFileName;
 }
 
-//____________________________________________________________
-void AliGenFLUKAsource::FlukaInit() 
+void AliGenFLUKAsource::AddFile(const Text_t *filename)
 {
-//
-// reset the existing file environment and open a new root file if
-// the pointer to the Fluka tree is null
+    fTreeChain->Add(filename);
     
-    TFile *File=0;
-    if (fTreeFluka==0) {
-        if (!File) {
-           File = new TFile(fFileName);
-           File->cd();
-           cout<<"I have opened "<<fFileName<<" file "<<endl;
-        }
-// get the tree address in the Fluka boundary source file
-       fTreeFluka = (TTree*)gDirectory->Get("h1");
-    } else {
-        File = fTreeFluka->GetCurrentFile();
-        File->cd();
-    }
+}
 
-    TTree *h2=fTreeFluka;
+
+//____________________________________________________________
+void AliGenFLUKAsource::FlukaInit() 
+{
+    TChain *h2=fTreeChain;
 //Set branch addresses
     h2->SetBranchAddress("Ip",&Ip);
     h2->SetBranchAddress("Ipp",&Ipp);
@@ -110,107 +128,122 @@ void AliGenFLUKAsource::FlukaInit()
 
 //____________________________________________________________
 void AliGenFLUKAsource::Generate()
-{
-
-  AliMC* pMC = AliMC::GetMC();
+{ 
+    AliMC* gMC = AliMC::GetMC();
 
-  const Int_t ifluge[28]={14, 15, 3, 2, 4, 4, 1, 13, 25, 5, 6, 10, 8, 9,
-                         11, 12,18, 26, 16, 21, 19, 20, 7, 16,16,0,0,0};
-  Float_t polar[3]= {0,0,0};
+    const Int_t ifluge[28]={kProton, kProtonBar, kElectron, kPositron,
+                         kNuE, kNuEBar, kGamma, kNeutron, kNeutronBar,
+                         kMuonPlus, kMuonMinus, kK0Long , kPiPlus, kPiMinus,
+                         kKPlus, kKMinus, kLambda0, kLambda0Bar, kK0Short,
+                         kSigmaMinus, kSigmaPlus, kSigma0, kPi0, kK0, kK0Bar,
+                         0,kNuMu,kNuMuBar};
+    Float_t polar[3]= {0,0,0};
   //
-  Float_t origin[3];
-  Float_t p[3];
-  Float_t prwn;
-  Float_t wgt, fwgt;
-  Float_t phi;
-  char name[100];
-  Float_t amass, charge, tlife;
-  Int_t itrtyp;
-  Int_t iwgt;
-  Int_t i, j, part, nt;
-  static Int_t irwn;
-  //
-  Float_t random[2];
+    Float_t origin[3];
+    Float_t p[3];
+    Float_t prwn;
+    Float_t wgt, fwgt;
+    Float_t phi;
+    char name[100];
+    Float_t amass, charge, tlife;
+    Int_t itrtyp;
+    Int_t iwgt;
+    Int_t i, j, part, nt;
+    static Int_t irwn;
+    //
+    Float_t random[2];
   //
-  FlukaInit();
-  TTree *h2=fTreeFluka;
-  Int_t nentries = (Int_t) h2->GetEntries();
-  if (fNpart == -1) fNpart=nentries;
+    FlukaInit();
+    TChain *h2=fTreeChain;
+    Int_t nentries = (Int_t) h2->GetEntries();
+    if (fNpart == -1) fNpart=Int_t(nentries*fFrac);
   
   // loop over number of particles
-  Int_t nb=0;
-  for (i=0; i<fNpart;i++) {
-    Int_t ev=pMC->CurrentEvent();
-    Int_t entry=fNpart*(ev-1)+i; 
-    nb = (Int_t)h2->GetEvent(entry); 
-    if (irwn > nentries) {
-      printf("No more entries in the FLUKA boundary source file\n");
-      return;
-    }   
-    if (Ip > 28 || Ip < 0) {
-      irwn++;
-      continue;
-    }
-    if ((Ip != fIkine && fIkine != 6 && fIkine != 9) || Age > fAgeMax){
-       irwn++;
-       continue;
-    } else if (fIkine == 9) {
-       if (Ip == 7 || Ip == 8 || Age > fAgeMax) { 
-          irwn++;
-          continue;
-       }
-    }
+    Int_t nb=0;
+    Int_t ev=gMC->CurrentEvent();
+    for (i=0; i<fNpart;i++) {
+       Int_t entry=fNpart*(ev-1)+i; 
+       nb = (Int_t)h2->GetEvent(entry); 
+       if (irwn > nentries) {
+           printf("No more entries in the FLUKA boundary source file\n");
+           TFile *File=0;
+           // Get AliRun object or create it 
+           if (!gAlice) {
+               gAlice = (AliRun*)File->Get("gAlice");
+               if (gAlice) printf("AliRun object found on file\n");
+               if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
+           }
+           TTree *fAli=gAlice->TreeK();
+           if (fAli) File =fAli->GetCurrentFile();
+           File->cd();
+           printf("Generate - I'm out \n");
+           return;
+       }   
+       if (Ip > 28 || Ip < 0) {
+           irwn++;
+           continue;
+       }
+       
+       if ((Ip != fIkine && fIkine != 6 && fIkine != 9 && fIkine != 10) || Age > fAgeMax){
+           irwn++;
+           continue;
+       } else if (fIkine == 9) {
+           if (Ip == 7 || Ip == 8 || Age > fAgeMax) { 
+               irwn++;
+               continue;
+           }
+       } else if (fIkine == 10) {
+           if (Ip == 8 || Age > fAgeMax) { 
+               irwn++;
+               continue;
+           }
+       }
     
 
-    irwn++;
-    printf("\n Particle type: %f \n \n ", Ip);
-    
-    if (Ip ==7){
-       prwn=Ekin;
-       part=1;
-    } else if (Ip == 8) {
-       prwn=sqrt(Ekin*Ekin + 2.*0.93956563);
-       part=13;
-    } else {
-       part=ifluge[int(Ip)-1];
-       pMC->Gfpart(part, name, itrtyp,  
-                  amass, charge, tlife); 
-       prwn=sqrt(Ekin*Ekin + 2.*amass);
-    }
-    origin[0]=Xi;
-    origin[1]=Yi;
-    origin[2]=Zi;
-
-    p[0]=Px*prwn;
-    p[1]=Py*prwn;
-    p[2]=Pz*prwn;
-
-    //handle particle weight correctly
-    wgt = (part == 13) ? Wgt*fAddWeight : Wgt;
-    iwgt=Int_t(wgt);
-    fwgt=wgt-Float_t(iwgt);
-    pMC->Rndm(random,2);
-    if (random[0] < fwgt) iwgt++;
-    if (part==1 && iwgt>100) iwgt=100;
-    Int_t nstack=0;
-    for (j=0; j<iwgt; j++) {
-       gAlice->SetTrack(1,-1,part,p,origin,polar,0,"Primary",nt);
-       pMC->Rndm(random,2);
-       phi=2*random[1]*TMath::Pi();
-       Float_t pn1=p[0]*TMath::Sin(phi) - p[1]*TMath::Cos(phi);
-       Float_t pn2=p[0]*TMath::Cos(phi) + p[1]*TMath::Sin(phi);
-       p[0]=pn1;
-       p[1]=pn2;
-       Float_t on1=origin[0]*TMath::Sin(phi)-origin[1]*TMath::Cos(phi);
-       Float_t on2=origin[0]*TMath::Cos(phi)+origin[1]*TMath::Sin(phi);
-       origin[0]=on1;
-       origin[1]=on2;
-       nstack++;
+       irwn++;
+//
+// PDG code from FLUKA particle type (Ip)
+       part=ifluge[int(Ip)-1]; 
+//
+// Calculate momentum from kinetic energy and mass of the particle
+       gMC->Gfpart(part, name, itrtyp,  
+                   amass, charge, tlife); 
+       prwn=Ekin*sqrt(1. + 2.*amass/Ekin);
+
+
+       origin[0]=Yi;
+       origin[1]=Xi;
+       origin[2]=Zi;
+       
+       p[0]=Py*prwn;
+       p[1]=Px*prwn;
+       p[2]=Pz*prwn;
+
+       //handle particle weight correctly
+       wgt = (part == 13) ? Wgt*fAddWeight : Wgt;
+       iwgt=Int_t(wgt);
+       fwgt=wgt-Float_t(iwgt);
+       gMC->Rndm(random,2);
+       if (random[0] < fwgt) iwgt++;
+       if (part==1 && iwgt>100) iwgt=100;
+       Int_t nstack=0;
+       for (j=0; j<iwgt; j++) {
+           gAlice->SetTrack(fTrackIt,-1,part,p,origin,polar,Age,"Primary",nt);
+           gMC->Rndm(random,2);
+           phi=2*random[1]*TMath::Pi();
+           Float_t pn1=p[0]*TMath::Sin(phi) - p[1]*TMath::Cos(phi);
+           Float_t pn2=p[0]*TMath::Cos(phi) + p[1]*TMath::Sin(phi);
+           p[0]=pn1;
+           p[1]=pn2;
+           Float_t on1=origin[0]*TMath::Sin(phi)-origin[1]*TMath::Cos(phi);
+           Float_t on2=origin[0]*TMath::Cos(phi)+origin[1]*TMath::Sin(phi);
+           origin[0]=on1;
+           origin[1]=on2;
+           nstack++;
+       }
+       if (nstack == 0) continue;
     }
-    if (nstack == 0) continue;
-  }
+    
     TFile *File=0;
 // Get AliRun object or create it 
     if (!gAlice) {