Skip to content

Commit

Permalink
Code for spectators from data driven fragmentation
Browse files Browse the repository at this point in the history 
- AliGenHijingEventHeader where I added data members for the number of free spectators calculated in HIJING and two other data members to keep track of HIJING settings about fragmentation (whether it is applied or not) and spectators (whether they are put in the stack or not)
- AliGenHijing here I added the code needed to calculate for any generated b value the number of free spectators from parametrized 2011 data
- AliZDC and AliZDCDigitizer to check if the fragmentation has already been considered in HIJING or not and then if the kinematics one would like to use for the generated spectators is from HIJING or from our parametrization (inclduing Fremi motion, beam crossing and beam divergence)

Chiara Oppedisano
  • Loading branch information
@amorsch
amorsch committed Nov 12, 2015
1 parent 207bbbc commit 85d5ca3
Show file tree
Hide file tree
Showing 8 changed files with 302 additions and 68 deletions.
9 changes: 7 additions & 2 deletions STEER/STEERBase/AliGenHijingEventHeader.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -25,7 +25,9 @@ AliGenHijingEventHeader::AliGenHijingEventHeader():
fJet1(0., 0., 0., 0.),
fJet2(0., 0., 0., 0.),
fJetFsr1(0., 0., 0., 0.),
fJetFsr2(0., 0., 0., 0.)
fJetFsr2(0., 0., 0., 0.),
fAreSpectatorsInTheStack(kFALSE),
fIsDataFragmentationSet(kFALSE)
{
// Constructor
}
Expand All @@ -38,7 +40,10 @@ AliGenHijingEventHeader::AliGenHijingEventHeader(const char* name):
fJet1(0., 0., 0., 0.),
fJet2(0., 0., 0., 0.),
fJetFsr1(0., 0., 0., 0.),
fJetFsr2(0., 0., 0., 0.)
fJetFsr2(0., 0., 0., 0.),
fAreSpectatorsInTheStack(kFALSE),
fIsDataFragmentationSet(kFALSE)

{
// Copy Constructor
}
19 changes: 18 additions & 1 deletion STEER/STEERBase/AliGenHijingEventHeader.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,12 @@ class AliGenHijingEventHeader : public AliGenEventHeader, public AliCollisionGeo
Float_t TotalEnergy() const {return fTotalEnergy;}
Int_t Trials() const {return fTrials;}
Int_t GetTrueNPart() const {return fNPart;}
Bool_t GetSpectatorsInTheStack() const {return fAreSpectatorsInTheStack;}
Bool_t GetFragmentationFromData() const {return fIsDataFragmentationSet;}
Int_t GetFreeProjSpecn() const {return fFreeProjSpecn;}
Int_t GetFreeProjSpecp() const {return fFreeProjSpecp;}
Int_t GetFreeTargSpecn() const {return fFreeTargSpecn;}
Int_t GetFreeTargSpecp() const {return fFreeTargSpecp;}

// Setters
void SetTotalEnergy(Float_t energy) {fTotalEnergy=energy;}
Expand All @@ -32,6 +38,10 @@ class AliGenHijingEventHeader : public AliGenEventHeader, public AliCollisionGeo
{jet1 = fJet1; jet2 = fJet2; jet3 = fJetFsr1; jet4 = fJetFsr2;}
void SetTrials(Int_t trials) {fTrials = trials;}
void SetTrueNPart(Int_t npart) {fNPart = npart;}
void SetSpectatorsInTheStack(Bool_t what) {fAreSpectatorsInTheStack=what;}
void SetDataFromFragmentation(Bool_t what) {fIsDataFragmentationSet=what;}
void SetFreeSpectators(Int_t specnproj, Int_t specpproj, Int_t specntarg, Int_t specptarg)
{fFreeProjSpecn=specnproj; fFreeProjSpecp=specpproj; fFreeTargSpecn=specntarg; fFreeTargSpecp=specptarg;}

protected:
Float_t fTotalEnergy; // Total energy of produced particles
Expand All @@ -41,8 +51,15 @@ class AliGenHijingEventHeader : public AliGenEventHeader, public AliCollisionGeo
TLorentzVector fJet2; // 4-Momentum-Vector of second triggered jet
TLorentzVector fJetFsr1; // 4-Momentum-Vector of first triggered jet
TLorentzVector fJetFsr2; // 4-Momentum-Vector of second triggered jet
// Added by Chiara O. for spectator generation
Bool_t fAreSpectatorsInTheStack; // check whether spectators are in the stack
Bool_t fIsDataFragmentationSet; // check if the data driven correction is switched on
Int_t fFreeProjSpecn; // Num. of spectator neutrons from projectile nucleus
Int_t fFreeProjSpecp; // Num. of spectator protons from projectile nucleus
Int_t fFreeTargSpecn; // Num. of spectator neutrons from target nucleus
Int_t fFreeTargSpecp; // Num. of spectator protons from target nucleus

ClassDef(AliGenHijingEventHeader,5) // Event header for hijing event
ClassDef(AliGenHijingEventHeader,6) // Event header for hijing event
};

#endif
95 changes: 88 additions & 7 deletions THijing/AliGenHijing.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -28,6 +28,8 @@
#include <TLorentzVector.h>
#include <TPDGCode.h>
#include <TParticle.h>
#include <TF1.h>
#include <TFile.h>

#include "AliGenHijing.h"
#include "AliGenHijingEventHeader.h"
Expand Down Expand Up @@ -76,7 +78,14 @@ AliGenHijing::AliGenHijing()
fNoHeavyQuarks(kFALSE),
fHeader(AliGenHijingEventHeader("Hijing")),
fSigmaNN(-1),
fNoElas(0)
fNoElas(0),
fDataFragmentation(kFALSE),
fFreeProjSpecn(0),
fFreeProjSpecp(0),
fFreeTargSpecn(0),
fFreeTargSpecp(0),
fFragmNeutrons(0x0),
fFragmProtons(0x0)
{
// Constructor
fEnergyCMS = 5500.;
Expand Down Expand Up @@ -122,7 +131,15 @@ AliGenHijing::AliGenHijing(Int_t npart)
fNoHeavyQuarks(kFALSE),
fHeader(AliGenHijingEventHeader("Hijing")),
fSigmaNN(-1),
fNoElas(0)
fNoElas(0),
fDataFragmentation(kFALSE),
fFreeProjSpecn(0),
fFreeProjSpecp(0),
fFreeTargSpecn(0),
fFreeTargSpecp(0),
fFragmNeutrons(0x0),
fFragmProtons(0x0)

{
// Default PbPb collisions at 5. 5 TeV
//
Expand Down Expand Up @@ -211,6 +228,14 @@ void AliGenHijing::Init()
fHijing->SetIHPR2(49, 0);
}

if(fDataFragmentation){
TFile *file = TFile::Open("$ALICE_ROOT/ZDC/fragmSpecDataDriven.root","READ");
if(!file->IsOpen()){
AliError("Could not open file $ALICE_ROOT/ZDC/fragmSpecDataDriven.root");
}
fFragmNeutrons = dynamic_cast<TF1*> (file->Get("funcorrn"));
fFragmProtons = dynamic_cast<TF1*> (file->Get("funcorrp"));
}

AliGenMC::Init();

Expand Down Expand Up @@ -244,7 +269,7 @@ void AliGenHijing::Generate()
//
Int_t nt = 0;
Int_t jev = 0;
Int_t j, kf, ks, ksp, imo;
Int_t j=0, kf=0, ks=0, ksp=0, imo=0;
kf = 0;


Expand Down Expand Up @@ -318,8 +343,8 @@ void AliGenHijing::Generate()
ks = iparticle->GetStatusCode();
if (kf == 92) continue;

if (!fSelectAll) selected = KinematicSelection(iparticle, 0) &&
SelectFlavor(kf);
if (!fSelectAll) selected = KinematicSelection(iparticle, 0) &&
SelectFlavor(kf);
hasSelectedDaughters = DaughtersSelection(iparticle);
//
// Put particle on the stack if it is either selected or
Expand Down Expand Up @@ -369,10 +394,20 @@ void AliGenHijing::Generate()
pSelected[i] = 1;
} // selected
} // particle loop final state


if(fDataFragmentation){
Float_t impPar = fHijing->GetBB();
fFreeProjSpecn = FreeSpectatorsn(impPar, fProjectileSpecn);
fFreeProjSpecp = FreeSpectatorsp(impPar, fProjectileSpecp);
fFreeTargSpecn = FreeSpectatorsn(impPar, fTargetSpecn);
fFreeTargSpecp = FreeSpectatorsp(impPar, fTargetSpecp);
//
printf("\n b %f fm - SPECTATORS* HIJING: %d %d %d %d -> DATA: %d %d %d %d \n", impPar, fProjectileSpecn, fProjectileSpecp, fTargetSpecn, fTargetSpecp, fFreeProjSpecn, fFreeProjSpecp, fFreeTargSpecn, fFreeTargSpecp);
}
//
// Write particles to stack

Int_t countSpecPn=0, countSpecPp=0;
Int_t countSpecTn=0, countSpecTp=0;
for (i = 0; i<np; i++) {
iparticle = (TParticle *) fParticles.At(i);
Bool_t hasMother = (iparticle->GetFirstMother() >=0);
Expand All @@ -396,6 +431,26 @@ void AliGenHijing::Generate()
imo = (mother->GetPdgCode() != 92) ? newPos[imo] : -1;
} // if has mother
Bool_t tFlag = (fTrackIt && !hasDaughter);
if(fDataFragmentation && (ksp==0 || ksp==1) && kf == kNeutron){ // Projectile neutrons
countSpecPn++;
if(countSpecPn>fFreeProjSpecn) continue;
//printf(" Putting spec.n from proj. %d into the stack\n", countSpecPn);
}
if(fDataFragmentation && (ksp==0 || ksp==1) && kf == kProton){ // Projectile protons
countSpecPp++;
if(countSpecPp>fFreeProjSpecp) continue;
//printf(" Putting spec.p from proj. %d into the stack\n", countSpecPp);
}
if(fDataFragmentation && (ksp==10 || ksp==11) && kf == kNeutron){ // Target neutrons
countSpecTn++;
if(countSpecTn>fFreeTargSpecn) continue;
//printf(" Putting spec.n from targ. %d into the stack\n", countSpecTn);
}
if(fDataFragmentation && (ksp==10 || ksp==11) && kf == kProton){ // Target protons
countSpecTp++;
if(countSpecTp>fFreeTargSpecp) continue;
//printf(" Putting spec.p from targ. %d into the stack\n", countSpecTp);
}
PushTrack(tFlag,imo,kf,p,origin,polar,tof,kPNoProcess,nt, 1., ks);
fNprimaries++;
KeepTrack(nt);
Expand Down Expand Up @@ -626,6 +681,10 @@ void AliGenHijing::MakeHeader()
}
}
fHeader.SetNDiffractive(nsd1, nsd2, ndd);
fHeader.SetFreeSpectators(fFreeProjSpecn, fFreeProjSpecp, fFreeTargSpecn, fFreeTargSpecp);
fHeader.SetSpectatorsInTheStack(fSpectators);
fHeader.SetDataFromFragmentation(fDataFragmentation);

AddHeader(&fHeader);
fCollisionGeometry = &fHeader;
}
Expand Down Expand Up @@ -687,3 +746,25 @@ Bool_t AliGenHijing::CheckTrigger()
} // fTrigger == 2
return triggered;
}

Int_t AliGenHijing::FreeSpectatorsn(Float_t b, Int_t nSpecn)
{
// Select no. of spectator neutrons to put on the stack
Float_t corr = fFragmNeutrons->Eval(b);

Float_t nSpecInStack = corr*nSpecn;
//printf(" FreeSpectatorsn: b %f corr. %f nspecn_corr %f \n",b, corr, nSpecInStack);

return (int) nSpecInStack;
}

Int_t AliGenHijing::FreeSpectatorsp(Float_t b, Int_t nSpecp)
{
// Select no. of spectator protons to put on the stack
Double_t corr = fFragmProtons->Eval(b);

Float_t nSpecInStack = corr*nSpecp;
//printf(" FreeSpectatorsp: corr. %f nspecn_corr %f \n",corr, nSpecInStack);

return (int) nSpecInStack;
}
22 changes: 20 additions & 2 deletions THijing/AliGenHijing.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,7 @@ class THijing;
class TParticle;
class TClonesArray;
class TGraph;
class TF1;


class AliGenHijing : public AliGenMC
Expand Down Expand Up @@ -61,7 +62,8 @@ class AliGenHijing : public AliGenMC
virtual void SetRandomPz(Bool_t flag = 0) {fRandomPz = flag;}
virtual void SwitchOffHeavyQuarks(Bool_t flag = kTRUE) {fNoHeavyQuarks = flag;}
virtual void SetSigmaNN(Float_t val) {fSigmaNN = val;}
virtual void SetNoElas(Bool_t b) {fNoElas = b; }
virtual void SetNoElas(Bool_t b) {fNoElas = b; }
virtual void SetDataDrivenSpectators() {fDataFragmentation = kTRUE;}

// Getters
virtual TString GetReferenceFrame() const {return fFrame;}
Expand All @@ -82,6 +84,11 @@ class AliGenHijing : public AliGenMC
{phimin = fPhiMinJet*180./TMath::Pi(); phimax = fPhiMaxJet*180./TMath::Pi();}
THijing *GetTHijing() const {return fHijing;}
virtual Float_t GetSigmaNN() const {return fSigmaNN;}
virtual Bool_t HasDataDrivenSpectators() const {return fDataFragmentation;}
virtual Int_t GetFreeProjSpecn() const {return fFreeProjSpecn;}
virtual Int_t GetFreeProjSpecp() const {return fFreeProjSpecp;}
virtual Int_t GetFreeTargSpecn() const {return fFreeTargSpecn;}
virtual Int_t GetFreeTargSpecp() const {return fFreeTargSpecp;}

// Physics Routines
virtual Bool_t ProvidesCollisionGeometry() const {return kTRUE;}
Expand Down Expand Up @@ -133,6 +140,14 @@ class AliGenHijing : public AliGenMC
AliGenHijingEventHeader fHeader; // MC Header
Float_t fSigmaNN; // If not -1 set sigmaNN (HIPR1)
Bool_t fNoElas; // If true switch off elastic scattering
// Added by Chiara to consider fragments production from data
Bool_t fDataFragmentation; // Spectators w. data driven correction
Int_t fFreeProjSpecn;// Num. of spectator neutrons from projectile nucleus
Int_t fFreeProjSpecp;// Num. of spectator protons from projectile nucleus
Int_t fFreeTargSpecn; // Num. of spectator neutrons from target nucleus
Int_t fFreeTargSpecp; // Num. of spectator protons from target nucleus
TF1 *fFragmNeutrons; // data driven correction for nuclear fragment formation
TF1 *fFragmProtons; // data driven correction for nuclear fragment formation

private:
AliGenHijing(const AliGenHijing &Hijing);
Expand All @@ -144,7 +159,10 @@ class AliGenHijing : public AliGenMC
Bool_t DaughtersSelection(const TParticle* iparticle);
// check if stable
Bool_t Stable(const TParticle* particle) const;
// calculate no. of free spectators from data driven parametrization
Int_t FreeSpectatorsn(Float_t b, Int_t nSpecn);
Int_t FreeSpectatorsp(Float_t b, Int_t nSpecp);

ClassDef(AliGenHijing, 9) // AliGenerator interface to Hijing
ClassDef(AliGenHijing, 10) // AliGenerator interface to Hijing
};
#endif
12 changes: 9 additions & 3 deletions ZDC/ZDCsim/AliZDC.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -68,7 +68,8 @@ AliZDC::AliZDC() :
fIspASystem(kFALSE),
fIsRELDISgen(kFALSE),
fOnlyZEM(kFALSE),
fFindMother(kFALSE)
fFindMother(kFALSE),
fSpectatorParam(1)
{
//
// Default constructor for the Zero Degree Calorimeter base class
Expand All @@ -95,7 +96,9 @@ AliZDC::AliZDC(const char *name, const char *title) :
fIspASystem(kFALSE),
fIsRELDISgen(kFALSE),
fOnlyZEM(kFALSE),
fFindMother(kFALSE)
fFindMother(kFALSE),
fSpectatorParam(1)

{
//
// Standard constructor for the Zero Degree Calorimeter base class
Expand Down Expand Up @@ -140,7 +143,9 @@ fBeamEnergy(ZDC.fBeamEnergy),
fIspASystem(ZDC.fIspASystem),
fIsRELDISgen(ZDC.fIsRELDISgen),
fOnlyZEM(ZDC.fOnlyZEM),
fFindMother(ZDC.fFindMother)
fFindMother(ZDC.fFindMother),
fSpectatorParam(ZDC.fSpectatorParam)

{
// copy constructor
}
Expand Down Expand Up @@ -446,6 +451,7 @@ void AliZDC::Hits2SDigits()
AliDigitizer* AliZDC::CreateDigitizer(AliDigitizationInput* digInput) const{
// Create the digitizer for ZDC
AliZDCDigitizer *zdcDigitizer = new AliZDCDigitizer(digInput);
if(fSpectatorParam!=1) zdcDigitizer->SetSpectatorParam(fSpectatorParam);
if(fSpectatorTracked==0) zdcDigitizer->SetSpectators2Track();
if(fBeamEnergy>0.01) zdcDigitizer->SetBeamEnergy(fBeamEnergy);
if(fIspASystem==kTRUE) zdcDigitizer->SetpAsystem();
Expand Down
13 changes: 7 additions & 6 deletions ZDC/ZDCsim/AliZDC.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -70,7 +70,7 @@ class AliZDC : public AliDetector {
virtual void SetYZPA(Float_t /*yZPA*/) {}

virtual void SetSwitchOnTrackreferences() {}

//Calibration methods
void SetZDCCalibFName(const char *name);
char* GetZDCCalibFName() const {return (char*)fZDCCalibFName.Data();}
Expand All @@ -84,9 +84,9 @@ class AliZDC : public AliDetector {
virtual AliTriggerDetector* CreateTriggerDetector() const
{return new AliZDCTrigger();}


void SetSpectatorParam(int imod=1) {fSpectatorParam=imod;}
void SetSpectatorsTrack() {fSpectatorTracked=0;}
Int_t SpectatorsTracked() const {return fSpectatorTracked;}
Int_t AreSpectatorsTracked() const {return fSpectatorTracked;}
void SetBeamEnergy(Float_t beamEnergy) {fBeamEnergy = beamEnergy;}
void SetpAsystem() {fIspASystem = kTRUE;}
void SetRELDISGenerator() {fIsRELDISgen = kTRUE;}
Expand Down Expand Up @@ -114,10 +114,11 @@ class AliZDC : public AliDetector {
Bool_t fIspASystem; // Configuring pA collisions (MC only)
Bool_t fIsRELDISgen; // Is RELDIS used as generator

Bool_t fOnlyZEM; // build only ZEM (no had. calorimeters!)
Bool_t fFindMother; // look for particle mothers in the stack in StepManager
Bool_t fOnlyZEM; // build only ZEM (no had. calorimeters!)
Bool_t fFindMother; // look for particle mothers in the stack in StepManager
Int_t fSpectatorParam; // kinematic model fro spectators (1=from AliGenZDC(DEFAULT), 2=from HIJING)

ClassDef(AliZDC,14) // Zero Degree Calorimeter base class
ClassDef(AliZDC,15) // Zero Degree Calorimeter base class
};

// Calibration
Expand Down
Loading

0 comments on commit 85d5ca3

Please sign in to comment.