Pythia8Generator Class Reference

#include <Pythia8Generator.h>

Inheritance diagram for Pythia8Generator:

Collaboration diagram for Pythia8Generator:

Public Member Functions
	Pythia8Generator ()
	~Pythia8Generator () override
Bool_t	ReadEvent (FairPrimaryGenerator *) override
void	SetParameters (char *)
void	Print ()
Bool_t	Init () override
Bool_t	Init (const char *inFile) override
Bool_t	Init (const char *inFile, int startEvent) override
void	SetMom (Double_t mom)
void	SetId (Double_t id)
void	UseRandom1 ()
void	UseRandom3 ()
void	SetfFDs (Double_t z)
void	SetTarget (TString s, Double_t x, Double_t y)
void	SetTargetCoordinates (Double_t start_z, Double_t end_z)
Int_t	nrOfRetries ()
Public Member Functions inherited from SHiP::Generator
	Generator ()
virtual	~Generator ()
virtual Bool_t	Init (const char *, int)=0
virtual Bool_t	Init (const char *)=0
virtual Bool_t	Init (const std::vector< std::string > &inFiles, int startNumber)
virtual Bool_t	Init (const std::vector< std::string > &inFiles)
virtual void	UseExternalFile (std::string x, Int_t i)
virtual void	UseExternalFile (std::vector< std::string > &inFiles, Int_t i)

Protected Attributes
Double_t	fMom
Int_t	fId
Bool_t	fUseRandom1
Bool_t	fUseRandom3
Float_t	hpx [1]
Float_t	hpy [1]
Float_t	hpz [1]
Float_t	hE [1]
Float_t	hM [1]
Float_t	mpx [1]
Float_t	mpy [1]
Float_t	mpz [1]
Float_t	mE [1]
Float_t	hid [1]
Float_t	mid [1]
Float_t	ck [1]
Float_t	ancestors [16]
Float_t	subprocCodes [16]
Int_t	fNevents
Int_t	fn
Int_t	fShipEventNr
TFile *	fInputFile
FairLogger *	fLogger
	pointer to a file
TTree *	fTree
	don't make it persistent, magic ROOT command
Pythia8::Pythia *	fPythia
Double_t	fFDs
Int_t	fnRetries
GenieGenerator *	fMaterialInvestigator
TString	targetName
Double_t	xOff
Double_t	yOff
Double_t	start [3]
Double_t	end [3]
Double_t	bparam
Double_t	mparam [10]
Double_t	startZ
Double_t	endZ
Bool_t	targetFromGeometry
Double_t	maxCrossSection
Protected Attributes inherited from SHiP::Generator
std::optional< std::string >	fextFile
Int_t	firstEvent = 0

Private Attributes
std::shared_ptr< Pythia8::RndmEngine >	fRandomEngine

Detailed Description

Definition at line 20 of file Pythia8Generator.h.

Constructor & Destructor Documentation

Pythia8Generator()

Pythia8Generator::Pythia8Generator

(

)

default constructor

Definition at line 27 of file Pythia8Generator.cxx.

                                   {
  fUseRandom1 = kFALSE;
  fUseRandom3 = kTRUE;
  fId = 2212;         // proton
  fMom = 400;         // proton
  fFDs = 7.7 / 10.4;  // correction for Pythia6 to match measured Ds production
  fInputFile = nullptr;
  targetName = "";
  xOff = 0;
  yOff = 0;
  targetFromGeometry = kFALSE;
  fPythia = new Pythia8::Pythia();
}

~Pythia8Generator()

Pythia8Generator::~Pythia8Generator ( )

override

destructor

Definition at line 151 of file Pythia8Generator.cxx.

{}

Member Function Documentation

Init() [1/3]

Bool_t Pythia8Generator::Init ( )

override

Definition at line 43 of file Pythia8Generator.cxx.

                              {
  if (fUseRandom1) fRandomEngine = std::make_shared<PyTr1Rng>();
  if (fUseRandom3) fRandomEngine = std::make_shared<PyTr3Rng>();
  if (fextFile) {
    fInputFile = TFile::Open(fextFile->c_str());
    LOG(info) << "Open external file with charm or beauty hadrons: "
              << *fextFile;
    if (!fInputFile) {
      LOG(fatal) << "Error opening input file.";
      return kFALSE;
    }
    fTree = fInputFile->Get<TTree>("pythia6");
    fNevents = fTree->GetEntries();
    fn = firstEvent;
    fTree->SetBranchAddress("id", &hid);  // particle id
    fTree->SetBranchAddress("px", &hpx);  // momentum
    fTree->SetBranchAddress("py", &hpy);
    fTree->SetBranchAddress("pz", &hpz);
    fTree->SetBranchAddress("E", &hE);
    fTree->SetBranchAddress("M", &hM);
    fTree->SetBranchAddress("mid", &mid);  // mother
    fTree->SetBranchAddress("mpx", &mpx);  // momentum
    fTree->SetBranchAddress("mpy", &mpy);
    fTree->SetBranchAddress("mpz", &mpz);
    fTree->SetBranchAddress("mE", &mE);
    if (fTree->GetBranch("k")) {
      fTree->SetBranchAddress("k", &ck);
      if (fTree->GetBranch("a0")) {
        for (Int_t i = 0; i < 16; i++) {
          TString na = "a";
          na += i;
          fTree->SetBranchAddress(na, &(ancestors[i]));
          TString ns = "s";
          ns += i;
          fTree->SetBranchAddress(ns, &(subprocCodes[i]));
        }
      }
    } else {
      ck[0] = 1;
      subprocCodes[0] = 88;
      ancestors[0] = 2212;
    }
    fPythia->readString("ProcessLevel:all = off");
    // find all long lived particles in pythia
    Int_t n = 1;
    while (n != 0) {
      n = fPythia->particleData.nextId(n);
      std::shared_ptr<Pythia8::ParticleDataEntry> p =
          fPythia->particleData.particleDataEntryPtr(n);
      if (p->tau0() > 1) {
        std::string particle = std::to_string(n) + ":mayDecay = false";
        fPythia->readString(particle);
        LOGF(info, "Made %s stable for Pythia, should decay in Geant4",
             p->name().c_str());
      }
    }
  } else {
    fPythia->setRndmEnginePtr(fRandomEngine);
    fPythia->settings.mode("Beams:idA", fId);
    fPythia->settings.mode("Beams:idB", 2212);
    fPythia->settings.mode("Beams:frameType", 2);
    fPythia->settings.parm("Beams:eA", fMom);  // codespell:ignore parm
    fPythia->settings.parm("Beams:eB", 0.);    // codespell:ignore parm
  }
  fPythia->init();
  if (targetFromGeometry) {
    // Use geometry-based coordinates passed from run_simScript.py
    fMaterialInvestigator = new GenieGenerator();
    start[0] = xOff;
    start[1] = yOff;
    start[2] = startZ;
    end[0] = xOff;
    end[1] = yOff;
    end[2] = endZ;
    LOG(info)
        << "Pythia8Generator: Using geometry-based target coordinates startZ="
        << startZ << " endZ=" << endZ;
    // find maximum interaction length
    bparam = shipgen::MeanMaterialBudget(start, end, mparam);
    maxCrossSection = mparam[9];
  } else if (targetName != "") {
    // Fallback to fragile TGeo navigation for backward compatibility
    fMaterialInvestigator = new GenieGenerator();
    TGeoVolume* top = gGeoManager->GetTopVolume();
    TGeoNode* target = top->FindNode(targetName);
    if (!target) {
      LOGF(error, "target not found, %s, program will crash",
           targetName.Data());
    }
    Double_t z_middle = target->GetMatrix()->GetTranslation()[2];
    TGeoBBox* sha = dynamic_cast<TGeoBBox*>(target->GetVolume()->GetShape());
    startZ = z_middle - sha->GetDZ();
    endZ = z_middle + sha->GetDZ();
    start[0] = xOff;
    start[1] = yOff;
    start[2] = startZ;
    end[0] = xOff;
    end[1] = yOff;
    end[2] = endZ;
    // find maximum interaction length
    bparam = shipgen::MeanMaterialBudget(start, end, mparam);
    maxCrossSection = mparam[9];
  }
  return kTRUE;
}

Init() [2/3]

Bool_t Pythia8Generator::Init ( const char *  inFile )

inlineoverridevirtual

Implements SHiP::Generator.

Definition at line 35 of file Pythia8Generator.h.

{ return Init(inFile, 0); };

Init() [3/3]

Bool_t Pythia8Generator::Init ( const char *  inFile, int  startEvent  )

inlineoverridevirtual

Implements SHiP::Generator.

Definition at line 37 of file Pythia8Generator.h.

                                                           {
    LOG(warning) << "Init with files not implemented for Pythia8Generator. "
                    "Using default Init() instead";
    return Init();
  };

nrOfRetries()

Int_t Pythia8Generator::nrOfRetries ( )

inline

Definition at line 65 of file Pythia8Generator.h.

{ return fnRetries; };

Print()

void Pythia8Generator::Print

(

)

ReadEvent()

Bool_t Pythia8Generator::ReadEvent ( FairPrimaryGenerator *  cpg )

override

public method ReadEvent

Definition at line 155 of file Pythia8Generator.cxx.

                                                            {
  Double_t x, y, z, px, py, pz, dl, e, tof;
  Int_t im, id, key;
  fnRetries = 0;
  // take charm hadrons from external file
  // correct eventually for too much primary Ds produced by pythia6
  key = 0;
  bool l = true;
  while (l) {
    if (fn == fNevents) {
      LOG(warning) << "End of input file. Rewind.";
    }
    fTree->GetEntry((fn + 1) % fNevents);
    // check that this and next entry is charm, otherwise continue reading
    l = false;
    if (static_cast<int>(mE[0]) == 0) {
      l = true;
    }
    bool isDs = false;
    if (static_cast<int>(fabs(hid[0])) == 431) {
      isDs = true;
    }
    fTree->GetEntry(fn % fNevents);
    if (static_cast<int>(mE[0]) == 0) {
      l = true;
    }
    fn++;
    if (static_cast<int>(fabs(hid[0])) == 431 || isDs) {
      Double_t rnr = gRandom->Uniform(0, 1);
      if (rnr > fFDs) {
        l = true;
        fn++;
      }
    }
  }
  Double_t zinter = 0;
  if (targetName != "") {
    // calculate primary proton interaction point:
    // loop over trajectory between start and end to pick an interaction point,
    // copied from GenieGenerator and adapted to hadrons
    Double_t prob2int = -1.;
    Double_t rndm = 0.;
    Double_t sigma;
    Double_t zinterStart = start[2];
    // simulate more downstream interaction points for interactions down in the
    // cascade
    Int_t nInter = ck[0];
    if (nInter > 16) {
      nInter = 16;
    }
    for (Int_t nI = 0; nI < nInter; nI++) {
      // if (!subprocCodes[nI]<90){continue;}  //if process is not inelastic, go
      // to next. Changed by taking now collision length
      prob2int = -1.;
      Double_t intLengthFactor = 1;  // for nucleons
      if (TMath::Abs(ancestors[nI]) < 1000) {
        intLengthFactor = 1.16;
      }  // for mesons
      // Fe: nuclear /\ 16.77 cm pion 20.42 cm  f=1.22
      // W:  nuclear /\ 9.946 cm pion 11.33 cm  f=1.14
      // Mo: nuclear /\ 15.25 cm pion 17.98 cm  f=1.18
      while (prob2int < rndm) {
        // place x,y,z uniform along path
        zinter = gRandom->Uniform(zinterStart, end[2]);
        Double_t point[3] = {xOff, yOff, zinter};
        bparam = shipgen::MeanMaterialBudget(start, point, mparam);
        Double_t interLength =
            mparam[8] * intLengthFactor *
            1.7;  // 1.7 = interaction length / collision length from PDG Tables
        TGeoNode* node = gGeoManager->FindNode(point[0], point[1], point[2]);
        TGeoMaterial* mat = nullptr;
        if (node && !gGeoManager->IsOutside()) {
          mat = node->GetVolume()->GetMaterial();
          Double_t n = mat->GetDensity() / mat->GetA();
          sigma = 1. / (n * mat->GetIntLen()) /
                  mbarn;  // no need to multiply with intLengthFactor, will
                          // cancel in next equation.
          prob2int = TMath::Exp(-interLength) * sigma / maxCrossSection;
        } else {
          prob2int = 0.;
        }
        rndm = gRandom->Uniform(0., 1.);
      }
      zinterStart = zinter;
    }
    zinter = zinter * cm;
  }
  for (Int_t c = 0; c < 2; c++) {
    if (c > 0) {
      fTree->GetEntry(fn % fNevents);
      fn++;
    }
    fPythia->event.reset();
    id = static_cast<Int_t>(hid[0]);
    fPythia->event.append(id, 1, 0, 0, hpx[0], hpy[0], hpz[0], hE[0], hM[0], 0.,
                          9.);
    // simulate displaced vertex, Pythia8 will not do it
    Double_t tau0 = fPythia->particleData.tau0(id);  // ctau in mm
    dl = gRandom->Exp(tau0) / hM[0];                 // mm/GeV
    fPythia->next();
    Int_t addedParticles = 0;
    if (c == 0) {
      // original particle responsible for interaction, "mother of charm" from
      // external file
      px = mpx[0];
      py = mpy[0];
      pz = mpz[0];
      x = 0.;
      y = 0.;
      z = zinter;
      id = mid[0];
      cpg->AddTrack(id, px, py, pz, x / cm, y / cm, z / cm, -1, false);
      addedParticles += 1;
    }
    for (Int_t ii = 1; ii < fPythia->event.size(); ii++) {
      id = fPythia->event[ii].id();
      Bool_t wanttracking = false;
      if (fPythia->event[ii].isFinal()) {
        wanttracking = true;
      }
      if (ii > 1) {
        z = fPythia->event[ii].zProd() + dl * fPythia->event[1].pz() + zinter;
        x = fPythia->event[ii].xProd() + dl * fPythia->event[1].px();
        y = fPythia->event[ii].yProd() + dl * fPythia->event[1].py();
        tof = fPythia->event[ii].tProd() / (10 * c_light) +
              dl * fPythia->event[1].e() / cm / c_light;
      } else {
        z = fPythia->event[ii].zProd() + zinter;
        x = fPythia->event[ii].xProd();
        y = fPythia->event[ii].yProd();
        tof =
            fPythia->event[ii].tProd() / (10 * c_light);  // to go from mm to s
      }
      pz = fPythia->event[ii].pz();
      px = fPythia->event[ii].px();
      py = fPythia->event[ii].py();
      e = fPythia->event[ii].e();
      im = fPythia->event[ii].mother1() + key;
 
      if (ii == 1) {
        im = 0;
      }
      cpg->AddTrack(id, px, py, pz, x / cm, y / cm, z / cm, im, wanttracking, e,
                    tof, 1.);
      addedParticles += 1;
    }
    key += addedParticles - 1;  // pythia counts from 1
  }
  counter += 1;
  // now the underlying event
  bool lx = true;
  while (lx) {
    fTree->GetEntry(fn % fNevents);
    lx = false;
    if (mE[0] == 0) {
      lx = true;
      fn++;
      cpg->AddTrack((Int_t)hid[0], hpx[0], hpy[0], hpz[0],
                    (mpx[0] + fPythia->event[0].xProd()) / cm,
                    (mpy[0] + fPythia->event[0].yProd()) / cm,
                    (mpz[0] + fPythia->event[0].zProd()) / cm + zinter / cm, -1,
                    true);
      // mpx,mpy,mpz are the vertex coordinates with respect to charm hadron,
      // first particle in Pythia is (system)
    }
  }
 
  return kTRUE;
}

SetfFDs()

void Pythia8Generator::SetfFDs ( Double_t  z )

inline

Definition at line 54 of file Pythia8Generator.h.

{ fFDs = z; };

SetId()

void Pythia8Generator::SetId ( Double_t  id )

inline

Definition at line 44 of file Pythia8Generator.h.

{ fId = id; };

SetMom()

void Pythia8Generator::SetMom ( Double_t  mom )

inline

Definition at line 43 of file Pythia8Generator.h.

{ fMom = mom; };

SetParameters()

void Pythia8Generator::SetParameters

(

char *

)

SetTarget()

void Pythia8Generator::SetTarget ( TString  s, Double_t  x, Double_t  y  )

inline

Definition at line 55 of file Pythia8Generator.h.

                                                    {
    targetName = s;
    xOff = x;
    yOff = y;
  };

SetTargetCoordinates()

void Pythia8Generator::SetTargetCoordinates ( Double_t  start_z, Double_t  end_z  )

inline

Definition at line 60 of file Pythia8Generator.h.

                                                              {
    startZ = start_z;
    endZ = end_z;
    targetFromGeometry = true;
  };

UseRandom1()

void Pythia8Generator::UseRandom1 ( )

inline

Definition at line 45 of file Pythia8Generator.h.

                    {
    fUseRandom1 = kTRUE;
    fUseRandom3 = kFALSE;
  };

UseRandom3()

void Pythia8Generator::UseRandom3 ( )

inline

Definition at line 49 of file Pythia8Generator.h.

                    {
    fUseRandom1 = kFALSE;
    fUseRandom3 = kTRUE;
  };

Member Data Documentation

ancestors

Float_t Pythia8Generator::ancestors[16]

protected

Definition at line 77 of file Pythia8Generator.h.

bparam

Double_t Pythia8Generator::bparam

protected

Definition at line 91 of file Pythia8Generator.h.

ck

Float_t Pythia8Generator::ck[1]

protected

Definition at line 76 of file Pythia8Generator.h.

end

Double_t Pythia8Generator::end[3]

protected

Definition at line 90 of file Pythia8Generator.h.

endZ

Double_t Pythia8Generator::endZ

protected

Definition at line 94 of file Pythia8Generator.h.

fFDs

Double_t Pythia8Generator::fFDs

protected

Definition at line 83 of file Pythia8Generator.h.

fId

Int_t Pythia8Generator::fId

protected

Definition at line 72 of file Pythia8Generator.h.

fInputFile

TFile* Pythia8Generator::fInputFile

protected

Definition at line 79 of file Pythia8Generator.h.

fLogger

FairLogger* Pythia8Generator::fLogger

protected

pointer to a file

Definition at line 80 of file Pythia8Generator.h.

fMaterialInvestigator

GenieGenerator* Pythia8Generator::fMaterialInvestigator

protected

Definition at line 85 of file Pythia8Generator.h.

fMom

Double_t Pythia8Generator::fMom

protected

Definition at line 71 of file Pythia8Generator.h.

fn

Int_t Pythia8Generator::fn

protected

Definition at line 78 of file Pythia8Generator.h.

fNevents

Int_t Pythia8Generator::fNevents

protected

Definition at line 78 of file Pythia8Generator.h.

fnRetries

Int_t Pythia8Generator::fnRetries

protected

Definition at line 84 of file Pythia8Generator.h.

fPythia

Pythia8::Pythia* Pythia8Generator::fPythia

protected

Definition at line 82 of file Pythia8Generator.h.

fRandomEngine

std::shared_ptr<Pythia8::RndmEngine> Pythia8Generator::fRandomEngine

private

Definition at line 68 of file Pythia8Generator.h.

fShipEventNr

Int_t Pythia8Generator::fShipEventNr

protected

Definition at line 78 of file Pythia8Generator.h.

fTree

TTree* Pythia8Generator::fTree

protected

don't make it persistent, magic ROOT command

Definition at line 81 of file Pythia8Generator.h.

fUseRandom1

Bool_t Pythia8Generator::fUseRandom1

protected

Definition at line 73 of file Pythia8Generator.h.

fUseRandom3

Bool_t Pythia8Generator::fUseRandom3

protected

Definition at line 74 of file Pythia8Generator.h.

hE

Float_t Pythia8Generator::hE[1]

protected

Definition at line 75 of file Pythia8Generator.h.

hid

Float_t Pythia8Generator::hid[1]

protected

Definition at line 76 of file Pythia8Generator.h.

hM

Float_t Pythia8Generator::hM[1]

protected

Definition at line 75 of file Pythia8Generator.h.

hpx

Float_t Pythia8Generator::hpx[1]

protected

Definition at line 75 of file Pythia8Generator.h.

hpy

Float_t Pythia8Generator::hpy[1]

protected

Definition at line 75 of file Pythia8Generator.h.

hpz

Float_t Pythia8Generator::hpz[1]

protected

Definition at line 75 of file Pythia8Generator.h.

maxCrossSection

Double_t Pythia8Generator::maxCrossSection

protected

Definition at line 96 of file Pythia8Generator.h.

mE

Float_t Pythia8Generator::mE[1]

protected

Definition at line 75 of file Pythia8Generator.h.

mid

Float_t Pythia8Generator::mid[1]

protected

Definition at line 76 of file Pythia8Generator.h.

mparam

Double_t Pythia8Generator::mparam[10]

protected

Definition at line 92 of file Pythia8Generator.h.

mpx

Float_t Pythia8Generator::mpx[1]

protected

Definition at line 75 of file Pythia8Generator.h.

mpy

Float_t Pythia8Generator::mpy[1]

protected

Definition at line 75 of file Pythia8Generator.h.

mpz

Float_t Pythia8Generator::mpz[1]

protected

Definition at line 75 of file Pythia8Generator.h.

start

Double_t Pythia8Generator::start[3]

protected

Definition at line 89 of file Pythia8Generator.h.

startZ

Double_t Pythia8Generator::startZ

protected

Definition at line 93 of file Pythia8Generator.h.

subprocCodes

Float_t Pythia8Generator::subprocCodes[16]

protected

Definition at line 77 of file Pythia8Generator.h.

targetFromGeometry

Bool_t Pythia8Generator::targetFromGeometry

protected

Definition at line 95 of file Pythia8Generator.h.

targetName

TString Pythia8Generator::targetName

protected

Definition at line 86 of file Pythia8Generator.h.

xOff

Double_t Pythia8Generator::xOff

protected

Definition at line 87 of file Pythia8Generator.h.

yOff

Double_t Pythia8Generator::yOff

protected

Definition at line 88 of file Pythia8Generator.h.

---

The documentation for this class was generated from the following files:

• shipgen/Pythia8Generator.h
• shipgen/Pythia8Generator.cxx