diff --git a/Generators/SFGen_i/CMakeLists.txt b/Generators/SFGen_i/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..0afcfdfd563c7c14813f1d555696f3f80f8a4b22
--- /dev/null
+++ b/Generators/SFGen_i/CMakeLists.txt
@@ -0,0 +1,19 @@
+# Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+
+# Declare the package name:
+atlas_subdir( SFGen_i )
+
+# External dependencies:
+
+find_package( Apfel )
+find_package( SFGen )
+
+# Install files from the package:
+atlas_install_python_modules( python/*.py )
+atlas_install_joboptions( share/*.py share/common/*.py )
+
+
+# Set SFGen specific environment variable(s).
+set( SFGenEnvironment_DIR ${CMAKE_CURRENT_SOURCE_DIR}
+ CACHE PATH "Location of SFGenEnvironment.cmake" )
+find_package( SFGenEnvironment )
diff --git a/Generators/SFGen_i/README.md b/Generators/SFGen_i/README.md
new file mode 100644
index 0000000000000000000000000000000000000000..c9c9721e60df499e760cef78b38297e805fbc79d
--- /dev/null
+++ b/Generators/SFGen_i/README.md
@@ -0,0 +1,72 @@
+[[_TOC_]]
+
+The *Athena SFGen interface* described here was developed by *Lorenzo Primomo*, who works for *ATLAS*. For interface questions please contact *Giancarlo Panizzo* or *Dominic Hirschbuehl*. For generator-specific questions, you might want to contact the *SFGen* author *Lucian Harland-Lang*.
+
+List of responsibles:
+
+* Lorenzo Primomo (lprimomo@cern.ch);
+* Giancarlo Panizzo (giancarlo.panizzo@uniud.it);
+* Dominic Hirschbuehl (dominic.hirschbuehl@cern.ch);
+* Lucian Harland-Lang (l.harland-lang@ucl.ac.uk));
+
+# Introduction
+
+*SFGen1.03* ([ref](https://sfgen.hepforge.org/) and [manual](https://sfgen.hepforge.org/SFGen1.03.pdf)) is a Monte Carlo generator which provides a publicly available tool for lepton pair production and lepton-lepton scattering within the *Structure Function approach*, including initial-state γ, Z and mixed Z/+q contributions. Arbitrary distributions and unweighted events can be generated, and errors due to the experiment uncertainty on the structure function (the equivalent of *PDF* uncertainties in the photon *PDF* framework) can be calculated on-the-fly. The LO collinear photon-initiated predictions can also be calculated for comparison, for arbitrary factorization and renormalization scales.
+
+The standalone executables of a more recent parched version can be found at:
+
+* `/cvmfs/sft.cern.ch/lcg/releases/LCG_88b/MCGenerators/SFGen/1.03p0/x86_64-centos7-gcc62-opt`;
+* `/cvmfs/sft.cern.ch/lcg/releases/LCG_104d_ATLAS_4/MCGenerators/SFGen/1.03.atlas2/`;
+* `/cvmfs/sft.cern.ch/lcg/releases/LCG_104d_ATLAS_3/MCGenerators/SFGen/1.03.atlas2/`;
+
+It will be available in *ATLAS cvsmfs* soon.
+
+## Intro to standalone SFGen
+
+When running independently, *SFGen* requires:
+1) Folders named `outputs/` and `evrecs/` in the directory that you want to run from;
+2) An input configuration file, named *input.DAT*. The format must follow that found in the directories with the standalone executables referenced above;
+3) The packages *APFEL* and *LHAPDF* must also be in the `ROOT_INCLUDE_PATH` and `LD_LIBRARY_PATH` of your environment. This is taken care of if you source a recent AthGeneration 21.6.106 release or the more recent 23.6.24 release;
+
+-----
+
+Standalone SFGen is run in only one step:
+
+The main MC code may be run with "./SFGen < input.DAT", using the same *input.DAT* file. This step generates an *LHE/HepMC/HEP-EVT* in the `evrecs/` folder that contains all of the generated events' information. This *LHE* file can then be fed into *Pythia* for showering.
+
+-----
+
+# SFGen Interface
+
+The new *Athena SFGen interface* proceeds primarily as standalone *SFGen* would. It proceeds fairly linearly:
+1) The necessary `outputs/` and `evrecs/` folders are created if they don't already exist;
+2) An *input.DAT* file is generated using arguments that are given to "Init". The class type for "Init" is defined here: `/athena/Generators/SFGen_i/python/SFGenUtils.py`. You can see all of the various input variables that can be given. The file is written with the "writeInputDAT" function in that python macro;
+3) "SFGenRun" also runs the SFGen executable, creating an *LHE* file. It is important that the output format be *LHE*, because the next step, which creates an *EVNT* file, is only compatible with the *LHE* format;
+4) The *EVNT* file is filled using the `/athena/Generators/SFGen_i/python/EventFiller.py` macro. To use a format other than *LHE*, a different version of this macro would have to be written;
+
+An example *jobOptions* file with *Pythia* showering can be found at `/athena/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_SDA.py`. The key difference here is of course that *Pythia8* is added to the "genSeq", and that we have to be careful about what *Pythia* command options are included. The `/athena/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_DD.py`, `/athena/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_SDA.py`, `/athena/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_SDB.py`, and `/athena/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_El.py` contain the latest recommended options for the case of photon-induced di-lepton for the cases of double-dissociated production, single dissociative production on one side and the other, and elastic production, respectively.
+
+## Some Notes
+
+* The "Init.diff" should be either 'el', 'sda', 'sdb', or 'dd'. Be careful to make sure that the "BeamRemnants:unresolvedHadron = 0" component of the *Pythia* commands is adjusted accordingly below (el -> 3, sda -> 2, sdb ->1, and dd -> 0).
+* The list of processes for "Init.proc" can be found in *Section 7* of the [manual](https://sfgen.hepforge.org/SFGen1.03.pdf).
+* Output kinematics can of course be changed from what is in the samples. Two body decays use "a" and "b" particles (as in "pta" and "ptb").
+ * Cuts that have "x" in them refer to the hard process as a whole, so e.g. the di-lepton system in a di-lepton event.
+* When running without showering with *Pythia*, an important drawback of using *LHE* files should be noted: they do not use high enough precision. What I mean is the following: In `EventFiller.py`, you can see that you enter in a four vector that contains every particle's px, py, pz, and E. The mass of the particles is then computed based on those values. Unfortunately, for highly boosted objects (i.e. those where E >> m) the calculation of four-vector mass with the precision given in the *LHE* files is slightly off. So in the output *EVNT* files, you will see protons with masses in the range of 900 - 1000 MeV. I have set the "generated_mass" of each particle to the right mass, but if you retrieve "m_m" for a truth-level particle, it will be slightly incorrect.
+
+# Running locally
+
+To run locally:
+1) Make a folder to hold the *Athena release*, go inside it, and run: asetup 21.6,latest,AthGeneration (21.6.106 or any later release should do);
+2) (Optional) Make a run folder and go inside it;
+3) mkdir 999999;
+4) Put a *jobOptions file* in `999999/`;
+5) Run with e.g.:
+
+ ```
+ Gen_tf.py --ecmEnergy=13000.0 --maxEvents=1000 --firstEvent=1 --randomSeed=13 --outputEVNTFile=test.pool.root --jobConfig=999999
+ ```
+
+# Releases
+
+The *SFGen* Athena interface will be soon incorporated into *AthGeneration release 23.6.25*.
diff --git a/Generators/SFGen_i/SFGenEnvironmentConfig.cmake b/Generators/SFGen_i/SFGenEnvironmentConfig.cmake
new file mode 100644
index 0000000000000000000000000000000000000000..3f9d0f31a284d9146de0435a5a54bc1dd89291d8
--- /dev/null
+++ b/Generators/SFGen_i/SFGenEnvironmentConfig.cmake
@@ -0,0 +1,24 @@
+# Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
+#
+# This module is used to set up the environment for SFGen
+#
+#
+
+# Set the environment variable(s):
+
+find_package( SFGen )
+find_package( Lhapdf )
+
+if( SFGEN_FOUND )
+
+ set( SFGENENVIRONMENT_ENVIRONMENT
+ FORCESET SFGENVER ${SFGEN_LCGVERSION}
+ FORCESET SFGENPATH ${SFGEN_LCGROOT}
+ APPEND LHAPATH ${LHAPDF_DATA_PATH}
+ APPEND LHAPDF_DATA_PATH ${LHAPDF_DATA_PATH} )
+
+endif()
+
+
+# Silently declare the module found:
+set( SFGENENVIRONMENT_FOUND TRUE )
diff --git a/Generators/SFGen_i/doc/SFGen.md b/Generators/SFGen_i/doc/SFGen.md
new file mode 100644
index 0000000000000000000000000000000000000000..3b089e038896b77fd8a7c93f7588ecf4b163b96f
--- /dev/null
+++ b/Generators/SFGen_i/doc/SFGen.md
@@ -0,0 +1,102 @@
+[[_TOC_]]
+
+# SFGen_i: An interface between SFGen Monte Carlo event generator for Athena version in release 21.6
+
+Author:
+
+* Lorenzo Primomo (lprimomo@cern.ch).
+
+# Introduction
+
+SFGen1.03 Monte Carlo event generator for central exclusive production.
+
+This interface package runs SFGen Monte Carlo Generator (v1.03) within
+Athena framework, and stores the events output into the transient store
+in HepMC format.
+
+SFGen Monte Carlo event generator is dedicated for lepton pair
+production and lepton-lepton scattering within the Structure Function
+approach, including initial-state $\gamma$, $Z$ and mixed $Z/\gamma + q$
+contributions. This documentation gives you some details about setting
+up the SFGen interface (which is prepared using the latest available
+SFGen version 1.03) and activating the module using the JobOptions.
+
+The code makes use of structure function grids in the LHAPDF format, as
+calculated using APFEL. These structure function grids are provided with
+the MC, and are calculated with MSHT20qe_nnlo PDFs in the high $Q^2$
+region.
+
+# Usage of Job Options
+
+The example job option file can be found under the following path.
+
+**SFGen_i/share/mc.SFGen.py**
+
+Before running this script, first setup the Athena environment:
+
+**asetup 23.6,latest,AthGeneration,slc6**
+
+Above command sets up the latest AthGeneration cache.
+
+The SFGen input parameters can be set from the job options service. The
+default parameters initializes SFGen for proton-proton beams at
+center-of-mass energy of 13 TeV for the scattering process
+\[$\gamma \gamma \rightarrow \mu^+ \mu^-$ (process n. 2)\].
+
+**All the parameters passed to SFGen are in the units specified in the
+SFGen manual https://sfgen.hepforge.org/SFGen1.03.pdf**
+
+The default mc.SFGen.py file can be copied to your test run directory.
+The initialization parameters can be changed via the following line in
+the python file.
+
+ SFGen.Initialize = ["parameter_1 value_1", "parameter_2 value_2"]
+
+Each quoted string sets one parameter value in the fortran variable
+format. You can set all the input parameters separated by commas,
+however, the important ones are listed below.
+
+**parameter_1:** must be one of the following variable names, an error
+message is returned if the specified variable is not in the input
+parameter list.
+
+**value_1:** is the value of the input parameter.\
+JO Example:\
+The following command generates 10 events for proton-proton collisions
+at 13 TeV center-of-mass energy along with important input parameters
+for process 2 i.e. $\gamma \gamma \rightarrow \mu^+ \mu^-$.\
+Running the Job Option to produce events:
+
+ Gen_tf.py --ecmEnergy=13000.0 --maxEvents=10 --firstEvent=1 --randomSeed=13 --outputEVNTFile=test.pool.root --jobConfig=mc.SFGen.py
+
+Example of the initialization of the input parameters in the JO is shown
+below.
+
+ genSeq.SFGen.Initialize = \
+ ["rts 13d3", # set the COM collision energy (in fortran syntax)
+ "diff 'tot'", # elastic ('el'), single/double dissociation ('sd'/'dd')
+ "PDFname 'SF_MSHT20qed_nnlo'", # PDF set name
+ "PDFmember 0", # PDF member
+ "proc 2", # Process number (1:gg->ee,2:gg->mumu,3:gg->tautau)
+ "outtag 'out'", # for output file name
+ "SFerror .false." # Include error from SF input - increases run time
+ ]
+
+
+One can also add the kinematical cuts to the outgoing particles or the
+system, details can be found in the manual.
+
+# Running SFGen in Standalone way
+
+One can directly carry out following commands in the fresh terminal to
+produce SFGen events in standalone way.
+
+ source /cvmfs/sft.cern.ch/lcg/releases/LCG_88b/MCGenerators/SFGen/1.03p0/x86_64-centos7-gcc62-opt/sfgenenv-genser.sh
+ cp -rf /cvmfs/sft.cern.ch/lcg/releases/LCG_88b/MCGenerators/SFGen/1.03p0/x86_64-centos7-gcc62-opt/bin/tmp/cern_username/
+ cd /tmp/cern_username/bin/
+ export LHAPATH=/cvmfs/sft.cern.ch/lcg/external/lhapdfsets/current/
+ ./SFGen < input.DAT
+
+You can change the input.DAT file to change the C.O.M collision energy
+for the particular process, process no., number of events etc. Above
+commands will produce output stored under the directory **'evrecs'** .
diff --git a/Generators/SFGen_i/python/EventFiller.py b/Generators/SFGen_i/python/EventFiller.py
new file mode 100644
index 0000000000000000000000000000000000000000..ce1c5f7ac48026866c4ef84a173ec815d9e96466
--- /dev/null
+++ b/Generators/SFGen_i/python/EventFiller.py
@@ -0,0 +1,134 @@
+# Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
+
+# This module fills the EVNT file using the information from the LHE file generated by SFGen in the 'evrecs/' directory
+
+from GeneratorModules.EvgenAlg import EvgenAlg
+from AthenaPython.PyAthena import StatusCode
+try:
+ from AthenaPython.PyAthena import HepMC3 as HepMC
+except ImportError:
+ from AthenaPython.PyAthena import HepMC as HepMC
+
+import os
+import math
+import ROOT
+
+# Getting pz in MeV unit
+
+def init_pz(energy):
+ return math.sqrt(energy**2 - 938.272046**2)
+
+class LheEVNTFiller(EvgenAlg):
+
+# Init definition
+
+ def __init__(self, ecmEnergy, name="LheEVNTFiller"):
+ super(LheEVNTFiller, self).__init__(name=name)
+ self.beamEnergyMeV = ecmEnergy * 1000. / 2.
+
+ fileName = "evrecs/evrecout.dat"
+ outputFileName = "outputs/outputout.dat"
+ eventsProcessed = 0
+
+# Checking out evrec and output directories
+
+ def initialize(self):
+
+ if(os.path.isfile(self.fileName) and os.path.isfile(self.outputFileName)):
+ print(self.fileName)
+ return StatusCode.Success
+ else:
+ return StatusCode.Failure
+
+# Filling the EVNT file
+
+ def fillEvent(self, evt):
+
+ eventsSeen = 0
+ firstLine = True
+
+# Reading cross section in the outputout.dat file in nb unit
+
+ with open(self.outputFileName, 'r') as inputOutputFile:
+ if self.eventsProcessed == 0:
+ for line in inputOutputFile:
+ if 'Cross section =' in line:
+ splitLine = line.split()
+ factor = 1.
+ if(splitLine[-1] == "pb"):
+ factor = 0.001
+ if(splitLine[-1] == "fb"):
+ factor = 0.000001
+ if(splitLine[-1] == "ub"):
+ factor = 1000.
+ if(splitLine[-1] == "mb"):
+ factor = 1000000.
+ print("MetaData: cross-section (nb)= "+str(float(splitLine[3])*factor))
+
+ break
+
+ with open(self.fileName, 'r') as inputfile:
+ event = False
+ for line in inputfile:
+ if not event and '<event>' not in line:
+ continue
+ # Check if we are just starting an event
+ if not event and '<event>' in line and eventsSeen == self.eventsProcessed:
+ event = True
+ continue
+ # Check if we have finished and are doing something else
+ if '<' in line or '>' in line:
+ event = False
+ eventsSeen += 1
+ firstLine = True
+ continue
+ if event and firstLine:
+ firstLine = False
+ evt.weights().push_back(float(line.split()[2]))
+
+ # Add the initial state protons
+ pos = HepMC.FourVector(0.0, 0.0, 0.0, 0.0)
+ gv = HepMC.GenVertex(pos)
+ ROOT.SetOwnership(gv, False)
+ evt.add_vertex(gv)
+ mom = HepMC.FourVector(0., 0., init_pz(self.beamEnergyMeV), self.beamEnergyMeV)
+ gp = HepMC.GenParticle()
+ gp.set_status( 2 )
+ gp.set_pdg_id( 2212 )
+ gp.set_momentum(mom)
+ gp.set_generated_mass(938.272046)
+ ROOT.SetOwnership(gp, False)
+ gv.add_particle_out(gp)
+
+ pos = HepMC.FourVector(0.0, 0.0, 0.0, 0.0)
+ gv = HepMC.GenVertex(pos)
+ ROOT.SetOwnership(gv, False)
+ evt.add_vertex(gv)
+ mom = HepMC.FourVector(0., 0., -init_pz(self.beamEnergyMeV), self.beamEnergyMeV)
+ gp = HepMC.GenParticle()
+ gp.set_status( 2 )
+ gp.set_pdg_id( 2212 )
+ gp.set_momentum(mom)
+ gp.set_generated_mass(938.272046)
+ ROOT.SetOwnership(gp, False)
+ gv.add_particle_out(gp)
+
+ continue
+ if event:
+ pos = HepMC.FourVector(0.0, 0.0, 0.0, 0.0)
+ gv = HepMC.GenVertex(pos)
+ ROOT.SetOwnership(gv, False)
+ evt.add_vertex(gv)
+ mom = HepMC.FourVector( float(line.split()[6])*1000. , float(line.split()[7])*1000. , float(line.split()[8])*1000. , float(line.split()[9])*1000. )
+ gp = HepMC.GenParticle()
+ gp.set_status(int(line.split()[1]) )
+ gp.set_pdg_id(int(line.split()[0]) )
+ gp.set_momentum(mom)
+ gp.set_generated_mass(float(line.split()[10]) * 1000.)
+ ROOT.SetOwnership(gp, False)
+ gv.add_particle_out(gp)
+
+ self.eventsProcessed += 1
+
+ return StatusCode.Success
+
diff --git a/Generators/SFGen_i/python/LheConverter.py b/Generators/SFGen_i/python/LheConverter.py
new file mode 100644
index 0000000000000000000000000000000000000000..aca5dcc853185e9e7ae6c047f53d442a922681ac
--- /dev/null
+++ b/Generators/SFGen_i/python/LheConverter.py
@@ -0,0 +1,102 @@
+# Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
+
+import os
+import re
+import math
+import xml.dom.minidom
+
+from GeneratorModules.EvgenAlg import EvgenAlg
+from AthenaPython.PyAthena import StatusCode
+
+
+class LheConverter(EvgenAlg):
+ '''
+ Class for converting output LHE file from SFGen to LHE recognised by Pythia.
+ '''
+
+# Init definition
+
+ def __init__(self, name='LheConverter'):
+ super(LheConverter, self).__init__(name=name)
+
+ inFileName = 'evrecs/evrecout.dat'
+ outFileName = 'events.lhe'
+ done = False
+ leptons = ['11', '-11', '13', '-13', '15', '-15'] # e, µ, τ
+
+# Checking out evrec directory
+
+ def initialize(self):
+ if(os.path.isfile(self.inFileName)):
+ print(self.fileName)
+
+ return self.convert()
+ else:
+ return StatusCode.Failure
+
+ def convert(self):
+ '''
+ Converts `evrecs/evrecout.dat` input file to `events.lhe` recognised by the Pythia:
+ - Replaces init block
+ - Changes photon px and py to zero and energy to pz
+ - Recalculates energy scales as an average of lepton pair transverse momentum
+ '''
+
+ if not self.done:
+ DOMTree = xml.dom.minidom.parse(self.inFileName)
+ collection = DOMTree.documentElement
+
+ # Replace init block
+ init = collection.getElementsByTagName('init')
+ init_repl = r'''
+ 13 -13 2.510000e+03 2.510000e+03 0 0 0 0 3 1
+ 1.000000e+00 0.000000e+00 1.000000e+00 9999
+ '''
+ init[0].firstChild.data = init_repl
+
+ # 1 -- energy scale
+ event_header = r'^(\s*\S*\s*\S*\s*\S*\s*)(\S*)(.*)$'
+ # pdgid px py pz e
+ event_particle = r'^(\s*)([0-9-]+)(\s+[0-9-]+\s+[0-9-]+\s+[0-9-]+\s+[0-9-]+\s+[0-9-]+\s+)(\S+)(\s+)(\S+)(\s+)(\S+)(\s+)(\S+)(.*)$'
+
+ events = collection.getElementsByTagName('event')
+ for i, event in enumerate(events):
+ new_particles = []
+ lepton_mom = []
+
+ particles = re.findall(event_particle, event.firstChild.data, re.MULTILINE)
+ for particle in particles:
+ particle = list(particle)
+ if particle[1] == '22': # photon
+ # Zero photon transverse momentum, change energy to be equal pz
+ particle[3] = '0.'
+ particle[5] = '0.'
+ particle[9] = particle[7]
+ new_particles.append(''.join(particle))
+ elif particle[1] in self.leptons:
+ # Read leptons px and py
+ lepton_mom.append(float(particle[3]))
+ lepton_mom.append(float(particle[5]))
+ new_particles.append(''.join(particle))
+ else:
+ new_particles.append(''.join(particle))
+
+ # Calculate new energy scale
+ l1_px, l1_py, l2_px, l2_py = lepton_mom
+ l1_pt = math.sqrt(l1_px**2 + l1_py**2)
+ l2_pt = math.sqrt(l2_px**2 + l2_py**2)
+ energy_scale = '{expression:.9E}'.format(expression=(l1_pt + l2_pt)/2.)
+
+
+ header = re.search(event_header, event.firstChild.data, re.MULTILINE)
+ header = list(header.groups())
+ header[1] = energy_scale
+
+ event.firstChild.data = '\n'.join([''.join(header)] + new_particles) + '\n '
+
+ with open(self.outFileName, 'w') as output:
+ output.write(DOMTree.toxml().replace('<?xml version="1.0" ?>', ' '))
+
+ self.done = True
+
+ return StatusCode.Success
diff --git a/Generators/SFGen_i/python/SFGenUtils.py b/Generators/SFGen_i/python/SFGenUtils.py
new file mode 100644
index 0000000000000000000000000000000000000000..7077a37afe66b3d49d5d852d58566ebe5a68b849
--- /dev/null
+++ b/Generators/SFGen_i/python/SFGenUtils.py
@@ -0,0 +1,257 @@
+# Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
+
+
+import subprocess, os, shlex, re
+
+from AthenaCommon import Logging
+
+## Get handle to Athena logging
+logger = Logging.logging.getLogger("SFGen_i")
+
+# This module defines the SFGen Utilities in the Athena framework
+
+class SFGenConfig:
+
+# Init definition
+
+ def __init__(self, runArgs):
+ self.sfgenpath = os.environ['SFGENPATH']
+
+ #SFGen specific variables for input.DAT, see writeInputDAT function for more elaboration
+
+ self.rts = 13000 #collision energy (GeV)
+ if hasattr(runArgs,"ecmEnergy"):
+ self.rts = runArgs.ecmEnergy
+
+ self.coll = False
+ self.PDFname = 'SF_MSHT20qed_nnlo'
+ self.PDFmember = 0
+ self.proc = 1 #Please consult SFGen Manual https://sfgen.hepforge.org/
+ self.outtag = 'out'
+ self.diff = 'tot' #Interaction: elastic ('el'), single ('sd','sda','sdb') and double ('dd') dissociation.
+ self.SFerror = False
+ self.mixed = False
+ self.Zinit = False
+ self.subt = False
+ self.lep1 = 'mu+'
+ self.lep2 = 'mu+'
+ self.kmu = 1.
+ self.ncall = 50000
+ self.itmx = 10
+ self.prec = 1
+ self.ncall1 = 10000
+ self.inccall = 10000
+ self.itend = 1000
+
+ self.iseed = 13
+ if hasattr(runArgs,"randomSeed"):
+ self.iseed = runArgs.randomSeed
+
+ self.genunw = False
+
+ self.nev = "10"
+ if hasattr(runArgs,"maxEvents"):
+ self.nev = runArgs.maxEvents
+
+ self.erec = 'lhe' #output file type
+ self.ymin = -2.4 # Minimum dilepton rapidity
+ self.ymax = 2.4 # Maximum dilepton rapidity
+ self.mmin = 10. # Minimum dilepton mass
+ self.mmax = 1000. # Maximum dilepton mass
+
+ self.gencuts = True # Generate cuts below
+
+ self.ptamin = 20. # Minimum pT of outgoing object a
+ self.ptbmin = 20. # Minimum pT of outgoing object b
+ self.etaamin = -2.4 # Minimum eta of outgoing object a
+ self.etaamax = 2.4 # Maximum eta of outgoing object a
+ self.etabmin = -2.4 # Minimum eta of outgoing object b
+ self.etabmax = 2.4 # Maximum eta of outgoing object b
+ self.ptllmin = 0. # Minimun pT of dilepton pair
+
+# Writing of the input.DAT file
+
+ def toFortran(self):
+
+ def fortDouble(x):
+ return str(x)+"d0"
+ def fortInt(x):
+ return str(x)
+ def fortBool(x):
+ return '.true.' if x else '.false.'
+ def fortStr(x):
+ return "'{}'".format(x)
+
+ conf = ""
+ conf+="***********************************************************************************\n"
+ conf+="***********************************************************************************\n"
+ conf+="***********************************************************************************\n"
+ conf+=fortDouble(self.rts) + " ! [rts] : CMS collision energy (GeV) \n"
+ conf+="***********************************************************************************\n"
+ conf+="***********************************************************************************\n"
+ conf+=fortInt(self.proc) + " ! [proc] : Process number \n"
+ conf+=fortStr(self.outtag) + " ! [outtag] : for output file \n"
+ conf+=fortStr(self.diff) + " ! [diff] : elastic ('el'), single/double dissociation ('sd'/'dd') \n"
+ conf+=fortBool(self.SFerror) + " ! [SFerror] : Include error from SF input - increases run time \n"
+ conf+=fortBool(self.mixed) + " ! [mixed] : include mixed gam/Z + q diagrams \n"
+ conf+=fortBool(self.Zinit) + " ! [Zinit] : include Z bosons in initial state \n"
+ conf+=fortBool(self.subt) + " ! [subt] : calculate *only* (positive) subtraction term \n"
+ conf+=fortStr(self.lep1) + " ! [lep1] : for lepton-lepton scattering (proc=4) \n"
+ conf+=fortStr(self.lep2) + " ! [lep2] : for lepton-lepton scattering (proc=4) \n"
+ conf+=fortDouble(self.kmu) + " ! [kmu] : = mu(f,r)/mu0 \n"
+ conf+="***********************************************************************************\n"
+ conf+="************************** To run in collinear mode *****************************\n"
+ conf+="***********************************************************************************\n"
+ conf+=fortBool(self.coll) + " ! [coll] : use collinear approach + lhapdf \n"
+ conf+=fortStr(self.PDFname) + " ! [PDFname] : PDF set \n"
+ conf+=fortInt(self.PDFmember) + " ! [PDFmember] : PDF member \n"
+ conf+="***********************************************************************************\n"
+ conf+="*************Integration parameters************************************************\n"
+ conf+="***********************************************************************************\n"
+ conf+=fortInt(self.ncall) + " ! [ncall] : Number of calls for preconditioning \n"
+ conf+=fortInt(self.itmx) + " ! [itmx] : Number of iterations for preconditioning \n"
+ conf+=fortDouble(self.prec) + " ! [prec] : Relative accuracy (in %) in main run \n"
+ conf+=fortInt(self.ncall1) + " ! [ncall1] : Number of calls in first iteration \n"
+ conf+=fortInt(self.inccall) + " ! [inccall] : Number of increase calls per iteration \n"
+ conf+=fortInt(self.itend) + " ! [itend] : Maximum number of iterations \n"
+ conf+=fortInt(self.iseed) + " ! [iseed] : Random number seed (integer > 0) \n"
+ conf+="***********************************************************************************\n"
+ conf+="********************Unweighted events**********************************************\n"
+ conf+="***********************************************************************************\n"
+ conf+=fortBool(self.genunw) + " ! [genunw] : Generate unweighted events \n"
+ conf+=fortInt(int(self.nev)) + " ! [nev] : Number of events (preferably controlled by maxEvents option in Gen_tf command) \n"
+ conf+=fortStr(self.erec) + " ! [erec] : Event record format ('hepmc','lhe','hepevt') \n"
+ conf+="***********************************************************************************\n"
+ conf+="******************* general cuts ************************************************\n"
+ conf+="***********************************************************************************\n"
+ conf+=fortDouble(self.ymin) + " ! [ymin] : Minimum dilepton rapidity \n"
+ conf+=fortDouble(self.ymax) + " ! [ymax] : Maximum dilepton rapidity \n"
+ conf+=fortDouble(self.mmin) + " ! [mmin] : Minimum dilepton mass \n"
+ conf+=fortDouble(self.mmax) + " ! [mmax] : Maximum dilepton mass \n"
+ conf+=fortBool(self.gencuts) + " ! [gencuts] : Generate cuts below \n"
+ conf+="***********************************************************************************\n"
+ conf+="********** 2 body final states : p(a) + p(b) **************************************\n"
+ conf+="***********************************************************************************\n"
+ conf+=fortDouble(self.ptamin) + " ! [ptamin] \n"
+ conf+=fortDouble(self.ptbmin) + " ! [ptbmin] \n"
+ conf+=fortDouble(self.etaamin) + " ! [etaamin] \n"
+ conf+=fortDouble(self.etaamax) + " ! [etaamax] \n"
+ conf+=fortDouble(self.etabmin) + " ! [etabmin] \n"
+ conf+=fortDouble(self.etabmax) + " ! [etabmax] \n"
+ conf+=fortDouble(self.ptllmin) + " ! [ptllmin] \n"
+ conf+="***********************************************************************************\n"
+ conf+="***********************************************************************************\n"
+
+ return conf
+
+ def outputLHEFile(self):
+ return "evrecs/evrec"+self.outtag+".dat"
+
+
+def writeInputDAT(Init):
+
+ with open("input.DAT", "w") as outF:
+ outF.write(Init.toFortran())
+
+ return
+
+
+def run_command(command, stdin = None ):
+ """
+ Run a command and print output continuously
+ """
+ process = subprocess.Popen(shlex.split(command), stdout=subprocess.PIPE, stdin=stdin)
+ while True:
+ output = process.stdout.readline().decode("utf-8")
+ if output == '' and process.poll() is not None:
+ break
+ if output:
+ # remove ANSI escape formatting characters
+ reaesc = re.compile(r'(\x9B|\x1B\[)[0-?]*[ -\/]*[@-~]')
+ text = reaesc.sub('', output.strip())
+ logger.info(text)
+
+ rc = process.poll()
+ return rc
+
+# SFGen initialization
+
+def SFGenInitialize(Init, stdin=None):
+
+ logger.info("Starting SFGen Initialization")
+
+ if not os.path.exists('evrecs'):
+ os.makedirs('evrecs')
+ if not os.path.exists('outputs'):
+ os.makedirs('outputs')
+
+
+ try:
+ inputDAT = open('input.DAT')
+
+ except IOError:
+ raise Exception("Problem with file IO; potentially input.DAT not created correctly")
+ else:
+
+ try:
+ rc = run_command(Init.sfgenpath+"/bin/SFGen", inputDAT)
+
+ except OSError:
+ raise Exception("File not found")
+
+ except Exception:
+ raise Exception("Non-OSError or IOError in execution block")
+
+ if rc:
+ raise Exception('Unexpected error in sfgen execution in SFGenInitialize')
+
+ return
+
+# SFGen execution
+
+def SFGenExecute(Init):
+
+ logger.info("Starting SFGen Itself")
+
+ if not os.path.exists('evrecs'):
+ os.makedirs('evrecs')
+ if not os.path.exists('outputs'):
+ os.makedirs('outputs')
+
+
+ try:
+ inputDAT = open('input.DAT')
+ except IOError:
+ raise Exception ("Problem with IO; potentially input.DAT not created correctly")
+ else:
+
+ try:
+
+ rc = run_command(Init.sfgenpath+'/bin/SFGen', stdin=inputDAT)
+
+ except OSError:
+ raise Exception("SFGen executable or file not found")
+
+ except Exception:
+ raise Exception("Non-OSError or IOError in SFGen execution block")
+
+ if rc:
+ raise Exception('Unexpected error in sfgen execution in SFGenExecute')
+
+ return
+
+# SFGen running
+
+def SFGenRun(Init, genSeq):
+
+ # dump the job configuration for fortran code
+ print(Init.toFortran())
+
+ # attach config to genSequence for later using in showering
+ genSeq.SFGenConfig = Init
+
+ writeInputDAT(Init)
+ SFGenInitialize(Init)
+ SFGenExecute(Init)
+
+ return
diff --git a/Generators/SFGen_i/share/common/LheEventFiller_Common.py b/Generators/SFGen_i/share/common/LheEventFiller_Common.py
new file mode 100644
index 0000000000000000000000000000000000000000..9922aa9472aa2a053f9dc9112ce0195bd5c6f950
--- /dev/null
+++ b/Generators/SFGen_i/share/common/LheEventFiller_Common.py
@@ -0,0 +1,8 @@
+"""
+Configuration for transfering SFGen output from Lhe to ENVT.
+"""
+
+import SFGen_i.EventFiller as EF
+
+ef = EF.LheEVNTFiller(runArgs.ecmEnergy)
+genSeq += ef
diff --git a/Generators/SFGen_i/share/common/Pythia8_Base_Common.py b/Generators/SFGen_i/share/common/Pythia8_Base_Common.py
new file mode 100644
index 0000000000000000000000000000000000000000..a021ac50926525eee0f3dfed6cca336b411f2bea
--- /dev/null
+++ b/Generators/SFGen_i/share/common/Pythia8_Base_Common.py
@@ -0,0 +1,19 @@
+from AthenaCommon import Logging
+logger = Logging.logging.getLogger("SFGen_i")
+
+## Shower config for Pythia8 with SFGen for elastic production
+
+from Pythia8_i.Pythia8_iConf import Pythia8_i
+genSeq += Pythia8_i("Pythia8")
+evgenConfig.generators += ["Pythia8"]
+
+## Control storing LHE in the HepMC record
+if "StoreLHE" in genSeq.Pythia8.__slots__.keys():
+ print("Pythia8_Base_Fragment.py: DISABLING storage of LHE record in HepMC by default. Please re-enable storage if desired")
+ genSeq.Pythia8.StoreLHE = False
+
+genSeq.Pythia8.LHEFile = genSeq.SFGenConfig.outputLHEFile()
+genSeq.Pythia8.CollisionEnergy = int(runArgs.ecmEnergy)
+
+testSeq.TestHepMC.MaxTransVtxDisp = 1000000
+testSeq.TestHepMC.MaxVtxDisp = 1000000000
diff --git a/Generators/SFGen_i/share/common/Pythia8_DD_Common.py b/Generators/SFGen_i/share/common/Pythia8_DD_Common.py
new file mode 100644
index 0000000000000000000000000000000000000000..f9c6c5d0b3eec3c08fee47f827141e4bf7f52c8f
--- /dev/null
+++ b/Generators/SFGen_i/share/common/Pythia8_DD_Common.py
@@ -0,0 +1,20 @@
+include('SFGen_i/Pythia8_Base_Common.py')
+
+# Pythia shower configuration flags in the double dissociation case
+
+if genSeq.SFGenConfig.diff != 'dd':
+ raise Exception("DD Pythia8 shower configuration can only be used with diff='dd'")
+
+genSeq.Pythia8.Commands += [
+ "PartonLevel:MPI = on",
+ "PartonLevel:Remnants = off",
+ "Check:event = off",
+ "BeamRemnants:primordialKT = off",
+ "LesHouches:matchInOut = off"
+ "PartonLevel:FSR = on",
+ "SpaceShower:dipoleRecoil = on",
+ "SpaceShower:pTmaxMatch = 2",
+ "SpaceShower:QEDshowerByQ = off",
+ "SpaceShower:pTdampMatch=1",
+ "BeamRemnants:unresolvedHadron = 0",
+]
diff --git a/Generators/SFGen_i/share/common/Pythia8_EL_Common.py b/Generators/SFGen_i/share/common/Pythia8_EL_Common.py
new file mode 100644
index 0000000000000000000000000000000000000000..0f5b086174dcd21019d6e4e7ae13bf1e6478bd22
--- /dev/null
+++ b/Generators/SFGen_i/share/common/Pythia8_EL_Common.py
@@ -0,0 +1,20 @@
+include('SFGen_i/Pythia8_Base_Common.py')
+
+# Pythia shower configuration flags in the elastic case
+
+if genSeq.SFGenConfig.diff != 'el':
+ raise Exception("EL Pythia8 shower configuration can only be used with diff='el'")
+
+genSeq.Pythia8.Commands += [
+ "PartonLevel:MPI = off",
+ "PartonLevel:Remnants = off",
+ "Check:event = off",
+ "BeamRemnants:primordialKT = off",
+ "LesHouches:matchInOut = off"
+ "PartonLevel:FSR = on",
+ "SpaceShower:dipoleRecoil = on",
+ "SpaceShower:pTmaxMatch = 2",
+ "SpaceShower:QEDshowerByQ = off",
+ "SpaceShower:pTdampMatch=1",
+ "BeamRemnants:unresolvedHadron = 3",
+]
diff --git a/Generators/SFGen_i/share/common/Pythia8_SD_Common.py b/Generators/SFGen_i/share/common/Pythia8_SD_Common.py
new file mode 100644
index 0000000000000000000000000000000000000000..d154b7f03b95ed0f2d24ae6cb4777838aca6a920
--- /dev/null
+++ b/Generators/SFGen_i/share/common/Pythia8_SD_Common.py
@@ -0,0 +1,26 @@
+include('SFGen_i/Pythia8_Base_Common.py')
+
+# Pythia shower configuration flags in the single dissociation case
+
+if not (genSeq.SFGenConfig.diff == 'sda' or genSeq.SFGenConfig.diff == 'sdb'):
+ raise Exception("SD Pythia8 shower configuration can only be used with diff='sda' or 'sdb'")
+
+unresolvedHadron = -1
+if genSeq.SFGenConfig.diff=='sda':
+ unresolvedHadron = 2
+elif genSeq.SFGenConfig.diff=='sdb':
+ unresolvedHadron = 1
+
+genSeq.Pythia8.Commands += [
+ "PartonLevel:MPI = off",
+ "PartonLevel:Remnants = off",
+ "Check:event = off",
+ "BeamRemnants:primordialKT = off",
+ "LesHouches:matchInOut = off"
+ "PartonLevel:FSR = on",
+ "SpaceShower:dipoleRecoil = on",
+ "SpaceShower:pTmaxMatch = 2",
+ "SpaceShower:QEDshowerByQ = off",
+ "SpaceShower:pTdampMatch=1",
+ "BeamRemnants:unresolvedHadron = {}".format(unresolvedHadron),
+]
diff --git a/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_DD.py b/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_DD.py
new file mode 100644
index 0000000000000000000000000000000000000000..c976a8bfb6412e71c84a1ee890f5467d23f13841
--- /dev/null
+++ b/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_DD.py
@@ -0,0 +1,31 @@
+evgenConfig.description = "SFGen MC gamma + gamma pp collisions at 13000 GeV to 2 muons with dissociation on both side"
+evgenConfig.keywords = ["2photon","2muon","dissociation"]
+evgenConfig.contact = ["lucian.harland-lang@physics.ox.ac.uk"]
+
+from SFGen_i.SFGenUtils import SFGenConfig, SFGenRun
+
+#class with the sfgen initialization parameters. Please see SFGenUtils for a complete list of tunable parameters.
+scConfig = SFGenConfig(runArgs)
+
+scConfig.PDFname = 'SF_MSHT20qed_nnlo' # PDF set name
+scConfig.PDFmember = 0 # PDF member
+scConfig.proc = 2 # Process number, see SFGen Manual for labelling at https://sfgen.hepforge.org/
+
+scConfig.diff = 'dd' # interaction: elastic ('el'), single ('sd','sda','sdb') and double ('dd') dissociation.
+scConfig.genunw = True
+scConfig.ymin = -5.0 # Minimum dilepton rapidity
+scConfig.ymax = 5.0 # Maximum dilepton rapidity
+scConfig.mmin = 20 # Minimum dilepton mass
+scConfig.mmax = 2000 # Maximum dilepton mass
+scConfig.gencuts = True # Generate cuts below
+
+scConfig.ptamin = 10.0 # Minimum pT of outgoing object a
+scConfig.ptbmin = 10. # Minimum pT of outgoing object b
+scConfig.etaamin = -2.5 # Minimum eta of outgoing object a
+scConfig.etaamax = 2.5 # Maximum eta of outgoing object a
+scConfig.etabmin = -2.5 # Minimum eta of outgoing object b
+scConfig.etabmax = 2.5 # Maximum eta of outgoing object b
+
+
+SFGenRun(scConfig, genSeq)
+include('SFGen_i/Pythia8_DD_Common.py')
diff --git a/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_El.py b/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_El.py
new file mode 100644
index 0000000000000000000000000000000000000000..ffd894ba4ef5ef1b3bde38f3220495be14446663
--- /dev/null
+++ b/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_El.py
@@ -0,0 +1,31 @@
+evgenConfig.description = "SFGen MC gamma + gamma pp collisions at 13000 GeV to 2 muons without dissociation"
+evgenConfig.keywords = ["2photon","2muon","dissociation"]
+evgenConfig.contact = ["lucian.harland-lang@physics.ox.ac.uk"]
+
+from SFGen_i.SFGenUtils import SFGenConfig, SFGenRun
+
+#class with the sfgen initialization parameters. Please see SFGenUtils for a complete list of tunable parameters.
+scConfig = SFGenConfig(runArgs)
+
+scConfig.PDFname = 'SF_MSHT20qed_nnlo' # PDF set name
+scConfig.PDFmember = 0 # PDF member
+scConfig.proc = 2 # Process number, see SFGen Manual for labelling at https://sfgen.hepforge.org/
+
+scConfig.diff = 'el' # interaction: elastic ('el'), single ('sd','sda','sdb') and double ('dd') dissociation.
+scConfig.genunw = True
+scConfig.ymin = -5.0 # Minimum dilepton rapidity
+scConfig.ymax = 5.0 # Maximum dilepton rapidity
+scConfig.mmin = 20 # Minimum dilepton mass
+scConfig.mmax = 2000 # Maximum dilepton mass
+scConfig.gencuts = True # Generate cuts below
+
+scConfig.ptamin = 10.0 # Minimum pT of outgoing object a
+scConfig.ptbmin = 10. # Minimum pT of outgoing object b
+scConfig.etaamin = -2.5 # Minimum eta of outgoing object a
+scConfig.etaamax = 2.5 # Maximum eta of outgoing object a
+scConfig.etabmin = -2.5 # Minimum eta of outgoing object b
+scConfig.etabmax = 2.5 # Maximum eta of outgoing object b
+
+
+SFGenRun(scConfig, genSeq)
+include('SFGen_i/Pythia8_EL_Common.py')
diff --git a/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_SDA.py b/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_SDA.py
new file mode 100644
index 0000000000000000000000000000000000000000..791a9099c804b8938bb83cecc18c7a2cf796550b
--- /dev/null
+++ b/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_SDA.py
@@ -0,0 +1,31 @@
+evgenConfig.description = "SFGen MC gamma + gamma pp collisions at 13000 GeV to 2 muons with dissociation on single side"
+evgenConfig.keywords = ["2photon","2muon","dissociation"]
+evgenConfig.contact = ["lucian.harland-lang@physics.ox.ac.uk"]
+
+from SFGen_i.SFGenUtils import SFGenConfig, SFGenRun
+
+#class with the sfgen initialization parameters. Please see SFGenUtils for a complete list of tunable parameters.
+scConfig = SFGenConfig(runArgs)
+
+scConfig.PDFname = 'SF_MSHT20qed_nnlo' # PDF set name
+scConfig.PDFmember = 0 # PDF member
+scConfig.proc = 2 # Process number, see SFGen Manual for labelling at https://sfgen.hepforge.org/
+
+scConfig.diff = 'sda' # interaction: elastic ('el'), single ('sd','sda','sdb') and double ('dd') dissociation.
+scConfig.genunw = True
+scConfig.ymin = -5.0 # Minimum dilepton rapidity
+scConfig.ymax = 5.0 # Maximum dilepton rapidity
+scConfig.mmin = 20 # Minimum dilepton mass
+scConfig.mmax = 2000 # Maximum dilepton mass
+scConfig.gencuts = True # Generate cuts below
+
+scConfig.ptamin = 10.0 # Minimum pT of outgoing object a
+scConfig.ptbmin = 10. # Minimum pT of outgoing object b
+scConfig.etaamin = -2.5 # Minimum eta of outgoing object a
+scConfig.etaamax = 2.5 # Maximum eta of outgoing object a
+scConfig.etabmin = -2.5 # Minimum eta of outgoing object b
+scConfig.etabmax = 2.5 # Maximum eta of outgoing object b
+
+
+SFGenRun(scConfig, genSeq)
+include('SFGen_i/Pythia8_SD_Common.py')
diff --git a/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_SDB.py b/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_SDB.py
new file mode 100644
index 0000000000000000000000000000000000000000..8d465c8b553d1dcf4c82ab8c70f7e21b28263692
--- /dev/null
+++ b/Generators/SFGen_i/share/mc.SFGenPy8_MuMu_SDB.py
@@ -0,0 +1,31 @@
+evgenConfig.description = "SFGen MC gamma + gamma pp collisions at 13000 GeV to 2 muons with dissociation on single side"
+evgenConfig.keywords = ["2photon","2muon","dissociation"]
+evgenConfig.contact = ["lucian.harland-lang@physics.ox.ac.uk"]
+
+from SFGen_i.SFGenUtils import SFGenConfig, SFGenRun
+
+#class with the sfgen initialization parameters. Please see SFGenUtils for a complete list of tunable parameters.
+scConfig = SFGenConfig(runArgs)
+
+scConfig.PDFname = 'SF_MSHT20qed_nnlo' # PDF set name
+scConfig.PDFmember = 0 # PDF member
+scConfig.proc = 2 # Process number, see SFGen Manual for labelling at https://sfgen.hepforge.org/
+
+scConfig.diff = 'sdb' # interaction: elastic ('el'), single ('sd','sda','sdb') and double ('dd') dissociation.
+scConfig.genunw = True
+scConfig.ymin = -5.0 # Minimum dilepton rapidity
+scConfig.ymax = 5.0 # Maximum dilepton rapidity
+scConfig.mmin = 20 # Minimum dilepton mass
+scConfig.mmax = 2000 # Maximum dilepton mass
+scConfig.gencuts = True # Generate cuts below
+
+scConfig.ptamin = 10.0 # Minimum pT of outgoing object a
+scConfig.ptbmin = 10. # Minimum pT of outgoing object b
+scConfig.etaamin = -2.5 # Minimum eta of outgoing object a
+scConfig.etaamax = 2.5 # Maximum eta of outgoing object a
+scConfig.etabmin = -2.5 # Minimum eta of outgoing object b
+scConfig.etabmax = 2.5 # Maximum eta of outgoing object b
+
+
+SFGenRun(scConfig, genSeq)
+include('SFGen_i/Pythia8_SD_Common.py')