main201
Back to index.
// main201.cc is a part of the PYTHIA event generator.
// Copyright (C) 2021 Torbjorn Sjostrand.
// PYTHIA is licenced under the GNU GPL v2 or later, see COPYING for details.
// Please respect the MCnet Guidelines, see GUIDELINES for details.
// This is a simple test program to compare Pythia and Vincia on
// inclusive jet rates at the LHC, for a sample with pThat > 100 GeV.
// Also illustrates simple use of OpenMP (if enabled) to run two instances
// of Pythia in parallel, here initialised for Pythia and Vincia shower
// models respectively.
// Authors:
// Peter Skands
// Keywords:
// Vincia
// Dire
// OpenMP
#include "Pythia8/Pythia.h"
using namespace Pythia8;
int main() {
// Common parameters used for both runs
const int nEvent = 1000;
const int nListJets = 5;
//************************************************************************
// Histograms.
Hist nJetsModel1("Model1 number of jets", 20, -0.5, 19.5);
Hist eTjetsModel1("Model1 pT for jets", 50, 0., 500.);
Hist yJetsModel1("Model1 y for jets", 20, -4., 4.);
Hist phiJetsModel1("Model1 phi for jets", 25, -M_PI, M_PI);
Hist distJetsModel1("Model1 R distance between jets", 50, 0., 10.);
Hist nJetsModel2("Model2 number of jets", 20, -0.5, 19.5);
Hist eTjetsModel2("Model2 pT for jets", 50, 0., 500.);
Hist yJetsModel2("Model2 y for jets", 20, -4., 4.);
Hist phiJetsModel2("Model2 phi for jets", 25, -M_PI, M_PI);
Hist distJetsModel2("Model2 R distance between jets", 50, 0., 10.);
Hist nJetsRatio("Model2/Model1 number of jets", 20, -0.5, 19.5);
Hist eTjetsRatio("Model2/Model1 pT for jets", 50, 0., 500.);
Hist yJetsRatio("Model2/Model1 y for jets", 20, -4., 4.);
Hist phiJetsRatio("Model2/Model1 phi for jets", 25, -M_PI, M_PI);
Hist distJetsRatio("Model2/Model1 R distance between jets", 50, 0., 10.);
// Loop over generators. Use OpenMP parallelisation if enabled.
#ifdef OPENMP
#pragma omp parallel for
#endif
for (int iRun=1; iRun<=2; ++iRun) {
Pythia pythia;
// Settings common to both runs
pythia.readString("Beams:eCM = 7000.");
pythia.readString("HardQCD:all = on");
pythia.readString("PhaseSpace:pTHatMin = 100.");
pythia.readString("Next:numberShowInfo = 0");
pythia.readString("Next:numberShowProcess = 1");
pythia.readString("Next:numberShowEvent = 1");
pythia.readString("PartonLevel:MPI = on");
pythia.readString("HadronLevel:all = on");
// Settings specific to second run
if (iRun == 2) {
// Switch to VINCIA shower model
pythia.readString("PartonShowers:Model = 2");
}
// Initialise generator for this run
if(!pythia.init()) {continue;}
// Set histogram pointers
Hist* nJetsPtr = &nJetsModel1;
Hist* eTjetsPtr = &eTjetsModel1;
Hist* yJetsPtr = &yJetsModel1;
Hist* phiJetsPtr = &phiJetsModel1;
Hist* distJetsPtr = &distJetsModel1;
// Switch to Model2 histograms for second
if (iRun == 2) {
nJetsPtr = &nJetsModel2;
eTjetsPtr = &eTjetsModel2;
yJetsPtr = &yJetsModel2;
phiJetsPtr = &phiJetsModel2;
distJetsPtr = &distJetsModel2;
}
// Set up SlowJet jet finder, with anti-kT clustering, R = 0.7,
// pT > 20 GeV, |eta| < 4, and pion mass assumed for non-photons
double etaMax = 4.;
double radius = 0.7;
double pTjetMin = 20.;
// Exclude neutrinos (and other invisible) from study.
int nSel = 2;
SlowJet slowJet( -1, radius, pTjetMin, etaMax, nSel, 1);
// Begin event loop.
double sumWeights = 0.;
for (int iEvent = 0; iEvent < nEvent; ++iEvent) {
// Generate event.
if (!pythia.next()) {
cout << " Event generation aborted prematurely, owing to error!\n";
cout << " Run number was : "<