#include <TDatabasePDG.h>
#include <TFile.h>
#include "AliConst.h"
+#include "AliLog.h"
#include "AliGenMUONLMR.h"
#include "AliMC.h"
#include "AliRun.h"
ClassImp(AliGenMUONLMR)
-AliGenMUONLMR::AliGenMUONLMR () : AliGenMC(), fNMuMin(2), fGenSingleProc(-1) {
+
+AliGenMUONLMR::AliGenMUONLMR () : AliGenMC(), fNMuMin(2), fGenSingleProc(-1), fCosTheta(0), fRhoLineShape(0), fHMultMu(0), fHNProc(0) {
//
// default constructor
//
fScaleMult[kPionLMR] = 0; // set pion multiplicity to zero
fScaleMult[kKaonLMR] = 0; // set kaon multiplicity to zero
Int_t pdg[7] = {211, 321, 221, 113, 223, 333, 331};
- char* fptname[7] = {"fPtPion","fPtKaon","fPtEta","fPtRho","fPtOmega","fPtPhi","fPtEtaPrime"};
- char* fyname[7] = {"fYPion","fYKaon","fYEta","fYRho","fYOmega","fYPhi","fYEtaPrime"};
- char* fnname[7] = {"fMultPion","fMultKaon","fMultEta","fMultRho","fMultOmega","fMultPhi","fMultEtaPrime"};
- char* fdname[2] = {"fDecPion","fDecKaon"};
+ const char* fptname[7] = {"fPtPion","fPtKaon","fPtEta","fPtRho","fPtOmega","fPtPhi","fPtEtaPrime"};
+ const char* fyname[7] = {"fYPion","fYKaon","fYEta","fYRho","fYOmega","fYPhi","fYEtaPrime"};
+ const char* fnname[7] = {"fMultPion","fMultKaon","fMultEta","fMultRho","fMultOmega","fMultPhi","fMultEtaPrime"};
+ const char* fdname[2] = {"fDecPion","fDecKaon"};
Double_t ptparam[7][3] = { {1,0.427,2.52}, // pions from Pythia
{1,0.58,2.57}, // kaons from Pythia
{1,0.641,2.62}, // eta from Pythia
}
Double_t delta = md3 * md3 / (mres * mres);
Double_t epsilon = mumass * mumass / (mres * mres);
- Int_t nbins = fDalitz[index]->GetNbinsX();
- Double_t xmin = fDalitz[index]->GetXaxis()->GetXmin();
+ Int_t nbins0 = fDalitz[index]->GetNbinsX();
+ Double_t xmin0 = fDalitz[index]->GetXaxis()->GetXmin();
Double_t deltax = fDalitz[index]->GetBinWidth(1);
- Double_t xd = xmin - deltax/2.;
- for (Int_t ibin = 0; ibin< nbins; ibin++) {
+ Double_t xd = xmin0 - deltax/2.;
+ for (Int_t ibin = 0; ibin< nbins0; ibin++) {
Double_t dalval = 0;
xd += deltax;
if (xd > 4. *epsilon) {
//-----------------------------------------------------------
-Double_t AliGenMUONLMR::YDistr(Double_t *px, Double_t *par){
+Double_t AliGenMUONLMR::YDistr(const Double_t *px, const Double_t *par){
// function for rapidity distribution: plateau at par[0] +
// gaussian tails centered at par[1] and with par[2]=sigma
Double_t x = TMath::Abs(px[0]);
//-----------------------------------------------------------
-Double_t AliGenMUONLMR::PtDistr(Double_t *px, Double_t *par){
+Double_t AliGenMUONLMR::PtDistr(const Double_t *px, const Double_t *par){
// pt distribution: power law
Double_t x = px[0];
Double_t func = par[0] * x / TMath::Power((1+(x/par[1])*(x/par[1])),par[2]);
Int_t nmuons = -1, npartPushed = 0, pdgPushed[100];
Double_t polar[3]= {0,0,0}; // Polarisation of the parent particle (for GEANT tracking)
Double_t origin0[3]; // Origin of the generated parent particle (for GEANT tracking)
+ Double_t time0; // Time0 of the generated parent particle
// Calculating vertex position per event
for (Int_t j=0;j<3;j++) origin0[j]=fOrigin[j];
+ time0 = fTimeOrigin;
if(fVertexSmear==kPerEvent) {
Vertex();
for (Int_t j=0;j<3;j++) origin0[j]=fVertex[j];
+ time0 = fTime;
}
printf ("In Generate()\n");
if (TMath::Abs(pdgPushed[ipart]) != 13) { // particle is not a muon, hence it's a mother
PushTrack(0,-1,pdgPushed[ipart],
pxPushed[ipart],pyPushed[ipart],pzPushed[ipart],ePushed[ipart],
- origin0[0],origin0[1],origin0[2],0.,
+ origin0[0],origin0[1],origin0[2],time0,
polar[0],polar[1],polar[2],
kPPrimary,ntmother,1,11);
KeepTrack(ntmother);
else {
PushTrack(1,ntmother,pdgPushed[ipart],
pxPushed[ipart],pyPushed[ipart],pzPushed[ipart],ePushed[ipart],
- origin0[0],origin0[1],origin0[2],0.,
+ origin0[0],origin0[1],origin0[2],time0,
polar[0],polar[1],polar[2],
kPDecay,ntchild,1,1);
KeepTrack(ntchild);
SetHighWaterMark(ntchild);
AliGenEventHeader* header = new AliGenEventHeader("LMR");
header->SetPrimaryVertex(fVertex);
+ header->SetInteractionTime(fTime);
header->SetNProduced(fNprimaries);
AddHeader(header);
}
//------------------------------------------------------------------
-void AliGenMUONLMR::Decay2Body(TParticle *mother){
+void AliGenMUONLMR::Decay2Body(const TParticle *mother){
// performs decay in two muons of the low mass resonances
Double_t md1 = fMu[0]->GetMass();
Int_t pdg = mother->GetPdgCode();
if (mother->GetPdgCode()>0) fMu[0]->SetMomentum(v1.Px(),v1.Py(),v1.Pz(),v1.E());
else fMu[1]->SetMomentum(v1.Px(),v1.Py(),v1.Pz(),v1.E());
- Int_t idmother = -1;
- if (TMath::Abs(mother->GetPdgCode())== 211) idmother = 0;
- if (TMath::Abs(mother->GetPdgCode())== 321) idmother = 1;
+ Int_t idmother = 0;
+ if (TMath::Abs(mother->GetPdgCode())== 211) {
+ idmother = 0;
+ } else if (TMath::Abs(mother->GetPdgCode())== 321) {
+ idmother = 1;
+ } else {
+ AliWarning("Mother particle is not a pion or kaon \n");
+ }
+
Double_t gammaRes = mother->Energy()/mres;
Double_t zResCM = fDecay[idmother]->GetRandom();
Double_t zResLab = gammaRes*zResCM;
//-------------------------------------------------------------------
-void AliGenMUONLMR::DalitzDecay(TParticle *mother){
+void AliGenMUONLMR::DalitzDecay(const TParticle *mother){
//
// perform dalitz decays of eta, omega and etaprime
//
Double_t mumass = fMu[0]->GetMass();
Double_t md3 = 0; // unless differently specified, third particle is a photon
if (mother->GetPdgCode() == 223) md3 = 0.134977; // if mother is an omega, third particle is a pi0
- Int_t index = -1;
- if (mother->GetPdgCode() == 221) index = 0; // eta
- else if (mother->GetPdgCode() == 223) index = 1; // omega
- else if (mother->GetPdgCode() == 331) index = 2; // etaPrime
+
+ Int_t index = 0;
+ if (TMath::Abs(mother->GetPdgCode())== 221) {
+ // eta
+ index = 0;
+ } else if (TMath::Abs(mother->GetPdgCode())== 223) {
+ // omega
+ index = 1;
+ } else if (mother->GetPdgCode() == 331) {
+ // eta'
+ index = 2;
+ } else {
+ AliWarning("Mother particle is not a eta, eta' or omega \n");
+ }
+
Int_t flag = 0;
Double_t xd=0, mvirt2=0;
Double_t countIt = 0;
//____________________________________________________________
-Double_t AliGenMUONLMR::RhoLineShapeNew(Double_t *x, Double_t *para){
+Double_t AliGenMUONLMR::RhoLineShapeNew(const Double_t *x, const Double_t* /*para*/){
//new parameterization implemented by Hiroyuki Sako (GSI)
Double_t mass = *x;
double r, GammaTot;