Enable usage of Graniitti MC generator in O2 simulations

MC/config/PWGUD/external/generator/GeneratorGraniitti.C

@@ -0,0 +1,104 @@
+namespace o2 {
+namespace eventgen {
+class GeneratorGraniitti_class : public Generator {
+public:
+  GeneratorGraniitti_class() { };
+  ~GeneratorGraniitti_class() = default;
+  bool setJsonFile(std::string fname) {
+    jsonFile = fname;
+    // check if jsonFile exists
+    if (gSystem->AccessPathName(jsonFile.c_str())) {
+      return false;
+    }
+
+    return setHepMCFile();
+  }
+
+  bool setHepMCFile() {
+    // item "OUTPUT" defines the hepmcFile
+    // find
+    //  "OUTPUT" : hepmcFile
+    std::string cmd = "grep \"OUTPUT\" "+jsonFile;
+    auto res = gSystem->GetFromPipe(cmd.c_str());
+    auto lines = res.Tokenize("\n");
+    if (lines->GetEntries() != 1) {
+      return false;
+    }
+
+    auto parts = ((TObjString*)lines->At(0))->GetString().Tokenize(":");
+    if (parts->GetEntries() != 2) {
+      return false;
+    }
+
+    auto fname = ((TObjString*)parts->At(1))->GetString();
+    hepmcFile =(std::string)fname.ReplaceAll("\"", "").ReplaceAll(",", "") + ".hepmc3";
+    return true;
+  }
+
+  bool createHepMCFile() {
+    // run graniitti with
+    //  jsonFile as input
+    auto cmd = "$Graniitti_ROOT/bin/gr -i " + jsonFile;
+    std::cout << "Generating events ...";
+    auto res = gSystem->GetFromPipe(cmd.c_str());
+    std::cout << "done!\n";
+
+    // check res to be ok
+
+    return true;
+  }
+
+  void openHepMCFile() {
+    reader = new o2::eventgen::GeneratorHepMC();
+    reader->setFileNames(hepmcFile);
+    if (!reader->Init())
+    {
+      std::cout << "GenerateEvent: Graniitti class/object not properly initialized"
+                << std::endl;
+    }
+  };
+
+  bool generateEvent() override {
+    return reader->generateEvent();
+  };
+
+  bool importParticles() override {
+    mParticles.clear();
+    if (!reader->importParticles()) {
+      return false;
+    }
+    printParticles();
+
+    return true;
+  };
+
+  void printParticles()
+  {
+    std::cout << "\n\n";
+    for (auto& particle : reader->getParticles())
+      particle.Print();
+  }
+
+ private:
+  o2::eventgen::GeneratorHepMC *reader = 0x0;
+  std::string jsonFile;
+  std::string hepmcFile;
+
+};
+
+} // namespace eventgen
+} // namespace o2
+
+FairGenerator* GeneratorGraniitti(std::string jsonFile) {
+
+  // create generator
+  auto gen = new o2::eventgen::GeneratorGraniitti_class();
+  if (gen->setJsonFile(jsonFile)) {
+    if (gen->createHepMCFile()) {
+      // preparing reader
+      gen->openHepMCFile();
+    }
+  }
+
+  return gen;
+}

MC/config/PWGUD/external/generator/GeneratorUpcgen.C

@@ -1,5 +1,5 @@
 R__LOAD_LIBRARY(libUpcgenlib.so)
-R__ADD_INCLUDE_PATH($upcgen_ROOT/include)
+R__ADD_INCLUDE_PATH($Upcgen_ROOT/include)
 
 #include "UpcGenerator.h"
 

MC/config/PWGUD/ini/makeGraniittiConfig.py

@@ -0,0 +1,132 @@
+#! /usr/bin/env python3
+
+### @author: Paul Buehler
+### @email: [email protected]
+
+import argparse
+import os
+import subprocess
+
+def createJson(args):
+  templateFile = os.getenv("O2DPG_ROOT")+"/MC/config/PWGUD/templates/ALICE_Graniitti.temp"
+  jsonFile = "ALICE_Graniitti.json"
+  processes = {
+    "kCon_pipi" : {
+      "OUTPUT" : "ALICE_Con_pipi",
+      "ENERGY" : 13600,
+      "PROCESS" : "PP[RES+CON]<C> -> pi+ pi-",
+      "RES" : ""
+    },
+    "kConRes_pipi" : {
+      "OUTPUT" : "ALICE_Con_pipi",
+      "ENERGY" : 13600,
+      "PROCESS" : "PP[RES+CON]<C> -> pi+ pi-",
+      "RES" : '["f0_500", "rho_770", "f0_980", "phi_1020", "f2_1270", "f0_1500", "f2_1525", "f0_1710", "f2_2150"]'
+    },
+    "kCon_KK" : {
+      "OUTPUT" : "ALICE_Con_pipi",
+      "ENERGY" : 13600,
+      "PROCESS" : "PP[RES+CON]<C> -> pi+ pi-",
+      "RES" : ""
+    },
+    "kConRes_KK" : {
+      "OUTPUT" : "ALICE_Con_pipi",
+      "ENERGY" : 13600,
+      "PROCESS" : "PP[RES+CON]<C> -> pi+ pi-",
+      "RES" : '["f0_500", "rho_770", "f0_980", "phi_1020", "f2_1270", "f0_1500", "f2_1525", "f0_1710", "f2_2150"]'
+    }
+  }
+
+  # is process defined?
+  if not args.process in processes.keys():
+    print("FATAL ERROR: ")
+    print("  Process ", args.process)
+    print("  is not defined!")
+    exit()
+  procdefs = processes[args.process]
+
+  # copy templateFile to jsonFile
+  cmd = "cp "+templateFile+" "+jsonFile
+  if subprocess.call(cmd, shell=True) > 0:
+    print("FATAL ERROR: ")
+    print("  ", templateFile)
+    print("  can not be copied to")
+    print("  ", jsonFile)
+    exit()
+
+  # update jsonFile
+  stat = 0
+  # OUTPUT
+  nl = '    "OUTPUT"     : "' + procdefs["OUTPUT"] + '",'
+  cmd = "sed -i '/\"OUTPUT\"/c\\" + nl + "' " + jsonFile
+  stat = stat + subprocess.call(cmd, shell=True)
+  # NEVENTS
+  nl = '    "NEVENTS"    : ' + args.nEvents + ','
+  cmd = "sed -i '/\"NEVENTS\"/c\\" + nl + "' " + jsonFile
+  stat = stat + subprocess.call(cmd, shell=True)
+  # ENERGY
+  beamEne = str(int(args.eCM)/2)
+  nl = '    "ENERGY"  : [' + beamEne + ', ' + beamEne + '],'
+  cmd = "sed -i '/\"ENERGY\"/c\\" + nl + "' " + jsonFile
+  stat = stat + subprocess.call(cmd, shell=True)
+  # PROCESS
+  nl = '    "PROCESS" : "' + procdefs["PROCESS"] + '",'
+  cmd = "sed -i '/\"PROCESS\"/c\\" + nl + "' " + jsonFile
+  stat = stat + subprocess.call(cmd, shell=True)
+  # RES
+  if procdefs["RES"] == "":
+    nl = '    "//RES"   :'
+  else:
+    nl = '    "RES"     : ' + procdefs["RES"] + ','
+  cmd = "sed -i '/\"RES\"/c\\" + nl + "' " + jsonFile
+  stat = stat + subprocess.call(cmd, shell=True)
+
+  return jsonFile
+
+# main
+
+parser = argparse.ArgumentParser(description='Make Graniitti configuration',
+                                 formatter_class=argparse.ArgumentDefaultsHelpFormatter)
+
+parser.add_argument('--process',default=None, choices=['kCon_pipi', 'kConRes_pipi', 'kCon_KK', 'kConRes_KK'],
+                    help='Process to switch on')
+
+parser.add_argument('--nEvents', default='100',
+                    help='Number of events to generate per TF')
+
+parser.add_argument('--eCM', type=float, default='13600',
+                    help='Centre-of-mass energy')
+
+parser.add_argument('--rapidity', default='cent', choices=['cent_eta', 'muon_eta'],
+                    help='Rapidity to select')
+
+parser.add_argument('--output', default='GenGraniitti.ini',
+                    help='Where to write the configuration')
+
+args = parser.parse_args()
+
+### prepare the json configuration file for graniitti
+jsonFile = createJson(args)
+
+### open output file
+fout = open(args.output, 'w')
+
+### Generator
+fout.write('[GeneratorExternal] \n')
+fout.write('fileName = ${O2DPG_ROOT}/MC/config/PWGUD/external/generator/GeneratorGraniitti.C \n')
+fout.write('funcName = GeneratorGraniitti("%s")  \n' % ("../"+jsonFile))
+
+###Trigger
+fout.write('[TriggerExternal] \n')
+fout.write('fileName = ${O2DPG_ROOT}/MC/config/PWGUD/trigger/selectParticlesInAcceptance.C \n')
+if args.rapidity == 'cent_rap':
+    fout.write('funcName = selectMotherPartInAcc(-0.9,0.9) \n')
+if args.rapidity == 'muon_rap':
+    fout.write('funcName = selectMotherPartInAcc(-4.0,-2.5) \n')
+if args.rapidity == 'cent_eta':
+    fout.write('funcName = selectDaughterPartInAcc(-0.95,0.95) \n')
+if args.rapidity == 'muon_eta':
+    fout.write('funcName = selectDaughterPartInAcc(-4.05,-2.45) \n')
+
+### close outout file
+fout.close()

MC/config/PWGUD/templates/ALICE_Graniitti.temp

@@ -0,0 +1,123 @@
+// JSON(5) + single line/multiline comments allowed, trailing comma allowed
+//
+// The format is strictly upper/lower case sensitive.
+//
+// [email protected], 2021
+{
+
+  // ----------------------------------------------------------------------
+  // GENERAL setup
+
+  "GENERALPARAM" : {
+
+    "OUTPUT"     :
+    "FORMAT"     : "hepmc3",       // hepmc3, hepmc2, hepevt
+    "CORES"      : 0,              // Number of CPU threads (0 for automatic)
+    "NEVENTS"    :
+    "INTEGRATOR" : "VEGAS",        // Integrator (VEGAS, FLAT)
+    "WEIGHTED"   : false,          // Weighted events (default false)
+    "MODELPARAM" : "TUNE0"         // General model tune
+  },
+
+  // ----------------------------------------------------------------------
+  // Process setup
+
+  "PROCESSPARAM" : {
+
+    "BEAM"    : ["p+","p+"],     // Beam PDG-ID / Name
+    "ENERGY"  :
+    "PROCESS" :
+    "RES"     :
+    "POMLOOP" : false,           // Eikonal Pomeron loop screening
+    "NSTARS"  : 0,               // N* excitation (0 = elastic, 1 = single, 2 = double)
+    "LHAPDF"  : "CT10nlo",       // LHAPDF parton distribution set
+    "HIST"    : 0,               // On-the-flight histograms (0 = off,1,2)
+    "RNDSEED" : 0,               // Random seed (uint)
+  },
+
+  // ----------------------------------------------------------------------
+  // Monte Carlo integrator setup
+
+  "INTEGRALPARAM" : {
+
+    // Loop integration
+    "POMLOOP" : {
+      "ND" : 0,               // Integral discretization [larger -> more discrete] (0 gives default discretization) (int)
+    },
+
+    // VEGAS integrator
+    "VEGAS" : {
+
+      "BINS" : 512,         // Maximum number of bins per dimension (NOTE: EVEN NUMBER)
+      "LAMBDA" : 1.25,      // Regularization parameter
+      "NCALL" : 40000,      // Number of function calls per iteration (5000 .. 50000 or more)
+      "ITER"  : 50,         // Number of initialization iteration (10 or more)
+      "CHI2MAX" : 10.0,     // Maximum Chi^2 in initialization
+      "PRECISION"  : 0.05,  // Integral relative precision target
+
+      "DEBUG" : -1          // Debug output (default -1)
+    },
+
+    // FLAT integrator
+    "FLAT" : {
+
+      "PRECISION"  : 0.01,      // Integral relative precision
+      "MIN_EVENTS" : 1000000    // Minimum number of events to be sampled
+    }
+  },
+
+
+  // ----------------------------------------------------------------------
+  // Generation cuts
+
+  "GENCUTS" : {
+
+    "<C>" : {
+      "Rap" : [-1.0, 1.0] // Rapidity boundaries of the final states (<C> class)
+    },
+    "<F>" : {
+      "Rap" : [-2.0, 2.0], // Rapidity boundaries of the system (<F> class)
+      "M"   : [ 0.0, 5.0]
+    },
+    "<Q>" : {
+      "Xi"  : [0.0, 0.05] // Invariant scale M^2/s (forward excitation)
+    }
+  },
+
+
+  // ----------------------------------------------------------------------
+  // Central system fiducial cuts
+  // All central system particles need to fullfill the conditions
+
+  "FIDCUTS" : {
+
+    "active" : false,
+
+     // Central system final states
+    "PARTICLE" : {
+
+      "Eta"  : [-1.0, 1.0],       // Pseudorapidity
+      "Rap"  : [-10.0, 10.0],     // Rapidity
+      "Pt"   : [0.1, 100000.0],   // Transverse momentum
+      "Et"   : [0.0, 100000.0]    // Transverse energy
+    },
+
+    // Central system
+    "SYSTEM" : {
+
+      "M"    : [0.0, 100000.0], // Mass
+      "Rap"  : [-10.0, 10.0],   // Rapidity
+      "Pt"   : [0.0, 100000.0]  // Transverse momentum
+    },
+
+    // Forward system
+    "FORWARD" : {
+      "M"    : [0.0, 100000.0], // Mass
+      "t"    : [0.0, 100000.0], // Forward protons/N* Mandelstam t (absolute value) (GeV^2)
+    },
+
+    // Custom user cuts (default false, otherwise identifier ID)
+    "USERCUTS" : false
+  }
+
+}

MC/run/PWGUD/runGraniittiANCHOR.sh

@@ -0,0 +1,30 @@
+# Run as:  ${O2DPG_ROOT}/GRID/utils/grid_submit.sh --script ./runGraniittiANCHOR.sh --jobname SLtest --outputspec "*.log@disk=1","*.root@disk=2" --packagespec "VO_ALICE@O2sim::v20240626-1" --wait --fetch-output --asuser pbuhler --local
+
+export ALIEN_JDL_LPMANCHORPASSNAME=apass4
+export ALIEN_JDL_MCANCHOR=apass4
+export ALIEN_JDL_COLLISIONSYSTEM=pp
+export ALIEN_JDL_CPULIMIT=8
+export ALIEN_JDL_LPMPASSNAME=apass4
+export ALIEN_JDL_LPMRUNNUMBER=535084
+export ALIEN_JDL_LPMPRODUCTIONTYPE=MC
+export ALIEN_JDL_LPMINTERACTIONTYPE=PbPb
+export ALIEN_JDL_LPMPRODUCTIONTAG=MyPass2Test
+export ALIEN_JDL_LPMANCHORRUN=535084
+export ALIEN_JDL_LPMANCHORPRODUCTION=LHC23f
+export ALIEN_JDL_LPMANCHORYEAR=2023
+
+export NTIMEFRAMES=2
+export NSIGEVENTS=100
+export NBKGEVENTS=1
+export SPLITID=20
+export PRODSPLIT=100
+export CYCLE=30
+export ALIEN_PROC_ID=2963436952
+
+
+#export ALIEN_JDL_ANCHOR_SIM_OPTIONS="-gen external -ini ${PWD}/GenGraniitti.ini --embedding -nb ${NBKGEVENTS} -colBkg PbPb -genBkg pythia8 -procBkg heavy_ion"
+
+${O2DPG_ROOT}/MC/config/PWGUD/ini/makeGraniittiConfig.py --process kConRes_pipi --eCM 13600 --nEvents 300 --rapidity cent_eta
+export ALIEN_JDL_ANCHOR_SIM_OPTIONS="-gen external -ini ${PWD}/GenGraniitti.ini"
+
+${O2DPG_ROOT}/MC/run/ANCHOR/anchorMC.sh