// Uses the TDPMjet implementation of TGenerator.
#include <TDPMjet.h>
+#include "DPMcommon.h"
#include <TRandom.h>
#include <TArrayI.h>
#include <TParticle.h>
#include "AliGenDPMjetEventHeader.h"
#include "AliRun.h"
#include "AliDpmJetRndm.h"
+#include "AliIonPDGCodes.h"
#include "AliHeader.h"
#include "AliStack.h"
#include "AliMC.h"
//______________________________________________________________________________
AliGenDPMjet::AliGenDPMjet()
:AliGenMC(),
- fBeamEn(2750.),
+ fBeamEn(0.),
fMinImpactParam(0.),
- fMaxImpactParam(5.),
+ fMaxImpactParam(100.),
fICentr(0),
fSelectAll(0),
fFlavor(0),
fTrials(0),
+ fNprimaries(0),
fSpectators(1),
fSpecn(0),
fSpecp(0),
fDecayAll(0),
fGenImpPar(0.),
fProcess(kDpmMb),
+ fTriggerParticle(0),
+ fTriggerEta(0.9),
+ fTriggerMinPt(-1),
+ fTriggerMaxPt(1000),
fTriggerMultiplicity(0),
fTriggerMultiplicityEta(0),
fTriggerMultiplicityPtMin(0),
fkTuneForDiff(0),
- fProcDiff(0)
+ fProcDiff(0),
+ fFragmentation(kFALSE),
+ fHeader(AliGenDPMjetEventHeader("DPMJET"))
+
{
// Constructor
fEnergyCMS = 5500.;
//______________________________________________________________________________
AliGenDPMjet::AliGenDPMjet(Int_t npart)
:AliGenMC(npart),
- fBeamEn(2750.),
+ fBeamEn(0.),
fMinImpactParam(0.),
- fMaxImpactParam(5.),
+ fMaxImpactParam(100.),
fICentr(0),
fSelectAll(0),
fFlavor(0),
fTrials(0),
+ fNprimaries(0),
fSpectators(1),
fSpecn(0),
fSpecp(0),
fDecayAll(0),
fGenImpPar(0.),
fProcess(kDpmMb),
+ fTriggerParticle(0),
+ fTriggerEta(0.9),
+ fTriggerMinPt(-1),
+ fTriggerMaxPt(1000),
fTriggerMultiplicity(0),
fTriggerMultiplicityEta(0),
fTriggerMultiplicityPtMin(0),
fkTuneForDiff(0),
- fProcDiff(0)
+ fProcDiff(0),
+ fFragmentation(kFALSE),
+ fHeader(AliGenDPMjetEventHeader("DPMJET"))
+
+
{
// Default PbPb collisions at 5. 5 TeV
//
AliGenDPMjet::AliGenDPMjet(const AliGenDPMjet &/*Dpmjet*/)
:AliGenMC(),
- fBeamEn(2750.),
+ fBeamEn(0.),
fMinImpactParam(0.),
- fMaxImpactParam(5.),
+ fMaxImpactParam(100.),
fICentr(0),
fSelectAll(0),
fFlavor(0),
fTrials(0),
+ fNprimaries(0),
fSpectators(1),
fSpecn(0),
fSpecp(0),
fDecayAll(0),
fGenImpPar(0.),
fProcess(kDpmMb),
+ fTriggerParticle(0),
+ fTriggerEta(0.9),
+ fTriggerMinPt(-1),
+ fTriggerMaxPt(1000),
fTriggerMultiplicity(0),
fTriggerMultiplicityEta(0),
fTriggerMultiplicityPtMin(0),
fkTuneForDiff(0),
- fProcDiff(0)
+ fProcDiff(0),
+ fFragmentation(kFALSE),
+ fHeader(0x0)
+
+
{
// Dummy copy constructor
fEnergyCMS = 5500.;
{
// Initialization
+ if(fEnergyCMS>0. && fBeamEn<0.1) fBeamEn = fEnergyCMS/2;
SetMC(new TDPMjet(fProcess, fAProjectile, fZProjectile, fATarget, fZTarget,
fBeamEn,fEnergyCMS));
fDPMjet=(TDPMjet*) fMCEvGen;
- //
+ if (fAProjectile == 1 && TMath::Abs(fZProjectile == 1)) fDPMjet->SetfIdp(1);
+
// **** Flag to force central production
// fICentr=1. central production forced
// fICentr<0 && fICentr>-100 -> bmin = fMinImpactParam, bmax = fMaxImpactParam
// fICentr<-99 -> fraction of x-sec. = XSFRAC
// fICentr=-1. -> evaporation/fzc suppressed
- // fICentr<-1. -> evaporation/fzc allowed
- if (fAProjectile == 1 && TMath::Abs(fZProjectile == 1)) fDPMjet->SetfIdp(1);
-
+ // fICentr<-1. -> evaporation/fzc allowed
fDPMjet->SetfFCentr(fICentr);
+
fDPMjet->SetbRange(fMinImpactParam, fMaxImpactParam);
fDPMjet->SetPi0Decay(fPi0Decay);
fDPMjet->SetDecayAll(fDecayAll);
-
- AliGenMC::Init();
+ fDPMjet->SetFragmentProd(fFragmentation);
//
// Initialize DPMjet
Int_t kf, ks, imo;
kf = 0;
fTrials = 0;
+ fNprimaries = 0;
// Set collision vertex position
if (fVertexSmear == kPerEvent) Vertex();
Printf("Triggered on event with multiplicity of %d >= %d", multiplicity, fTriggerMultiplicity);
}
+ //Trigger on the presence of a given particle in some phase space
+ if (fTriggerParticle) {
+ Bool_t triggered = kFALSE;
+ for (Long_t i = 0; i < np; i++) {
+ TParticle * iparticle = (TParticle *) fParticles.At(i);
+ kf = CheckPDGCode(iparticle->GetPdgCode());
+ if (kf != fTriggerParticle) continue;
+ if (iparticle->Pt() == 0.) continue;
+ if (TMath::Abs(iparticle->Eta()) > fTriggerEta) continue;
+ if ( iparticle->Pt() > fTriggerMaxPt || iparticle->Pt() < fTriggerMinPt ) continue;
+ triggered = kTRUE;
+ break;
+ }
+ if (!triggered) continue;
+ }
if(fkTuneForDiff && ( (TMath::Abs(fEnergyCMS - 900) < 1) || (TMath::Abs(fEnergyCMS - 2760) < 1) || (TMath::Abs(fEnergyCMS - 7000) < 1)) ) {
if(!CheckDiffraction() ) continue;
- Bool_t tFlag = (fTrackIt && (ks == 1));
+ Bool_t tFlag = (fTrackIt && (ks==1 || ks==-1));
+ //printf(" AliGemDPMJet->PushTrack: kf %d ks %d flag %d\n",kf,ks,tFlag);
+ if(kf>10000 && (ks==-1 || ks==1000 || ks==1001)) kf += 1000000000;
PushTrack(tFlag, imo, kf,
p[0], p[1], p[2], p[3],
origin[0], origin[1], origin[2], tof,
polar[0], polar[1], polar[2],
kPNoProcess, nt, 1., ks);
KeepTrack(nt);
+ fNprimaries++;
newPos[i] = nt;
} // if selected
} // particle loop
{
// Return true for a stable particle
//
-
-// if (particle->GetFirstDaughter() < 0 ) return kTRUE;
- if (particle->GetStatusCode() == 1) return kTRUE;
+ int st = particle->GetStatusCode();
+ if(st == 1 || st == -1) return kTRUE;
else return kFALSE;
}
//______________________________________________________________________________
void AliGenDPMjet::MakeHeader()
{
-// printf("MakeHeader %13.3f \n", fDPMjet->GetBImpac());
// Builds the event header, to be called after each event
- AliGenEventHeader* header = new AliGenDPMjetEventHeader("DPMJET");
- ((AliGenDPMjetEventHeader*) header)->SetNProduced(fDPMjet->GetNumStablePc());
- ((AliGenDPMjetEventHeader*) header)->SetImpactParameter(fDPMjet->GetBImpac());
- ((AliGenDPMjetEventHeader*) header)->SetTotalEnergy(fDPMjet->GetTotEnergy());
- ((AliGenDPMjetEventHeader*) header)->SetParticipants(fDPMjet->GetProjParticipants(),
- fDPMjet->GetTargParticipants());
-
- if(fProcDiff>0){
- ((AliGenDPMjetEventHeader*) header)->SetProcessType(fProcDiff);
- }
- else
- ((AliGenDPMjetEventHeader*) header)->SetProcessType(fDPMjet->GetProcessCode());
+ fHeader.SetNProduced(fNprimaries);
+ fHeader.SetImpactParameter(fDPMjet->GetBImpac());
+ fHeader.SetTotalEnergy(fDPMjet->GetTotEnergy());
+ fHeader.SetParticipants(fDPMjet->GetProjParticipants(),
+ fDPMjet->GetTargParticipants());
+
+ fHeader.SetCollisions(DTGLCP.ncp, DTGLCP.nct,
+ fDPMjet->GetProjWounded(),fDPMjet->GetTargWounded());
+
+ if(fProcDiff>0) fHeader.SetProcessType(fProcDiff);
+ else fHeader.SetProcessType(fDPMjet->GetProcessCode());
// Bookkeeping for kinematic bias
- ((AliGenDPMjetEventHeader*) header)->SetTrials(fTrials);
+ fHeader.SetTrials(fTrials);
// Event Vertex
- header->SetPrimaryVertex(fVertex);
- header->SetInteractionTime(fTime);
- gAlice->SetGenEventHeader(header);
- AddHeader(header);
+ fHeader.SetPrimaryVertex(fVertex);
+ fHeader.SetInteractionTime(fTime);
+ fHeader.SetNDiffractive(POEVT1.nsd1, POEVT1.nsd2, POEVT1.ndd);
+// gAlice->SetGenEventHeader(fHeader);
+ AddHeader(&fHeader);
+ fCollisionGeometry = &fHeader;
}
-void AliGenDPMjet::AddHeader(AliGenEventHeader* header)
+//______________________________________________________________________________
+/*void AliGenDPMjet::AddHeader(AliGenEventHeader* fHeader)
{
- // Add header to container or runloader
+ // Add fHeader to container or runloader
if (fContainer) {
- fContainer->AddHeader(header);
+ fContainer->AddHeader(fHeader);
} else {
- AliRunLoader::Instance()->GetHeader()->SetGenEventHeader(header);
+ AliRunLoader::Instance()->GetHeader()->SetGenEventHeader(fHeader);
}
-}
+}*/
//______________________________________________________________________________
TParticle * part = (TParticle *) fParticles.At(iPart);
Double_t E= part->Energy();
Double_t P= part->P();
- M= TMath::Sqrt((fEnergyCMS-E-P)*(fEnergyCMS-E+P));
+ Double_t M2 = (fEnergyCMS-E-P)*(fEnergyCMS-E+P);
+ if(M2<0) return kFALSE;
+ M= TMath::Sqrt(M2);
}
Double_t Mmin, Mmax, wSD, wDD, wND;
return kTRUE;
}
+// -------------------------------------------------------
+void AliGenDPMjet::SetIonPDGCodes()
+{
+ // Defining PDG codes for the ions
+ AliIonPDGCodes *pdgcodes = new AliIonPDGCodes();
+ pdgcodes->AddParticlesToPdgDataBase();
+}
-
+// -------------------------------------------------------
Bool_t AliGenDPMjet::GetWeightsDiffraction(Double_t M, Double_t &Mmin, Double_t &Mmax,
Double_t &wSD, Double_t &wDD, Double_t &wND)
{