* provided "as is" without express or implied warranty. *
**************************************************************************/
-/*
-$Log$
-Revision 1.6 1999/09/29 09:24:12 fca
-Introduction of the Copyright and cvs Log
+/* $Id$ */
-*/
+// Read background particles from a FLUKA boundary source file
+// This is a very special generator that works for background studies for the muon-spectrometer.
+// The input files come from FLUKA simulations.
+// Files can be chained.
+// Author: andreas.morsch@cern.ch
#include "AliGenFLUKAsource.h"
-#include "AliGenMUONlib.h"
-#include "AliMC.h"
-#include "AliRun.h"
-#include "AliPDG.h"
-#include <TDirectory.h>
+#include <stdlib.h>
+
+#include <TDatabasePDG.h>
+#include <TPDGCode.h>
+#include <RVersion.h>
+#include <TChain.h>
#include <TFile.h>
#include <TTree.h>
-#include <stdlib.h>
- ClassImp(AliGenFLUKAsource)
- AliGenFLUKAsource::AliGenFLUKAsource()
- :AliGenerator(-1)
+#include <TVirtualMC.h>
+
+#include "AliRun.h"
+
+ClassImp(AliGenFLUKAsource)
+
+AliGenFLUKAsource::AliGenFLUKAsource()
+ :AliGenerator(-1),
+ fIkine(6),
+ fAgeMax(1.e-5),
+ fAddWeight(1.),
+ fZshift(0.),
+ fFrac(0.),
+ fSourceId(-1),
+ fFileName(0),
+ fTreeChain(0),
+ fTreeFluka(0),
+ fIp(0.),
+ fIpp(0.),
+ fXi(0.),
+ fYi(0.),
+ fZi(0.),
+ fPx(0.),
+ fPy(0.),
+ fPz(0.),
+ fEkin(0.),
+ fZv(0.),
+ fRv(0.),
+ fItra(0.),
+ fIgas(0.),
+ fWgt(0.),
+ fEtag(0.),
+ fPtg(0.),
+ fAge(0.)
{
- //
+ // Constructor
fName="FLUKA";
fTitle="FLUKA Boundary Source";
// Read in all particle types by default
- fIkine=6;
// Set maximum admitted age of particles to 1.e-05 by default
- fAgeMax=1.e-05;
// Do not add weight
- fAddWeight=1.;
// Shift the z-coordinate of the impact point by 4.5 cm only if it reads
// from specific boundary source to the chamber (fZshift=4.5;),else there
// is no need to shift as it reads boundary source for the whole volume of
// the Muon Arm; the default value corresponds to boundary source for the
// whole volume of the MUON Arm
- fZshift=0;
// Set the default file
- fFileName="flukasource.root";
-
- fTreeFluka=0;
fTreeChain = new TChain("h1");
//
// Read all particles
fNpart=-1;
-
}
AliGenFLUKAsource::AliGenFLUKAsource(Int_t npart)
- :AliGenerator(npart)
+ :AliGenerator(npart),
+ fIkine(6),
+ fAgeMax(1.e-5),
+ fAddWeight(1.),
+ fZshift(0.),
+ fFrac(0.),
+ fSourceId(-1),
+ fFileName(""),
+ fTreeChain(new TChain("h1")),
+ fTreeFluka(0),
+ fIp(0.),
+ fIpp(0.),
+ fXi(0.),
+ fYi(0.),
+ fZi(0.),
+ fPx(0.),
+ fPy(0.),
+ fPz(0.),
+ fEkin(0.),
+ fZv(0.),
+ fRv(0.),
+ fItra(0.),
+ fIgas(0.),
+ fWgt(0.),
+ fEtag(0.),
+ fPtg(0.),
+ fAge(0.)
{
- //
- fName="FLUKA";
- fTitle="FLUKA Boundary Source";
- // Read in all particle types by default
- fIkine=6;
- // Set maximum admitted age of particles to 1.e-05 by default
- fAgeMax=1.e-05;
- // Do not add weight
- fAddWeight=1.;
- // Shift the z-coordinate of the impact point by 4.5 cm only if it reads
- // from specific boundary source to the chamber (fZshift=4.5;),else there
- // is no need to shift as it reads boundary source for the whole volume of
- // the Muon Arm; the default value corresponds to boundary source for the
- // whole volume of the MUON Arm
- fZshift=0;
- // Set the default file
- fFileName="flukasource.root";
+ // Constructor
+ fName = "FLUKA";
+ fTitle = "FLUKA Boundary Source";
+}
- fTreeFluka=0;
- fTreeChain = new TChain("h1");
+AliGenFLUKAsource::AliGenFLUKAsource(const AliGenFLUKAsource & FLUKAsource):
+ AliGenerator(FLUKAsource),
+ fIkine(6),
+ fAgeMax(1.e-5),
+ fAddWeight(1.),
+ fZshift(0.),
+ fFrac(0.),
+ fSourceId(-1),
+ fFileName(0),
+ fTreeChain(0),
+ fTreeFluka(0),
+ fIp(0.),
+ fIpp(0.),
+ fXi(0.),
+ fYi(0.),
+ fZi(0.),
+ fPx(0.),
+ fPy(0.),
+ fPz(0.),
+ fEkin(0.),
+ fZv(0.),
+ fRv(0.),
+ fItra(0.),
+ fIgas(0.),
+ fWgt(0.),
+ fEtag(0.),
+ fPtg(0.),
+ fAge(0.)
+{
+// Copy constructor
+ FLUKAsource.Copy(*this);
}
+
//____________________________________________________________
AliGenFLUKAsource::~AliGenFLUKAsource()
{
+// Destructor
if (fTreeFluka) delete fTreeFluka;
if (fTreeChain) delete fTreeChain;
- // if (fFileName) delete fFileName;
}
void AliGenFLUKAsource::AddFile(const Text_t *filename)
{
+// Add a file to the chain
fTreeChain->Add(filename);
}
//____________________________________________________________
void AliGenFLUKAsource::FlukaInit()
{
+// Set branch addresses of data entries
TChain *h2=fTreeChain;
//Set branch addresses
- h2->SetBranchAddress("Ip",&Ip);
- h2->SetBranchAddress("Ipp",&Ipp);
- h2->SetBranchAddress("Xi",&Xi);
- h2->SetBranchAddress("Yi",&Yi);
- h2->SetBranchAddress("Zi",&Zi);
- h2->SetBranchAddress("Px",&Px);
- h2->SetBranchAddress("Py",&Py);
- h2->SetBranchAddress("Pz",&Pz);
- h2->SetBranchAddress("Ekin",&Ekin);
- h2->SetBranchAddress("Zv",&Zv);
- h2->SetBranchAddress("Rv",&Rv);
- h2->SetBranchAddress("Itra",&Itra);
- h2->SetBranchAddress("Igas",&Igas);
- h2->SetBranchAddress("Wgt",&Wgt);
- h2->SetBranchAddress("Etag",&Etag);
- h2->SetBranchAddress("Ptg",&Ptg);
- h2->SetBranchAddress("Age",&Age);
+ h2->SetBranchAddress("Ip",&fIp);
+ h2->SetBranchAddress("Ipp",&fIpp);
+ h2->SetBranchAddress("Xi",&fXi);
+ h2->SetBranchAddress("Yi",&fYi);
+ h2->SetBranchAddress("Zi",&fZi);
+ h2->SetBranchAddress("Px",&fPx);
+ h2->SetBranchAddress("Py",&fPy);
+ h2->SetBranchAddress("Pz",&fPz);
+ h2->SetBranchAddress("Ekin",&fEkin);
+ h2->SetBranchAddress("Zv",&fZv);
+ h2->SetBranchAddress("Rv",&fRv);
+ h2->SetBranchAddress("Itra",&fItra);
+ h2->SetBranchAddress("Igas",&fIgas);
+ h2->SetBranchAddress("Wgt",&fWgt);
+ h2->SetBranchAddress("Etag",&fEtag);
+ h2->SetBranchAddress("Ptg",&fPtg);
+ h2->SetBranchAddress("Age",&fAge);
}
//____________________________________________________________
void AliGenFLUKAsource::Generate()
-{
- AliMC* gMC = AliMC::GetMC();
+{
+// Generate one event
- const Int_t ifluge[28]={kProton, kProtonBar, kElectron, kPositron,
+ const Int_t kIfluge[28]={kProton, kProtonBar, kElectron, kPositron,
kNuE, kNuEBar, kGamma, kNeutron, kNeutronBar,
kMuonPlus, kMuonMinus, kK0Long , kPiPlus, kPiMinus,
kKPlus, kKMinus, kLambda0, kLambda0Bar, kK0Short,
Float_t prwn;
Float_t wgt, fwgt;
Float_t phi;
- char name[100];
- Float_t amass, charge, tlife;
- Int_t itrtyp;
+ Float_t amass;
Int_t iwgt;
Int_t i, j, part, nt;
- static Int_t irwn;
+ static Int_t irwn=0;
//
Float_t random[2];
//
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;
Int_t ev=gMC->CurrentEvent();
nb = (Int_t)h2->GetEvent(entry);
if (irwn > nentries) {
printf("No more entries in the FLUKA boundary source file\n");
- TFile *File=0;
+ TFile *pFile=0;
// Get AliRun object or create it
if (!gAlice) {
- gAlice = (AliRun*)File->Get("gAlice");
+ gAlice = (AliRun*)pFile->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();
+ if (fAli) pFile =fAli->GetCurrentFile();
+ pFile->cd();
printf("Generate - I'm out \n");
return;
}
- if (Ip > 28 || Ip < 0) {
+
+ Int_t ifip = Int_t(fIp);
+
+
+ if (fSourceId != -1 && fIgas !=fSourceId) {
irwn++;
continue;
}
- if ((Ip != fIkine && fIkine != 6 && fIkine != 9 && fIkine != 10) || Age > fAgeMax){
+ if (ifip > 28 || ifip < 0) {
irwn++;
continue;
- } else if (fIkine == 9) {
- if (Ip == 7 || Ip == 8 || Age > fAgeMax) {
+ }
+
+ if ((ifip != fIkine && fIkine != kAll && fIkine != kCharged
+ && fIkine != 10) || fAge > fAgeMax){
+ irwn++;
+ continue;
+ } else if (fIkine == kCharged) {
+ if (ifip == 7 || ifip == 8 || fAge > fAgeMax) {
irwn++;
continue;
}
- } else if (fIkine == 10) {
- if (Ip == 8 || Age > fAgeMax) {
+ } else if (fIkine == kNoNeutron) {
+ if (ifip == 8 || fAge > fAgeMax) {
irwn++;
continue;
}
irwn++;
- printf("\n Particle type: %f \n \n ", Ip);
-
- if (Ip ==7){
- prwn=Ekin;
- part=1;
- } else if (Ip == 8) {
- prwn=Ekin*sqrt(1. + 2.*0.93956563/Ekin);
- part=13;
- } else {
- part=ifluge[int(Ip)-1];
- gMC->Gfpart(part, name, itrtyp,
- amass, charge, tlife);
- prwn=Ekin*sqrt(1. + 2.*amass/Ekin);
- }
- origin[0]=Yi;
- origin[1]=Xi;
- origin[2]=Zi;
+//
+// PDG code from FLUKA particle type (ifip)
+ part=kIfluge[int(ifip)-1];
+//
+// Calculate momentum from kinetic energy and mass of the particle
+#if ROOT_VERSION_CODE > 197895
+ amass = gMC->ParticleMass(part);
+#else
+ amass = (TDatabasePDG::Instance())->GetParticle(part)->Mass();
+#endif
+ prwn=fEkin*sqrt(1. + 2.*amass/fEkin);
+
+
+ origin[0]=fYi;
+ origin[1]=fXi;
+ origin[2]=fZi;
- p[0]=Py*prwn;
- p[1]=Px*prwn;
- p[2]=Pz*prwn;
+ p[0]=fPy*prwn;
+ p[1]=fPx*prwn;
+ p[2]=fPz*prwn;
//handle particle weight correctly
- wgt = (part == 13) ? Wgt*fAddWeight : Wgt;
+ wgt = (part == 13) ? fWgt*fAddWeight : fWgt;
iwgt=Int_t(wgt);
fwgt=wgt-Float_t(iwgt);
- gMC->Rndm(random,2);
+ 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);
+ PushTrack(fTrackIt,-1,part,p,origin,polar,fAge,kPPrimary,nt);
+ 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);
if (nstack == 0) continue;
}
- TFile *File=0;
+ TFile *pFile=0;
// Get AliRun object or create it
if (!gAlice) {
- gAlice = (AliRun*)File->Get("gAlice");
+ gAlice = (AliRun*)pFile->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();
+ if (fAli) pFile =fAli->GetCurrentFile();
+ pFile->cd();
}
+AliGenFLUKAsource& AliGenFLUKAsource::operator=(const AliGenFLUKAsource& rhs)
+{
+// Assignment operator
+ rhs.Copy(*this);
+ return (*this);
+}
+void AliGenFLUKAsource::Copy(TObject &) const
+{
+ Fatal("Copy","Not implemented!\n");
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff