diff --git a/.github/workflows/clean-ci.yml b/.github/workflows/clean-ci.yml new file mode 100644 index 000000000..aa5364983 --- /dev/null +++ b/.github/workflows/clean-ci.yml @@ -0,0 +1,39 @@ +--- +name: Clean PR builds + +'on': + workflow_dispatch: + inputs: + pr: + description: PR number in this repo to be cleaned + type: string # can't use number here + required: true + + # Warning: GitHub limits the total number of inputs to 10, so a maximum of + # 9 checks is allowed here! + # Warning: the check_* keys are magic and must consist of the string + # "check_" followed by the applicable check name exactly. The + # "description" field is only the human-readable label for the input. + 'check_build/O2DPG/sim/o2': + description: build/O2DPG/sim/o2 + type: boolean + default: true + 'check_build/O2DPG/O2fst/o2': + description: build/O2DPG/O2fst/o2 + type: boolean + default: true + +permissions: {} + +jobs: + clean: + name: Clean PR checks + uses: alisw/ali-bot/.github/workflows/clean-pr-checks.yml@master + with: + owner: ${{ github.event.repository.owner.login }} + repo: ${{ github.event.repository.name }} + pr: ${{ github.event.inputs.pr }} + checks: ${{ toJSON(github.event.inputs) }} + permissions: + pull-requests: read # to get last commit for pr (octokit/graphql-action) + statuses: write # for set-github-status diff --git a/DATA/production/configurations/asyncReco/setenv_extra.sh b/DATA/production/configurations/asyncReco/setenv_extra.sh index 6ccb976f7..56a473d53 100644 --- a/DATA/production/configurations/asyncReco/setenv_extra.sh +++ b/DATA/production/configurations/asyncReco/setenv_extra.sh @@ -420,7 +420,7 @@ elif [[ $ALIGNLEVEL == 1 ]]; then export TPC_CORR_SCALING+=" --enable-M-shape-correction " fi - if [[ $ALIEN_JDL_LPMANCHORYEAR == "2023" ]] && [[ $BEAMTYPE == "PbPb" ]] ; then + if [[ $ALIEN_JDL_LPMANCHORYEAR -ge 2023 ]] && [[ $BEAMTYPE == "PbPb" ]] ; then # adding additional cluster errors # the values below should be squared, but the validation of those values (0.01 and 0.0225) is ongoing TPCEXTRAERR=";GPU_rec_tpc.clusterError2AdditionalYSeeding=0.1;GPU_rec_tpc.clusterError2AdditionalZSeeding=0.15;" diff --git a/DATA/production/qc-async/emc.json b/DATA/production/qc-async/emc.json index d726cdc96..1eb272482 100644 --- a/DATA/production/qc-async/emc.json +++ b/DATA/production/qc-async/emc.json @@ -79,7 +79,7 @@ "cycleDurationSeconds": "60", "dataSource": { "type": "direct", - "query": "emcal-triggers:EMC/CELLSTRGR/0;ctp-digits:CTP/DIGITS;ctp-config:CTP/CONFIG/0?lifetime=condition&ccdb-path=CTP/Config/Config&ccdb-run-dependent=true" + "query": "emcal-triggers:EMC/CELLSTRGR/0;ctp-digits:CTP/DIGITS;ctp-config:CTP/CONFIG/0?lifetime=condition&ccdb-path=CTP/Config/Config&ccdb-run-dependent=1" }, "taskParameters": { } diff --git a/DATA/production/qc-async/emc_PbPb.json b/DATA/production/qc-async/emc_PbPb.json index 0380c2283..e9328d890 100644 --- a/DATA/production/qc-async/emc_PbPb.json +++ b/DATA/production/qc-async/emc_PbPb.json @@ -93,7 +93,7 @@ "cycleDurationSeconds": "60", "dataSource": { "type": "direct", - "query": "emcal-triggers:EMC/CELLSTRGR/0;ctp-digits:CTP/DIGITS;ctp-config:CTP/CONFIG/0?lifetime=condition&ccdb-path=CTP/Config/Config&ccdb-run-dependent=true" + "query": "emcal-triggers:EMC/CELLSTRGR/0;ctp-digits:CTP/DIGITS;ctp-config:CTP/CONFIG/0?lifetime=condition&ccdb-path=CTP/Config/Config&ccdb-run-dependent=1" }, "taskParameters": { "AliasMB" : "CMTVXTSC,C0TVXTSC,C0TVXTCE" diff --git a/DATA/production/qc-async/tpc.json b/DATA/production/qc-async/tpc.json index a8c7f2ac9..32548958e 100644 --- a/DATA/production/qc-async/tpc.json +++ b/DATA/production/qc-async/tpc.json @@ -67,7 +67,8 @@ "cutMaxpTPC": "20.", "cutMinpTPCMIPs": "0.45", "cutMaxpTPCMIPs": "0.55", - "turnOffHistosForAsync": "true" + "turnOffHistosForAsync": "true", + "getdEdxVspHypoHist": "true" } }, "TPCTrackClusters": { diff --git a/DATA/tools/epn/gen_topo_o2dpg.sh b/DATA/tools/epn/gen_topo_o2dpg.sh index 03a236f20..39bdfc98a 100755 --- a/DATA/tools/epn/gen_topo_o2dpg.sh +++ b/DATA/tools/epn/gen_topo_o2dpg.sh @@ -82,15 +82,21 @@ while true; do break fi fi - if [[ ! -d O2DPG ]]; then git clone https://github.com/AliceO2Group/O2DPG.git 1>&2 || { echo O2DPG checkout failed 1>&2; exit 1; }; fi - cd O2DPG + for CHECKOUTATTEMPT in 1 2; do + if [[ ! -d O2DPG ]]; then git clone https://github.com/AliceO2Group/O2DPG.git 1>&2 || { echo O2DPG checkout failed 1>&2; exit 1; }; fi + cd O2DPG + rm -f DATA/core_dump_* + git reset --hard HEAD &> /dev/null && git clean -d -f &> /dev/null && break + [[ $CHECKOUTATTEMPT -eq 2 ]] && { echo git reset error 1>&2; exit 1; } + echo "Clean-up of O2DPG repository failed. Removing repository and cloning it from scratch" 1>&2 + cd $GEN_TOPO_WORKDIR || { echo Cannot enter work dir 1>&2; exit 1; } + rm -rf O2DPG + done git checkout $GEN_TOPO_SOURCE &> /dev/null if [[ $? != 0 ]]; then git fetch --tags origin 1>&2 || { echo Repository update failed 1>&2; exit 1; } git checkout $GEN_TOPO_SOURCE &> /dev/null || { echo commit does not exist 1>&2; exit 1; } fi - git reset --hard $GEN_TOPO_SOURCE &> /dev/null || { echo git reset error 1>&2; exit 1; } - rm -f DATA/core_dump_* # At a tag, or a detached non-dirty commit, but not on a branch if ! git describe --exact-match --tags HEAD &> /dev/null && ( git symbolic-ref -q HEAD &> /dev/null || ! git diff-index --quiet HEAD &> /dev/null ); then unset GEN_TOPO_CACHEABLE diff --git a/GRID/utils/grid_submit.sh b/GRID/utils/grid_submit.sh index d8fb95397..84cc73fcc 100755 --- a/GRID/utils/grid_submit.sh +++ b/GRID/utils/grid_submit.sh @@ -194,7 +194,7 @@ ONGRID=0 JOBTTL=82000 CPUCORES=8 -PRODSPLIT=1 +PRODSPLIT=${PRODSPLIT:-1} # this tells us to continue an existing job --> in this case we don't create a new workdir while [ $# -gt 0 ] ; do case $1 in diff --git a/MC/bin/o2_hybrid_gen.py b/MC/bin/o2_hybrid_gen.py new file mode 100755 index 000000000..12b5b7c73 --- /dev/null +++ b/MC/bin/o2_hybrid_gen.py @@ -0,0 +1,118 @@ +#!/usr/bin/env python3 + +# This script creates a template JSON file for the configuration of the hybrid generator +# The generators to be used are passed as a list of strings via the --gen flag, while the +# output filename is set using --output. All the generators available in O2sim are supported. +# Example: +# $O2DPG_ROOT/MC/bin/o2_hybrid_gen.py --gen pythia8 boxgen external extkinO2 hepmc pythia8hf \ +# --output config.json + +import argparse +import json + +def main(): + parser = argparse.ArgumentParser(description='Create a JSON file from command line flags.') + parser.add_argument('--gen', type=str, nargs='+', required=True, help='List of generators to be used') + parser.add_argument('--output', type=str, required=True, help='Output JSON file path') + + args = parser.parse_args() + + # put in a list all the elementes in the gen flag + noConfGen = ["pythia8pp", "pythia8hf", "pythia8hi", "pythia8powheg"] + gens = [] + for gen in args.gen: + if gen == "pythia8": + gens.append({ + 'name': 'pythia8', + 'config': { + "config": "$O2_ROOT/share/Generators/egconfig/pythia8_inel.cfg", + "hooksFileName": "", + "hooksFuncName": "", + "includePartonEvent": False, + "particleFilter": "", + "verbose": 0 + } + }) + elif gen == "external": + gens.append({ + 'name': 'external', + 'config': { + "fileName": "${O2DPG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorParamPromptJpsiToElectronEvtGen_pp13TeV.C", + "funcName": "GeneratorParamPromptJpsiToElectronEvtGen_pp13TeV()" + } + }) + elif gen == "extkinO2": + gens.append({ + 'name': 'extkinO2', + 'config': { + "skipNonTrackable": True, + "continueMode": False, + "roundRobin": False, + "randomize": False, + "rngseed": 0, + "randomphi": False, + "fileName": "/path/to/filename.root" + } + }) + elif gen == "hepmc": + gens.append({ + "name": "hepmc", + "config": { + "configcmd": { + "fileNames": "", + "cmd": "" + }, + "confighepmc": { + "version": 2, + "eventsToSkip": 0, + "fileName": "/path/to/filename.hepmc", + "prune": False + } + } + }) + elif gen == "boxgen": + gens.append({ + "name": "boxgen", + "config": { + "pdg": 13, + "number": 1, + "eta": [ + -4, + -2.5 + ], + "prange": [ + 0.1, + 5 + ], + "phirange": [ + 0, + 360 + ] + } + }) + elif gen in noConfGen: + gens.append({ + "name": gen, + "config": "" + }) + else: + print(f"Generator {gen} not found in the list of available generators") + exit(1) + + # fill fractions with 1 for each generator + fractions = [1] * len(gens) + + # Put gens and fractions in the data dictionary + data = { + "generators": gens, + "fractions": fractions + } + + # Write the data dictionary to a JSON file + with open(args.output, 'w') as f: + json.dump(data, f, indent=2) + + print(f"JSON file created at {args.output}") + +if __name__ == "__main__": + main() diff --git a/MC/bin/o2dpg_sim_workflow.py b/MC/bin/o2dpg_sim_workflow.py index 6febf07da..09e9ecd0b 100755 --- a/MC/bin/o2dpg_sim_workflow.py +++ b/MC/bin/o2dpg_sim_workflow.py @@ -21,10 +21,11 @@ import importlib.util import argparse from os import environ, mkdir -from os.path import join, dirname, isdir, isabs +from os.path import join, dirname, isdir, isabs, isfile import random import json import itertools +import math import requests, re pandas_available = True try: @@ -61,6 +62,7 @@ parser.add_argument('-ini',help='generator init parameters file (full paths required), for example: ${O2DPG_ROOT}/MC/config/PWGHF/ini/GeneratorHF.ini', default='') parser.add_argument('-confKey',help='generator or trigger configuration key values, for example: "GeneratorPythia8.config=pythia8.cfg;A.x=y"', default='') parser.add_argument('--readoutDets',help='comma separated string of detectors readout (does not modify material budget - only hit creation)', default='all') +parser.add_argument('--make-evtpool', help='Generate workflow for event pool creation.', action='store_true') parser.add_argument('-interactionRate',help='Interaction rate, used in digitization', default=-1) parser.add_argument('-bcPatternFile',help='Bunch crossing pattern file, used in digitization (a file name or "ccdb")', default='') @@ -116,9 +118,9 @@ parser.add_argument('--no-tpc-digitchunking', action='store_true', help=argparse.SUPPRESS) parser.add_argument('--no-strangeness-tracking', action='store_true', default=False, help="Disable strangeness tracking") parser.add_argument('--combine-tpc-clusterization', action='store_true', help=argparse.SUPPRESS) #<--- useful for small productions (pp, low interaction rate, small number of events) -parser.add_argument('--first-orbit', default=0, type=int, help=argparse.SUPPRESS) # to set the first orbit number of the run for HBFUtils (only used when anchoring) +parser.add_argument('--first-orbit', default=256, type=int, help=argparse.SUPPRESS) # to set the first orbit number of the run for HBFUtils (only used when anchoring); default 256 for convenience to allow for some orbits-early # (consider doing this rather in O2 digitization code directly) -parser.add_argument('--orbits-early', default=0, type=int, help=argparse.SUPPRESS) # number of orbits to start simulating earlier +parser.add_argument('--orbits-early', default=0, type=float, help=argparse.SUPPRESS) # number of orbits to start simulating earlier # to reduce start of timeframe effects in MC --> affects collision context parser.add_argument('--sor', default=-1, type=int, help=argparse.SUPPRESS) # may pass start of run with this (otherwise it is autodetermined from run number) parser.add_argument('--run-anchored', action='store_true', help=argparse.SUPPRESS) @@ -242,21 +244,44 @@ def retrieve_sor(run_number): """ retrieves start of run (sor) from the RCT/Info/RunInformation table with a simple http request - in case of problems, 0 will be returned + in case of problems, 0 will be returned. Simple http request has advantage + of not needing to initialize a Ccdb object. """ + url="http://alice-ccdb.cern.ch/browse/RCT/Info/RunInformation/"+str(run_number) ansobject=requests.get(url) tokens=ansobject.text.split("\n") + # determine start of run, earlier values take precedence (see also implementation in BasicCCDBManager::getRunDuration) + STF=0 + # extract SOR by pattern matching + for t in tokens: + match_object=re.match(r"\s*(STF\s*=\s*)([0-9]*)\s*", t) + if match_object != None: + STF=int(match_object[2]) + break + if STF > 0: + return STF + + SOX=0 + # extract SOX by pattern matching + for t in tokens: + match_object=re.match(r"\s*(STF\s*=\s*)([0-9]*)\s*", t) + if match_object != None: + SOX=int(match_object[2]) + break + if SOX > 0: + return SOX + SOR=0 # extract SOR by pattern matching for t in tokens: match_object=re.match(r"\s*(SOR\s*=\s*)([0-9]*)\s*", t) if match_object != None: - SOR=match_object[2] + SOR=int(match_object[2]) break - - return int(SOR) + + return SOR # check and sanitize config-key values (extract and remove diamond vertex arguments into finalDiamondDict) @@ -626,7 +651,7 @@ def getDPL_global_options(bigshm=False, ccdbbackend=True): PreCollContextTask=createTask(name='precollcontext_' + str(tf), needs=precollneeds, tf=tf, cwd=timeframeworkdir, cpu='1') PreCollContextTask['cmd']='${O2_ROOT}/bin/o2-steer-colcontexttool -i ' + signalprefix + ',' + str(INTRATE) + ',' + str(args.ns) + ':' + str(args.ns) \ + ' --show-context ' + ' --timeframeID ' + str(tf-1 + int(args.production_offset)*NTIMEFRAMES) \ - + ' --orbitsPerTF ' + str(orbitsPerTF) + ' --orbits ' + str(orbitsPerTF + args.orbits_early) \ + + ' --orbitsPerTF ' + str(orbitsPerTF) + ' --orbits ' + str(orbitsPerTF + math.ceil(args.orbits_early)) \ + ' --seed ' + str(TFSEED) + ' --noEmptyTF --first-orbit ' + str(args.first_orbit - args.orbits_early) PreCollContextTask['cmd'] += ' --bcPatternFile ccdb' # <--- the object should have been set in (local) CCDB if includeQED: @@ -699,9 +724,14 @@ def getDPL_global_options(bigshm=False, ccdbbackend=True): # possible generator) workflow['stages'].append(SGN_CONFIG_task) + + # default flags for extkinO2 signal simulation (no transport) + extkinO2Config = '' + if GENERATOR == 'extkinO2': + extkinO2Config = ';GeneratorFromO2Kine.randomize=true;GeneratorFromO2Kine.rngseed=' + str(TFSEED) # determine final conf key for signal simulation - CONFKEY = constructConfigKeyArg(create_geant_config(args, args.confKey)) + CONFKEY = constructConfigKeyArg(create_geant_config(args, args.confKey + extkinO2Config)) # ----------------- # transport signals # ----------------- @@ -741,7 +771,11 @@ def getDPL_global_options(bigshm=False, ccdbbackend=True): # determine the skip number cmd = 'export HEPMCEVENTSKIP=$(${O2DPG_ROOT}/UTILS/ReadHepMCEventSkip.sh ../HepMCEventSkip.json ' + str(tf) + ');' SGNGENtask['cmd'] = cmd - SGNGENtask['cmd'] +='${O2_ROOT}/bin/o2-sim --noGeant -j 1 --field ccdb --vertexMode kCCDB' \ + + # No vertexing for event pool generation + vtxmode = 'kNoVertex' if args.make_evtpool else 'kCCDB' + + SGNGENtask['cmd'] +='${O2_ROOT}/bin/o2-sim --noGeant -j 1 --field ccdb --vertexMode ' + vtxmode \ + ' --run ' + str(args.run) + ' ' + str(CONFKEY) + str(TRIGGER) \ + ' -g ' + str(GENERATOR) + ' ' + str(INIFILE) + ' -o genevents ' + embeddinto \ + ('', ' --timestamp ' + str(args.timestamp))[args.timestamp!=-1] \ @@ -754,6 +788,11 @@ def getDPL_global_options(bigshm=False, ccdbbackend=True): if sep_event_mode == True: workflow['stages'].append(SGNGENtask) signalneeds = signalneeds + [SGNGENtask['name']] + if args.make_evtpool: + continue + + # GeneratorFromO2Kine parameters are needed only before the transport + CONFKEY = re.sub(r'GeneratorFromO2Kine.*?;', '', CONFKEY) sgnmem = 6000 if COLTYPE == 'PbPb' else 4000 SGNtask=createTask(name='sgnsim_'+str(tf), needs=signalneeds, tf=tf, cwd='tf'+str(tf), lab=["GEANT"], @@ -1520,26 +1559,32 @@ def addQCPerTF(taskName, needs, readerCommand, configFilePath, objectsFile=''): TFcleanup['cmd'] += 'rm *cluster*.root' workflow['stages'].append(TFcleanup) -# AOD merging as one global final step -aodmergerneeds = ['aod_' + str(tf) for tf in range(1, NTIMEFRAMES + 1)] -AOD_merge_task = createTask(name='aodmerge', needs = aodmergerneeds, lab=["AOD"], mem='2000', cpu='1') -AOD_merge_task['cmd'] = ' set -e ; [ -f aodmerge_input.txt ] && rm aodmerge_input.txt; ' -AOD_merge_task['cmd'] += ' for i in `seq 1 ' + str(NTIMEFRAMES) + '`; do echo "tf${i}/AO2D.root" >> aodmerge_input.txt; done; ' -AOD_merge_task['cmd'] += ' o2-aod-merger --input aodmerge_input.txt --output AO2D.root' -# produce MonaLisa event stat file -AOD_merge_task['cmd'] += ' ; ${O2DPG_ROOT}/MC/bin/o2dpg_determine_eventstat.py' -workflow['stages'].append(AOD_merge_task) - -job_merging = False -if includeFullQC: - workflow['stages'].extend(include_all_QC_finalization(ntimeframes=NTIMEFRAMES, standalone=False, run=args.run, productionTag=args.productionTag, conditionDB=args.conditionDB, qcdbHost=args.qcdbHost)) - - -if includeAnalysis: - # include analyses and potentially final QC upload tasks - add_analysis_tasks(workflow["stages"], needs=[AOD_merge_task["name"]], is_mc=True, collision_system=COLTYPE) - if QUALITYCONTROL_ROOT: - add_analysis_qc_upload_tasks(workflow["stages"], args.productionTag, args.run, "passMC") +if not args.make_evtpool: + # AOD merging as one global final step + aodmergerneeds = ['aod_' + str(tf) for tf in range(1, NTIMEFRAMES + 1)] + AOD_merge_task = createTask(name='aodmerge', needs = aodmergerneeds, lab=["AOD"], mem='2000', cpu='1') + AOD_merge_task['cmd'] = ' set -e ; [ -f aodmerge_input.txt ] && rm aodmerge_input.txt; ' + AOD_merge_task['cmd'] += ' for i in `seq 1 ' + str(NTIMEFRAMES) + '`; do echo "tf${i}/AO2D.root" >> aodmerge_input.txt; done; ' + AOD_merge_task['cmd'] += ' o2-aod-merger --input aodmerge_input.txt --output AO2D.root' + # produce MonaLisa event stat file + AOD_merge_task['cmd'] += ' ; ${O2DPG_ROOT}/MC/bin/o2dpg_determine_eventstat.py' + workflow['stages'].append(AOD_merge_task) + + job_merging = False + if includeFullQC: + workflow['stages'].extend(include_all_QC_finalization(ntimeframes=NTIMEFRAMES, standalone=False, run=args.run, productionTag=args.productionTag, conditionDB=args.conditionDB, qcdbHost=args.qcdbHost)) + + if includeAnalysis: + # include analyses and potentially final QC upload tasks + add_analysis_tasks(workflow["stages"], needs=[AOD_merge_task["name"]], is_mc=True, collision_system=COLTYPE) + if QUALITYCONTROL_ROOT: + add_analysis_qc_upload_tasks(workflow["stages"], args.productionTag, args.run, "passMC") +else: + wfneeds=['sgngen_' + str(tf) for tf in range(1, NTIMEFRAMES + 1)] + tfpool=['tf' + str(tf) + '/genevents_Kine.root' for tf in range(1, NTIMEFRAMES + 1)] + POOL_merge_task = createTask(name='poolmerge', needs=wfneeds, lab=["POOL"], mem='2000', cpu='1') + POOL_merge_task['cmd'] = '${O2DPG_ROOT}/UTILS/root_merger.py -o evtpool.root -i ' + ','.join(tfpool) + workflow['stages'].append(POOL_merge_task) # adjust for alternate (RECO) software environments adjust_RECO_environment(workflow, args.alternative_reco_software) diff --git a/MC/config/PWGDQ/EvtGen/DecayTablesEvtgen/PSITOJPSIPIPI.DEC b/MC/config/PWGDQ/EvtGen/DecayTablesEvtgen/PSITOJPSIPIPI.DEC new file mode 100644 index 000000000..70bc72517 --- /dev/null +++ b/MC/config/PWGDQ/EvtGen/DecayTablesEvtgen/PSITOJPSIPIPI.DEC @@ -0,0 +1,8 @@ +Decay J/psi +### from DECAY.DEC +1.000 e+ e- PHOTOS VLL; +Enddecay +Decay psi(2S) +1.000 J/psi pi+ pi- VVPIPI; +Enddecay +End diff --git a/MC/config/PWGDQ/EvtGen/DecayTablesEvtgen/X3872TOJPSIPIPI.DEC b/MC/config/PWGDQ/EvtGen/DecayTablesEvtgen/X3872TOJPSIPIPI.DEC new file mode 100644 index 000000000..8451b0777 --- /dev/null +++ b/MC/config/PWGDQ/EvtGen/DecayTablesEvtgen/X3872TOJPSIPIPI.DEC @@ -0,0 +1,10 @@ +Decay X_1(3872) +#### Breit-Wigner function for the pion distribution (S-Wave approximation) +1.000 J/psi pi+ pi- XJPSIRHO0PIPI_S; +Enddecay +Decay J/psi +### from DECAY.DEC +1.000 e+ e- PHOTOS VLL; +Enddecay + +End diff --git a/MC/config/PWGDQ/external/generator/GeneratorBottomonia.C b/MC/config/PWGDQ/external/generator/GeneratorBottomonia.C new file mode 100644 index 000000000..25665ede2 --- /dev/null +++ b/MC/config/PWGDQ/external/generator/GeneratorBottomonia.C @@ -0,0 +1,280 @@ +// +// generators for bottomonia considering at midrapidity and forward rapidity +// + +R__ADD_INCLUDE_PATH($O2DPG_ROOT/MC/config/PWGDQ/EvtGen) +R__LOAD_LIBRARY(libpythia6) +#include "GeneratorCocktail.C" +#include "GeneratorEvtGen.C" + +namespace o2 +{ +namespace eventgen +{ + +///////////////////////////////////////////////////////////////////////////// +class O2_GeneratorParamUpsilon1SFwdY : public GeneratorTGenerator +{ + + public: + O2_GeneratorParamUpsilon1SFwdY() : GeneratorTGenerator("ParamUpsilon1S") + { + paramUpsilon1S = new GeneratorParam(1, -1, PtUpsilon1Spp13TeV, YUpsilon1Spp13TeV, V2Upsilon1Spp13TeV, IpUpsilon1Spp13TeV); + paramUpsilon1S->SetMomentumRange(0., 1.e6); + paramUpsilon1S->SetPtRange(0, 999.); + paramUpsilon1S->SetYRange(-4.2, -2.3); + paramUpsilon1S->SetPhiRange(0., 360.); + paramUpsilon1S->SetDecayer(new TPythia6Decayer()); + paramUpsilon1S->SetForceDecay(kNoDecay); // particle left undecayed + // - - paramUpsilon1S->SetTrackingFlag(1); // (from AliGenParam) -> check this + setTGenerator(paramUpsilon1S); + }; + + ~O2_GeneratorParamUpsilon1SFwdY() + { + delete paramUpsilon1S; + }; + + Bool_t Init() override + { + GeneratorTGenerator::Init(); + paramUpsilon1S->Init(); + return true; + } + + void SetNSignalPerEvent(Int_t nsig) { paramUpsilon1S->SetNumberParticles(nsig); } + + //-------------------------------------------------------------------------// + static Double_t PtUpsilon1Spp13TeV(const Double_t* px, const Double_t* /*dummy*/) + { + // Upsilon1S pt shape from LHCb pp@13TeV arXiv:1804.09214 + Double_t x = *px; + Float_t p0, p1, p2, p3; + + p0 = 4.11195e+02; + p1 = 1.03097e+01; + p2 = 1.62309e+00; + p3 = 4.84709e+00; + + return (p0 * x / TMath::Power(1. + TMath::Power(x / p1, p2), p3)); + } + + //-------------------------------------------------------------------------// + static Double_t YUpsilon1Spp13TeV(const Double_t* py, const Double_t* /*dummy*/) + { + // Upsilon1S y shape from LHCb pp@13TeV arXiv:1804.09214 + Double_t x = *py; + Float_t p0, p1; + + p0 = 3.07931e+03; + p1 = -3.53102e-02; + + return (p0 * (1. + p1 * x * x)); + } + + //-------------------------------------------------------------------------// + static Double_t V2Upsilon1Spp13TeV(const Double_t* /*dummy*/, const Double_t* /*dummy*/) + { + // Upsilon(1S) v2 + return 0.; + } + + //-------------------------------------------------------------------------// + static Int_t IpUpsilon1Spp13TeV(TRandom*) + { + return 553; + } + + private: + GeneratorParam* paramUpsilon1S = nullptr; +}; + +///////////////////////////////////////////////////////////////////////////// +class O2_GeneratorParamUpsilon2SFwdY : public GeneratorTGenerator +{ + + public: + O2_GeneratorParamUpsilon2SFwdY() : GeneratorTGenerator("ParamUpsilon2S") + { + paramUpsilon2S = new GeneratorParam(1, -1, PtUpsilon2Spp13TeV, YUpsilon2Spp13TeV, V2Upsilon2Spp13TeV, IpUpsilon2Spp13TeV); + paramUpsilon2S->SetMomentumRange(0., 1.e6); + paramUpsilon2S->SetPtRange(0, 999.); + paramUpsilon2S->SetYRange(-4.2, -2.3); + paramUpsilon2S->SetPhiRange(0., 360.); + paramUpsilon2S->SetDecayer(new TPythia6Decayer()); + paramUpsilon2S->SetForceDecay(kNoDecay); // particle left undecayed + // - - paramUpsilon2S->SetTrackingFlag(1); // (from AliGenParam) -> check this + setTGenerator(paramUpsilon2S); + }; + + ~O2_GeneratorParamUpsilon2SFwdY() + { + delete paramUpsilon2S; + }; + + Bool_t Init() override + { + GeneratorTGenerator::Init(); + paramUpsilon2S->Init(); + return true; + } + + void SetNSignalPerEvent(Int_t nsig) { paramUpsilon2S->SetNumberParticles(nsig); } + + //-------------------------------------------------------------------------// + static Double_t PtUpsilon2Spp13TeV(const Double_t* px, const Double_t* /*dummy*/) + { + // Upsilon2S pt shape from LHCb pp@13TeV arXiv:1804.09214 + Double_t x = *px; + Float_t p0, p1, p2, p3; + + p0 = 8.15699e+01; + p1 = 1.48060e+01; + p2 = 1.50018e+00; + p3 = 6.34208e+00; + + return (p0 * x / TMath::Power(1. + TMath::Power(x / p1, p2), p3)); + } + + //-------------------------------------------------------------------------// + static Double_t YUpsilon2Spp13TeV(const Double_t* py, const Double_t* /*dummy*/) + { + // Upsilon2s y shape from LHCb pp@13TeV arXiv:1804.09214 + Double_t x = *py; + Float_t p0, p1; + + p0 = 7.50409e+02; + p1 = -3.57039e-02; + + return (p0 * (1. + p1 * x * x)); + } + + //-------------------------------------------------------------------------// + static Double_t V2Upsilon2Spp13TeV(const Double_t* /*dummy*/, const Double_t* /*dummy*/) + { + // Upsilon(2S) v2 + return 0.; + } + + //-------------------------------------------------------------------------// + static Int_t IpUpsilon2Spp13TeV(TRandom*) + { + return 100553; + } + + private: + GeneratorParam* paramUpsilon2S = nullptr; +}; + +///////////////////////////////////////////////////////////////////////////// +class O2_GeneratorParamUpsilon3SFwdY : public GeneratorTGenerator +{ + + public: + O2_GeneratorParamUpsilon3SFwdY() : GeneratorTGenerator("ParamUpsilon3S") + { + paramUpsilon3S = new GeneratorParam(1, -1, PtUpsilon3Spp13TeV, YUpsilon3Spp13TeV, V2Upsilon3Spp13TeV, IpUpsilon3Spp13TeV); + paramUpsilon3S->SetMomentumRange(0., 1.e6); + paramUpsilon3S->SetPtRange(0, 999.); + paramUpsilon3S->SetYRange(-4.2, -2.3); + paramUpsilon3S->SetPhiRange(0., 360.); + paramUpsilon3S->SetDecayer(new TPythia6Decayer()); + paramUpsilon3S->SetForceDecay(kNoDecay); // particle left undecayed + // - - paramUpsilon3S->SetTrackingFlag(1); // (from AliGenParam) -> check this + setTGenerator(paramUpsilon3S); + }; + + ~O2_GeneratorParamUpsilon3SFwdY() + { + delete paramUpsilon3S; + }; + + Bool_t Init() override + { + GeneratorTGenerator::Init(); + paramUpsilon3S->Init(); + return true; + } + + void SetNSignalPerEvent(Int_t nsig) { paramUpsilon3S->SetNumberParticles(nsig); } + + //-------------------------------------------------------------------------// + static Double_t PtUpsilon3Spp13TeV(const Double_t* px, const Double_t* /*dummy*/) + { + // Upsilon3S pt shape from LHCb pp@13TeV arXiv:1804.09214 + Double_t x = *px; + Float_t p0, p1, p2, p3; + + p0 = 3.51590e+01; + p1 = 2.30813e+01; + p2 = 1.40822e+00; + p3 = 9.38026e+00; + + return (p0 * x / TMath::Power(1. + TMath::Power(x / p1, p2), p3)); + } + + //-------------------------------------------------------------------------// + static Double_t YUpsilon3Spp13TeV(const Double_t* py, const Double_t* /*dummy*/) + { + // Upsilon3s y shape from LHCb pp@13TeV arXiv:1804.09214 + Double_t x = *py; + Float_t p0, p1; + + p0 = 3.69961e+02; + p1 = -3.54650e-02; + + return (p0 * (1. + p1 * x * x)); + } + + //-------------------------------------------------------------------------// + static Double_t V2Upsilon3Spp13TeV(const Double_t* /*dummy*/, const Double_t* /*dummy*/) + { + // Upsilon(3S) v2 + return 0.; + } + + //-------------------------------------------------------------------------// + static Int_t IpUpsilon3Spp13TeV(TRandom*) + { + return 200553; + } + + private: + GeneratorParam* paramUpsilon3S = nullptr; +}; + + +} // namespace eventgen +} // namespace o2 + +FairGenerator* GeneratorCocktailBottomoniaToMuonEvtGen_pp13TeV() +{ + + auto genCocktailEvtGen = new o2::eventgen::GeneratorEvtGen(); + + auto genUpsilon1S = new o2::eventgen::O2_GeneratorParamUpsilon1SFwdY; + genUpsilon1S->SetNSignalPerEvent(1); // 1 Upsilon(1S) generated per event by GeneratorParam + + auto genUpsilon2S = new o2::eventgen::O2_GeneratorParamUpsilon2SFwdY; + genUpsilon2S->SetNSignalPerEvent(1); // 1 Upsilon(2S) generated per event by GeneratorParam + + auto genUpsilon3S = new o2::eventgen::O2_GeneratorParamUpsilon3SFwdY; + genUpsilon3S->SetNSignalPerEvent(1); // 1 Upsilon(3S) generated per event by GeneratorParam + + genCocktailEvtGen->AddGenerator(genUpsilon1S, 1); // add Upsilon(1S) generator + genCocktailEvtGen->AddGenerator(genUpsilon2S, 1); // add Upsilon(2S) generator + genCocktailEvtGen->AddGenerator(genUpsilon3S, 1); // add Upsilon(3S) generator + + TString pdgs = "553;100553;200553"; + std::string spdg; + TObjArray* obj = pdgs.Tokenize(";"); + genCocktailEvtGen->SetSizePdg(obj->GetEntriesFast()); + for (int i = 0; i < obj->GetEntriesFast(); i++) { + spdg = obj->At(i)->GetName(); + genCocktailEvtGen->AddPdg(std::stoi(spdg), i); + printf("PDG %d \n", std::stoi(spdg)); + } + genCocktailEvtGen->SetForceDecay(kEvtDiMuon); + + return genCocktailEvtGen; +} \ No newline at end of file diff --git a/MC/config/PWGDQ/external/generator/GeneratorPromptCharmonia.C b/MC/config/PWGDQ/external/generator/GeneratorPromptCharmonia.C index 2be7d7d8b..61bbf2c26 100644 --- a/MC/config/PWGDQ/external/generator/GeneratorPromptCharmonia.C +++ b/MC/config/PWGDQ/external/generator/GeneratorPromptCharmonia.C @@ -638,6 +638,75 @@ class O2_GeneratorParamPsiPbPb5TeV : public GeneratorTGenerator GeneratorParam* paramPsi = nullptr; }; + +class O2_GeneratorParamX3872MidY : public GeneratorTGenerator +{ + + public: + O2_GeneratorParamX3872MidY() : GeneratorTGenerator("ParamX3872MidY") + { + paramX3872 = new GeneratorParam(1, -1, PtX3872pp13TeV, YX3872pp13TeV, V2X3872pp13TeV, IpX3872pp13TeV); + paramX3872->SetMomentumRange(0., 1.e6); + paramX3872->SetPtRange(0., 1000.); + paramX3872->SetYRange(-1.0, 1.0); + paramX3872->SetPhiRange(0., 360.); + paramX3872->SetDecayer(new TPythia6Decayer()); // Pythia + paramX3872->SetForceDecay(kNoDecay); // particle left undecayed + setTGenerator(paramX3872); + }; + + ~O2_GeneratorParamX3872MidY() + { + delete paramX3872; + }; + + Bool_t Init() override + { + GeneratorTGenerator::Init(); + paramX3872->Init(); + return true; + } + + void SetNSignalPerEvent(Int_t nsig) { paramX3872->SetNumberParticles(nsig); } + + //-------------------------------------------------------------------------// + static Double_t PtX3872pp13TeV(const Double_t* px, const Double_t* /*dummy*/) + { + // prompt X3872 pT + // pp, 13TeV (tuned LHCb pp 13 TeV) + // + const Double_t kC = 7.64519e+00 ; + const Double_t kpt0 = 5.30628e+00; + const Double_t kn = 3.30887e+00; + Double_t pt = px[0]; + return kC * pt / TMath::Power((1. + (pt / kpt0) * (pt / kpt0)), kn); + } + + //-------------------------------------------------------------------------// + static Double_t YX3872pp13TeV(const Double_t* py, const Double_t* /*dummy*/) + { + // flat rapidity distribution assumed at midrapidity + return 1.; + } + + //-------------------------------------------------------------------------// + static Double_t V2X3872pp13TeV(const Double_t* /*dummy*/, const Double_t* /*dummy*/) + { + // X3872 v2 + return 0.; + } + + //-------------------------------------------------------------------------// + static Int_t IpX3872pp13TeV(TRandom*) + { + return 9920443; + } + + private: + GeneratorParam* paramX3872 = nullptr; +}; + + } // namespace eventgen } // namespace o2 @@ -741,6 +810,31 @@ FairGenerator* GeneratorCocktailPromptCharmoniaToMuonEvtGen_pp13TeV() return genCocktailEvtGen; } + +FairGenerator* + GeneratorParamPromptPsiToJpsiPiPiEvtGen_pp13TeV(TString pdgs = "100443") +{ + auto gen = new o2::eventgen::GeneratorEvtGen(); + gen->SetNSignalPerEvent(1); // number of jpsis per event + + std::string spdg; + TObjArray* obj = pdgs.Tokenize(";"); + gen->SetSizePdg(obj->GetEntriesFast()); + for (int i = 0; i < obj->GetEntriesFast(); i++) { + spdg = obj->At(i)->GetName(); + gen->AddPdg(std::stoi(spdg), i); + printf("PDG %d \n", std::stoi(spdg)); + } + TString pathO2 = gSystem->ExpandPathName("$O2DPG_ROOT/MC/config/PWGDQ/EvtGen/DecayTablesEvtgen"); + gen->SetDecayTable(Form("%s/PSITOJPSIPIPI.DEC", pathO2.Data())); + + // print debug + gen->PrintDebug(); + + return gen; +} + + FairGenerator* GeneratorParamPromptJpsiToMuonEvtGen_pp13TeV(TString pdgs = "443") { @@ -844,6 +938,29 @@ FairGenerator* } +FairGenerator* + GeneratorParamX3872ToJpsiEvtGen_pp13TeV(TString pdgs = "9920443") +{ + auto gen = new o2::eventgen::GeneratorEvtGen(); + gen->SetNSignalPerEvent(1); // number of jpsis per event + + std::string spdg; + TObjArray* obj = pdgs.Tokenize(";"); + gen->SetSizePdg(obj->GetEntriesFast()); + for (int i = 0; i < obj->GetEntriesFast(); i++) { + spdg = obj->At(i)->GetName(); + gen->AddPdg(std::stoi(spdg), i); + printf("PDG %d \n", std::stoi(spdg)); + } + TString pathO2 = gSystem->ExpandPathName("$O2DPG_ROOT/MC/config/PWGDQ/EvtGen/DecayTablesEvtgen"); + gen->SetDecayTable(Form("%s/X3872TOJPSIPIPI.DEC", pathO2.Data())); + + // print debug + gen->PrintDebug(); + + return gen; +} + diff --git a/MC/config/PWGDQ/external/generator/generator_pythia8_NonPromptSignals_gaptriggered_dq.C b/MC/config/PWGDQ/external/generator/generator_pythia8_NonPromptSignals_gaptriggered_dq.C index eac21670c..8a9f5e7c0 100644 --- a/MC/config/PWGDQ/external/generator/generator_pythia8_NonPromptSignals_gaptriggered_dq.C +++ b/MC/config/PWGDQ/external/generator/generator_pythia8_NonPromptSignals_gaptriggered_dq.C @@ -9,7 +9,6 @@ R__ADD_INCLUDE_PATH($O2DPG_ROOT/MC/config/PWGDQ/EvtGen) #include using namespace o2::eventgen; -using namespace Pythia8; namespace o2 { @@ -140,7 +139,7 @@ Pythia8::Event mOutputEvent; // Control gap-triggering unsigned long long mGeneratedEvents; int mInverseTriggerRatio; - Pythia pythiaMBgen; // minimum bias event + Pythia8::Pythia pythiaMBgen; // minimum bias event TString mConfigMBdecays; int mPDG; std::vector mHadronsPDGs; diff --git a/MC/config/PWGDQ/external/generator/generator_pythia8_withInjectedBottomoniaSignals_gaptriggered_dq.C b/MC/config/PWGDQ/external/generator/generator_pythia8_withInjectedBottomoniaSignals_gaptriggered_dq.C new file mode 100644 index 000000000..de4ab2b7a --- /dev/null +++ b/MC/config/PWGDQ/external/generator/generator_pythia8_withInjectedBottomoniaSignals_gaptriggered_dq.C @@ -0,0 +1,76 @@ +#include "FairGenerator.h" +#include "Generators/GeneratorPythia8.h" +#include "Pythia8/Pythia.h" +#include "TRandom.h" +#include "GeneratorBottomonia.C" +#include + +using namespace o2::eventgen; +using namespace Pythia8; + +class GeneratorPythia8BottomoniaInjectedGapTriggeredDQ : public o2::eventgen::GeneratorPythia8 { +public: + + /// default constructor + GeneratorPythia8BottomoniaInjectedGapTriggeredDQ() = default; + + /// constructor + GeneratorPythia8BottomoniaInjectedGapTriggeredDQ(int inputTriggerRatio = 5, int gentype = 0) { + + mGeneratedEvents = 0; + mGeneratorParam = 0x0; + mInverseTriggerRatio = inputTriggerRatio; + switch (gentype) { + case 0: // generate bottomonia cocktail at forward rapidity + mGeneratorParam = (Generator*)GeneratorCocktailBottomoniaToMuonEvtGen_pp13TeV(); + break; + } + mGeneratorParam->Init(); + } + + /// Destructor + ~GeneratorPythia8BottomoniaInjectedGapTriggeredDQ() = default; + +protected: + +Bool_t importParticles() override + { + GeneratorPythia8::importParticles(); + bool genOk = false; + if (mGeneratedEvents % mInverseTriggerRatio == 0) { // add injected prompt signals to the stack + bool genOk = false; + while (!genOk) { + genOk = (mGeneratorParam->generateEvent() && mGeneratorParam->importParticles()) ? true : false ; + } + int originalSize = mParticles.size(); + for (int ipart=0; ipart < mGeneratorParam->getParticles().size(); ipart++) { + TParticle part = TParticle(mGeneratorParam->getParticles().at(ipart)); + if(part.GetFirstMother() >= 0) part.SetFirstMother(part.GetFirstMother() + originalSize); + if(part.GetFirstDaughter() >= 0) part.SetFirstDaughter(part.GetFirstDaughter() + originalSize); + if(part.GetLastDaughter() >= 0) part.SetLastDaughter(part.GetLastDaughter() + originalSize); + mParticles.push_back(part); + // encodeParticleStatusAndTracking method already called in GeneratorEvtGen.C + } + mGeneratorParam->clearParticles(); + } + + mGeneratedEvents++; + return true; + } + + +private: + Generator* mGeneratorParam; + // Control gap-triggering + unsigned long long mGeneratedEvents; + int mInverseTriggerRatio; +}; + +// Predefined generators: +FairGenerator *GeneratorPythia8InjectedBottomoniaGapTriggered(int inputTriggerRatio, int gentype) { + auto myGen = new GeneratorPythia8BottomoniaInjectedGapTriggeredDQ(inputTriggerRatio,gentype); + auto seed = (gRandom->TRandom::GetSeed() % 900000000); + myGen->readString("Random:setSeed on"); + myGen->readString("Random:seed " + std::to_string(seed)); + return myGen; +} diff --git a/MC/config/PWGDQ/external/generator/generator_pythia8_withInjectedPromptSignals_gaptriggered_dq.C b/MC/config/PWGDQ/external/generator/generator_pythia8_withInjectedPromptSignals_gaptriggered_dq.C index f0910dcad..ce6c7a7a5 100644 --- a/MC/config/PWGDQ/external/generator/generator_pythia8_withInjectedPromptSignals_gaptriggered_dq.C +++ b/MC/config/PWGDQ/external/generator/generator_pythia8_withInjectedPromptSignals_gaptriggered_dq.C @@ -45,7 +45,12 @@ public: case 7: // generate prompt charmonia cocktail at forward rapidity at 5TeV mGeneratorParam = (Generator*)GeneratorCocktailPromptCharmoniaToMuonEvtGen_PbPb5TeV(); break; - + case 8: // generate prompt X_1(3872) to Jpsi pi pi at midrapidity + mGeneratorParam = (Generator*)GeneratorParamX3872ToJpsiEvtGen_pp13TeV("9920443"); + break; + case 9: // generate prompt psi2S to Jpsi pi pi at midrapidity + mGeneratorParam = (Generator*)GeneratorParamPromptPsiToJpsiPiPiEvtGen_pp13TeV("100443"); + break; } mGeneratorParam->Init(); } diff --git a/MC/config/PWGDQ/ini/GeneratorHF_bbbarToBplus_midy_triggerGap.ini b/MC/config/PWGDQ/ini/GeneratorHF_bbbarToBplus_midy_triggerGap.ini index b1e836ff4..32e12f000 100644 --- a/MC/config/PWGDQ/ini/GeneratorHF_bbbarToBplus_midy_triggerGap.ini +++ b/MC/config/PWGDQ/ini/GeneratorHF_bbbarToBplus_midy_triggerGap.ini @@ -12,9 +12,3 @@ funcName = GeneratorBplusToJpsiKaon_EvtGen() config = ${O2DPG_MC_CONFIG_ROOT}/MC/config/common/pythia8/generator/pythia8_hf.cfg hooksFileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/hooks/pythia8_userhooks_qqbar.C hooksFuncName = pythia8_userhooks_bbbar(-1.5,1.5) - -### The setup inhibits transport of primary particles which are produce at forward rapidity. -### The settings below only transports particles in the barrel, which is currently defined as |eta| < 2 - -[Stack] -transportPrimary = barrel diff --git a/MC/config/PWGDQ/ini/GeneratorHF_bbbar_PsiAndJpsi_midy_triggerGap.ini b/MC/config/PWGDQ/ini/GeneratorHF_bbbar_PsiAndJpsi_midy_triggerGap.ini index 5073a28b0..2ba063903 100644 --- a/MC/config/PWGDQ/ini/GeneratorHF_bbbar_PsiAndJpsi_midy_triggerGap.ini +++ b/MC/config/PWGDQ/ini/GeneratorHF_bbbar_PsiAndJpsi_midy_triggerGap.ini @@ -12,9 +12,3 @@ funcName = GeneratorBeautyToPsiAndJpsi_EvtGenMidY() config = ${O2DPG_MC_CONFIG_ROOT}/MC/config/common/pythia8/generator/pythia8_hf.cfg hooksFileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/hooks/pythia8_userhooks_qqbar.C hooksFuncName = pythia8_userhooks_bbbar(-1.5,1.5) - -### The setup inhibits transport of primary particles which are produce at forward rapidity. -### The settings below only transports particles in the barrel, which is currently defined as |eta| < 2 - -[Stack] -transportPrimary = barrel diff --git a/MC/config/PWGDQ/ini/Generator_InjectedBottomoniaFwdy_TriggerGap.ini b/MC/config/PWGDQ/ini/Generator_InjectedBottomoniaFwdy_TriggerGap.ini new file mode 100755 index 000000000..410ec7738 --- /dev/null +++ b/MC/config/PWGDQ/ini/Generator_InjectedBottomoniaFwdy_TriggerGap.ini @@ -0,0 +1,7 @@ +### The external generator derives from GeneratorPythia8. +[GeneratorExternal] +fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/external/generator/generator_pythia8_withInjectedBottomoniaSignals_gaptriggered_dq.C +funcName=GeneratorPythia8InjectedBottomoniaGapTriggered(5,0) + +[GeneratorPythia8] +config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_triggerGap.cfg diff --git a/MC/config/PWGDQ/ini/Generator_InjectedPromptPsi2SToJpsiPiPiMidy_TriggerGap.ini b/MC/config/PWGDQ/ini/Generator_InjectedPromptPsi2SToJpsiPiPiMidy_TriggerGap.ini new file mode 100644 index 000000000..30082ec0f --- /dev/null +++ b/MC/config/PWGDQ/ini/Generator_InjectedPromptPsi2SToJpsiPiPiMidy_TriggerGap.ini @@ -0,0 +1,6 @@ +### The external generator derives from GeneratorPythia8. +[GeneratorExternal] +fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/external/generator/generator_pythia8_withInjectedPromptSignals_gaptriggered_dq.C +funcName=GeneratorPythia8InjectedPromptCharmoniaGapTriggered(5,9) +[GeneratorPythia8] +config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_triggerGap.cfg diff --git a/MC/config/PWGDQ/ini/Generator_InjectedPromptX3782Midy_TriggerGap.ini b/MC/config/PWGDQ/ini/Generator_InjectedPromptX3782Midy_TriggerGap.ini new file mode 100644 index 000000000..097d4377d --- /dev/null +++ b/MC/config/PWGDQ/ini/Generator_InjectedPromptX3782Midy_TriggerGap.ini @@ -0,0 +1,6 @@ +### The external generator derives from GeneratorPythia8. +[GeneratorExternal] +fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/external/generator/generator_pythia8_withInjectedPromptSignals_gaptriggered_dq.C +funcName=GeneratorPythia8InjectedPromptCharmoniaGapTriggered(5,8) +[GeneratorPythia8] +config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_triggerGap.cfg diff --git a/MC/config/PWGDQ/ini/tests/Generator_InjectedBottomoniaFwdy_TriggerGap.C b/MC/config/PWGDQ/ini/tests/Generator_InjectedBottomoniaFwdy_TriggerGap.C new file mode 100644 index 000000000..b7a166069 --- /dev/null +++ b/MC/config/PWGDQ/ini/tests/Generator_InjectedBottomoniaFwdy_TriggerGap.C @@ -0,0 +1,103 @@ +int External() +{ + int checkPdgSignal[] = {553,100553,200553}; + int checkPdgDecay = 13; + double rapiditymin = -4.3; double rapiditymax = -2.3; + std::string path{"o2sim_Kine.root"}; + std::cout << "Check for\nsignal PDG " << checkPdgSignal << "\ndecay PDG " << checkPdgDecay << "\n"; + TFile file(path.c_str(), "READ"); + if (file.IsZombie()) { + std::cerr << "Cannot open ROOT file " << path << "\n"; + return 1; + } + + auto tree = (TTree*)file.Get("o2sim"); + std::vector* tracks{}; + tree->SetBranchAddress("MCTrack", &tracks); + + int nLeptons{}; + int nAntileptons{}; + int nLeptonPairs{}; + int nLeptonPairsToBeDone{}; + int nSignalUpsilon1S{}; + int nSignalUpsilon2S{}; + int nSignalUpsilon3S{}; + int nSignalUpsilon1SWithinAcc{}; + int nSignalUpsilon2SWithinAcc{}; + int nSignalUpsilon3SWithinAcc{}; + auto nEvents = tree->GetEntries(); + o2::steer::MCKinematicsReader mcreader("o2sim", o2::steer::MCKinematicsReader::Mode::kMCKine); + Bool_t isInjected = kFALSE; + + for (int i = 0; i < nEvents; i++) { + tree->GetEntry(i); + for (auto& track : *tracks) { + auto pdg = track.GetPdgCode(); + auto rapidity = track.GetRapidity(); + auto idMoth = track.getMotherTrackId(); + if (pdg == checkPdgDecay) { + // count leptons + nLeptons++; + } else if(pdg == -checkPdgDecay) { + // count anti-leptons + nAntileptons++; + } else if (pdg == checkPdgSignal[0] || pdg == checkPdgSignal[1] || pdg == checkPdgSignal[2]) { + if(idMoth < 0){ + // count signal PDG + if (pdg == checkPdgSignal[0]) { + nSignalUpsilon1S++; + } else if (pdg == checkPdgSignal[1]) { + nSignalUpsilon2S++; + } else { + nSignalUpsilon3S++; + } + + // count signal PDG within acceptance + if (rapidity > rapiditymin && rapidity < rapiditymax) { + if (pdg == checkPdgSignal[0]) { + nSignalUpsilon1SWithinAcc++; + } else if (pdg == checkPdgSignal[1]) { + nSignalUpsilon2SWithinAcc++; + } else { + nSignalUpsilon3SWithinAcc++; + } + } + } + auto child0 = o2::mcutils::MCTrackNavigator::getDaughter0(track, *tracks); + auto child1 = o2::mcutils::MCTrackNavigator::getDaughter1(track, *tracks); + if (child0 != nullptr && child1 != nullptr) { + // check for parent-child relations + auto pdg0 = child0->GetPdgCode(); + auto pdg1 = child1->GetPdgCode(); + std::cout << "First and last children of parent " << checkPdgSignal << " are PDG0: " << pdg0 << " PDG1: " << pdg1 << "\n"; + if (std::abs(pdg0) == checkPdgDecay && std::abs(pdg1) == checkPdgDecay && pdg0 == -pdg1) { + nLeptonPairs++; + if (child0->getToBeDone() && child1->getToBeDone()) { + nLeptonPairsToBeDone++; + } + } + } + } + } + } + std::cout << "#events: " << nEvents << "\n" + << "#leptons: " << nLeptons << "\n" + << "#antileptons: " << nAntileptons << "\n" + << "#signal (Upsilon(1S)): " << nSignalUpsilon1S << "; within acceptance " << rapiditymin << " < y < " << rapiditymax << " : " << nSignalUpsilon1SWithinAcc << "\n" + << "#signal (Upsilon(2S)): " << nSignalUpsilon2S << "; within acceptance " << rapiditymin << " < y < " << rapiditymax << " : " << nSignalUpsilon2SWithinAcc << "\n" + << "#signal (Upsilon(2S)): " << nSignalUpsilon3S << "; within acceptance " << rapiditymin << " < y < " << rapiditymax << " : " << nSignalUpsilon3SWithinAcc << "\n" + << "#lepton pairs: " << nLeptonPairs << "\n" + << "#lepton pairs to be done: " << nLeptonPairs << "\n"; + + + if (nLeptonPairs == 0 || nLeptons == 0 || nAntileptons == 0) { + std::cerr << "Number of leptons, number of anti-leptons as well as number of lepton pairs should all be greater than 1.\n"; + return 1; + } + if (nLeptonPairs != nLeptonPairsToBeDone) { + std::cerr << "The number of lepton pairs should be the same as the number of lepton pairs which should be transported.\n"; + return 1; + } + + return 0; +} \ No newline at end of file diff --git a/MC/config/PWGDQ/ini/tests/Generator_InjectedPromptPsi2SToJpsiPiPiMidy_TriggerGap.C b/MC/config/PWGDQ/ini/tests/Generator_InjectedPromptPsi2SToJpsiPiPiMidy_TriggerGap.C new file mode 100644 index 000000000..9254fdab1 --- /dev/null +++ b/MC/config/PWGDQ/ini/tests/Generator_InjectedPromptPsi2SToJpsiPiPiMidy_TriggerGap.C @@ -0,0 +1,99 @@ +int External() +{ + int checkPdgSignal[] = {100443}; + int checkPdgDecay[] = {443, 211, -211}; + int leptonPdg = 11; + Double_t rapidityWindow = 1.0; + std::string path{"o2sim_Kine.root"}; + std::cout << "Check for\nsignal PDG " << checkPdgSignal[0] << "\n decay PDG " << checkPdgDecay[0] << ", " << checkPdgDecay[1] << ", " << checkPdgDecay[2] << "\n"; + TFile file(path.c_str(), "READ"); + if (file.IsZombie()) { + std::cerr << "Cannot open ROOT file " << path << "\n"; + return 1; + } + auto tree = (TTree*)file.Get("o2sim"); + std::vector* tracks{}; + tree->SetBranchAddress("MCTrack", &tracks); + + int nLeptons{}; + int nAntileptons{}; + int nLeptonPairs{}; + int nLeptonPairsToBeDone{}; + int nSignalJpsi{}; + int nSignalPionsPos{}; + int nSignalPionsNeg{}; + int nSignalPsi2S{}; + int nSignalJpsiWithinAcc{}; + int nSignalPionsPosWithinAcc{}; + int nSignalPionsNegWithinAcc{}; + auto nEvents = tree->GetEntries(); + o2::steer::MCKinematicsReader mcreader("o2sim", o2::steer::MCKinematicsReader::Mode::kMCKine); + Bool_t hasPsi2SMoth = kFALSE; + + for (int i = 0; i < nEvents; i++) { + tree->GetEntry(i); + for (auto& track : *tracks) { + auto pdg = track.GetPdgCode(); + auto rapidity = track.GetRapidity(); + auto idMoth = track.getMotherTrackId(); + int idJpsi = -1; int IdChild0 = -1; int IdChild1 = -1; + if (pdg == leptonPdg) { + // count leptons + nLeptons++; + } else if(pdg == -leptonPdg) { + // count anti-leptons + nAntileptons++; + } else if (pdg == checkPdgSignal[0]) { + // check daughters + std::cout << "Signal PDG: " << pdg << "\n"; + for (int j{track.getFirstDaughterTrackId()}; j <= track.getLastDaughterTrackId(); ++j) { + auto pdgDau = tracks->at(j).GetPdgCode(); + std::cout << "Daughter " << j << " is: " << pdgDau << "\n"; + if(TMath::Abs(pdgDau) == checkPdgDecay[0] ) { nSignalJpsi++; if( std::abs(track.GetRapidity()) < rapidityWindow) nSignalJpsiWithinAcc++; idJpsi = j; } + if(pdgDau == checkPdgDecay[1] ) { nSignalPionsPos++; if( std::abs(track.GetRapidity()) < rapidityWindow) nSignalPionsPosWithinAcc++; } + if(pdgDau == checkPdgDecay[2] ) { nSignalPionsNeg++; if( std::abs(track.GetRapidity()) < rapidityWindow) nSignalPionsNegWithinAcc++; } + } + + auto trackJpsi = tracks->at(idJpsi); + for (int j{trackJpsi.getFirstDaughterTrackId()}; j <= trackJpsi.getLastDaughterTrackId(); ++j) { + auto pdgDau = tracks->at(j).GetPdgCode(); + if(pdgDau == leptonPdg ) IdChild0 = j; + if(pdgDau == -leptonPdg ) IdChild1 = j; + } + auto child0 = tracks->at(IdChild0); + auto child1 = tracks->at(IdChild1); + // check for parent-child relations + auto pdg0 = child0.GetPdgCode(); + auto pdg1 = child1.GetPdgCode(); + std::cout << "Lepton daughter particles of mother " << trackJpsi.GetPdgCode() << " are PDG0: " << pdg0 << " PDG1: " << pdg1 << "\n"; + if (std::abs(pdg0) == leptonPdg && std::abs(pdg1) == leptonPdg && pdg0 == -pdg1) { + nLeptonPairs++; + if (child0.getToBeDone() && child1.getToBeDone()) { + nLeptonPairsToBeDone++; + } + } + //} + } + } + } + std::cout << "#events: " << nEvents << "\n" + << "#leptons: " << nLeptons << "\n" + << "#antileptons: " << nAntileptons << "\n" + << "#signal (jpsi <- psi2S): " << nSignalJpsi << "; within acceptance (|y| < " << rapidityWindow << "): " << nSignalJpsiWithinAcc << "\n" + << "#signal (pi+ <- psi2S): " << nSignalPionsPos << "; within acceptance (|y| < " << rapidityWindow << "): " << nSignalPionsPosWithinAcc << "\n" + << "#signal (pi- <- psi2S): " << nSignalPionsNeg << "; within acceptance (|y| < " << rapidityWindow << "): " << nSignalPionsNegWithinAcc << "\n" + << "#lepton pairs: " << nLeptonPairs << "\n" + << "#lepton pairs to be done: " << nLeptonPairs << "\n"; + + + if (nLeptonPairs == 0 || nLeptons == 0 || nAntileptons == 0) { + std::cerr << "Number of leptons, number of anti-leptons as well as number of lepton pairs should all be greater than 1.\n"; + return 1; + } + if (nLeptonPairs != nLeptonPairsToBeDone) { + std::cerr << "The number of lepton pairs should be the same as the number of lepton pairs which should be transported.\n"; + return 1; + } + + return 0; +} diff --git a/MC/config/PWGDQ/ini/tests/Generator_InjectedPromptX3782Midy_TriggerGap.C b/MC/config/PWGDQ/ini/tests/Generator_InjectedPromptX3782Midy_TriggerGap.C new file mode 100644 index 000000000..836c8df30 --- /dev/null +++ b/MC/config/PWGDQ/ini/tests/Generator_InjectedPromptX3782Midy_TriggerGap.C @@ -0,0 +1,99 @@ +int External() +{ + int checkPdgSignal[] = {9920443}; // pdg code X3872 + int checkPdgDecay[] = {443, 211, -211}; + int leptonPdg = 11; + Double_t rapidityWindow = 1.0; + std::string path{"o2sim_Kine.root"}; + std::cout << "Check for\nsignal PDG " << checkPdgSignal[0] << "\n decay PDG " << checkPdgDecay[0] << ", " << checkPdgDecay[1] << ", " << checkPdgDecay[2] << "\n"; + TFile file(path.c_str(), "READ"); + if (file.IsZombie()) { + std::cerr << "Cannot open ROOT file " << path << "\n"; + return 1; + } + auto tree = (TTree*)file.Get("o2sim"); + std::vector* tracks{}; + tree->SetBranchAddress("MCTrack", &tracks); + + int nLeptons{}; + int nAntileptons{}; + int nLeptonPairs{}; + int nLeptonPairsToBeDone{}; + int nSignalX3872{}; + int nSignalPionsPos{}; + int nSignalPionsNeg{}; + int nSignalPsi2S{}; + int nSignalX3872WithinAcc{}; + int nSignalPionsPosWithinAcc{}; + int nSignalPionsNegWithinAcc{}; + auto nEvents = tree->GetEntries(); + o2::steer::MCKinematicsReader mcreader("o2sim", o2::steer::MCKinematicsReader::Mode::kMCKine); + Bool_t hasPsi2SMoth = kFALSE; + + for (int i = 0; i < nEvents; i++) { + tree->GetEntry(i); + for (auto& track : *tracks) { + auto pdg = track.GetPdgCode(); + auto rapidity = track.GetRapidity(); + auto idMoth = track.getMotherTrackId(); + int idX3872 = -1; int IdChild0 = -1; int IdChild1 = -1; + if (pdg == leptonPdg) { + // count leptons + nLeptons++; + } else if(pdg == -leptonPdg) { + // count anti-leptons + nAntileptons++; + } else if (pdg == checkPdgSignal[0]) { + // check daughters + std::cout << "Signal PDG: " << pdg << "\n"; + for (int j{track.getFirstDaughterTrackId()}; j <= track.getLastDaughterTrackId(); ++j) { + auto pdgDau = tracks->at(j).GetPdgCode(); + std::cout << "Daughter " << j << " is: " << pdgDau << "\n"; + if(TMath::Abs(pdgDau) == checkPdgDecay[0] ) { nSignalX3872++; if( std::abs(track.GetRapidity()) < rapidityWindow) nSignalX3872WithinAcc++; idX3872 = j; } + if(pdgDau == checkPdgDecay[1] ) { nSignalPionsPos++; if( std::abs(track.GetRapidity()) < rapidityWindow) nSignalPionsPosWithinAcc++; } + if(pdgDau == checkPdgDecay[2] ) { nSignalPionsNeg++; if( std::abs(track.GetRapidity()) < rapidityWindow) nSignalPionsNegWithinAcc++; } + } + + auto trackX3872 = tracks->at(idX3872); + for (int j{trackX3872.getFirstDaughterTrackId()}; j <= trackX3872.getLastDaughterTrackId(); ++j) { + auto pdgDau = tracks->at(j).GetPdgCode(); + if(pdgDau == leptonPdg ) IdChild0 = j; + if(pdgDau == -leptonPdg ) IdChild1 = j; + } + auto child0 = tracks->at(IdChild0); + auto child1 = tracks->at(IdChild1); + // check for parent-child relations + auto pdg0 = child0.GetPdgCode(); + auto pdg1 = child1.GetPdgCode(); + std::cout << "Lepton daughter particles of mother " << trackX3872.GetPdgCode() << " are PDG0: " << pdg0 << " PDG1: " << pdg1 << "\n"; + if (std::abs(pdg0) == leptonPdg && std::abs(pdg1) == leptonPdg && pdg0 == -pdg1) { + nLeptonPairs++; + if (child0.getToBeDone() && child1.getToBeDone()) { + nLeptonPairsToBeDone++; + } + } + //} + } + } + } + std::cout << "#events: " << nEvents << "\n" + << "#leptons: " << nLeptons << "\n" + << "#antileptons: " << nAntileptons << "\n" + << "#signal (jpsi <- X3872): " << nSignalX3872 << "; within acceptance (|y| < " << rapidityWindow << "): " << nSignalX3872WithinAcc << "\n" + << "#signal (pi+ <- X3872): " << nSignalPionsPos << "; within acceptance (|y| < " << rapidityWindow << "): " << nSignalPionsPosWithinAcc << "\n" + << "#signal (pi- <- X3872): " << nSignalPionsNeg << "; within acceptance (|y| < " << rapidityWindow << "): " << nSignalPionsNegWithinAcc << "\n" + << "#lepton pairs: " << nLeptonPairs << "\n" + << "#lepton pairs to be done: " << nLeptonPairs << "\n"; + + + if (nLeptonPairs == 0 || nLeptons == 0 || nAntileptons == 0) { + std::cerr << "Number of leptons, number of anti-leptons as well as number of lepton pairs should all be greater than 1.\n"; + return 1; + } + if (nLeptonPairs != nLeptonPairsToBeDone) { + std::cerr << "The number of lepton pairs should be the same as the number of lepton pairs which should be transported.\n"; + return 1; + } + + return 0; +} diff --git a/MC/config/PWGEM/decaytables/decaytable_LMee.dat b/MC/config/PWGEM/decaytables/decaytable_LMee.dat index f35b10657..3dbbc5a76 100644 --- a/MC/config/PWGEM/decaytables/decaytable_LMee.dat +++ b/MC/config/PWGEM/decaytables/decaytable_LMee.dat @@ -277,7 +277,8 @@ 111 pi0 0 0 0 0.13498 0.000000 0.00000 2.56000E-05 0 1 1 0 0.988230 22 22 0 0 0 1 2 0.011740 22 -11 11 0 0 - 1 0 0.000033 -11 11 -11 11 0 + 1 0 0.000000 -11 11 -11 11 0 + 1 0 0.000000 -11 11 0 0 0 221 eta 0 0 0 0.54786 0.000001 0.00000 1.50886E-07 0 1 1 0 0.394074 22 22 0 0 0 1 0 0.326777 111 111 111 0 0 @@ -287,8 +288,8 @@ 1 2 0.006899 22 11 -11 0 0 1 0 0.000309 22 13 -13 0 0 1 0 0.000006 -13 13 0 0 0 - 1 0 0.000024 -11 11 -11 11 0 - 1 0 0.000268 211 -211 11 -11 0 + 1 0 0.000000 -11 11 -11 11 0 + 1 0 0.000000 211 -211 11 -11 0 331 eta' 0 0 0 0.95778 0.000197 0.002 1.00336E-09 0 1 1 0 0.428744 211 -211 221 0 0 1 0 0.290822 22 113 0 0 0 @@ -301,8 +302,9 @@ 1 0 0.003794 211 -211 111 0 0 1 0 0.000001 211 -211 211 -211 0 1 0 0.000180 211 -211 111 111 0 - 1 0 0.002384 211 -211 -11 11 0 + 1 0 0.000000 211 -211 -11 11 0 1 2 0.000473 22 11 -11 0 0 + 1 0 0.000000 -11 11 -11 11 0 113 rho0 0 0 0 0.77526 0.149100 0.40000 1.33000E-12 0 1 1 3 0.988927 211 -211 0 0 0 1 0 0.009900 211 -211 22 0 0 diff --git a/MC/config/PWGEM/external/generator/GeneratorEMCocktailV2.C b/MC/config/PWGEM/external/generator/GeneratorEMCocktailV2.C index 13770644d..0511ecaf7 100644 --- a/MC/config/PWGEM/external/generator/GeneratorEMCocktailV2.C +++ b/MC/config/PWGEM/external/generator/GeneratorEMCocktailV2.C @@ -435,11 +435,11 @@ GenerateEMCocktail(Int_t collisionsSystem = GeneratorParamEMlibV2::kpp7TeV, Double_t yGenRange = 0.1, TString useLMeeDecaytable = "", Int_t weightingMode = 1) { - TString O2DPG_ROOT = TString(getenv("O2DPG_ROOT")); - paramFile=paramFile.ReplaceAll("$O2DPG_ROOT",O2DPG_ROOT); - paramFile=paramFile.ReplaceAll("${O2DPG_ROOT}",O2DPG_ROOT); - useLMeeDecaytable=useLMeeDecaytable.ReplaceAll("$O2DPG_ROOT",O2DPG_ROOT); - useLMeeDecaytable=useLMeeDecaytable.ReplaceAll("${O2DPG_ROOT}",O2DPG_ROOT); + TString O2DPG_MC_CONFIG_ROOT = TString(getenv("O2DPG_MC_CONFIG_ROOT")); + paramFile=paramFile.ReplaceAll("$O2DPG_MC_CONFIG_ROOT",O2DPG_MC_CONFIG_ROOT); + paramFile=paramFile.ReplaceAll("${O2DPG_MC_CONFIG_ROOT}",O2DPG_MC_CONFIG_ROOT); + useLMeeDecaytable=useLMeeDecaytable.ReplaceAll("$O2DPG_MC_CONFIG_ROOT",O2DPG_MC_CONFIG_ROOT); + useLMeeDecaytable=useLMeeDecaytable.ReplaceAll("${O2DPG_MC_CONFIG_ROOT}",O2DPG_MC_CONFIG_ROOT); if (paramFile.BeginsWith("alien://")){ TGrid::Connect("alien://"); } diff --git a/MC/config/PWGEM/ini/GeneratorEMCocktailwithFlow3040.ini b/MC/config/PWGEM/ini/GeneratorEMCocktailwithFlow3040.ini new file mode 100644 index 000000000..d91e2e9ea --- /dev/null +++ b/MC/config/PWGEM/ini/GeneratorEMCocktailwithFlow3040.ini @@ -0,0 +1,6 @@ +### The setup uses an external event generator +### This part sets the path of the file and the function call to retrieve it + +[GeneratorExternal] +fileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGEM/external/generator/GeneratorEMCocktailV2.C +funcName=GenerateEMCocktail(400,0,3,63,"${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGEM/parametrizations/PbPb5TeV_central.json","5TeV_3040_wRatio",350,0.0,30.0,10000,1,1,0,0,"5TeV_3040_wRatio",0,1.1,"${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGEM/decaytables/decaytable_LMee.dat",1) diff --git a/MC/config/PWGEM/ini/GeneratorEMCocktailwithFlow4050.ini b/MC/config/PWGEM/ini/GeneratorEMCocktailwithFlow4050.ini new file mode 100644 index 000000000..f97517bab --- /dev/null +++ b/MC/config/PWGEM/ini/GeneratorEMCocktailwithFlow4050.ini @@ -0,0 +1,6 @@ +### The setup uses an external event generator +### This part sets the path of the file and the function call to retrieve it + +[GeneratorExternal] +fileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGEM/external/generator/GeneratorEMCocktailV2.C +funcName=GenerateEMCocktail(400,0,3,63,"${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGEM/parametrizations/PbPb5TeV_central.json","5TeV_4050_wRatio",350,0.0,30.0,10000,1,1,0,0,"5TeV_4050_wRatio",0,1.1,"${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGEM/decaytables/decaytable_LMee.dat",1) diff --git a/MC/config/PWGEM/ini/Pythia8_Beauty_Cocktail.ini b/MC/config/PWGEM/ini/Pythia8_Beauty_Cocktail.ini index 3fcdf4660..c40d93ef6 100644 --- a/MC/config/PWGEM/ini/Pythia8_Beauty_Cocktail.ini +++ b/MC/config/PWGEM/ini/Pythia8_Beauty_Cocktail.ini @@ -1,11 +1,11 @@ [GeneratorExternal] fileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGEM/external/generator/Generator_pythia8_forcedDecays.C -funcName=GeneratePythia8ForcedDecays("411;421;431;4122;4232;4132;4332;") +funcName=GeneratePythia8ForcedDecays("411;421;431;4122;4232;4132;4332") [GeneratorPythia8] config = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGEM/pythia8/generator/pythia8_hf_cocktail.cfg hooksFileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/hooks/pythia8_userhooks_qqbar.C -hooksFuncName = pythia8_userhooks_bbbar(-9999.,9999.) +hooksFuncName = pythia8_userhooks_bbbar(-12.,12.) includePartonEvent=true [DecayerPythia8] diff --git a/MC/config/PWGEM/ini/Pythia8_Charm_Cocktail.ini b/MC/config/PWGEM/ini/Pythia8_Charm_Cocktail.ini index d0b05b8c9..ad7392644 100644 --- a/MC/config/PWGEM/ini/Pythia8_Charm_Cocktail.ini +++ b/MC/config/PWGEM/ini/Pythia8_Charm_Cocktail.ini @@ -5,7 +5,7 @@ funcName=GeneratePythia8ForcedDecays("411;421;431;4122;4232;4132;4332") [GeneratorPythia8] config = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGEM/pythia8/generator/pythia8_hf_cocktail.cfg hooksFileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/hooks/pythia8_userhooks_qqbar.C -hooksFuncName = pythia8_userhooks_ccbar(-9999.,9999.) +hooksFuncName = pythia8_userhooks_ccbar(-12.,12.) includePartonEvent=true [DecayerPythia8] diff --git a/MC/config/PWGEM/parametrizations/PbPb5TeV_central.json b/MC/config/PWGEM/parametrizations/PbPb5TeV_central.json index 85e0fec17..6cc9cb0bd 100644 --- a/MC/config/PWGEM/parametrizations/PbPb5TeV_central.json +++ b/MC/config/PWGEM/parametrizations/PbPb5TeV_central.json @@ -27,7 +27,10 @@ }, "111_pt": "TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841))", "221_pt": "(TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841)))*((0.317347/(1+exp(-(x-0.475129)/0.245111))+-0.265641*TMath::Gaus(x,0.276828,1.17016)+-0.0289536*TMath::Gaus(x,5.40691,2.58911)+0.249861+((0.550454*exp((((0.995785*x)-sqrt((x*x)+(0.5478530000*0.5478530000)))/sqrt(1-(0.995785*0.995785)))/1.71099))+(0.0106936*(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/(exp((((0.995785*x)-sqrt((x*x)+(0.1349770000*0.1349770000)))/sqrt(1-(0.995785*0.995785)))/1.71099)+(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/2.)", - "333_pt": "(TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841)))*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071)" + "333_pt": "(TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841)))*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071)", + "111_v2_def": "x*x*(0.29722*exp(-x/0.318227)+0.560172*pow(1+x/0.401715,-2.04434)+0.193576*exp(-x/1.2557))", + "310_v2_def": "x*x*(0.0280321*exp(-x/0.49621)+0.328197*pow(1+x/0.0128527,-0.853155)+0.214485*exp(-x/1.29581))", + "333_v2_def": "x*x*(1.2367*exp(-x/0.629473)+-3.9919*pow(1+x/0.886275,-3.11899)+0.274776*exp(-x/1.82705))" }, "5TeV_4050_wRatio": { "histoMtScaleFactor": { @@ -37,7 +40,10 @@ }, "111_pt": "TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313))", "221_pt": "(TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313)))*((0.298649/(1+exp(-(x-0.494749)/0.23628))+-0.255064*TMath::Gaus(x,0.170043,1.31558)+-0.024999*TMath::Gaus(x,8.22033,2.71373)+0.251236+((0.550454*exp((((0.995785*x)-sqrt((x*x)+(0.5478530000*0.5478530000)))/sqrt(1-(0.995785*0.995785)))/1.71099))+(0.0106936*(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/(exp((((0.995785*x)-sqrt((x*x)+(0.1349770000*0.1349770000)))/sqrt(1-(0.995785*0.995785)))/1.71099)+(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/2.)", - "333_pt": "(TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313)))*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985)" + "333_pt": "(TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313)))*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985)", + "111_v2_def": "x*x*(0.171282*exp(-x/0.353318)+1.13749*pow(1+x/0.290699,-2.02699)+0.20604*exp(-x/1.19178))", + "310_v2_def": "x*x*(-0.524086*exp(-x/0.402378)+0.891134*pow(1+x/1.15413,-3.07832)+0.0812177*exp(-x/1.22691))", + "333_v2_def": "x*x*(-2.78876*exp(-x/0.581864)+3.6669*pow(1+x/3.01301,-6.60505)+0.0436791*exp(-x/1.23285))" }, "5TeV_4050_wRatio_pi0up": { "histoMtScaleFactor": { @@ -47,7 +53,10 @@ }, "111_pt": "TMath::TwoPi()*x*(10.1537*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.423934)+708.793*pow(exp(-0.595237*x-0.468542*x*x)+x/0.502946,-5.71377))", "221_pt": "(TMath::TwoPi()*x*(10.1537*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.423934)+708.793*pow(exp(-0.595237*x-0.468542*x*x)+x/0.502946,-5.71377)))*((0.298649/(1+exp(-(x-0.494749)/0.23628))+-0.255064*TMath::Gaus(x,0.170043,1.31558)+-0.024999*TMath::Gaus(x,8.22033,2.71373)+0.251236+((0.550454*exp((((0.995785*x)-sqrt((x*x)+(0.5478530000*0.5478530000)))/sqrt(1-(0.995785*0.995785)))/1.71099))+(0.0106936*(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/(exp((((0.995785*x)-sqrt((x*x)+(0.1349770000*0.1349770000)))/sqrt(1-(0.995785*0.995785)))/1.71099)+(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/2.)", - "333_pt": "(TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313)))*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985)" + "333_pt": "(TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313)))*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985)", + "111_v2_def": "x*x*(0.171282*exp(-x/0.353318)+1.13749*pow(1+x/0.290699,-2.02699)+0.20604*exp(-x/1.19178))", + "310_v2_def": "x*x*(-0.524086*exp(-x/0.402378)+0.891134*pow(1+x/1.15413,-3.07832)+0.0812177*exp(-x/1.22691))", + "333_v2_def": "x*x*(-2.78876*exp(-x/0.581864)+3.6669*pow(1+x/3.01301,-6.60505)+0.0436791*exp(-x/1.23285))" }, "5TeV_3040_wRatio_pi0up": { "histoMtScaleFactor": { @@ -57,7 +66,10 @@ }, "111_pt": "TMath::TwoPi()*x*(26.8208*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.408437)+1019.42*pow(exp(-0.662498*x-0.498822*x*x)+x/0.495404,-5.68121))", "221_pt": "(TMath::TwoPi()*x*(26.8208*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.408437)+1019.42*pow(exp(-0.662498*x-0.498822*x*x)+x/0.495404,-5.68121)))*((0.317347/(1+exp(-(x-0.475129)/0.245111))+-0.265641*TMath::Gaus(x,0.276828,1.17016)+-0.0289536*TMath::Gaus(x,5.40691,2.58911)+0.249861+((0.550454*exp((((0.995785*x)-sqrt((x*x)+(0.5478530000*0.5478530000)))/sqrt(1-(0.995785*0.995785)))/1.71099))+(0.0106936*(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/(exp((((0.995785*x)-sqrt((x*x)+(0.1349770000*0.1349770000)))/sqrt(1-(0.995785*0.995785)))/1.71099)+(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/2.)", - "333_pt": "(TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841)))*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071)" + "333_pt": "(TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841)))*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071)", + "111_v2_def": "x*x*(0.29722*exp(-x/0.318227)+0.560172*pow(1+x/0.401715,-2.04434)+0.193576*exp(-x/1.2557))", + "310_v2_def": "x*x*(0.0280321*exp(-x/0.49621)+0.328197*pow(1+x/0.0128527,-0.853155)+0.214485*exp(-x/1.29581))", + "333_v2_def": "x*x*(1.2367*exp(-x/0.629473)+-3.9919*pow(1+x/0.886275,-3.11899)+0.274776*exp(-x/1.82705))" }, "5TeV_0510_wRatio_pi0up": { "histoMtScaleFactor": { @@ -107,7 +119,10 @@ }, "111_pt": "TMath::TwoPi()*x*(22.2048*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.41549)+802.356*pow(exp(-0.814916*x-0.343908*x*x)+x/0.515668,-5.69624))", "221_pt": "(TMath::TwoPi()*x*(22.2048*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.41549)+802.356*pow(exp(-0.814916*x-0.343908*x*x)+x/0.515668,-5.69624)))*((0.317347/(1+exp(-(x-0.475129)/0.245111))+-0.265641*TMath::Gaus(x,0.276828,1.17016)+-0.0289536*TMath::Gaus(x,5.40691,2.58911)+0.249861+((0.550454*exp((((0.995785*x)-sqrt((x*x)+(0.5478530000*0.5478530000)))/sqrt(1-(0.995785*0.995785)))/1.71099))+(0.0106936*(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/(exp((((0.995785*x)-sqrt((x*x)+(0.1349770000*0.1349770000)))/sqrt(1-(0.995785*0.995785)))/1.71099)+(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/2.)", - "333_pt": "(TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841)))*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071)" + "333_pt": "(TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841)))*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071)", + "111_v2_def": "x*x*(0.29722*exp(-x/0.318227)+0.560172*pow(1+x/0.401715,-2.04434)+0.193576*exp(-x/1.2557))", + "310_v2_def": "x*x*(0.0280321*exp(-x/0.49621)+0.328197*pow(1+x/0.0128527,-0.853155)+0.214485*exp(-x/1.29581))", + "333_v2_def": "x*x*(1.2367*exp(-x/0.629473)+-3.9919*pow(1+x/0.886275,-3.11899)+0.274776*exp(-x/1.82705))" }, "5TeV_4050_wRatio_pi0down": { "histoMtScaleFactor": { @@ -117,7 +132,10 @@ }, "111_pt": "TMath::TwoPi()*x*(8.04609*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.432276)+566.92*pow(exp(-0.727741*x-0.316795*x*x)+x/0.523032,-5.73361))", "221_pt": "(TMath::TwoPi()*x*(8.04609*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.432276)+566.92*pow(exp(-0.727741*x-0.316795*x*x)+x/0.523032,-5.73361)))*((0.298649/(1+exp(-(x-0.494749)/0.23628))+-0.255064*TMath::Gaus(x,0.170043,1.31558)+-0.024999*TMath::Gaus(x,8.22033,2.71373)+0.251236+((0.550454*exp((((0.995785*x)-sqrt((x*x)+(0.5478530000*0.5478530000)))/sqrt(1-(0.995785*0.995785)))/1.71099))+(0.0106936*(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/(exp((((0.995785*x)-sqrt((x*x)+(0.1349770000*0.1349770000)))/sqrt(1-(0.995785*0.995785)))/1.71099)+(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/2.)", - "333_pt": "(TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313)))*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985)" + "333_pt": "(TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313)))*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985)", + "111_v2_def": "x*x*(0.171282*exp(-x/0.353318)+1.13749*pow(1+x/0.290699,-2.02699)+0.20604*exp(-x/1.19178))", + "310_v2_def": "x*x*(-0.524086*exp(-x/0.402378)+0.891134*pow(1+x/1.15413,-3.07832)+0.0812177*exp(-x/1.22691))", + "333_v2_def": "x*x*(-2.78876*exp(-x/0.581864)+3.6669*pow(1+x/3.01301,-6.60505)+0.0436791*exp(-x/1.23285))" }, "5TeV_4050_wRatio_ratiosdown": { "histoMtScaleFactor": { @@ -127,7 +145,10 @@ }, "111_pt": "TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313))", "221_pt": "(TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313)))*(((0.550454*exp((((0.995785*x)-sqrt((x*x)+(0.5478530000*0.5478530000)))/sqrt(1-(0.995785*0.995785)))/1.71099))+(0.0106936*(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/(exp((((0.995785*x)-sqrt((x*x)+(0.1349770000*0.1349770000)))/sqrt(1-(0.995785*0.995785)))/1.71099)+(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))", - "333_pt": "(TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313)))*((x<=0.5)*(0.3/0.5*x+0.5)*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985)+(x>0.5)*0.8*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985))" + "333_pt": "(TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313)))*((x<=0.5)*(0.3/0.5*x+0.5)*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985)+(x>0.5)*0.8*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985))", + "111_v2_def": "x*x*(0.171282*exp(-x/0.353318)+1.13749*pow(1+x/0.290699,-2.02699)+0.20604*exp(-x/1.19178))", + "310_v2_def": "x*x*(-0.524086*exp(-x/0.402378)+0.891134*pow(1+x/1.15413,-3.07832)+0.0812177*exp(-x/1.22691))", + "333_v2_def": "x*x*(-2.78876*exp(-x/0.581864)+3.6669*pow(1+x/3.01301,-6.60505)+0.0436791*exp(-x/1.23285))" }, "5TeV_3040_wRatio_ratiosdown": { "histoMtScaleFactor": { @@ -137,7 +158,10 @@ }, "111_pt": "TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841))", "221_pt": "(TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841)))*(((0.550454*exp((((0.995785*x)-sqrt((x*x)+(0.5478530000*0.5478530000)))/sqrt(1-(0.995785*0.995785)))/1.71099))+(0.0106936*(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))/(exp((((0.995785*x)-sqrt((x*x)+(0.1349770000*0.1349770000)))/sqrt(1-(0.995785*0.995785)))/1.71099)+(0.620714*TMath::Power(1+((x/29.189)*(x/29.189)),-(506.895)))))", - "333_pt": "(TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841)))*((x<=0.5)*(0.3/0.5*x+0.5)*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071)+(x>0.5)*0.8*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071))" + "333_pt": "(TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841)))*((x<=0.5)*(0.3/0.5*x+0.5)*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071)+(x>0.5)*0.8*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071))", + "111_v2_def": "x*x*(0.29722*exp(-x/0.318227)+0.560172*pow(1+x/0.401715,-2.04434)+0.193576*exp(-x/1.2557))", + "310_v2_def": "x*x*(0.0280321*exp(-x/0.49621)+0.328197*pow(1+x/0.0128527,-0.853155)+0.214485*exp(-x/1.29581))", + "333_v2_def": "x*x*(1.2367*exp(-x/0.629473)+-3.9919*pow(1+x/0.886275,-3.11899)+0.274776*exp(-x/1.82705))" }, "5TeV_0510_wRatio_ratiosdown": { "histoMtScaleFactor": { @@ -187,7 +211,10 @@ }, "111_pt": "TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841))", "221_pt": "(TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841)))*(0.317347/(1+exp(-(x-0.475129)/0.245111))+-0.265641*TMath::Gaus(x,0.276828,1.17016)+-0.0289536*TMath::Gaus(x,5.40691,2.58911)+0.249861)", - "333_pt": "(TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841)))*((x<=0.5)*(-0.3/0.5*x+1.5)*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071)+(x>0.5)*1.2*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071))" + "333_pt": "(TMath::TwoPi()*x*(24.3786*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.412015)+908.266*pow(exp(-0.735883*x-0.425837*x*x)+x/0.505017,-5.68841)))*((x<=0.5)*(-0.3/0.5*x+1.5)*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071)+(x>0.5)*1.2*(0.106572/(1+exp(-(x-1.71328)/1.68645))+-0.54012*TMath::Gaus(x,0.754082,1.88026)+0.452996*TMath::Gaus(x,1.49283,2.4567)+0.0943071))", + "111_v2_def": "x*x*(0.29722*exp(-x/0.318227)+0.560172*pow(1+x/0.401715,-2.04434)+0.193576*exp(-x/1.2557))", + "310_v2_def": "x*x*(0.0280321*exp(-x/0.49621)+0.328197*pow(1+x/0.0128527,-0.853155)+0.214485*exp(-x/1.29581))", + "333_v2_def": "x*x*(1.2367*exp(-x/0.629473)+-3.9919*pow(1+x/0.886275,-3.11899)+0.274776*exp(-x/1.82705))" }, "5TeV_4050_wRatio_ratiosup": { "histoMtScaleFactor": { @@ -197,7 +224,10 @@ }, "111_pt": "TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313))", "221_pt": "(TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313)))*(0.298649/(1+exp(-(x-0.494749)/0.23628))+-0.255064*TMath::Gaus(x,0.170043,1.31558)+-0.024999*TMath::Gaus(x,8.22033,2.71373)+0.251236)", - "333_pt": "(TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313)))*((x<=0.5)*(-0.3/0.5*x+1.5)*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985)+(x>0.5)*1.2*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985))" + "333_pt": "(TMath::TwoPi()*x*(9.00063*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.428379)+636.492*pow(exp(-0.659016*x-0.397367*x*x)+x/0.512428,-5.72313)))*((x<=0.5)*(-0.3/0.5*x+1.5)*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985)+(x>0.5)*1.2*(0.0979137/(1+exp(-(x-1.53259)/0.370136))+-0.537275*TMath::Gaus(x,0.706186,2.03543)+0.45201*TMath::Gaus(x,1.18876,2.40911)+0.104985))", + "111_v2_def": "x*x*(0.171282*exp(-x/0.353318)+1.13749*pow(1+x/0.290699,-2.02699)+0.20604*exp(-x/1.19178))", + "310_v2_def": "x*x*(-0.524086*exp(-x/0.402378)+0.891134*pow(1+x/1.15413,-3.07832)+0.0812177*exp(-x/1.22691))", + "333_v2_def": "x*x*(-2.78876*exp(-x/0.581864)+3.6669*pow(1+x/3.01301,-6.60505)+0.0436791*exp(-x/1.23285))" }, "5TeV_0005_wRatio_etatest": { "histoMtScaleFactor": { @@ -399,4 +429,4 @@ "221_pt": "(TMath::TwoPi()*x*(219.382*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.384303)+2471.6*pow(exp(-1.12373*x-0.752507*x*x)+x/0.438343,-5.45943)))*((((x<9.445400)*(TMath::Max(TMath::Max(1.68883e-09*TMath::Power((exp(- -1.16032*sqrt(x*x+0.547*0.547-0.139*0.139)-0.435153*(x*x+0.547*0.547-0.139*0.139))+TMath::Power(0.480215/1.68883e-09,-1./16.7077)*sqrt(x*x+0.547*0.547-0.139*0.139)/0.0356546),-16.7077)/TMath::Power((exp(- -1.16032*x-0.435153*x*x)+x/0.0356546),-16.7077)+1.75797*TMath::Landau(x,2.41484,1.01421),0.00469052*TMath::Power((exp(- -0.124633*sqrt(x*x+0.547*0.547-0.139*0.139)-0.310906*(x*x+0.547*0.547-0.139*0.139))+TMath::Power(0.479677/0.00469052,-1./0.657178)*sqrt(x*x+0.547*0.547-0.139*0.139)/0.0210313),-0.657178)/TMath::Power((exp(- -0.124633*x-0.310906*x*x)+x/0.0210313),-0.657178)+1.42452*TMath::Landau(x,2.46537,0.969902)),TMath::Max(4.88363e-10*TMath::Power((exp(- -2.82354*sqrt(x*x+0.547*0.547-0.139*0.139)-0.728315*(x*x+0.547*0.547-0.139*0.139))+TMath::Power(0.483053/4.88363e-10,-1./16.2073)*sqrt(x*x+0.547*0.547-0.139*0.139)/0.0021492),-16.2073)/TMath::Power((exp(- -2.82354*x-0.728315*x*x)+x/0.0021492),-16.2073)+1.57731*TMath::Landau(x,2.40466,1.00957),0.00469052*TMath::Power((exp(- -0.124633*sqrt(x*x+0.547*0.547-0.139*0.139)-0.310906*(x*x+0.547*0.547-0.139*0.139))+TMath::Power(0.479677/0.00469052,-1./0.657178)*sqrt(x*x+0.547*0.547-0.139*0.139)/0.0210313),-0.657178)/TMath::Power((exp(- -0.124633*x-0.310906*x*x)+x/0.0210313),-0.657178)+1.42452*TMath::Landau(x,2.46537,0.969902))))+(x>=9.445400)*1.06*(((0.52091*1.)+(159.301*(-0.00307058*TMath::Power(1.+((x/1.03888)*(x/1.03888)),-(0.567478)))))/(1.+(-0.00307058*TMath::Power(1.+((x/1.03888)*(x/1.03888)),-(0.567478))))))))", "333_pt": "(TMath::TwoPi()*x*(219.382*exp(-(sqrt(x*x+0.139571*0.139571)-0.139571)/0.384303)+2471.6*pow(exp(-1.12373*x-0.752507*x*x)+x/0.438343,-5.45943)))*(0.1191/(1+exp(-(x- -0.316437)/0.554431))+-0.507618*TMath::Gaus(x,1.69654,1.62649)+0.494136*TMath::Gaus(x,2.56805,1.6193)+0.0789659)" } -} \ No newline at end of file +} diff --git a/MC/config/PWGLF/ini/GeneratorLF_Resonances_PbPb.ini b/MC/config/PWGLF/ini/GeneratorLF_Resonances_PbPb.ini new file mode 100644 index 000000000..b66e66584 --- /dev/null +++ b/MC/config/PWGLF/ini/GeneratorLF_Resonances_PbPb.ini @@ -0,0 +1,10 @@ +[GeneratorExternal] +fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator_pythia8_LF.C +funcName=generateLF("${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator/resonancelistgun_pbpb.json", true, 4) + +# [GeneratorPythia8] +# config=${O2_ROOT}/share/Generators/egconfig/pythia8_hi.cfg + +[DecayerPythia8] # after for transport code! +config[0]=${O2DPG_MC_CONFIG_ROOT}/MC/config/common/pythia8/decayer/base.cfg +config[1]=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator/resonances.cfg \ No newline at end of file diff --git a/MC/config/PWGLF/ini/GeneratorLF_Resonances_PbPb_exotic.ini b/MC/config/PWGLF/ini/GeneratorLF_Resonances_PbPb_exotic.ini new file mode 100644 index 000000000..e2e575af0 --- /dev/null +++ b/MC/config/PWGLF/ini/GeneratorLF_Resonances_PbPb_exotic.ini @@ -0,0 +1,10 @@ +[GeneratorExternal] +fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator_pythia8_LF.C +funcName=generateLF("${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator/resonancelistgun_exotic_pbpb.json", true, 4) + +# [GeneratorPythia8] +# config=${O2_ROOT}/share/Generators/egconfig/pythia8_hi.cfg + +[DecayerPythia8] # after for transport code! +config[0]=${O2DPG_MC_CONFIG_ROOT}/MC/config/common/pythia8/decayer/base.cfg +config[1]=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator/resonances.cfg \ No newline at end of file diff --git a/MC/config/PWGLF/ini/GeneratorLF_Resonances_pp.ini b/MC/config/PWGLF/ini/GeneratorLF_Resonances_pp.ini new file mode 100644 index 000000000..8e91af28b --- /dev/null +++ b/MC/config/PWGLF/ini/GeneratorLF_Resonances_pp.ini @@ -0,0 +1,10 @@ +[GeneratorExternal] +fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator_pythia8_LF.C +funcName=generateLF("${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator/resonancelistgun.json", true, 4) + +# [GeneratorPythia8] # if triggered then this will be used as the background event +# config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator/pythia8_inel_136tev.cfg + +[DecayerPythia8] # after for transport code! +config[0]=${O2DPG_MC_CONFIG_ROOT}/MC/config/common/pythia8/decayer/base.cfg +config[1]=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator/resonances.cfg \ No newline at end of file diff --git a/MC/config/PWGLF/ini/GeneratorLF_Resonances_pp_exotic.ini b/MC/config/PWGLF/ini/GeneratorLF_Resonances_pp_exotic.ini new file mode 100644 index 000000000..a49804cc8 --- /dev/null +++ b/MC/config/PWGLF/ini/GeneratorLF_Resonances_pp_exotic.ini @@ -0,0 +1,10 @@ +[GeneratorExternal] +fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator_pythia8_LF.C +funcName=generateLF("${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator/resonancelistgun_exotic.json", true, 4) + +# [GeneratorPythia8] # if triggered then this will be used as the background event +# config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator/pythia8_inel_136tev.cfg + +[DecayerPythia8] # after for transport code! +config[0]=${O2DPG_MC_CONFIG_ROOT}/MC/config/common/pythia8/decayer/base.cfg +config[1]=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator/resonances.cfg \ No newline at end of file diff --git a/MC/config/PWGLF/ini/GeneratorLF_SyntheFlowXi.ini b/MC/config/PWGLF/ini/GeneratorLF_SyntheFlowXi.ini new file mode 100644 index 000000000..4752520bf --- /dev/null +++ b/MC/config/PWGLF/ini/GeneratorLF_SyntheFlowXi.ini @@ -0,0 +1,6 @@ +[GeneratorExternal] +fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator_pythia8_syntheFlowXi.C +funcName=generator_syntheFlowXi() + +[GeneratorPythia8] +config=${O2DPG_MC_CONFIG_ROOT}/MC/config/common/pythia8/generator/pythia8_hi.cfg diff --git a/MC/config/PWGLF/ini/tests/GeneratorLF_Resonances_PbPb.C b/MC/config/PWGLF/ini/tests/GeneratorLF_Resonances_PbPb.C new file mode 100644 index 000000000..fd735962c --- /dev/null +++ b/MC/config/PWGLF/ini/tests/GeneratorLF_Resonances_PbPb.C @@ -0,0 +1,182 @@ +int External() +{ + std::string path{"o2sim_Kine.root"}; + int numberOfInjectedSignalsPerEvent{1}; + int numberOfGapEvents{4}; + int numberOfEventsProcessed{0}; + int numberOfEventsProcessedWithoutInjection{0}; + std::vector injectedPDGs = { + 323, // K*+- + -323, // K*bar+- + 9010221, // f_0(980) + 113, // rho(770)0 + 213, // rho(770)+ + -213, // rho(770)bar- + 3224, // Sigma(1385)+ + -3224, // Sigma(1385)bar- + 3114, // Sigma(1385)- + -3114, // Sigma(1385)bar+ + 3124, // Lambda(1520)0 + -3124, // Lambda(1520)0bar + 3324, // Xi(1530)0 + -3324, // Xi(1530)0bar + 10323, // K1(1270)+ + -10323, // K1(1270)-bar + 2224, // Delta(1232)+ + -2224, // Delta(1232)bar- + 2114, // Delta(1232)0 + -2114, // Delta(1232)0bar + 123314, // Xi(1820)- + -123314, // Xi(1820)+ + 123324, // Xi(1820)0 + -123324 // Xi(1820)0bar + }; + std::vector> decayDaughters = { + {311, 211}, // K*+- + {-311, -211}, // K*bar+- + {211, -211}, // f_0(980) + {211, -211}, // rho(770)0 + {211, 111}, // rho(770)+ + {-211, 111}, // rho(770)bar- + {3122, 211}, // Sigma(1385)+ + {-3122, -211}, // Sigma(1385)bar- + {3122, -211}, // Sigma(1385)- + {-3122, 211}, // Sigma(1385)bar+ + {2212, -321}, // Lambda(1520)0 + {-2212, 321}, // Lambda(1520)0bar + {3312, 211}, // Xi(1530)0 + {-3312, -211}, // Xi(1530)0bar + {321, 211}, // K1(1270)+ + {-321, -211}, // K1(1270)-bar + {2212, 211}, // Delta(1232)+ + {-2212, -211}, // Delta(1232)bar- + {2212, -211}, // Delta(1232)- + {-2212, 211}, // Delta(1232)bar+ + {3122, -311}, // Xi(1820)- + {3122, 311}, // Xi(1820)+ + {3122, 310}, // Xi(1820)0 + {-3122, 310} // Xi(1820)0bar + }; + + auto nInjection = injectedPDGs.size(); + + TFile file(path.c_str(), "READ"); + if (file.IsZombie()) + { + std::cerr << "Cannot open ROOT file " << path << "\n"; + return 1; + } + + auto tree = (TTree *)file.Get("o2sim"); + if (!tree) + { + std::cerr << "Cannot find tree o2sim in file " << path << "\n"; + return 1; + } + std::vector *tracks{}; + tree->SetBranchAddress("MCTrack", &tracks); + + std::vector nSignal; + for (int i = 0; i < nInjection; i++) + { + nSignal.push_back(0); + } + std::vector> nDecays; + std::vector nNotDecayed; + for (int i = 0; i < nInjection; i++) + { + std::vector nDecay; + for (int j = 0; j < decayDaughters[i].size(); j++) + { + nDecay.push_back(0); + } + nDecays.push_back(nDecay); + nNotDecayed.push_back(0); + } + auto nEvents = tree->GetEntries(); + bool hasInjection = false; + for (int i = 0; i < nEvents; i++) + { + hasInjection = false; + numberOfEventsProcessed++; + auto check = tree->GetEntry(i); + for (int idxMCTrack = 0; idxMCTrack < tracks->size(); ++idxMCTrack) + { + auto track = tracks->at(idxMCTrack); + auto pdg = track.GetPdgCode(); + auto it = std::find(injectedPDGs.begin(), injectedPDGs.end(), pdg); + int index = std::distance(injectedPDGs.begin(), it); // index of injected PDG + if (it != injectedPDGs.end()) // found + { + // count signal PDG + nSignal[index]++; + if (track.getFirstDaughterTrackId() < 0) + { + nNotDecayed[index]++; + continue; + } + for (int j{track.getFirstDaughterTrackId()}; j <= track.getLastDaughterTrackId(); ++j) + { + auto pdgDau = tracks->at(j).GetPdgCode(); + bool foundDau = false; + // count decay PDGs + for (int idxDaughter = 0; idxDaughter < decayDaughters[index].size(); ++idxDaughter) + { + if (pdgDau == decayDaughters[index][idxDaughter]) + { + nDecays[index][idxDaughter]++; + foundDau = true; + hasInjection = true; + break; + } + } + if (!foundDau) + { + std::cerr << "Decay daughter not found: " << pdg << " -> " << pdgDau << "\n"; + } + } + } + } + if (!hasInjection) + { + numberOfEventsProcessedWithoutInjection++; + } + } + std::cout << "--------------------------------\n"; + std::cout << "# Events: " << nEvents << "\n"; + for (int i = 0; i < nInjection; i++) + { + std::cout << "# Mother \n"; + std::cout << injectedPDGs[i] << " generated: " << nSignal[i] << ", " << nNotDecayed[i] << " did not decay\n"; + if (nSignal[i] == 0) + { + std::cerr << "No generated: " << injectedPDGs[i] << "\n"; + // return 1; // At least one of the injected particles should be generated + } + for (int j = 0; j < decayDaughters[i].size(); j++) + { + std::cout << "# Daughter " << decayDaughters[i][j] << ": " << nDecays[i][j] << "\n"; + } + // if (nSignal[i] != nEvents * numberOfInjectedSignalsPerEvent) + // { + // std::cerr << "Number of generated: " << injectedPDGs[i] << ", lower than expected\n"; + // // return 1; // Don't need to return 1, since the number of generated particles is not the same for each event + // } + } + std::cout << "--------------------------------\n"; + std::cout << "Number of events processed: " << numberOfEventsProcessed << "\n"; + std::cout << "Number of input for the gap events: " << numberOfGapEvents << "\n"; + std::cout << "Number of events processed without injection: " << numberOfEventsProcessedWithoutInjection << "\n"; + // injected event + numberOfGapEvents*gap events + injected event + numberOfGapEvents*gap events + ... + // total fraction of the gap event: numberOfEventsProcessedWithoutInjection/numberOfEventsProcessed + float ratioOfNormalEvents = numberOfEventsProcessedWithoutInjection / numberOfEventsProcessed; + if (ratioOfNormalEvents > 0.75) + { + std::cout << "The number of injected event is loo low!!" << std::endl; + return 1; + } + + return 0; +} + +void GeneratorLF_Resonances_PbPb5360_injection() { External(); } diff --git a/MC/config/PWGLF/ini/tests/GeneratorLF_Resonances_PbPb_exotic.C b/MC/config/PWGLF/ini/tests/GeneratorLF_Resonances_PbPb_exotic.C new file mode 100644 index 000000000..becf5ce3e --- /dev/null +++ b/MC/config/PWGLF/ini/tests/GeneratorLF_Resonances_PbPb_exotic.C @@ -0,0 +1,161 @@ +int External() +{ + std::string path{"o2sim_Kine.root"}; + int numberOfInjectedSignalsPerEvent{1}; + int numberOfGapEvents{4}; + int numberOfEventsProcessed{0}; + int numberOfEventsProcessedWithoutInjection{0}; + std::vector injectedPDGs = { + 9010221, // f_0(980) + 10221, // f_0(1370) + 9030221, // f_0(1500) + 10331, // f_0(1710) + 20223, // f_1(1285) + 20333, // f_1(1420) + 335, // f_2(1525) + 10323, // K1(1270)+ + -10323, // K1(1270)-bar + 123314, // Xi(1820)- + -123314, // Xi(1820)+ + 123324, // Xi(1820)0 + -123324 // Xi(1820)0bar + }; + std::vector> decayDaughters = { + {211, -211}, // f_0(980) + {310, 310}, // f_0(1370) + {310, 310}, // f_0(1500) + {310, 310}, // f_0(1710) + {310, -321, 211}, // f_1(1285) + {310, -321, 211}, // f_1(1420) + {310, 310}, // f_2(1525) + {321, 211}, // K1(1270)+ + {-321, -211}, // K1(1270)-bar + {2212, 211}, // Delta(1232)+ + {3122, -311}, // Xi(1820)- + {3122, 311}, // Xi(1820)+ + {3122, 310}, // Xi(1820)0 + {-3122, 310} // Xi(1820)0bar + }; + + auto nInjection = injectedPDGs.size(); + + TFile file(path.c_str(), "READ"); + if (file.IsZombie()) + { + std::cerr << "Cannot open ROOT file " << path << "\n"; + return 1; + } + + auto tree = (TTree *)file.Get("o2sim"); + if (!tree) + { + std::cerr << "Cannot find tree o2sim in file " << path << "\n"; + return 1; + } + std::vector *tracks{}; + tree->SetBranchAddress("MCTrack", &tracks); + + std::vector nSignal; + for (int i = 0; i < nInjection; i++) + { + nSignal.push_back(0); + } + std::vector> nDecays; + std::vector nNotDecayed; + for (int i = 0; i < nInjection; i++) + { + std::vector nDecay; + for (int j = 0; j < decayDaughters[i].size(); j++) + { + nDecay.push_back(0); + } + nDecays.push_back(nDecay); + nNotDecayed.push_back(0); + } + auto nEvents = tree->GetEntries(); + bool hasInjection = false; + for (int i = 0; i < nEvents; i++) + { + hasInjection = false; + numberOfEventsProcessed++; + auto check = tree->GetEntry(i); + for (int idxMCTrack = 0; idxMCTrack < tracks->size(); ++idxMCTrack) + { + auto track = tracks->at(idxMCTrack); + auto pdg = track.GetPdgCode(); + auto it = std::find(injectedPDGs.begin(), injectedPDGs.end(), pdg); + int index = std::distance(injectedPDGs.begin(), it); // index of injected PDG + if (it != injectedPDGs.end()) // found + { + // count signal PDG + nSignal[index]++; + if (track.getFirstDaughterTrackId() < 0) + { + nNotDecayed[index]++; + continue; + } + for (int j{track.getFirstDaughterTrackId()}; j <= track.getLastDaughterTrackId(); ++j) + { + auto pdgDau = tracks->at(j).GetPdgCode(); + bool foundDau = false; + // count decay PDGs + for (int idxDaughter = 0; idxDaughter < decayDaughters[index].size(); ++idxDaughter) + { + if (pdgDau == decayDaughters[index][idxDaughter]) + { + nDecays[index][idxDaughter]++; + foundDau = true; + hasInjection = true; + break; + } + } + if (!foundDau) + { + std::cerr << "Decay daughter not found: " << pdg << " -> " << pdgDau << "\n"; + } + } + } + } + if (!hasInjection) + { + numberOfEventsProcessedWithoutInjection++; + } + } + std::cout << "--------------------------------\n"; + std::cout << "# Events: " << nEvents << "\n"; + for (int i = 0; i < nInjection; i++) + { + std::cout << "# Mother \n"; + std::cout << injectedPDGs[i] << " generated: " << nSignal[i] << ", " << nNotDecayed[i] << " did not decay\n"; + if (nSignal[i] == 0) + { + std::cerr << "No generated: " << injectedPDGs[i] << "\n"; + // return 1; // At least one of the injected particles should be generated + } + for (int j = 0; j < decayDaughters[i].size(); j++) + { + std::cout << "# Daughter " << decayDaughters[i][j] << ": " << nDecays[i][j] << "\n"; + } + // if (nSignal[i] != nEvents * numberOfInjectedSignalsPerEvent) + // { + // std::cerr << "Number of generated: " << injectedPDGs[i] << ", lower than expected\n"; + // // return 1; // Don't need to return 1, since the number of generated particles is not the same for each event + // } + } + std::cout << "--------------------------------\n"; + std::cout << "Number of events processed: " << numberOfEventsProcessed << "\n"; + std::cout << "Number of input for the gap events: " << numberOfGapEvents << "\n"; + std::cout << "Number of events processed without injection: " << numberOfEventsProcessedWithoutInjection << "\n"; + // injected event + numberOfGapEvents*gap events + injected event + numberOfGapEvents*gap events + ... + // total fraction of the gap event: numberOfEventsProcessedWithoutInjection/numberOfEventsProcessed + float ratioOfNormalEvents = numberOfEventsProcessedWithoutInjection / numberOfEventsProcessed; + if (ratioOfNormalEvents > 0.75) + { + std::cout << "The number of injected event is loo low!!" << std::endl; + return 1; + } + + return 0; +} + +void GeneratorLF_Resonances_PbPb5360_injection() { External(); } diff --git a/MC/config/PWGLF/ini/tests/GeneratorLF_Resonances_pp.C b/MC/config/PWGLF/ini/tests/GeneratorLF_Resonances_pp.C new file mode 100644 index 000000000..c8edc36d7 --- /dev/null +++ b/MC/config/PWGLF/ini/tests/GeneratorLF_Resonances_pp.C @@ -0,0 +1,182 @@ +int External() +{ + std::string path{"o2sim_Kine.root"}; + int numberOfInjectedSignalsPerEvent{1}; + int numberOfGapEvents{4}; + int numberOfEventsProcessed{0}; + int numberOfEventsProcessedWithoutInjection{0}; + std::vector injectedPDGs = { + 323, // K*+- + -323, // K*bar+- + 9010221, // f_0(980) + 113, // rho(770)0 + 213, // rho(770)+ + -213, // rho(770)bar- + 3224, // Sigma(1385)+ + -3224, // Sigma(1385)bar- + 3114, // Sigma(1385)- + -3114, // Sigma(1385)bar+ + 3124, // Lambda(1520)0 + -3124, // Lambda(1520)0bar + 3324, // Xi(1530)0 + -3324, // Xi(1530)0bar + 10323, // K1(1270)+ + -10323, // K1(1270)-bar + 2224, // Delta(1232)+ + -2224, // Delta(1232)bar- + 2114, // Delta(1232)0 + -2114, // Delta(1232)0bar + 123314, // Xi(1820)- + -123314, // Xi(1820)+ + 123324, // Xi(1820)0 + -123324 // Xi(1820)0bar + }; + std::vector> decayDaughters = { + {311, 211}, // K*+- + {-311, -211}, // K*bar+- + {211, -211}, // f_0(980) + {211, -211}, // rho(770)0 + {211, 111}, // rho(770)+ + {-211, 111}, // rho(770)bar- + {3122, 211}, // Sigma(1385)+ + {-3122, -211}, // Sigma(1385)bar- + {3122, -211}, // Sigma(1385)- + {-3122, 211}, // Sigma(1385)bar+ + {2212, -321}, // Lambda(1520)0 + {-2212, 321}, // Lambda(1520)0bar + {3312, 211}, // Xi(1530)0 + {-3312, -211}, // Xi(1530)0bar + {321, 211}, // K1(1270)+ + {-321, -211}, // K1(1270)-bar + {2212, 211}, // Delta(1232)+ + {-2212, -211}, // Delta(1232)bar- + {2212, -211}, // Delta(1232)- + {-2212, 211}, // Delta(1232)bar+ + {3122, -311}, // Xi(1820)- + {3122, 311}, // Xi(1820)+ + {3122, 310}, // Xi(1820)0 + {-3122, 310} // Xi(1820)0bar + }; + + auto nInjection = injectedPDGs.size(); + + TFile file(path.c_str(), "READ"); + if (file.IsZombie()) + { + std::cerr << "Cannot open ROOT file " << path << "\n"; + return 1; + } + + auto tree = (TTree *)file.Get("o2sim"); + if (!tree) + { + std::cerr << "Cannot find tree o2sim in file " << path << "\n"; + return 1; + } + std::vector *tracks{}; + tree->SetBranchAddress("MCTrack", &tracks); + + std::vector nSignal; + for (int i = 0; i < nInjection; i++) + { + nSignal.push_back(0); + } + std::vector> nDecays; + std::vector nNotDecayed; + for (int i = 0; i < nInjection; i++) + { + std::vector nDecay; + for (int j = 0; j < decayDaughters[i].size(); j++) + { + nDecay.push_back(0); + } + nDecays.push_back(nDecay); + nNotDecayed.push_back(0); + } + auto nEvents = tree->GetEntries(); + bool hasInjection = false; + for (int i = 0; i < nEvents; i++) + { + hasInjection = false; + numberOfEventsProcessed++; + auto check = tree->GetEntry(i); + for (int idxMCTrack = 0; idxMCTrack < tracks->size(); ++idxMCTrack) + { + auto track = tracks->at(idxMCTrack); + auto pdg = track.GetPdgCode(); + auto it = std::find(injectedPDGs.begin(), injectedPDGs.end(), pdg); + int index = std::distance(injectedPDGs.begin(), it); // index of injected PDG + if (it != injectedPDGs.end()) // found + { + // count signal PDG + nSignal[index]++; + if (track.getFirstDaughterTrackId() < 0) + { + nNotDecayed[index]++; + continue; + } + for (int j{track.getFirstDaughterTrackId()}; j <= track.getLastDaughterTrackId(); ++j) + { + auto pdgDau = tracks->at(j).GetPdgCode(); + bool foundDau = false; + // count decay PDGs + for (int idxDaughter = 0; idxDaughter < decayDaughters[index].size(); ++idxDaughter) + { + if (pdgDau == decayDaughters[index][idxDaughter]) + { + nDecays[index][idxDaughter]++; + foundDau = true; + hasInjection = true; + break; + } + } + if (!foundDau) + { + std::cerr << "Decay daughter not found: " << pdg << " -> " << pdgDau << "\n"; + } + } + } + } + if (!hasInjection) + { + numberOfEventsProcessedWithoutInjection++; + } + } + std::cout << "--------------------------------\n"; + std::cout << "# Events: " << nEvents << "\n"; + for (int i = 0; i < nInjection; i++) + { + std::cout << "# Mother \n"; + std::cout << injectedPDGs[i] << " generated: " << nSignal[i] << ", " << nNotDecayed[i] << " did not decay\n"; + if (nSignal[i] == 0) + { + std::cerr << "No generated: " << injectedPDGs[i] << "\n"; + // return 1; // At least one of the injected particles should be generated + } + for (int j = 0; j < decayDaughters[i].size(); j++) + { + std::cout << "# Daughter " << decayDaughters[i][j] << ": " << nDecays[i][j] << "\n"; + } + // if (nSignal[i] != nEvents * numberOfInjectedSignalsPerEvent) + // { + // std::cerr << "Number of generated: " << injectedPDGs[i] << ", lower than expected\n"; + // // return 1; // Don't need to return 1, since the number of generated particles is not the same for each event + // } + } + std::cout << "--------------------------------\n"; + std::cout << "Number of events processed: " << numberOfEventsProcessed << "\n"; + std::cout << "Number of input for the gap events: " << numberOfGapEvents << "\n"; + std::cout << "Number of events processed without injection: " << numberOfEventsProcessedWithoutInjection << "\n"; + // injected event + numberOfGapEvents*gap events + injected event + numberOfGapEvents*gap events + ... + // total fraction of the gap event: numberOfEventsProcessedWithoutInjection/numberOfEventsProcessed + float ratioOfNormalEvents = numberOfEventsProcessedWithoutInjection / numberOfEventsProcessed; + if (ratioOfNormalEvents > 0.75) + { + std::cout << "The number of injected event is loo low!!" << std::endl; + return 1; + } + + return 0; +} + +void GeneratorLF_Resonances_pp1360_injection() { External(); } diff --git a/MC/config/PWGLF/ini/tests/GeneratorLF_Resonances_pp_exotic.C b/MC/config/PWGLF/ini/tests/GeneratorLF_Resonances_pp_exotic.C new file mode 100644 index 000000000..61d739eab --- /dev/null +++ b/MC/config/PWGLF/ini/tests/GeneratorLF_Resonances_pp_exotic.C @@ -0,0 +1,163 @@ +int External() +{ + std::string path{"o2sim_Kine.root"}; + int numberOfInjectedSignalsPerEvent{1}; + int numberOfGapEvents{4}; + int numberOfEventsProcessed{0}; + int numberOfEventsProcessedWithoutInjection{0}; + std::vector injectedPDGs = { + 9010221, // f_0(980) + 10221, // f_0(1370) + 9030221, // f_0(1500) + 10331, // f_0(1710) + 20223, // f_1(1285) + 20333, // f_1(1420) + 335, // f_2(1525) + 115, // a_2(1230) + 10323, // K1(1270)+ + -10323, // K1(1270)-bar + 123314, // Xi(1820)- + -123314, // Xi(1820)+ + 123324, // Xi(1820)0 + -123324 // Xi(1820)0bar + }; + std::vector> decayDaughters = { + {211, -211}, // f_0(980) + {310, 310}, // f_0(1370) + {310, 310}, // f_0(1500) + {310, 310}, // f_0(1710) + {310, -321, 211}, // f_1(1285) + {310, -321, 211}, // f_1(1420) + {310, 310}, // f_2(1525) + {310, 310}, // a_2(1230) + {321, 211}, // K1(1270)+ + {-321, -211}, // K1(1270)-bar + {2212, 211}, // Delta(1232)+ + {3122, -311}, // Xi(1820)- + {3122, 311}, // Xi(1820)+ + {3122, 310}, // Xi(1820)0 + {-3122, 310} // Xi(1820)0bar + }; + + auto nInjection = injectedPDGs.size(); + + TFile file(path.c_str(), "READ"); + if (file.IsZombie()) + { + std::cerr << "Cannot open ROOT file " << path << "\n"; + return 1; + } + + auto tree = (TTree *)file.Get("o2sim"); + if (!tree) + { + std::cerr << "Cannot find tree o2sim in file " << path << "\n"; + return 1; + } + std::vector *tracks{}; + tree->SetBranchAddress("MCTrack", &tracks); + + std::vector nSignal; + for (int i = 0; i < nInjection; i++) + { + nSignal.push_back(0); + } + std::vector> nDecays; + std::vector nNotDecayed; + for (int i = 0; i < nInjection; i++) + { + std::vector nDecay; + for (int j = 0; j < decayDaughters[i].size(); j++) + { + nDecay.push_back(0); + } + nDecays.push_back(nDecay); + nNotDecayed.push_back(0); + } + auto nEvents = tree->GetEntries(); + bool hasInjection = false; + for (int i = 0; i < nEvents; i++) + { + hasInjection = false; + numberOfEventsProcessed++; + auto check = tree->GetEntry(i); + for (int idxMCTrack = 0; idxMCTrack < tracks->size(); ++idxMCTrack) + { + auto track = tracks->at(idxMCTrack); + auto pdg = track.GetPdgCode(); + auto it = std::find(injectedPDGs.begin(), injectedPDGs.end(), pdg); + int index = std::distance(injectedPDGs.begin(), it); // index of injected PDG + if (it != injectedPDGs.end()) // found + { + // count signal PDG + nSignal[index]++; + if (track.getFirstDaughterTrackId() < 0) + { + nNotDecayed[index]++; + continue; + } + for (int j{track.getFirstDaughterTrackId()}; j <= track.getLastDaughterTrackId(); ++j) + { + auto pdgDau = tracks->at(j).GetPdgCode(); + bool foundDau = false; + // count decay PDGs + for (int idxDaughter = 0; idxDaughter < decayDaughters[index].size(); ++idxDaughter) + { + if (pdgDau == decayDaughters[index][idxDaughter]) + { + nDecays[index][idxDaughter]++; + foundDau = true; + hasInjection = true; + break; + } + } + if (!foundDau) + { + std::cerr << "Decay daughter not found: " << pdg << " -> " << pdgDau << "\n"; + } + } + } + } + if (!hasInjection) + { + numberOfEventsProcessedWithoutInjection++; + } + } + std::cout << "--------------------------------\n"; + std::cout << "# Events: " << nEvents << "\n"; + for (int i = 0; i < nInjection; i++) + { + std::cout << "# Mother \n"; + std::cout << injectedPDGs[i] << " generated: " << nSignal[i] << ", " << nNotDecayed[i] << " did not decay\n"; + if (nSignal[i] == 0) + { + std::cerr << "No generated: " << injectedPDGs[i] << "\n"; + // return 1; // At least one of the injected particles should be generated + } + for (int j = 0; j < decayDaughters[i].size(); j++) + { + std::cout << "# Daughter " << decayDaughters[i][j] << ": " << nDecays[i][j] << "\n"; + } + // if (nSignal[i] != nEvents * numberOfInjectedSignalsPerEvent) + // { + // std::cerr << "Number of generated: " << injectedPDGs[i] << ", lower than expected\n"; + // // return 1; // Don't need to return 1, since the number of generated particles is not the same for each event + // } + } + std::cout << "--------------------------------\n"; + std::cout << "Number of events processed: " << numberOfEventsProcessed << "\n"; + std::cout << "Number of input for the gap events: " << numberOfGapEvents << "\n"; + std::cout << "Number of events processed without injection: " << numberOfEventsProcessedWithoutInjection << "\n"; + // injected event + numberOfGapEvents*gap events + injected event + numberOfGapEvents*gap events + ... + // total fraction of the gap event: numberOfEventsProcessedWithoutInjection/numberOfEventsProcessed + float ratioOfNormalEvents = numberOfEventsProcessedWithoutInjection / numberOfEventsProcessed; + if (ratioOfNormalEvents > 0.75) + { + std::cout << "The number of injected event is loo low!!" << std::endl; + return 1; + } + + return 0; +} + +void GeneratorLF_Resonances_pp1360_injection() { External(); } diff --git a/MC/config/PWGLF/ini/tests/GeneratorLF_SyntheFlowXi.C b/MC/config/PWGLF/ini/tests/GeneratorLF_SyntheFlowXi.C new file mode 100644 index 000000000..ab7eeb695 --- /dev/null +++ b/MC/config/PWGLF/ini/tests/GeneratorLF_SyntheFlowXi.C @@ -0,0 +1,28 @@ +int External() +{ + std::string path{"o2sim_Kine.root"}; + + TFile file(path.c_str(), "READ"); + if (file.IsZombie()) + { + std::cerr << "Cannot open ROOT file " << path << "\n"; + return 1; + } + + auto tree = (TTree *)file.Get("o2sim"); + if (!tree) + { + std::cerr << "Cannot find tree o2sim in file " << path << "\n"; + return 1; + } + std::vector *tracks{}; + tree->SetBranchAddress("MCTrack", &tracks); + + auto nEvents = tree->GetEntries(); + if (nEvents < 1) + { + std::cerr << "No events actually generated: not OK!"; + return 1; + } + return 0; +} diff --git a/MC/config/PWGLF/pythia8/generator/resonancelistgun_exotic.json b/MC/config/PWGLF/pythia8/generator/resonancelistgun_exotic.json new file mode 100644 index 000000000..1d1f07fcb --- /dev/null +++ b/MC/config/PWGLF/pythia8/generator/resonancelistgun_exotic.json @@ -0,0 +1,128 @@ +{ + "f_0(980)" : { + "pdg": 9010221, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "f_0(1370)" : { + "pdg": 10221, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "f_0(1500)" : { + "pdg": 9030221, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "f_0(1710)" : { + "pdg": 10331, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "f_1(1285)" : { + "pdg": 20223, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "f_1(1420)" : { + "pdg": 20333, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "f_2(1525)" : { + "pdg": 335, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "a_2(1230)" : { + "pdg": 115, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "K1(1270)+": { + "pdg": 10323, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "K1(1270)-": { + "pdg": -10323, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Xi(1820)0" : { + "pdg": 123314, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Anti-Xi(1820)0" : { + "pdg": -123314, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Xi(1820)-" : { + "pdg": 123324, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Xi(1820)+" : { + "pdg": -123324, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + } +} diff --git a/MC/config/PWGLF/pythia8/generator/resonancelistgun_exotic_pbpb.json b/MC/config/PWGLF/pythia8/generator/resonancelistgun_exotic_pbpb.json new file mode 100644 index 000000000..99b38064d --- /dev/null +++ b/MC/config/PWGLF/pythia8/generator/resonancelistgun_exotic_pbpb.json @@ -0,0 +1,128 @@ +{ + "f_0(980)" : { + "pdg": 9010221, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "f_0(1370)" : { + "pdg": 10221, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "f_0(1500)" : { + "pdg": 9030221, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "f_0(1710)" : { + "pdg": 10331, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "f_1(1285)" : { + "pdg": 20223, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "f_1(1420)" : { + "pdg": 20333, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "f_2(1525)" : { + "pdg": 335, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "a_2(1230)" : { + "pdg": 115, + "n": 1, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "K1(1270)+": { + "pdg": 10323, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "K1(1270)-": { + "pdg": -10323, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Xi(1820)0" : { + "pdg": 123314, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Anti-Xi(1820)0" : { + "pdg": -123314, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Xi(1820)-" : { + "pdg": 123324, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Xi(1820)+" : { + "pdg": -123324, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + } +} diff --git a/MC/config/PWGLF/pythia8/generator/resonancelistgun_pbpb.json b/MC/config/PWGLF/pythia8/generator/resonancelistgun_pbpb.json new file mode 100644 index 000000000..89413907a --- /dev/null +++ b/MC/config/PWGLF/pythia8/generator/resonancelistgun_pbpb.json @@ -0,0 +1,218 @@ +{ + "K(892)+" : { + "pdg": 323, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "K(892)-" : { + "pdg": -323, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "f_0(980)" : { + "pdg": 9010221, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "rho(770)0" : { + "pdg": 113, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "rho(770)+" : { + "pdg": 213, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "rho(770)-" : { + "pdg": -213, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Sigma(1385)-" : { + "pdg": 3114, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "anti-Sigma(1385)+" : { + "pdg": -3114, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Sigma(1385)+" : { + "pdg": 3224, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "anti-Sigma(1385)-" : { + "pdg": -3224, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Lambda(1520)0" : { + "pdg": 102134, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "anti-Lambda(1520)0" : { + "pdg": -102134, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Xi(1530)0" : { + "pdg": 3324, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "anti-Xi(1530)0" : { + "pdg": -3324, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "K1(1270)+": { + "pdg": 10323, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "K1(1270)-": { + "pdg": -10323, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Delta++" : { + "pdg": 2224, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "anti-Delta++" : { + "pdg": -2224, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Delta0" : { + "pdg": 2114, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "anti-Delta0" : { + "pdg": -2114, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Xi(1820)0" : { + "pdg": 123314, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Anti-Xi(1820)0" : { + "pdg": -123314, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Xi(1820)-" : { + "pdg": 123324, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + }, + "Xi(1820)+" : { + "pdg": -123324, + "n": 10, + "ptMin": 0.0, + "ptMax": 20, + "etaMin": -1.2, + "etaMax": 1.2, + "genDecayed": true + } +} diff --git a/MC/config/PWGLF/pythia8/generator_pythia8_extraStrangeness.C b/MC/config/PWGLF/pythia8/generator_pythia8_extraStrangeness.C index 75651ff00..c7a57f956 100644 --- a/MC/config/PWGLF/pythia8/generator_pythia8_extraStrangeness.C +++ b/MC/config/PWGLF/pythia8/generator_pythia8_extraStrangeness.C @@ -5,6 +5,7 @@ #include "FairPrimaryGenerator.h" #include "Generators/GeneratorPythia8.h" #include "TRandom3.h" +#include "TF1.h" #include "TParticlePDG.h" #include "TDatabasePDG.h" @@ -18,10 +19,10 @@ public: GeneratorPythia8ExtraStrangeness() { genMinPt=0.0; genMaxPt=20.0; - genminY=-1.5; - genmaxY=1.5; - genminEta=-1.5; - genmaxEta=1.5; + genminY=-1.0; + genmaxY=1.0; + genminEta=-1.0; + genmaxEta=1.0; pdg=0; E=0; @@ -37,6 +38,22 @@ public: xProd=0.; yProd=0.; zProd=0.; fLVHelper = std::make_unique(); + + fSpectrumXi = std::make_unique("fSpectrumXi", this, &GeneratorPythia8ExtraStrangeness::boltzPlusPower, 0., genMaxPt, 5, "GeneratorPythia8ExtraStrangeness", "boltzPlusPower"); + + fSpectrumXi->FixParameter(0, 1.32171); + fSpectrumXi->FixParameter(1, 4.84e-1); + fSpectrumXi->FixParameter(2, 111.9); + fSpectrumXi->FixParameter(3, -2.56511e+00); + fSpectrumXi->FixParameter(4, 1.14011e-04); + + fSpectrumOm = std::make_unique("fSpectrumOm", this, &GeneratorPythia8ExtraStrangeness::boltzPlusPower, 0., genMaxPt, 5, "GeneratorPythia8ExtraStrangeness", "boltzPlusPower"); + + fSpectrumOm->FixParameter(0, 1.67245e+00); + fSpectrumOm->FixParameter(1, 5.18174e-01); + fSpectrumOm->FixParameter(2, 1.73747e+01); + fSpectrumOm->FixParameter(3, -2.56681e+00); + fSpectrumOm->FixParameter(4, 1.87513e-04); } Double_t y2eta(Double_t pt, Double_t mass, Double_t y){ @@ -85,7 +102,7 @@ public: ranGenerator->SetSeed(0); // generate transverse momentum - const double gen_pT = ranGenerator->Uniform(0,5); + const double gen_pT = fSpectrumXi->GetRandom(genMinPt,genMaxPt); //Actually could be something else without loss of generality but okay const double gen_phi = ranGenerator->Uniform(0,2*TMath::Pi()); @@ -111,7 +128,7 @@ public: ranGenerator->SetSeed(0); // generate transverse momentum - const double gen_pT = ranGenerator->Uniform(0,5); + const double gen_pT = fSpectrumOm->GetRandom(genMinPt,genMaxPt); //Actually could be something else without loss of generality but okay const double gen_phi = ranGenerator->Uniform(0,2*TMath::Pi()); @@ -129,6 +146,26 @@ public: set4momentum(fLVHelper->Px(),fLVHelper->Py(),fLVHelper->Pz()); } +Double_t boltzPlusPower(const Double_t *x, const Double_t *p) +{ + // a plain parametrization. not meant to be physics worthy. + // adjusted to match preliminary 5 TeV shape. + + Double_t pt = x[0]; + Double_t mass = p[0]; + Double_t mt = TMath::Sqrt(pt * pt + mass * mass); + Double_t T = p[1]; + Double_t norm = p[2]; + + Double_t lowptpart = mt * TMath::Exp(-mt / T); + Double_t highptpart = p[4]*TMath::Power(x[0], p[3]); + + Double_t mixup = 1./(1.+TMath::Exp((x[0]-4.5)/.1)); + + //return pt * norm * (mixup * mt * TMath::Exp(-mt / T) + (1.-mixup)*highptpart) ; + return pt * norm * (mt * TMath::Exp(-mt / T) + (1.-mixup)*highptpart) ; +} + //__________________________________________________________________ Bool_t generateEvent() override { @@ -235,9 +272,15 @@ private: double zProd; /// z-coordinate position production vertex [cm] std::unique_ptr fLVHelper; + std::unique_ptr fSpectrumXi; + std::unique_ptr fSpectrumOm; }; FairGenerator *generator_extraStrangeness() { - return new GeneratorPythia8ExtraStrangeness(); + auto generator = new GeneratorPythia8ExtraStrangeness(); + gRandom->SetSeed(0); + generator->readString("Random:setSeed = on"); + generator->readString("Random:seed =" + std::to_string(gRandom->Integer(900000000 - 2) + 1)); + return generator; } diff --git a/MC/config/PWGLF/pythia8/generator_pythia8_syntheFlow.C b/MC/config/PWGLF/pythia8/generator_pythia8_syntheFlow.C index e96dae792..5800c0b31 100644 --- a/MC/config/PWGLF/pythia8/generator_pythia8_syntheFlow.C +++ b/MC/config/PWGLF/pythia8/generator_pythia8_syntheFlow.C @@ -7,6 +7,9 @@ #include "TRandom3.h" #include "TParticlePDG.h" #include "TDatabasePDG.h" +#include "CCDB/BasicCCDBManager.h" +#include "TH1F.h" +#include "TH1D.h" #include #include @@ -19,15 +22,25 @@ public: lutGen = new o2::eventgen::FlowMapper(); // -------- CONFIGURE SYNTHETIC FLOW ------------ - // specify a v2 vs pT here - TFile *filehep = new TFile("/Users/daviddc/Downloads/HEPData-ins1116150-v1-Table_1.root", "READ"); - TH1D *hv = (TH1D*) filehep->Get("Table 1/Hist1D_y6"); - - TFile *fileEcc = new TFile("/Users/daviddc/Downloads/eccentricityvsb.root", "READ"); - TH1D *hEccentricities = (TH1D*) fileEcc->Get("hEccentricities"); + // establish connection to ccdb + o2::ccdb::CcdbApi ccdb_api; + ccdb_api.init("https://alice-ccdb.cern.ch"); + + // config was placed at midpoint of run 544122, retrieve that + std::map metadataRCT, headers; + headers = ccdb_api.retrieveHeaders("RCT/Info/RunInformation/544122", metadataRCT, -1); + int64_t tsSOR = atol(headers["SOR"].c_str()); + int64_t tsEOR = atol(headers["EOR"].c_str()); + int64_t midRun = 0.5*tsSOR+0.5*tsEOR; + + map metadata; // can be empty + auto list = ccdb_api.retrieveFromTFileAny("Users/d/ddobrigk/syntheflow", metadata, midRun); + + TH1D *hv2vspT = (TH1D*) list->FindObject("hFlowVsPt_ins1116150_v1_Table_1"); + TH1D *heccvsb = (TH1D*) list->FindObject("hEccentricityVsB"); cout<<"Generating LUT for flow test"<CreateLUT(hv, hEccentricities); + lutGen->CreateLUT(hv2vspT, heccvsb); cout<<"Finished creating LUT!"<SetSeed(0); + generator->readString("Random:setSeed = on"); + generator->readString("Random:seed =" + std::to_string(gRandom->Integer(900000000 - 2) + 1)); + return generator; } diff --git a/MC/config/PWGLF/pythia8/generator_pythia8_syntheFlowXi.C b/MC/config/PWGLF/pythia8/generator_pythia8_syntheFlowXi.C new file mode 100644 index 000000000..8042d58ab --- /dev/null +++ b/MC/config/PWGLF/pythia8/generator_pythia8_syntheFlowXi.C @@ -0,0 +1,343 @@ + +#include "Pythia8/Pythia.h" +#include "Pythia8/HeavyIons.h" +#include "FairGenerator.h" +#include "FairPrimaryGenerator.h" +#include "Generators/GeneratorPythia8.h" +#include "TRandom3.h" +#include "TParticlePDG.h" +#include "TDatabasePDG.h" +#include "CCDB/BasicCCDBManager.h" +#include "TH1F.h" +#include "TH1D.h" + +#include +#include + +class GeneratorPythia8SyntheFlowXi : public o2::eventgen::GeneratorPythia8 +{ +public: + /// Constructor + GeneratorPythia8SyntheFlowXi() { + lutGen = new o2::eventgen::FlowMapper(); + + // -------- CONFIGURE SYNTHETIC FLOW ------------ + // establish connection to ccdb + o2::ccdb::CcdbApi ccdb_api; + ccdb_api.init("https://alice-ccdb.cern.ch"); + + // config was placed at midpoint of run 544122, retrieve that + std::map metadataRCT, headers; + headers = ccdb_api.retrieveHeaders("RCT/Info/RunInformation/544122", metadataRCT, -1); + int64_t tsSOR = atol(headers["SOR"].c_str()); + int64_t tsEOR = atol(headers["EOR"].c_str()); + int64_t midRun = 0.5*tsSOR+0.5*tsEOR; + + map metadata; // can be empty + auto list = ccdb_api.retrieveFromTFileAny("Users/d/ddobrigk/syntheflow", metadata, midRun); + + TH1D *hv2vspT = (TH1D*) list->FindObject("hFlowVsPt_ins1116150_v1_Table_1"); + TH1D *heccvsb = (TH1D*) list->FindObject("hEccentricityVsB"); + + cout<<"Generating LUT for flow test"<CreateLUT(hv2vspT, heccvsb); + cout<<"Finished creating LUT!"<(); + + fSpectrumXi = std::make_unique("fSpectrumXi", this, &GeneratorPythia8SyntheFlowXi::boltzPlusPower, 0., genMaxPt, 5, "GeneratorPythia8ExtraStrangeness", "boltzPlusPower"); + + fSpectrumXi->FixParameter(0, 1.32171); + fSpectrumXi->FixParameter(1, 4.84e-1); + fSpectrumXi->FixParameter(2, 111.9); + fSpectrumXi->FixParameter(3, -2.56511e+00); + fSpectrumXi->FixParameter(4, 1.14011e-04); + + fSpectrumOm = std::make_unique("fSpectrumOm", this, &GeneratorPythia8SyntheFlowXi::boltzPlusPower, 0., genMaxPt, 5, "GeneratorPythia8ExtraStrangeness", "boltzPlusPower"); + + fSpectrumOm->FixParameter(0, 1.67245e+00); + fSpectrumOm->FixParameter(1, 5.18174e-01); + fSpectrumOm->FixParameter(2, 1.73747e+01); + fSpectrumOm->FixParameter(3, -2.56681e+00); + fSpectrumOm->FixParameter(4, 1.87513e-04); + } + + /// Destructor + ~GeneratorPythia8SyntheFlowXi() = default; + + Double_t y2eta(Double_t pt, Double_t mass, Double_t y){ + Double_t mt = TMath::Sqrt(mass * mass + pt * pt); + return TMath::ASinH(mt / pt * TMath::SinH(y)); + } + + /// set 4-momentum + void set4momentum(double input_px, double input_py, double input_pz){ + px = input_px; + py = input_py; + pz = input_pz; + E = sqrt( m*m+px*px+py*py+pz*pz ); + fourMomentum.px(px); + fourMomentum.py(py); + fourMomentum.pz(pz); + fourMomentum.e(E); + p = sqrt( px*px+py*py+pz*pz ); + y = 0.5*log( (E+pz)/(E-pz) ); + eta = 0.5*log( (p+pz)/(p-pz) ); + } + + //__________________________________________________________________ + Pythia8::Particle createParticle(){ + //std::cout << "createParticle() mass " << m << " pdgCode " << pdg << std::endl; + Pythia8::Particle myparticle; + myparticle.id(pdg); + myparticle.status(11); + myparticle.px(px); + myparticle.py(py); + myparticle.pz(pz); + myparticle.e(E); + myparticle.m(m); + myparticle.xProd(xProd); + myparticle.yProd(yProd); + myparticle.zProd(zProd); + + return myparticle; + } + + //_________________________________________________________________________________ + /// generate uniform eta and uniform momentum + void genSpectraMomentumEtaXi(double minP, double maxP, double minY, double maxY){ + // random generator + std::unique_ptr ranGenerator { new TRandom3() }; + ranGenerator->SetSeed(0); + + // generate transverse momentum + const double gen_pT = fSpectrumXi->GetRandom(genMinPt,genMaxPt); + + //Actually could be something else without loss of generality but okay + const double gen_phi = ranGenerator->Uniform(0,2*TMath::Pi()); + + // sample flat in rapidity, calculate eta + Double_t gen_Y=10, gen_eta=10; + + while( gen_eta>genmaxEta || gen_etaUniform(minY,maxY); + //(Double_t pt, Double_t mass, Double_t y) + gen_eta = y2eta(gen_pT, m, gen_Y); + } + + fLVHelper->SetPtEtaPhiM(gen_pT, gen_eta, gen_phi, m); + set4momentum(fLVHelper->Px(),fLVHelper->Py(),fLVHelper->Pz()); + } + + //_________________________________________________________________________________ + /// generate uniform eta and uniform momentum + void genSpectraMomentumEtaOm(double minP, double maxP, double minY, double maxY){ + // random generator + std::unique_ptr ranGenerator { new TRandom3() }; + ranGenerator->SetSeed(0); + + // generate transverse momentum + const double gen_pT = fSpectrumOm->GetRandom(genMinPt,genMaxPt); + + //Actually could be something else without loss of generality but okay + const double gen_phi = ranGenerator->Uniform(0,2*TMath::Pi()); + + // sample flat in rapidity, calculate eta + Double_t gen_Y=10, gen_eta=10; + + while( gen_eta>genmaxEta || gen_etaUniform(minY,maxY); + //(Double_t pt, Double_t mass, Double_t y) + gen_eta = y2eta(gen_pT, m, gen_Y); + } + + fLVHelper->SetPtEtaPhiM(gen_pT, gen_eta, gen_phi, m); + set4momentum(fLVHelper->Px(),fLVHelper->Py(),fLVHelper->Pz()); + } + + //_________________________________________________________________________________ + /// shape function + Double_t boltzPlusPower(const Double_t *x, const Double_t *p) + { + // a plain parametrization. not meant to be physics worthy. + // adjusted to match preliminary 5 TeV shape. + + Double_t pt = x[0]; + Double_t mass = p[0]; + Double_t mt = TMath::Sqrt(pt * pt + mass * mass); + Double_t T = p[1]; + Double_t norm = p[2]; + + Double_t lowptpart = mt * TMath::Exp(-mt / T); + Double_t highptpart = p[4]*TMath::Power(x[0], p[3]); + + Double_t mixup = 1./(1.+TMath::Exp((x[0]-4.5)/.1)); + + //return pt * norm * (mixup * mt * TMath::Exp(-mt / T) + (1.-mixup)*highptpart) ; + return pt * norm * (mt * TMath::Exp(-mt / T) + (1.-mixup)*highptpart) ; + } + + //__________________________________________________________________ + Bool_t generateEvent() override { + + // Generate PYTHIA event + Bool_t lPythiaOK = kFALSE; + while (!lPythiaOK){ + lPythiaOK = mPythia.next(); + } + + //+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+ + // add extra xi/omega content + // characterise event + Long_t nParticles = mPythia.event.size(); + Long_t nChargedParticlesAtMidRap = 0; + Long_t nPionsAtMidRap = 0; + for ( Long_t j=0; j < nParticles; j++ ) { + Int_t pypid = mPythia.event[j].id(); + Float_t pyrap = mPythia.event[j].y(); + Float_t pyeta = mPythia.event[j].eta(); + + // final only + if (!mPythia.event[j].isFinal()) continue; + + if ( TMath::Abs(pyrap) < 0.5 && TMath::Abs(pypid)==211 ) nPionsAtMidRap++; + if ( TMath::Abs(pyeta) < 0.5 && TMath::Abs(mPythia.event[j].charge())>1e-5 ) nChargedParticlesAtMidRap++; + } + + // now we have the multiplicity + //+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+ + // XI ABUNDANCE FIX + //Adjust relative abundance of multi-strange particles by injecting some + Double_t lExpectedXiToPion = TMath::Max(4.74929e-03 - 4.08255e-03*TMath::Exp(-nChargedParticlesAtMidRap/4.76660e+00) - 0.00211334,0.); + Double_t lExpectedXi = 5.0*nPionsAtMidRap*lExpectedXiToPion; // extra rich, factor 5 + Int_t lXiYield = gRandom->Poisson(3*lExpectedXi); //factor 3: fix the rapidity acceptance + m = 1.32171; + pdg = 3312; + cout<<"Adding extra xi: "<Uniform()>0.5?+1:-1; + xProd=0.0; + yProd=0.0; + zProd=0.0; + genSpectraMomentumEtaXi(genMinPt,genMaxPt,genminY,genmaxY); + Pythia8::Particle lAddedParticle = createParticle(); + mPythia.event.append(lAddedParticle); + //lAddedParticles++; + } + //+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+ + + //+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+ + // OMEGA ABUNDANCE FIX + //Adjust relative abundance of multi-strange particles by injecting some + Double_t lExpectedOmegaToPion = TMath::Max(8.55057e-04 - 7.38732e-04*TMath::Exp(-nChargedParticlesAtMidRap/2.40545e+01) - 6.56785e-05,0.); + Double_t lExpectedOmega = 5.0*nPionsAtMidRap*lExpectedOmegaToPion; // extra rich, factor 5 + Int_t lOmegaYield = gRandom->Poisson(3*lExpectedOmega); //factor 3: fix the rapidity acceptance + m = 1.67245; + pdg = 3334; + cout<<"Adding extra omegas: "<Uniform()>0.5?+1:-1; + xProd=0.0; + yProd=0.0; + zProd=0.0; + genSpectraMomentumEtaOm(genMinPt,genMaxPt,genminY,genmaxY); + Pythia8::Particle lAddedParticle = createParticle(); + mPythia.event.append(lAddedParticle); + //lAddedParticles++; + } + //+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+ + + //+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+ + // loop over the entire event record and rotate all particles + // synthetic flow exercise + // first: get event plane + float eventPlaneAngle = mPythia.info.hiInfo->phi(); + float impactParameter = mPythia.info.hiInfo->b(); + + for ( Long_t j=0; j < mPythia.event.size(); j++ ) { + float pyphi = mPythia.event[j].phi(); + float pypT = mPythia.event[j].pT(); + + // calculate delta with EP + float deltaPhiEP = pyphi - eventPlaneAngle; + float shift = 0.0; + while(deltaPhiEP<0.0){ + deltaPhiEP += 2*TMath::Pi(); + shift += 2*TMath::Pi(); + } + while(deltaPhiEP>2*TMath::Pi()){ + deltaPhiEP -= 2*TMath::Pi(); + shift -= 2*TMath::Pi(); + } + float newDeltaPhiEP = lutGen->MapPhi(deltaPhiEP, impactParameter, pypT); + float pyphiNew = newDeltaPhiEP - shift + eventPlaneAngle; + + if(pyphiNew>TMath::Pi()) + pyphiNew -= 2.0*TMath::Pi(); + if(pyphiNew<-TMath::Pi()) + pyphiNew += 2.0*TMath::Pi(); + mPythia.event[j].rot(0.0, pyphiNew-pyphi); + } + //+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+ + + return true; + } + +private: + double genMinPt; /// minimum 3-momentum for generated particles + double genMaxPt; /// maximum 3-momentum for generated particles + double genminY; /// minimum pseudorapidity for generated particles + double genmaxY; /// maximum pseudorapidity for generated particles + double genminEta; + double genmaxEta; + + Pythia8::Vec4 fourMomentum; /// four-momentum (px,py,pz,E) + + double E; /// energy: sqrt( m*m+px*px+py*py+pz*pz ) [GeV/c] + double m; /// particle mass [GeV/c^2] + int pdg; /// particle pdg code + double px; /// x-component momentum [GeV/c] + double py; /// y-component momentum [GeV/c] + double pz; /// z-component momentum [GeV/c] + double p; /// momentum + double y; /// rapidity + double eta; /// pseudorapidity + double xProd; /// x-coordinate position production vertex [cm] + double yProd; /// y-coordinate position production vertex [cm] + double zProd; /// z-coordinate position production vertex [cm] + + std::unique_ptr fLVHelper; + std::unique_ptr fSpectrumXi; + std::unique_ptr fSpectrumOm; + o2::eventgen::FlowMapper *lutGen; // for mapping phi angles +}; + + FairGenerator *generator_syntheFlowXi() + { + auto generator = new GeneratorPythia8SyntheFlowXi(); + gRandom->SetSeed(0); + generator->readString("Random:setSeed = on"); + generator->readString("Random:seed =" + std::to_string(gRandom->Integer(900000000 - 2) + 1)); + return generator; + } diff --git a/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.EL3PI.DEC b/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.EL3PI.DEC new file mode 100644 index 000000000..86fed720d --- /dev/null +++ b/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.EL3PI.DEC @@ -0,0 +1,7 @@ +Decay tau- +1.0 e- anti-nu_e nu_tau PHOTOS TAULNUNU; #[Reconstructed PDG2011] +Enddecay +Decay tau+ +1.0 pi+ pi+ pi- anti-nu_tau TAUHADNU -0.108 0.775 0.149 1.364 0.400 1.23 0.4; #[Reconstructed PDG2011] +Enddecay +End diff --git a/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.ELMU.DEC b/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.ELMU.DEC new file mode 100644 index 000000000..cc0976928 --- /dev/null +++ b/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.ELMU.DEC @@ -0,0 +1,9 @@ +Decay tau- +0.5 e- anti-nu_e nu_tau PHOTOS TAULNUNU; #[Reconstructed PDG2011] +0.5 mu- anti-nu_mu nu_tau PHOTOS TAULNUNU; #[Reconstructed PDG2011] +Enddecay +Decay tau+ +0.5 e+ anti-nu_tau nu_e PHOTOS TAULNUNU; #[Reconstructed PDG2011] +0.5 mu+ anti-nu_tau nu_mu PHOTOS TAULNUNU; #[Reconstructed PDG2011] +Enddecay +End diff --git a/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.ELPI.DEC b/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.ELPI.DEC new file mode 100644 index 000000000..cabb74aa9 --- /dev/null +++ b/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.ELPI.DEC @@ -0,0 +1,10 @@ +Decay tau- +1.0 e- anti-nu_e nu_tau PHOTOS TAULNUNU; #[Reconstructed PDG2011] +Enddecay +Decay tau+ +0.25 pi+ anti-nu_tau TAUSCALARNU; #[Reconstructed PDG2011] +0.25 pi+ pi0 anti-nu_tau TAUHADNU -0.108 0.775 0.149 1.364 0.400; #[Reconstructed PDG2011] +0.25 pi0 pi0 pi+ anti-nu_tau TAUHADNU -0.108 0.775 0.149 1.364 0.400 1.23 0.4; #[Reconstructed PDG2011] +0.25 anti-nu_tau pi+ pi0 pi0 pi0 PYTHIA 41; #[Reconstructed PDG2011] +Enddecay +End diff --git a/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.PO3PI.DEC b/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.PO3PI.DEC new file mode 100644 index 000000000..b17413d5c --- /dev/null +++ b/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.PO3PI.DEC @@ -0,0 +1,7 @@ +Decay tau- +1.0 pi- pi- pi+ nu_tau TAUHADNU -0.108 0.775 0.149 1.364 0.400 1.23 0.4; #[Reconstructed PDG2011] +Enddecay +Decay tau+ +1.0 e+ anti-nu_tau nu_e PHOTOS TAULNUNU; #[Reconstructed PDG2011] +Enddecay +End diff --git a/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.POPI.DEC b/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.POPI.DEC new file mode 100644 index 000000000..e66ddc4be --- /dev/null +++ b/MC/config/PWGUD/external/generator/DecayTablesEvtGen/TAUTAU.POPI.DEC @@ -0,0 +1,10 @@ +Decay tau+ +1.0 e+ nu_e anti-nu_tau PHOTOS TAULNUNU; #[Reconstructed PDG2011] +Enddecay +Decay tau- +0.25 pi- nu_tau TAUSCALARNU; #[Reconstructed PDG2011] +0.25 pi- pi0 nu_tau TAUHADNU -0.108 0.775 0.149 1.364 0.400; #[Reconstructed PDG2011] +0.25 pi0 pi0 pi- nu_tau TAUHADNU -0.108 0.775 0.149 1.364 0.400 1.23 0.4; #[Reconstructed PDG2011] +0.25 nu_tau pi- pi0 pi0 pi0 PYTHIA 41; #[Reconstructed PDG2011] +Enddecay +End diff --git a/MC/config/PWGUD/external/generator/GeneratorStarlight.C b/MC/config/PWGUD/external/generator/GeneratorStarlight.C index 10495d485..132e57444 100644 --- a/MC/config/PWGUD/external/generator/GeneratorStarlight.C +++ b/MC/config/PWGUD/external/generator/GeneratorStarlight.C @@ -1,8 +1,11 @@ +//R__LOAD_LIBRARY(libDPMJET.so) +//R__LOAD_LIBRARY(libDpmJetLib.so) R__LOAD_LIBRARY(libStarlib.so) R__ADD_INCLUDE_PATH($STARlight_ROOT/include) #include "randomgenerator.h" #include "upcXevent.h" +#include "upcevent.h" #include "starlight.h" #include "inputParameters.h" @@ -18,6 +21,7 @@ class GeneratorStarlight_class : public Generator GeneratorStarlight_class(){}; ~GeneratorStarlight_class() = default; void selectConfiguration(std::string val) { mSelectedConfiguration = val; }; + void setExtraParams(std::string val) { mExtraParams = val; }; void setCollisionSystem(float energyCM, int beam1Z, int beam1A, int beam2Z, int beam2A) {eCM = energyCM; projZ=beam1Z; projA=beam1A; targZ=beam2Z; targA=beam2A;}; bool setParameter(std::string line) { if (not mInputParameters.setParameter(line)){ @@ -37,6 +41,7 @@ class GeneratorStarlight_class : public Generator float gamma1 = beam1energy/0.938272; float gamma2 = beam2energy/0.938272; float rapMax = 4.1 + 0.5*(TMath::ACosH(gamma2)-TMath::ACosH(gamma1)); + float dy = 0.01; const struct SLConfig { const char* name; @@ -45,58 +50,66 @@ class GeneratorStarlight_class : public Generator int nwbins; float wmin; float wmax; - float dy; int pdg_mother; bool decay_EvtGen; } slConfig[] = { - {"kTwoGammaToMuLow", 1, 13, 292, 0.4, 15.0, 0.01, -1, 0 }, // from 0.4 to 15 GeV - {"kTwoGammaToElLow", 1, 11, 292, 0.4, 15.0, 0.01, -1, 0 }, // from 0.4 to 15 GeV - {"kTwoGammaToMuMedium", 1, 13, 264, 1.8, 15.0, 0.01, -1, 0 }, // from 1.8 to 15 GeV - {"kTwoGammaToElMedium", 1, 11, 264, 1.8, 15.0, 0.01, -1, 0 }, // from 1.8 to 15 GeV - {"kTwoGammaToMuHigh", 1, 13, 220, 4.0, 15.0, 0.01, -1, 0 }, // from 4.0 to 15 GeV - {"kTwoGammaToElHigh", 1, 11, 220, 4.0, 15.0, 0.01, -1, 0 }, // from 4.0 to 15 GeV - {"kTwoGammaToRhoRho", 1, 33, 20, -1.0, -1.0, 0.01, -1, 0 }, // - {"kTwoGammaToF2", 1, 225, 20, -1.0, -1.0, 0.01, -1, 0 }, // - {"kCohRhoToPi", 3, 113, 1200, -1.0, -1.0, 0.02, 113, 0 }, // - {"kCohRhoToElEl", 3, 113011, 1200, -1.0, -1.0, 0.02, 113, 0 }, // - {"kCohRhoToMuMu", 3, 113013, 1200, -1.0, -1.0, 0.02, 113, 0 }, // - {"kCohRhoToPiWithCont", 3, 913, 1200, -1.0, -1.0, 0.02, -1, 0 }, // - {"kCohRhoToPiFlat", 3, 113, 1, -1.0, 2.5, 0.02, 113, 0 }, // - {"kCohPhiToKa", 2, 333, 20, -1.0, -1.0, 0.01, 333, 0 }, // - {"kDirectPhiToKaKa", 3, 933, 20, -1.0, -1.0, 0.01, 333, 0 }, // - {"kCohOmegaTo2Pi", 2, 223, 20, -1.0, -1.0, 0.01, 223, 0 }, // - {"kCohOmegaTo3Pi", 2, 223, 20, -1.0, -1.0, 0.01, 223, 0 }, // - {"kCohOmegaToPiPiPi", 2, 223211111, 20, -1.0, -1.0, 0.01, 333, 0 }, // - {"kCohJpsiToMu", 2, 443013, 20, -1.0, -1.0, 0.01, 443, 0 }, // - {"kCohJpsiToEl", 2, 443011, 20, -1.0, -1.0, 0.01, 443, 0 }, // - {"kCohJpsiToElRad", 2, 443011, 20, -1.0, -1.0, 0.01, 443, 1 }, // - {"kCohJpsiToProton", 2, 4432212, 20, -1.0, -1.0, 0.01, 443, 0 }, // - {"kCohPsi2sToMu", 2, 444013, 20, -1.0, -1.0, 0.01, 100443, 0 }, // - {"kCohPsi2sToEl", 2, 444011, 20, -1.0, -1.0, 0.01, 100443, 0 }, // - {"kCohPsi2sToMuPi", 2, 444013, 20, -1.0, -1.0, 0.01, 100443, 1 }, // - {"kCohPsi2sToElPi", 2, 444011, 20, -1.0, -1.0, 0.01, 100443, 1 }, // - {"kCohUpsilonToMu", 2, 553013, 20, -1.0, -1.0, 0.01, 553, 0 }, // - {"kCohUpsilonToEl", 2, 553011, 20, -1.0, -1.0, 0.01, 553, 0 }, // - {"kIncohRhoToPi", 4, 113, 1200, -1.0, -1.0, 0.02, 113, 0 }, // - {"kIncohRhoToElEl", 4, 113011, 1200, -1.0, -1.0, 0.02, 113, 0 }, // - {"kIncohRhoToMuMu", 4, 113013, 1200, -1.0, -1.0, 0.02, 113, 0 }, // - {"kIncohRhoToPiWithCont",4, 913, 1200, -1.0, -1.0, 0.02, -1, 0 }, // - {"kIncohRhoToPiFlat", 4, 113, 1, -1.0, 2.5, 0.02, 113, 0 }, // - {"kIncohPhiToKa", 4, 333, 20, -1.0, -1.0, 0.01, 333, 0 }, // - {"kIncohOmegaTo2Pi", 4, 223, 20, -1.0, -1.0, 0.01, 223, 0 }, // - {"kIncohOmegaTo3Pi", 4, 223, 20, -1.0, -1.0, 0.01, 223, 0 }, // - {"kIncohOmegaToPiPiPi", 4, 223211111, 20, -1.0, -1.0, 0.01, 223, 0 }, // - {"kIncohJpsiToMu", 4, 443013, 20, -1.0, -1.0, 0.01, 443, 0 }, // - {"kIncohJpsiToEl", 4, 443011, 20, -1.0, -1.0, 0.01, 443, 0 }, // - {"kIncohJpsiToElRad", 4, 443011, 20, -1.0, -1.0, 0.01, 443, 1 }, // - {"kIncohJpsiToProton", 4, 4432212, 20, -1.0, -1.0, 0.01, 443, 0 }, // - {"kIncohJpsiToLLbar", 4, 4433122, 20, -1.0, -1.0, 0.01, 443, 0 }, // - {"kIncohPsi2sToMu", 4, 444013, 20, -1.0, -1.0, 0.01, 100443, 0 }, // - {"kIncohPsi2sToEl", 4, 444011, 20, -1.0, -1.0, 0.01, 100443, 0 }, // - {"kIncohPsi2sToMuPi", 4, 444013, 20, -1.0, -1.0, 0.01, 100443, 1 }, // - {"kIncohPsi2sToElPi", 4, 444011, 20, -1.0, -1.0, 0.01, 100443, 1 }, // - {"kIncohUpsilonToMu", 4, 553013, 20, -1.0, -1.0, 0.01, 553, 0 }, // - {"kIncohUpsilonToEl", 4, 553011, 20, -1.0, -1.0, 0.01, 553, 0 }, // + {"kTwoGammaToMuLow", 1, 13, 876, 0.4, 15.0, -1, 0 }, // from 0.4 to 15 GeV + {"kTwoGammaToElLow", 1, 11, 876, 0.4, 15.0, -1, 0 }, // from 0.4 to 15 GeV + {"kTwoGammaToMuMedium", 1, 13, 792, 1.8, 15.0, -1, 0 }, // from 1.8 to 15 GeV + {"kTwoGammaToElMedium", 1, 11, 792, 1.8, 15.0, -1, 0 }, // from 1.8 to 15 GeV + {"kTwoGammaToMuHigh", 1, 13, 660, 4.0, 15.0, -1, 0 }, // from 4.0 to 15 GeV + {"kTwoGammaToElHigh", 1, 11, 660, 4.0, 15.0, -1, 0 }, // from 4.0 to 15 GeV + {"kTwoGammaToRhoRho", 1, 33, 20, -1.0, -1.0, -1, 0 }, // + {"kTwoGammaToF2", 1, 225, 20, -1.0, -1.0, -1, 0 }, // + {"kCohRhoToPi", 3, 113, 1200, -1.0, -1.0, 113, 0 }, // + {"kCohRhoToElEl", 3, 113011, 1200, -1.0, -1.0, 113, 0 }, // + {"kCohRhoToMuMu", 3, 113013, 1200, -1.0, -1.0, 113, 0 }, // + {"kCohRhoToPiWithCont", 3, 913, 1200, -1.0, -1.0, 113, 0 }, // + {"kCohRhoToPiFlat", 3, 113, 1, -1.0, 2.5, 113, 0 }, // + {"kCohPhiToKa", 2, 333, 20, -1.0, -1.0, 333, 0 }, // + {"kCohPhiToEl", 2, 333011, 20, -1.0, -1.0, 333, 0 }, // + {"kDirectPhiToKaKa", 3, 933, 20, -1.0, -1.0, 333, 0 }, // + {"kCohOmegaTo2Pi", 2, 223, 20, -1.0, -1.0, 223, 0 }, // + {"kCohOmegaTo3Pi", 2, 223, 20, -1.0, -1.0, 223, 1 }, // + {"kCohOmegaToPiPiPi", 2, 223211111, 20, -1.0, -1.0, 233, 0 }, // + {"kCohRhoPrimeTo4Pi", 3, 999, 1200, -1.0, 5.0, 30113, 0 }, // + {"kCohJpsiToMu", 2, 443013, 20, -1.0, -1.0, 443, 0 }, // + {"kCohJpsiToEl", 2, 443011, 20, -1.0, -1.0, 443, 0 }, // + {"kCohJpsiToElRad", 2, 443011, 20, -1.0, -1.0, 443, 1 }, // + {"kCohJpsiToProton", 2, 4432212, 20, -1.0, -1.0, 443, 0 }, // + {"kCohPsi2sToMu", 2, 444013, 20, -1.0, -1.0, 100443, 0 }, // + {"kCohPsi2sToEl", 2, 444011, 20, -1.0, -1.0, 100443, 0 }, // + {"kCohPsi2sToMuPi", 2, 444013, 20, -1.0, -1.0, 100443, 1 }, // + {"kCohPsi2sToElPi", 2, 444011, 20, -1.0, -1.0, 100443, 1 }, // + {"kCohUpsilonToMu", 2, 553013, 20, -1.0, -1.0, 553, 0 }, // + {"kCohUpsilonToEl", 2, 553011, 20, -1.0, -1.0, 553, 0 }, // + {"kIncohRhoToPi", 4, 113, 1200, -1.0, -1.0, 113, 0 }, // + {"kIncohRhoToElEl", 4, 113011, 1200, -1.0, -1.0, 113, 0 }, // + {"kIncohRhoToMuMu", 4, 113013, 1200, -1.0, -1.0, 113, 0 }, // + {"kIncohRhoToPiWithCont",4, 913, 1200, -1.0, -1.0, 113, 0 }, // + {"kIncohRhoToPiFlat", 4, 113, 1, -1.0, 2.5, 113, 0 }, // + {"kIncohPhiToKa", 4, 333, 20, -1.0, -1.0, 333, 0 }, // + {"kIncohOmegaTo2Pi", 4, 223, 20, -1.0, -1.0, 223, 0 }, // + {"kIncohOmegaTo3Pi", 4, 223, 20, -1.0, -1.0, 223, 1 }, // + {"kIncohOmegaToPiPiPi", 4, 223211111, 20, -1.0, -1.0, 223, 0 }, // + {"kIncohRhoPrimeTo4Pi", 4, 999, 1200, -1.0, 5.0, 30113, 0 }, // + {"kIncohJpsiToMu", 4, 443013, 20, -1.0, -1.0, 443, 0 }, // + {"kIncohJpsiToEl", 4, 443011, 20, -1.0, -1.0, 443, 0 }, // + {"kIncohJpsiToElRad", 4, 443011, 20, -1.0, -1.0, 443, 1 }, // + {"kIncohJpsiToProton", 4, 4432212, 20, -1.0, -1.0, 443, 0 }, // + {"kIncohJpsiToLLbar", 4, 4433122, 20, -1.0, -1.0, 443, 0 }, // + {"kIncohPsi2sToMu", 4, 444013, 20, -1.0, -1.0, 100443, 0 }, // + {"kIncohPsi2sToEl", 4, 444011, 20, -1.0, -1.0, 100443, 0 }, // + {"kIncohPsi2sToMuPi", 4, 444013, 20, -1.0, -1.0, 100443, 1 }, // + {"kIncohPsi2sToElPi", 4, 444011, 20, -1.0, -1.0, 100443, 1 }, // + {"kIncohUpsilonToMu", 4, 553013, 20, -1.0, -1.0, 553, 0 }, // + {"kIncohUpsilonToEl", 4, 553011, 20, -1.0, -1.0, 553, 0 }, // +// {"kDpmjetSingle", 5, 113, 20, -1.0, -1.0, -1, 0 }, // + {"kTauLowToEl3Pi", 1, 15, 990, 3.5, 20.0, -1, 1 }, // from 0.4 to 15 GeV + {"kTauLowToPo3Pi", 1, 15, 990, 3.5, 20.0, -1, 1 }, // from 0.4 to 15 GeV + {"kTauLowToElMu", 1, 15, 990, 3.5, 20.0, -1, 1 }, // from 0.4 to 15 GeV + {"kTauLowToElPiPi0", 1, 15, 990, 3.5, 20.0, -1, 1 }, // from 0.4 to 15 GeV + {"kTauLowToPoPiPi0", 1, 15, 990, 3.5, 20.0, -1, 1 }, // from 0.4 to 15 GeV }; const int nProcess = sizeof(slConfig)/sizeof(SLConfig); @@ -131,7 +144,7 @@ class GeneratorStarlight_class : public Generator setParameter(Form("W_MIN = %.1f #Min value of w",slConfig[idx].wmin)); setParameter(Form("W_N_BINS = %3i #Bins i w",slConfig[idx].nwbins)); setParameter(Form("RAP_MAX = %.2f #max y",rapMax)); - setParameter(Form("RAP_N_BINS = %.0f #Bins i y",rapMax*2./slConfig[idx].dy)); + setParameter(Form("RAP_N_BINS = %.0f #Bins i y",rapMax*2./dy)); setParameter("CUT_PT = 0 #Cut in pT? 0 = (no, 1 = yes)"); setParameter("PT_MIN = 0 #Minimum pT in GeV"); setParameter("PT_MAX = 10 #Maximum pT in GeV"); @@ -146,11 +159,22 @@ class GeneratorStarlight_class : public Generator setParameter("IF_STRENGTH = 1. #% of interfernce (0.0 - 0.1)"); setParameter("INT_PT_MAX = 0.24 #Maximum pt considered, when interference is turned on"); setParameter("INT_PT_N_BINS = 120 #Number of pt bins when interference is turned on"); - setParameter("XSEC_METHOD = 1 # Set to 0 to use old method for calculating gamma-gamma luminosity"); //CM - setParameter("BSLOPE_DEFINITION = 0"); // using default slope + setParameter("XSEC_METHOD = 0 # Set to 0 to use old method for calculating gamma-gamma luminosity"); //CM + setParameter("BSLOPE_DEFINITION = 2"); // using default slope setParameter("BSLOPE_VALUE = 4.0"); // default slope value setParameter("PRINT_VM = 0"); // print cross sections and fluxes vs rapidity in stdout for VM photoproduction processes - + + // Photonuclear specific options, energies in Lab frame. These values should be within the range of the values specified in the DPMJet input file (when DPMJet is used) + if(slConfig[idx].prod_mode == 5 || slConfig[idx].prod_mode == 6 || slConfig[idx].prod_mode == 7){ + setParameter("MIN_GAMMA_ENERGY = 1000.0"); + setParameter("MAX_GAMMA_ENERGY = 600000.0"); + } + + TString extraPars(mExtraParams); + TString token; + Ssiz_t from = 0; + while(extraPars.Tokenize(token, from, ";"))setParameter(token.Data()); + if (not mInputParameters.init()) { std::cout << "InitStarLight parameter initialization has failed" << std::endl; return false; @@ -171,6 +195,11 @@ class GeneratorStarlight_class : public Generator return false; } + if (mInputParameters.interactionType() >= 5) { + mUpcEvent = mStarLight->produceUpcEvent(); + mUpcEvent.boost(0.5*(TMath::ACosH(mInputParameters.beam1LorentzGamma()) - TMath::ACosH(mInputParameters.beam2LorentzGamma()))); + } + mEvent = mStarLight->produceEvent(); // boost event to the experiment CM frame mEvent.boost(0.5*(TMath::ACosH(mInputParameters.beam1LorentzGamma()) - TMath::ACosH(mInputParameters.beam2LorentzGamma()))); @@ -185,17 +214,28 @@ class GeneratorStarlight_class : public Generator { int nVtx(0); float vtx(0), vty(0), vtz(0), vtt(0); - const std::vector* slVtx(mEvent.getVertices()); + const std::vector* slVtx; + const std::vector* slPartArr; + int npart = 0; + + if (mInputParameters.interactionType() >= 5) { + slVtx = mUpcEvent.getVertices(); + slPartArr = mUpcEvent.getParticles(); + npart = mUpcEvent.getParticles()->size(); + } + else{ + slVtx = mEvent.getVertices(); + slPartArr = mEvent.getParticles(); + npart = mEvent.getParticles()->size(); + } + if (slVtx == 0) { // not vertex assume 0,0,0,0; vtx = vty = vtz = vtt = 0.0; - } else { // a vertex exits - slVtx = mEvent.getVertices(); - nVtx = slVtx->size(); - } // end if - - const std::vector* slPartArr(mEvent.getParticles()); - const int npart(mEvent.getParticles()->size()); - + } + else { // a vertex exits + nVtx = slVtx->size(); + } // end if + if(mPdgMother != -1){ //Reconstruct mother particle for VM processes TLorentzVector lmoth; TLorentzVector ldaug; @@ -219,7 +259,7 @@ class GeneratorStarlight_class : public Generator mParticles.push_back(particle); o2::mcutils::MCGenHelper::encodeParticleStatusAndTracking(mParticles.back(), false); } - if(!mDecayEvtGen){ // Don't import daughters in case of external decayer + if(!mDecayEvtGen || mPdgMother == -1){ // Don't import daughters in case of external decayer for(int ipart=0;ipartat(ipart))); if (nVtx < 1) { // No verticies @@ -232,10 +272,10 @@ class GeneratorStarlight_class : public Generator } // end if TParticle particle(slPart->getPdgCode(), 1, - 0, + (mPdgMother != -1 ? 0 :-1), + -1, + -1, -1, - slPart->getFirstDaughter(), - slPart->getLastDaughter(), slPart->GetPx(), slPart->GetPy(), slPart->GetPz(), @@ -245,7 +285,7 @@ class GeneratorStarlight_class : public Generator mParticles.push_back(particle); o2::mcutils::MCGenHelper::encodeParticleStatusAndTracking(mParticles.back(), true); } - } + } return true; } @@ -254,7 +294,9 @@ class GeneratorStarlight_class : public Generator inputParameters mInputParameters; // simulation input information. randomGenerator mRandomGenerator; // STARLIGHT's own random generator upcXEvent mEvent; // object holding STARlight simulated event. + upcEvent mUpcEvent; std::string mSelectedConfiguration = ""; + std::string mExtraParams = ""; int mPdgMother = -1; bool mDecayEvtGen = 0; float eCM = 5020; //CMS energy @@ -270,11 +312,12 @@ class GeneratorStarlight_class : public Generator FairGenerator* - GeneratorStarlight(std::string configuration = "empty",float energyCM = 5020, int beam1Z = 82, int beam1A = 208, int beam2Z = 82, int beam2A = 208) + GeneratorStarlight(std::string configuration = "empty",float energyCM = 5020, int beam1Z = 82, int beam1A = 208, int beam2Z = 82, int beam2A = 208, std::string extrapars = "") { auto gen = new o2::eventgen::GeneratorStarlight_class(); gen->selectConfiguration(configuration); gen->setCollisionSystem(energyCM, beam1Z, beam1A, beam2Z, beam2A); + gen->setExtraParams(extrapars); return gen; } diff --git a/MC/config/PWGUD/external/generator/GeneratorStarlightToEvtGen.C b/MC/config/PWGUD/external/generator/GeneratorStarlightToEvtGen.C index 3fe6e17c0..84b6a6e4e 100644 --- a/MC/config/PWGUD/external/generator/GeneratorStarlightToEvtGen.C +++ b/MC/config/PWGUD/external/generator/GeneratorStarlightToEvtGen.C @@ -4,30 +4,38 @@ // o2-sim -j 4 -n 10 -g external -t external -m "PIPE ITS TPC" -o sgn --configFile GeneratorHF_bbbar_midy.ini // // -R__ADD_INCLUDE_PATH($O2DPG_ROOT/MC/config/PWGDQ/EvtGen) -R__ADD_INCLUDE_PATH($O2DPG_ROOT/MC/config/PWGUD/external/generator) +R__ADD_INCLUDE_PATH($O2DPG_MC_CONFIG_ROOT/MC/config/PWGDQ/EvtGen) +R__ADD_INCLUDE_PATH($O2DPG_MC_CONFIG_ROOT/MC/config/PWGUD/external/generator) #include "GeneratorEvtGen.C" #include "GeneratorStarlight.C" FairGenerator* - GeneratorStarlightToEvtGen(std::string configuration = "empty",float energyCM = 5020, int beam1Z = 82, int beam1A = 208, int beam2Z = 82, int beam2A = 208) + GeneratorStarlightToEvtGen(std::string configuration = "empty",float energyCM = 5020, int beam1Z = 82, int beam1A = 208, int beam2Z = 82, int beam2A = 208, std::string extrapars = "") { auto gen = new o2::eventgen::GeneratorEvtGen(); gen->selectConfiguration(configuration); gen->setCollisionSystem(energyCM, beam1Z, beam1A, beam2Z, beam2A); + gen->setExtraParams(extrapars); - gen->SetSizePdg(3); + gen->SetSizePdg(5); gen->AddPdg(443,0); gen->AddPdg(100443,1); gen->AddPdg(223,2); - gen->SetPolarization(1); //Transversal + gen->AddPdg(15,3); + gen->AddPdg(-15,4); + if (configuration.find("kTau") == std::string::npos) gen->SetPolarization(1); //Transversal - TString pathO2 = gSystem->ExpandPathName("$O2DPG_ROOT/MC/config/PWGUD/external/generator/DecayTablesEvtGen"); + TString pathO2 = gSystem->ExpandPathName("$O2DPG_MC_CONFIG_ROOT/MC/config/PWGUD/external/generator/DecayTablesEvtGen"); if (configuration.find("Psi2sToMuPi") != std::string::npos) gen->SetDecayTable(Form("%s/PSI2S.MUMUPIPI.DEC",pathO2.Data())); else if (configuration.find("Psi2sToElPi") != std::string::npos) gen->SetDecayTable(Form("%s/PSI2S.EEPIPI.DEC",pathO2.Data())); else if (configuration.find("RhoPrime") != std::string::npos) gen->SetDecayTable(Form("%s/RHOPRIME.RHOPIPI.DEC",pathO2.Data())); else if (configuration.find("OmegaTo3Pi") != std::string::npos) gen->SetDecayTable(Form("%s/OMEGA.3PI.DEC",pathO2.Data())); else if (configuration.find("JpsiToElRad") != std::string::npos) gen->SetDecayTable(Form("%s/JPSI.EE.DEC",pathO2.Data())); + else if (configuration.find("ToEl3Pi") != std::string::npos) gen->SetDecayTable(Form("%s/TAUTAU.EL3PI.DEC",pathO2.Data())); + else if (configuration.find("ToPo3Pi") != std::string::npos) gen->SetDecayTable(Form("%s/TAUTAU.PO3PI.DEC",pathO2.Data())); + else if (configuration.find("ToElMu") != std::string::npos) gen->SetDecayTable(Form("%s/TAUTAU.ELMU.DEC",pathO2.Data())); + else if (configuration.find("ToElPiPi0") != std::string::npos) gen->SetDecayTable(Form("%s/TAUTAU.ELPI.DEC",pathO2.Data())); + else if (configuration.find("ToPoPiPi0") != std::string::npos) gen->SetDecayTable(Form("%s/TAUTAU.POPI.DEC",pathO2.Data())); return gen; -} \ No newline at end of file +} diff --git a/MC/config/PWGUD/ini/dpmjet_PbPb5360.txt b/MC/config/PWGUD/ini/dpmjet_PbPb5360.txt new file mode 100644 index 000000000..a7e25193e --- /dev/null +++ b/MC/config/PWGUD/ini/dpmjet_PbPb5360.txt @@ -0,0 +1,45 @@ +********************************************************************** +* Example for a DTUNUC input file. +* Uncomment the input-cards according to your requirements. +* +* Format: A10,6E10.0,A8 +* (except for the section enclosed by "PHOINPUT" and "ENDINPUT" +* which is format-free) +* lines starting with "*" are comment lines +********************************************************************** +* +* projectile / target / Energy +* ------------------- +* 1 2 3 4 5 6 7 +*23456789012345678901234567890123456789012345678901234567890123456789012345678 +PROJPAR 0.0 PHOTON +TARPAR 208.0 82.0 +ENERGY 1000.0 600000.0 +*ENERGY 100.0 +* Initialize the random number generator +RNDMINIT 55.0 101.0 15.0 73.0 +* +* +* PHOJET-specific input +* --------------------- +* The following lines control the event-generation with PHOJET for +* individual photon/nucleon-nucleon collisions. +* For details see the PHOJET-manual available at +* http://lepton.bartol.udel.edu/~eng/phojet.html +* Any options explained in the PHOJET-manual can be used in between +* the "PHOINPUT" and "ENDINPUT" cards. +PHOINPUT +PROCESS 1 0 1 1 1 1 1 1 +ENDINPUT +* +* Output +* ------ +* some default output (particle multiplicities etc.) +HISTOGRAM 101.0 102.0 +* +* Start of event generation +* ------------------------- +*START 5000.0 0.0 +START 100.0 0.0 +STOP +*...+....1....+....2....+....3....+....4....+....5....+....6....+....7... diff --git a/MC/config/PWGUD/ini/makeStarlightConfig.py b/MC/config/PWGUD/ini/makeStarlightConfig.py index 17d09a442..2cc92456f 100755 --- a/MC/config/PWGUD/ini/makeStarlightConfig.py +++ b/MC/config/PWGUD/ini/makeStarlightConfig.py @@ -17,13 +17,16 @@ parser.add_argument('--rapidity', default='cent', choices=['cent_rap', 'muon_rap', 'cent_eta', 'muon_eta'], help='Rapidity to select') -parser.add_argument('--process',default=None, choices=['kTwoGammaToMuLow', 'kTwoGammaToElLow', 'kTwoGammaToMuMedium', 'kTwoGammaToElMedium', 'kTwoGammaToMuHigh', 'kTwoGammaToElHigh', 'kTwoGammaToRhoRho', 'kTwoGammaToF2', 'kCohRhoToPi', 'kCohRhoToElEl', 'kCohRhoToMuMu', 'kCohRhoToPiWithCont', 'kCohRhoToPiFlat', 'kCohPhiToKa', 'kDirectPhiToKaKa','kCohOmegaTo2Pi', 'kCohOmegaTo3Pi', 'kCohOmegaToPiPiPi', 'kCohJpsiToMu', 'kCohJpsiToEl', 'kCohJpsiToElRad', 'kCohJpsiToProton', 'kCohPsi2sToMu','kCohPsi2sToEl', 'kCohPsi2sToMuPi', 'kCohPsi2sToElPi', 'kCohUpsilonToMu', 'kCohUpsilonToEl', 'kIncohRhoToPi', 'kIncohRhoToElEl', 'kIncohRhoToMuMu', 'kIncohRhoToPiWithCont', 'kIncohRhoToPiFlat', 'kIncohPhiToKa', 'kIncohOmegaTo2Pi', 'kIncohOmegaTo3Pi', 'kIncohOmegaToPiPiPi', 'kIncohJpsiToMu', 'kIncohJpsiToEl', 'kIncohJpsiToElRad', 'kIncohJpsiToProton', 'kIncohJpsiToLLbar', 'kIncohPsi2sToMu', 'kIncohPsi2sToEl', 'kIncohPsi2sToMuPi', 'kIncohPsi2sToElPi', 'kIncohUpsilonToMu', 'kIncohUpsilonToEl'], +parser.add_argument('--process',default=None, choices=['kTwoGammaToMuLow', 'kTwoGammaToElLow', 'kTwoGammaToMuMedium', 'kTwoGammaToElMedium', 'kTwoGammaToMuHigh', 'kTwoGammaToElHigh', 'kTwoGammaToRhoRho', 'kTwoGammaToF2', 'kCohRhoToPi', 'kCohRhoToElEl', 'kCohRhoToMuMu', 'kCohRhoToPiWithCont', 'kCohRhoToPiFlat', 'kCohPhiToKa', 'kDirectPhiToKaKa', 'kCohPhiToEl', 'kCohOmegaTo2Pi', 'kCohOmegaTo3Pi', 'kCohOmegaToPiPiPi', 'kCohRhoPrimeTo4Pi', 'kCohJpsiToMu', 'kCohJpsiToEl', 'kCohJpsiToElRad', 'kCohJpsiToProton', 'kCohPsi2sToMu','kCohPsi2sToEl', 'kCohPsi2sToMuPi', 'kCohPsi2sToElPi', 'kCohUpsilonToMu', 'kCohUpsilonToEl', 'kIncohRhoToPi', 'kIncohRhoToElEl', 'kIncohRhoToMuMu', 'kIncohRhoToPiWithCont', 'kIncohRhoToPiFlat', 'kIncohPhiToKa', 'kIncohOmegaTo2Pi', 'kIncohOmegaTo3Pi', 'kIncohOmegaToPiPiPi', 'kIncohRhoPrimeTo4Pi', 'kIncohJpsiToMu', 'kIncohJpsiToEl', 'kIncohJpsiToElRad', 'kIncohJpsiToProton', 'kIncohJpsiToLLbar', 'kIncohPsi2sToMu', 'kIncohPsi2sToEl', 'kIncohPsi2sToMuPi', 'kIncohPsi2sToElPi', 'kIncohUpsilonToMu', 'kIncohUpsilonToEl', 'kDpmjetSingle', 'kTauLowToEl3Pi', 'kTauLowToPo3Pi', 'kTauMediumToEl3Pi', 'kTauMediumToPo3Pi', 'kTauHighToEl3Pi', 'kTauHighToPo3Pi', 'kTauLowToElMu', 'kTauLowToElPiPi0', 'kTauLowToPoPiPi0'], help='Process to switch on') parser.add_argument('--output', default='GenStarlight.ini', help='Where to write the configuration') +parser.add_argument('--extraPars', default='', + help='Extra parameters for SL config') + args = parser.parse_args() @@ -68,17 +71,17 @@ ### Generator fout.write('[GeneratorExternal] \n') -if 'Psi2sToMuPi' in args.process or 'Psi2sToElPi' in args.process or 'RhoPrime' in args.process or 'OmegaTo3Pi' in args.process or 'JpsiToElRad' in args.process : - fout.write('fileName = ${O2DPG_ROOT}/MC/config/PWGUD/external/generator/GeneratorStarlightToEvtGen.C \n') - fout.write('funcName = GeneratorStarlightToEvtGen("%s", %f, %d, %d, %d, %d) \n' % (args.process,args.eCM ,pZ,pA,tZ,tA)) +if 'Psi2sToMuPi' in args.process or 'Psi2sToElPi' in args.process or 'OmegaTo3Pi' in args.process or 'JpsiToElRad' in args.process or 'kTau' in args.process: + fout.write('fileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGUD/external/generator/GeneratorStarlightToEvtGen.C \n') + fout.write('funcName = GeneratorStarlightToEvtGen("%s", %f, %d, %d, %d, %d, "%s") \n' % (args.process,args.eCM ,pZ,pA,tZ,tA,args.extraPars)) else: - fout.write('fileName = ${O2DPG_ROOT}/MC/config/PWGUD/external/generator/GeneratorStarlight.C \n') - fout.write('funcName = GeneratorStarlight("%s", %f, %d, %d, %d, %d) \n' % (args.process,args.eCM ,pZ,pA,tZ,tA)) + fout.write('fileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGUD/external/generator/GeneratorStarlight.C \n') + fout.write('funcName = GeneratorStarlight("%s", %f, %d, %d, %d, %d, "%s") \n' % (args.process,args.eCM ,pZ,pA,tZ,tA,args.extraPars)) ###Trigger fout.write('[TriggerExternal] \n') -fout.write('fileName = ${O2DPG_ROOT}/MC/config/PWGUD/trigger/selectParticlesInAcceptance.C \n') -if 'kTwoGamma' in args.process: +fout.write('fileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGUD/trigger/selectParticlesInAcceptance.C \n') +if 'kTwoGamma' in args.process or 'kTau' in args.process: if args.rapidity == 'cent_eta': fout.write('funcName = selectDirectPartInAcc(-0.9,0.9) \n') if args.rapidity == 'muon_eta': diff --git a/MC/config/PWGUD/trigger/selectParticlesInAcceptance.C b/MC/config/PWGUD/trigger/selectParticlesInAcceptance.C index 9a8d39407..513a93c13 100644 --- a/MC/config/PWGUD/trigger/selectParticlesInAcceptance.C +++ b/MC/config/PWGUD/trigger/selectParticlesInAcceptance.C @@ -24,13 +24,24 @@ o2::eventgen::Trigger selectDaughterPartInAcc(double etaMin = -1., double etaMax for (const auto& particle : particles) { if (particle.GetFirstMother() == -1) if ((particle.Y() < etaMin) || (particle.Y() > etaMax)) return kFALSE; - if (particle.GetFirstMother() != -1 && particle.GetFirstDaughter() == -1 && particle.GetPdgCode() != 22) + if (particle.GetFirstMother() != -1 && particle.GetFirstDaughter() == -1 && particle.GetPdgCode() != 22 && TMath::Abs(particle.GetPdgCode()) != 12 && TMath::Abs(particle.GetPdgCode()) != 14 && TMath::Abs(particle.GetPdgCode()) != 16) if ((particle.Eta() < etaMin) || (particle.Eta() > etaMax)) return kFALSE; } return kTRUE; }; } +o2::eventgen::Trigger selectDileptonsInAcc(double etaMin = -1., double etaMax = -1.) +{ + return [etaMin, etaMax](const std::vector& particles) -> bool { + for (const auto& particle : particles) { + if (particle.GetFirstMother() != -1 && particle.GetFirstDaughter() == -1 && (TMath::Abs(particle.GetPdgCode()) == 11 || TMath::Abs(particle.GetPdgCode() == 13))) + if ((particle.Eta() < etaMin) || (particle.Eta() > etaMax)) return kFALSE; + } + return kTRUE; + }; +} + o2::eventgen::Trigger selectDirectPartInAcc(double etaMin = -1., double etaMax = -1.) { return [etaMin, etaMax](const std::vector& particles) -> bool { @@ -41,4 +52,4 @@ o2::eventgen::Trigger selectDirectPartInAcc(double etaMin = -1., double etaMax = } return kTRUE; }; -} \ No newline at end of file +} diff --git a/MC/config/QC/json/emc-bc-task.json b/MC/config/QC/json/emc-bc-task.json index 7cb788561..38e5c26eb 100644 --- a/MC/config/QC/json/emc-bc-task.json +++ b/MC/config/QC/json/emc-bc-task.json @@ -36,7 +36,7 @@ "maxNumberCycles": "-1", "dataSource": { "type": "direct", - "query": "emcal-triggers:EMC/CELLSTRGR/0;ctp-digits:CTP/DIGITS;ctp-config:CTP/CONFIG/0?lifetime=condition&ccdb-path=CTP/Config/Config&ccdb-run-dependent=true" + "query": "emcal-triggers:EMC/CELLSTRGR/0;ctp-digits:CTP/DIGITS;ctp-config:CTP/CONFIG/0?lifetime=condition&ccdb-path=CTP/Config/Config&ccdb-run-dependent=1" }, "taskParameters": { "AliasMB": "MINBIAS_TVX_L0,C0TVX,CMTVXTSC,C0TVXTSC,C0TVXTCE" diff --git a/MC/config/examples/epos4/epos.sh b/MC/config/examples/epos4/epos.sh new file mode 100755 index 000000000..c416e03f7 --- /dev/null +++ b/MC/config/examples/epos4/epos.sh @@ -0,0 +1,41 @@ +#!/bin/sh +# Script based on CRMC example +# EPOS4 option files must contain ihepmc set to 2 to print HepMC +# data on stdout. -hepmc flag is not needed anymore, but -hepstd is fundamental +# in order not to print useless information on stdout (a z-*optns*.mtr file will be created) + +optns="example" +seed=$RANDOM + +while test $# -gt 0 ; do + case $1 in + -i|--input) optns=$2 ; shift ;; + -s|--seed) seed=$2 ; shift ;; + -h|--help) usage; exit 0 ;; + esac + shift +done + +if [ ! -f $optns.optns ]; then + echo "Error: Options file $optns.optns not found" + exit 1 +fi + +if [ $seed -eq 0 ]; then + echo "Seed can't be 0, random number will be used" + seed=$RANDOM +fi + +# Check if the environment variable EPO4 is set (mandatory with o2dpg-sim-tests on CI machines) +# If not, set all the EPOS4 related variables, most likely they are unset as well. +if [ -z "${EPO4}" ]; then + export EPO4=$EPOS4_ROOT/epos4/ + export LIBDIR=$EPOS4_ROOT/epos4/bin + export EPO4VSN=4.0.0 + export OPT=./ + export HTO=./ + export CHK=./ +fi + +# Or filters the stdout with only HepMC2 useful data +$EPOS4_ROOT/epos4/scripts/epos -hepstd -s $seed $optns | sed -n 's/^\(HepMC::\|[EAUWVP] \)/\1/p' diff --git a/MC/config/examples/epos4/generator/example.optns b/MC/config/examples/epos4/generator/example.optns new file mode 100644 index 000000000..553ecafa8 --- /dev/null +++ b/MC/config/examples/epos4/generator/example.optns @@ -0,0 +1,32 @@ +!-------------------------------------------------------------------- +! proton-proton collision no hydro no hadronic cascade +!-------------------------------------------------------------------- + +!--------------------------------------- +! Define run +!--------------------------------------- + +application hadron !hadron-hadron, hadron-nucleus, or nucleus-nucleus +set laproj 1 !projectile atomic number +set maproj 1 !projectile mass number +set latarg 1 !target atomic number +set matarg 1 !target mass number +set ecms 13000 !sqrt(s)_pp +set istmax 25 !max status considered for storage + +ftime on !string formation time non-zero +!suppressed decays: +nodecays + 110 20 2130 -2130 2230 -2230 1130 -1130 1330 -1330 2330 -2330 3331 -3331 +end + +set ninicon 1 !number of initial conditions used for hydro evolution +core full !core/corona activated +hydro hlle !hydro activated +eos x3ff !eos activated (epos standard EoS) +hacas full !hadronic cascade activated (UrQMD) +set nfreeze 1 !number of freeze out events per hydro event +set modsho 1 !printout every modsho events +set centrality 0 !0=min bias +set ihepmc 2 !HepMC output enabled on stdout +set nfull 10 diff --git a/MC/config/examples/epos4/generator_EPOS4.C b/MC/config/examples/epos4/generator_EPOS4.C new file mode 100644 index 000000000..6d08a7b19 --- /dev/null +++ b/MC/config/examples/epos4/generator_EPOS4.C @@ -0,0 +1,61 @@ +#include +#include +#include +#include + +class GeneratorEPOS4 : public o2::eventgen::GeneratorHepMC +{ + public: + GeneratorEPOS4() = default; + ~GeneratorEPOS4() = default; + +}; + +// Next function takes the optns file as argument and edits the maximum number of events to be generated. +// When used as an external generator it is important that the events passed with the -n flag are the same +// or lower of the events set in the optns file, otherwise the generation will crash. That is why the .ini +// example file contains the maximum integer available, assuming less events than that are generated in a real +// life scenario. Unfortunately a larger number cannot be used at the moment since EPOS4 has a signed integer +// type for the nfull parameter. Might be changed in the future. +// When running locally, or on the GRID (not in hyperloop), the default parameters provided in the .ini file of the +// external generation can be overwritten using the confKeyValues option (or similar depending on the tool used). +FairGenerator* generateEPOS4(const std::string &name, const int& nEvents) +{ + // check if the file exists + auto filename = gSystem->ExpandPathName(name.c_str()); + if (!std::filesystem::exists(filename)) + { + LOG(fatal) << "Options file " << filename << " does not exist!"; + exit(1); + } + // cache all the lines of the optns file and replace the number of events + std::ifstream file(filename); + std::string line; + bool found = false; + std::stringstream buffer; + while (std::getline(file, line)) + { + if (line.find("nfull") != std::string::npos){ + // replace the number of events + found = true; + line = "set nfull " + std::to_string(nEvents); + } + buffer << line << "\n"; + } + file.close(); + // Write the updated content back to a file in the current directory + std::ofstream outFile("cfg.optns"); + outFile << buffer.str(); + if (!found) + { + outFile << "set nfull " + std::to_string(nEvents); + } + outFile.close(); + auto gen = new GeneratorEPOS4(); + auto& param0 = o2::eventgen::GeneratorFileOrCmdParam::Instance(); + auto& param = o2::eventgen::GeneratorHepMCParam::Instance(); + auto& conf = o2::conf::SimConfig::Instance(); + // setup the HepMC generator to run with automatic FIFOs + gen->setup(param0, param, conf); + return gen; +} diff --git a/MC/config/examples/ini/GeneratorEPOS4.ini b/MC/config/examples/ini/GeneratorEPOS4.ini new file mode 100644 index 000000000..33935cbe2 --- /dev/null +++ b/MC/config/examples/ini/GeneratorEPOS4.ini @@ -0,0 +1,10 @@ +[GeneratorExternal] +fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/examples/epos4/generator_EPOS4.C +funcName=generateEPOS4("${O2DPG_MC_CONFIG_ROOT}/MC/config/examples/epos4/generator/example.optns", 2147483647) + +[GeneratorFileOrCmd] +cmd=${O2DPG_MC_CONFIG_ROOT}/MC/config/examples/epos4/epos.sh -i cfg +bMaxSwitch=none + +[HepMC] +version=2 \ No newline at end of file diff --git a/MC/config/examples/ini/tests/GeneratorEPOS4.C b/MC/config/examples/ini/tests/GeneratorEPOS4.C new file mode 100644 index 000000000..5f735eddf --- /dev/null +++ b/MC/config/examples/ini/tests/GeneratorEPOS4.C @@ -0,0 +1,64 @@ +int External() +{ + std::string path{"o2sim_Kine.root"}; + // Check that file exists, can be opened and has the correct tree + TFile file(path.c_str(), "READ"); + if (file.IsZombie()) + { + std::cerr << "Cannot open ROOT file " << path << "\n"; + return 1; + } + auto tree = (TTree *)file.Get("o2sim"); + if (!tree) + { + std::cerr << "Cannot find tree o2sim in file " << path << "\n"; + return 1; + } + std::vector *tracks{}; + tree->SetBranchAddress("MCTrack", &tracks); + + // Check if all events are filled + auto nEvents = tree->GetEntries(); + for (Long64_t i = 0; i < nEvents; ++i) + { + tree->GetEntry(i); + if (tracks->empty()) + { + std::cerr << "Empty entry found at event " << i << "\n"; + return 1; + } + } + // Check if there are 100 events, as simulated in the o2dpg-test + if (nEvents != 100) + { + std::cerr << "Expected 100 events, got " << nEvents << "\n"; + return 1; + } + // check if each event has two protons with 6800 GeV of energy + // exits if the particle is not a proton + for (int i = 0; i < nEvents; i++) + { + auto check = tree->GetEntry(i); + int count = 0; + for (int idxMCTrack = 0; idxMCTrack < tracks->size(); ++idxMCTrack) + { + auto track = tracks->at(idxMCTrack); + double energy = track.GetEnergy(); + // Check if track energy is approximately equal to 6800 GeV (a tolerance of 65 keV is considered, straight equality does not work due to floating point precision) + if (std::abs(energy - 6800) < 1e-4) + { + if (track.GetPdgCode() != 2212){ + std::cerr << "Found 6800 GeV particle with pdgID " << track.GetPdgCode() << "\n"; + return 1; + } + count++; + } + } + if (count < 2) + { + std::cerr << "Event " << i << " has less than 2 protons at 6800 GeV\n"; + return 1; + } + } + return 0; +} \ No newline at end of file diff --git a/MC/run/ANCHOR/anchorMC.sh b/MC/run/ANCHOR/anchorMC.sh index 3c7b52336..83144bdc7 100755 --- a/MC/run/ANCHOR/anchorMC.sh +++ b/MC/run/ANCHOR/anchorMC.sh @@ -235,20 +235,28 @@ TIMESTAMP=`grep "Determined timestamp to be" timestampsampling_${ALIEN_JDL_LPMRU echo_info "TIMESTAMP IS ${TIMESTAMP}" # -- Create aligned geometry using ITS ideal alignment to avoid overlaps in geant -CCDBOBJECTS_IDEAL_MC="ITS/Calib/Align" -TIMESTAMP_IDEAL_MC=1 -${O2_ROOT}/bin/o2-ccdb-downloadccdbfile --host http://alice-ccdb.cern.ch/ -p ${CCDBOBJECTS_IDEAL_MC} -d ${ALICEO2_CCDB_LOCALCACHE} --timestamp ${TIMESTAMP_IDEAL_MC} -CCDB_RC="${?}" -if [ ! "${CCDB_RC}" == "0" ]; then - echo_error "Problem during CCDB prefetching of ${CCDBOBJECTS_IDEAL_MC}. Exiting." - exit ${CCDB_RC} +if [ "${ENABLE_PARALLEL_WORLD}" == "0" ]; then + CCDBOBJECTS_IDEAL_MC="ITS/Calib/Align" + TIMESTAMP_IDEAL_MC=1 + ${O2_ROOT}/bin/o2-ccdb-downloadccdbfile --host http://alice-ccdb.cern.ch/ -p ${CCDBOBJECTS_IDEAL_MC} -d ${ALICEO2_CCDB_LOCALCACHE} --timestamp ${TIMESTAMP_IDEAL_MC} + CCDB_RC="${?}" + if [ ! "${CCDB_RC}" == "0" ]; then + echo_error "Problem during CCDB prefetching of ${CCDBOBJECTS_IDEAL_MC}. Exiting." + exit ${CCDB_RC} + fi fi # TODO This can potentially be removed or if needed, should be taken over by o2dpg_sim_workflow_anchored.py and O2_dpg_workflow_runner.py -echo "run with echo in pipe" | ${O2_ROOT}/bin/o2-create-aligned-geometry-workflow --configKeyValues "HBFUtils.startTime=${TIMESTAMP}" --condition-remap=file://${ALICEO2_CCDB_LOCALCACHE}=ITS/Calib/Align -b --run +if [ "${ENABLE_PARALLEL_WORLD}" == "0" ]; then + echo "run with echo in pipe" | ${O2_ROOT}/bin/o2-create-aligned-geometry-workflow --configKeyValues "HBFUtils.startTime=${TIMESTAMP}" --condition-remap=file://${ALICEO2_CCDB_LOCALCACHE}=ITS/Calib/Align -b --run +else + echo "run with echo in pipe" | ${O2_ROOT}/bin/o2-create-aligned-geometry-workflow --configKeyValues "HBFUtils.startTime=${TIMESTAMP}" -b --run +fi mkdir -p $ALICEO2_CCDB_LOCALCACHE/GLO/Config/GeometryAligned ln -s -f $PWD/o2sim_geometry-aligned.root $ALICEO2_CCDB_LOCALCACHE/GLO/Config/GeometryAligned/snapshot.root -[[ -f $PWD/its_GeometryTGeo.root ]] && mkdir -p $ALICEO2_CCDB_LOCALCACHE/ITS/Config/Geometry && ln -s -f $PWD/its_GeometryTGeo.root $ALICEO2_CCDB_LOCALCACHE/ITS/Config/Geometry/snapshot.root +if [ "${ENABLE_PARALLEL_WORLD}" == "0" ]; then + [[ -f $PWD/its_GeometryTGeo.root ]] && mkdir -p $ALICEO2_CCDB_LOCALCACHE/ITS/Config/Geometry && ln -s -f $PWD/its_GeometryTGeo.root $ALICEO2_CCDB_LOCALCACHE/ITS/Config/Geometry/snapshot.root +fi [[ -f $PWD/mft_GeometryTGeo.root ]] && mkdir -p $ALICEO2_CCDB_LOCALCACHE/MFT/Config/Geometry && ln -s -f $PWD/mft_GeometryTGeo.root $ALICEO2_CCDB_LOCALCACHE/MFT/Config/Geometry/snapshot.root # -- RUN THE MC WORKLOAD TO PRODUCE AOD -- diff --git a/MC/run/ANCHOR/tests/test_anchor_2023_apass2_pp.sh b/MC/run/ANCHOR/tests/test_anchor_2023_apass2_pp.sh index 9d8378044..a399658f5 100755 --- a/MC/run/ANCHOR/tests/test_anchor_2023_apass2_pp.sh +++ b/MC/run/ANCHOR/tests/test_anchor_2023_apass2_pp.sh @@ -30,6 +30,9 @@ export SEED=5 # for pp and 50 events per TF, we launch only 4 workers. export NWORKERS=2 +export ENABLE_PARALLEL_WORLD=1 +export ALIEN_JDL_ANCHOR_SIM_OPTIONS="-gen pythia8 -confKey \"GeometryManagerParam.useParallelWorld=1;GeometryManagerParam.usePwGeoBVH=1;GeometryManagerParam.usePwCaching=1\"" + # run the central anchor steering script; this includes # * derive timestamp # * derive interaction rate diff --git a/MC/run/PWGDQ/runBeautyToJpsi_fwdy_PbPb.sh b/MC/run/PWGDQ/runBeautyToJpsi_fwdy_PbPb.sh index 5db05ff3d..3eaa90881 100644 --- a/MC/run/PWGDQ/runBeautyToJpsi_fwdy_PbPb.sh +++ b/MC/run/PWGDQ/runBeautyToJpsi_fwdy_PbPb.sh @@ -15,7 +15,7 @@ NWORKERS=${NWORKERS:-8} NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 5020 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ - -trigger "external" -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_fwdy.ini \ + -trigger "external" -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_fwdy.ini -interactionRate 500000 \ -genBkg pythia8 -procBkg "heavy_ion" -colBkg PbPb --embedding -nb ${NBKGEVENTS} --mft-assessment-full --fwdmatching-assessment-full --fwdmatching-save-trainingdata # run workflow diff --git a/MC/run/PWGDQ/runBeautyToJpsi_fwdy_pp.sh b/MC/run/PWGDQ/runBeautyToJpsi_fwdy_pp.sh index 25f2249e6..76f886fec 100755 --- a/MC/run/PWGDQ/runBeautyToJpsi_fwdy_pp.sh +++ b/MC/run/PWGDQ/runBeautyToJpsi_fwdy_pp.sh @@ -15,7 +15,7 @@ NWORKERS=${NWORKERS:-8} NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 13600 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ - -trigger "external" -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_fwdy.ini \ + -trigger "external" -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_fwdy.ini -interactionRate 500000 \ -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} --fwdmatching-4-param --fwdmatching-cut-4-param # run workflow diff --git a/MC/run/PWGDQ/runBeautyToJpsi_midy_pp.sh b/MC/run/PWGDQ/runBeautyToJpsi_midy_pp.sh index 1a09db46a..f8ef40880 100755 --- a/MC/run/PWGDQ/runBeautyToJpsi_midy_pp.sh +++ b/MC/run/PWGDQ/runBeautyToJpsi_midy_pp.sh @@ -15,7 +15,7 @@ NWORKERS=${NWORKERS:-8} NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 900 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ - -trigger "external" -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_midy.ini \ + -trigger "external" -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_midy.ini -interactionRate 500000 \ -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} # run workflow diff --git a/MC/run/PWGDQ/runBeautyToPsiAndJpsi_fwdy_pp_triggerGap.sh b/MC/run/PWGDQ/runBeautyToPsiAndJpsi_fwdy_pp_triggerGap.sh index 642bed9df..8a344f66d 100755 --- a/MC/run/PWGDQ/runBeautyToPsiAndJpsi_fwdy_pp_triggerGap.sh +++ b/MC/run/PWGDQ/runBeautyToPsiAndJpsi_fwdy_pp_triggerGap.sh @@ -15,7 +15,7 @@ NWORKERS=${NWORKERS:-8} NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 13600 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ - -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_PsiAndJpsi_fwdy_triggerGap.ini --mft-assessment-full --fwdmatching-assessment-full + -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_PsiAndJpsi_fwdy_triggerGap.ini --mft-assessment-full --fwdmatching-assessment-full -interactionRate 500000 # run workflow ${O2DPG_ROOT}/MC/bin/o2_dpg_workflow_runner.py -f workflow.json diff --git a/MC/run/PWGDQ/runBeautyToPsiAndJpsi_midy_pp_triggerGap.sh b/MC/run/PWGDQ/runBeautyToPsiAndJpsi_midy_pp_triggerGap.sh index 2a2edcc42..d9387ea7f 100755 --- a/MC/run/PWGDQ/runBeautyToPsiAndJpsi_midy_pp_triggerGap.sh +++ b/MC/run/PWGDQ/runBeautyToPsiAndJpsi_midy_pp_triggerGap.sh @@ -15,7 +15,7 @@ NWORKERS=${NWORKERS:-8} NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 13600 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ - -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_PsiAndJpsi_midy_triggerGap.ini + -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_PsiAndJpsi_midy_triggerGap.ini -interactionRate 500000 # run workflow ${O2DPG_ROOT}/MC/bin/o2_dpg_workflow_runner.py -f workflow.json diff --git a/MC/run/PWGDQ/runBeautyToPsi_fwdy_PbPb.sh b/MC/run/PWGDQ/runBeautyToPsi_fwdy_PbPb.sh index e791a1502..febbd98ac 100755 --- a/MC/run/PWGDQ/runBeautyToPsi_fwdy_PbPb.sh +++ b/MC/run/PWGDQ/runBeautyToPsi_fwdy_PbPb.sh @@ -15,7 +15,7 @@ NWORKERS=${NWORKERS:-4} NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 5020 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ - -trigger "external" -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_Psi2S_fwdy.ini \ + -trigger "external" -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_Psi2S_fwdy.ini -interactionRate 500000 \ -genBkg pythia8 -procBkg "heavy_ion" -colBkg PbPb --embedding -nb ${NBKGEVENTS} --mft-assessment-full --fwdmatching-assessment-full # run workflow diff --git a/MC/run/PWGDQ/runBeautyToPsi_fwdy_pp.sh b/MC/run/PWGDQ/runBeautyToPsi_fwdy_pp.sh index f1f18ae86..02b9bd27d 100755 --- a/MC/run/PWGDQ/runBeautyToPsi_fwdy_pp.sh +++ b/MC/run/PWGDQ/runBeautyToPsi_fwdy_pp.sh @@ -15,7 +15,7 @@ NWORKERS=${NWORKERS:-8} NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 900 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ - -trigger "external" -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_Psi2S_fwdy.ini \ + -trigger "external" -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_Psi2S_fwdy.ini -interactionRate 500000 \ -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} # run workflow diff --git a/MC/run/PWGDQ/runBeautyToPsi_midy_pp.sh b/MC/run/PWGDQ/runBeautyToPsi_midy_pp.sh index 001abad3b..941db0167 100755 --- a/MC/run/PWGDQ/runBeautyToPsi_midy_pp.sh +++ b/MC/run/PWGDQ/runBeautyToPsi_midy_pp.sh @@ -16,7 +16,7 @@ NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 900 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ -trigger "external" -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbar_Psi2S_midy.ini \ - -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} + -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} -interactionRate 500000 # run workflow ${O2DPG_ROOT}/MC/bin/o2_dpg_workflow_runner.py -f workflow.json diff --git a/MC/run/PWGDQ/runBplusToJpsi_midy_pp.sh b/MC/run/PWGDQ/runBplusToJpsi_midy_pp.sh index 509cb814a..c559a75ac 100755 --- a/MC/run/PWGDQ/runBplusToJpsi_midy_pp.sh +++ b/MC/run/PWGDQ/runBplusToJpsi_midy_pp.sh @@ -15,7 +15,7 @@ NWORKERS=${NWORKERS:-8} NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 13600 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ - -trigger "external" -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbarToBplus_midy.ini \ + -trigger "external" -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbarToBplus_midy.ini -interactionRate 500000 \ -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} # run workflow diff --git a/MC/run/PWGDQ/runBplusToJpsi_midy_pp_triggerGap.sh b/MC/run/PWGDQ/runBplusToJpsi_midy_pp_triggerGap.sh index 04402b17b..449e67216 100755 --- a/MC/run/PWGDQ/runBplusToJpsi_midy_pp_triggerGap.sh +++ b/MC/run/PWGDQ/runBplusToJpsi_midy_pp_triggerGap.sh @@ -15,7 +15,7 @@ NWORKERS=${NWORKERS:-8} NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 13600 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ - -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbarToBplus_midy_triggerGap.ini + -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/GeneratorHF_bbbarToBplus_midy_triggerGap.ini -interactionRate 500000 # run workflow ${O2DPG_ROOT}/MC/bin/o2_dpg_workflow_runner.py -f workflow.json diff --git a/MC/run/PWGDQ/runFwdMuBoxGen_PbPb.sh b/MC/run/PWGDQ/runFwdMuBoxGen_PbPb.sh index 7fd5435ee..0e5e257a9 100755 --- a/MC/run/PWGDQ/runFwdMuBoxGen_PbPb.sh +++ b/MC/run/PWGDQ/runFwdMuBoxGen_PbPb.sh @@ -17,7 +17,7 @@ NTIMEFRAMES=${NTIMEFRAMES:-1} NBOXMUONS=${NBOXMUONS:-2} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 5020 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 \ - -confKey "GeneratorExternal.fileName=${O2DPG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorBoxFwd.C;GeneratorExternal.funcName=fwdMuBoxGen()" \ + -confKey "GeneratorExternal.fileName=${O2DPG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorBoxFwd.C;GeneratorExternal.funcName=fwdMuBoxGen()" -interactionRate 500000 \ -genBkg pythia8 -procBkg "heavy_ion" -colBkg PbPb --embedding -nb ${NBKGEVENTS} --mft-assessment-full --fwdmatching-assessment-full --fwdmatching-save-trainingdata # run workflow (MFT-related tasks) diff --git a/MC/run/PWGDQ/runFwdMuBoxGen_pp.sh b/MC/run/PWGDQ/runFwdMuBoxGen_pp.sh index 6d154e06e..4284a5a0c 100755 --- a/MC/run/PWGDQ/runFwdMuBoxGen_pp.sh +++ b/MC/run/PWGDQ/runFwdMuBoxGen_pp.sh @@ -17,7 +17,7 @@ NTIMEFRAMES=${NTIMEFRAMES:-1} NBOXMUONS=${NBOXMUONS:-2} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 13600 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ - -confKey "GeneratorExternal.fileName=${O2DPG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorBoxFwd.C;GeneratorExternal.funcName=fwdMuBoxGen()" \ + -confKey "GeneratorExternal.fileName=${O2DPG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorBoxFwd.C;GeneratorExternal.funcName=fwdMuBoxGen()" -interactionRate 500000 \ -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} --mft-assessment-full --fwdmatching-assessment-full --fwdmatching-save-trainingdata # run workflow (MFT-related tasks) diff --git a/MC/run/PWGDQ/runPromptCharmonia_fwdy_PbPb.sh b/MC/run/PWGDQ/runPromptCharmonia_fwdy_PbPb.sh index d305c150e..82d0b10fd 100755 --- a/MC/run/PWGDQ/runPromptCharmonia_fwdy_PbPb.sh +++ b/MC/run/PWGDQ/runPromptCharmonia_fwdy_PbPb.sh @@ -16,7 +16,7 @@ NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 5020 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ -confKey "GeneratorExternal.fileName=${O2DPG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorCocktailPromptCharmoniaToMuonEvtGen_PbPb5TeV.C;GeneratorExternal.funcName=GeneratorCocktailPromptCharmoniaToMuonEvtGen_PbPb5TeV()" \ - -genBkg pythia8 -procBkg "heavy_ion" -colBkg PbPb --embedding -nb ${NBKGEVENTS} --mft-assessment-full --fwdmatching-assessment-full --fwdmatching-save-trainingdata + -genBkg pythia8 -procBkg "heavy_ion" -colBkg PbPb --embedding -nb ${NBKGEVENTS} --mft-assessment-full --fwdmatching-assessment-full --fwdmatching-save-trainingdata -interactionRate 500000 # run workflow ${O2DPG_ROOT}/MC/bin/o2_dpg_workflow_runner.py -f workflow.json diff --git a/MC/run/PWGDQ/runPromptCharmonia_fwdy_pp.sh b/MC/run/PWGDQ/runPromptCharmonia_fwdy_pp.sh index dc0c28fae..24f2c1fb5 100755 --- a/MC/run/PWGDQ/runPromptCharmonia_fwdy_pp.sh +++ b/MC/run/PWGDQ/runPromptCharmonia_fwdy_pp.sh @@ -17,7 +17,7 @@ TARGETTASK=${TARGETTASK:+-tt ${TARGETTASK}} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 13600 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ -confKey "GeneratorExternal.fileName=${O2DPG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorParamPromptJpsiToMuonEvtGen_pp13TeV.C;GeneratorExternal.funcName=GeneratorParamPromptJpsiToMuonEvtGen_pp13TeV()" \ - -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} --fwdmatching-4-param --fwdmatching-cut-4-param + -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} --fwdmatching-4-param --fwdmatching-cut-4-param -interactionRate 500000 # run workflow diff --git a/MC/run/PWGDQ/runPromptCharmonia_fwdy_pp_assessment.sh b/MC/run/PWGDQ/runPromptCharmonia_fwdy_pp_assessment.sh index e4f80946e..348863a2f 100755 --- a/MC/run/PWGDQ/runPromptCharmonia_fwdy_pp_assessment.sh +++ b/MC/run/PWGDQ/runPromptCharmonia_fwdy_pp_assessment.sh @@ -16,7 +16,7 @@ NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 900 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ -confKey "GeneratorExternal.fileName=${O2DPG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorCocktailPromptCharmoniaToMuonEvtGen_pp13TeV.C;GeneratorExternal.funcName=GeneratorCocktailPromptCharmoniaToMuonEvtGen_pp13TeV()" \ - -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} --mft-assessment-full --fwdmatching-assessment-full + -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} --mft-assessment-full --fwdmatching-assessment-full -interactionRate 500000 # run workflow diff --git a/MC/run/PWGDQ/runPromptCharmonia_fwdy_pp_assessment_triggerGap.sh b/MC/run/PWGDQ/runPromptCharmonia_fwdy_pp_assessment_triggerGap.sh index b64e13ee3..5da085250 100755 --- a/MC/run/PWGDQ/runPromptCharmonia_fwdy_pp_assessment_triggerGap.sh +++ b/MC/run/PWGDQ/runPromptCharmonia_fwdy_pp_assessment_triggerGap.sh @@ -15,7 +15,7 @@ NWORKERS=${NWORKERS:-8} NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 13600 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ - -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/Generator_InjectedPromptCharmoniaFwdy_TriggerGap.ini --mft-assessment-full --fwdmatching-assessment-full + -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/Generator_InjectedPromptCharmoniaFwdy_TriggerGap.ini --mft-assessment-full --fwdmatching-assessment-full -interactionRate 500000 # run workflow ${O2DPG_ROOT}/MC/bin/o2_dpg_workflow_runner.py -f workflow.json diff --git a/MC/run/PWGDQ/runPromptCharmonia_midy_pp.sh b/MC/run/PWGDQ/runPromptCharmonia_midy_pp.sh index 30aa3a865..5a1a503ad 100755 --- a/MC/run/PWGDQ/runPromptCharmonia_midy_pp.sh +++ b/MC/run/PWGDQ/runPromptCharmonia_midy_pp.sh @@ -16,7 +16,7 @@ NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 900 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ -confKey "GeneratorExternal.fileName=${O2DPG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorCocktailPromptCharmoniaToElectronEvtGen_pp13TeV.C;GeneratorExternal.funcName=GeneratorCocktailPromptCharmoniaToElectronEvtGen_pp13TeV()" \ - -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} + -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} -interactionRate 500000 # run workflow diff --git a/MC/run/PWGDQ/runPromptCharmonia_midy_pp_triggerGap.sh b/MC/run/PWGDQ/runPromptCharmonia_midy_pp_triggerGap.sh index d4ad2179a..401d6ced2 100755 --- a/MC/run/PWGDQ/runPromptCharmonia_midy_pp_triggerGap.sh +++ b/MC/run/PWGDQ/runPromptCharmonia_midy_pp_triggerGap.sh @@ -14,7 +14,7 @@ NBKGEVENTS=${NBKGEVENTS:-1} NWORKERS=${NWORKERS:-8} NTIMEFRAMES=${NTIMEFRAMES:-1} -${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 13600 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ +${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 13600 -interactionRate 500000 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/Generator_InjectedPromptCharmoniaMidy_TriggerGap.ini # run workflow diff --git a/MC/run/PWGDQ/runPromptJpsi_fwdy_pp_assessment_triggerGap.sh b/MC/run/PWGDQ/runPromptJpsi_fwdy_pp_assessment_triggerGap.sh index 5273b9398..a784b4737 100755 --- a/MC/run/PWGDQ/runPromptJpsi_fwdy_pp_assessment_triggerGap.sh +++ b/MC/run/PWGDQ/runPromptJpsi_fwdy_pp_assessment_triggerGap.sh @@ -15,7 +15,7 @@ NWORKERS=${NWORKERS:-8} NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 13600 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ - -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/Generator_InjectedPromptJpsiFwdy_TriggerGap.ini --mft-assessment-full --fwdmatching-assessment-full + -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/Generator_InjectedPromptJpsiFwdy_TriggerGap.ini --mft-assessment-full --fwdmatching-assessment-full -interactionRate 500000 # run workflow ${O2DPG_ROOT}/MC/bin/o2_dpg_workflow_runner.py -f workflow.json diff --git a/MC/run/PWGDQ/runPromptJpsi_midy_pp.sh b/MC/run/PWGDQ/runPromptJpsi_midy_pp.sh index eed212575..2398cd456 100755 --- a/MC/run/PWGDQ/runPromptJpsi_midy_pp.sh +++ b/MC/run/PWGDQ/runPromptJpsi_midy_pp.sh @@ -16,7 +16,7 @@ NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 900 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ -confKey "GeneratorExternal.fileName=${O2DPG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorParamPromptJpsiToElectronEvtGen_pp13TeV.C;GeneratorExternal.funcName=GeneratorParamPromptJpsiToElectronEvtGen_pp13TeV()" \ - -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} + -genBkg pythia8 -procBkg inel -colBkg pp --embedding -nb ${NBKGEVENTS} -interactionRate 500000 # run workflow ${O2DPG_ROOT}/MC/bin/o2_dpg_workflow_runner.py -f workflow.json diff --git a/MC/run/PWGDQ/runPromptPsi2S_fwdy_pp_assessment_triggerGap.sh b/MC/run/PWGDQ/runPromptPsi2S_fwdy_pp_assessment_triggerGap.sh index 26f6cc930..5baa5f223 100755 --- a/MC/run/PWGDQ/runPromptPsi2S_fwdy_pp_assessment_triggerGap.sh +++ b/MC/run/PWGDQ/runPromptPsi2S_fwdy_pp_assessment_triggerGap.sh @@ -15,7 +15,7 @@ NWORKERS=${NWORKERS:-8} NTIMEFRAMES=${NTIMEFRAMES:-1} ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 13600 -gen external -j ${NWORKERS} -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -e TGeant4 -mod "--skipModules ZDC" \ - -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/Generator_InjectedPromptPsi2SFwdy_TriggerGap.ini --mft-assessment-full --fwdmatching-assessment-full + -ini $O2DPG_ROOT/MC/config/PWGDQ/ini/Generator_InjectedPromptPsi2SFwdy_TriggerGap.ini --mft-assessment-full --fwdmatching-assessment-full -interactionRate 500000 # run workflow ${O2DPG_ROOT}/MC/bin/o2_dpg_workflow_runner.py -f workflow.json diff --git a/MC/run/PWGGAJE/run_jets.sh b/MC/run/PWGGAJE/run_jets.sh index 314450710..a15530f3b 100755 --- a/MC/run/PWGGAJE/run_jets.sh +++ b/MC/run/PWGGAJE/run_jets.sh @@ -43,7 +43,7 @@ ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM ${CONFIG_ENERGY} -col pp -gen py -ptHatMin ${PTHATMIN} -ptHatMax ${PTHATMAX} \ -tf ${NTIMEFRAMES} -ns ${NSIGEVENTS} -e ${SIMENGINE} \ -j ${NWORKERS} -mod "--skipModules ZDC" \ - -weightPow ${WEIGHTPOW} + -weightPow ${WEIGHTPOW} -interactionRate 500000 # run workflow ${O2DPG_ROOT}/MC/bin/o2_dpg_workflow_runner.py -f workflow.json diff --git a/MC/run/PWGGAJE/run_jets_embedding.sh b/MC/run/PWGGAJE/run_jets_embedding.sh index 88d69fd32..103f986b2 100755 --- a/MC/run/PWGGAJE/run_jets_embedding.sh +++ b/MC/run/PWGGAJE/run_jets_embedding.sh @@ -48,6 +48,6 @@ ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM ${CONFIG_ENERGY} \ -ptHatMin ${PTHATMIN} -ptHatMax ${PTHATMAX} \ -tf ${NTIMEFRAMES} -ns ${NSIGEVENTS} -e ${SIMENGINE} \ -j ${NWORKERS} -mod "--skipModules ZDC" \ - -weightPow ${WEIGHTPOW} + -weightPow ${WEIGHTPOW} -interactionRate 500000 # run workflow ${O2DPG_ROOT}/MC/bin/o2_dpg_workflow_runner.py -f workflow.json diff --git a/MC/run/PWGHF/analysis_benchmark.sh b/MC/run/PWGHF/analysis_benchmark.sh index 5d26d3d35..9f57cdb57 100755 --- a/MC/run/PWGHF/analysis_benchmark.sh +++ b/MC/run/PWGHF/analysis_benchmark.sh @@ -25,16 +25,13 @@ PYPROCESS=${PYPROCESS:-ccbar} #ccbar, bbar, ... # create simulation workflow ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 5020 -col pp -gen pythia8 -proc ${PYPROCESS} -tf ${NTIMEFRAMES} -nb ${NBKGEVENTS} \ - -ns ${NSIGEVENTS} -e ${SIMENGINE} \ - -j ${NWORKERS} --embedding -interactionRate 50000 + -ns ${NSIGEVENTS} -e ${SIMENGINE} -interactionRate 500000 \ + -j ${NWORKERS} -genBkg pythia8 --embedding # Simulating a user who extends this workflow by an analysis task -${O2DPG_ROOT}/MC/bin/o2dpg-workflow-tools.py create workflow_ana --workflow --add-task mchist -needs="" -for i in $(seq 1 $NTIMEFRAMES) -do - needs+="aodmerge_$i " -done +${O2DPG_ROOT}/MC/bin/o2dpg-workflow-tools.py create workflow_ana --add-task mchist +needs="aodmerge" + # Comments: # 1. The output AOD name is the one created by the simulation workflow # 2. base name of AOD merge tasks is known as well to be aodmerge diff --git a/MC/run/PWGLF/run_SyntheticFlowXi_PbPb.sh b/MC/run/PWGLF/run_SyntheticFlowXi_PbPb.sh new file mode 100755 index 000000000..088d145af --- /dev/null +++ b/MC/run/PWGLF/run_SyntheticFlowXi_PbPb.sh @@ -0,0 +1,44 @@ +#!/bin/bash + +# +# A example workflow MC->RECO->AOD for a simple pp min bias +# production, targetting test beam conditions. + +# make sure O2DPG + O2 is loaded +[ ! "${O2DPG_ROOT}" ] && echo "Error: This needs O2DPG loaded" && exit 1 +[ ! "${O2_ROOT}" ] && echo "Error: This needs O2 loaded" && exit 1 + +# ----------- CONFIGURE -------------------------- +export IGNORE_VALIDITYCHECK_OF_CCDB_LOCALCACHE=1 +#export ALICEO2_CCDB_LOCALCACHE=.ccdb + +# ----------- LOAD UTILITY FUNCTIONS -------------------------- +. ${O2_ROOT}/share/scripts/jobutils.sh + +# ----------- START ACTUAL JOB ----------------------------- + +NWORKERS=${NWORKERS:-8} +SIMENGINE=${SIMENGINE:-TGeant4} +NSIGEVENTS=${NSIGEVENTS:-1} +NTIMEFRAMES=${NTIMEFRAMES:-1} +INTRATE=${INTRATE:-50000} +SYSTEM=${SYSTEM:-PbPb} +ENERGY=${ENERGY:-5360} +CFGINIFILE=${CFGINIFILE:-"${O2DPG_ROOT}/MC/config/PWGLF/ini/GeneratorLF_SyntheFlowXi.ini"} +[[ ${SPLITID} != "" ]] && SEED="-seed ${SPLITID}" || SEED="" + +echo "NWORKERS = $NWORKERS" + +# create workflow +O2_SIM_WORKFLOW=${O2_SIM_WORKFLOW:-"${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py"} +$O2_SIM_WORKFLOW -eCM ${ENERGY} -col ${SYSTEM} -gen external \ + -j ${NWORKERS} \ + -ns ${NSIGEVENTS} -tf ${NTIMEFRAMES} -interactionRate ${INTRATE} \ + -confKey "Diamond.width[2]=6." \ + ${SEED} \ + -e ${SIMENGINE} \ + -ini $CFGINIFILE + +# run workflow +O2_SIM_WORKFLOW_RUNNER=${O2_SIM_WORKFLOW_RUNNER:-"${O2DPG_ROOT}/MC/bin/o2_dpg_workflow_runner.py"} +$O2_SIM_WORKFLOW_RUNNER -f workflow.json -tt aod --cpu-limit $NWORKERS diff --git a/MC/run/PWGUD/runSTARlightANCHOR.sh b/MC/run/PWGUD/runSTARlightANCHOR.sh index cfabe1c44..43a7b32b8 100755 --- a/MC/run/PWGUD/runSTARlightANCHOR.sh +++ b/MC/run/PWGUD/runSTARlightANCHOR.sh @@ -24,4 +24,4 @@ export ALIEN_PROC_ID=2963436952 export ALIEN_JDL_ANCHOR_SIM_OPTIONS="-gen external -ini ${PWD}/GenStarlight.ini --embedding -nb ${NBKGEVENTS} -colBkg PbPb -genBkg pythia8 -procBkg heavy_ion" ${O2DPG_ROOT}/MC/config/PWGUD/ini/makeStarlightConfig.py --process kCohPsi2sToMuPi --collType PbPb --eCM 5360 --rapidity cent_eta -${O2DPG_ROOT}/MC/run/ANCHOR/anchorMC.sh \ No newline at end of file +${O2DPG_ROOT}/MC/run/ANCHOR/anchorMC.sh diff --git a/MC/run/examples/event_pool.sh b/MC/run/examples/event_pool.sh new file mode 100644 index 000000000..04da7faa7 --- /dev/null +++ b/MC/run/examples/event_pool.sh @@ -0,0 +1,55 @@ +#!/usr/bin/env bash + +# Example on how to produce an event pool and how to feed it +# to the O2DPG simulation workflow + +# make sure O2DPG + O2 is loaded +[ ! "${O2DPG_ROOT}" ] && echo "Error: This needs O2DPG loaded" && exit 1 +[ ! "${O2_ROOT}" ] && echo "Error: This needs O2 loaded" && exit 1 + +# Parse arguments +MAKE=false +INPUT="" + +help() { + echo "Usage: $0 [--make] [-i|--input ]" + echo " --make: Create the event pool" + echo " -i|--input: Input event pool file to be used in the simulation workflow. Alien paths are supported." + echo " A full path must be provided (use of environment variables allowed), otherwise generation will fail." + echo " -h|--help: Display this help message" + exit 0 +} + +while [[ "$#" -gt 0 ]]; do + case $1 in + --make) MAKE=true ;; + -i|--input) INPUT="$2"; shift ;; + -h|--help) help ;; + *) echo "Unknown operation requested: $1"; help ;; + esac + shift +done + +if $MAKE; then + echo "Started generation of event pool" + # Workflow creation. All the parameters are used as examples + # No transport will be executed. The workflow will stop at the event generation and will conclude with the merging of all the + # kinematic root files of the timeframes in a file called evtpool.root in the current working directory + ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 14000 -col pp -gen pythia8 -proc cdiff -tf 2 -ns 5000 --make-evtpool -seed 546 -interactionRate 500000 -productionTag "evtpoolcreation" -o evtpool + # Workflow runner + ${O2DPG_ROOT}/MC/bin/o2dpg_workflow_runner.py -f evtpool.json -tt pool +elif [[ -n "$INPUT" ]]; then + echo "Input file provided: $INPUT" + if [[ -f "$INPUT" && -s "$INPUT" ]] || [[ "$INPUT" == alien://* ]]; then + # Workflow creation. Phi Rotation is set manually, while the event randomisation of the pool is set by default + ${O2DPG_ROOT}/MC/bin/o2dpg_sim_workflow.py -eCM 14000 -confKey "GeneratorFromO2Kine.randomphi=true;GeneratorFromO2Kine.fileName=$INPUT" -gen extkinO2 -tf 2 -ns 10 -e TGeant4 -j 4 -interactionRate 500000 -seed 546 -productionTag "evtpooltest" + # Workflow runner. The rerun option is set in case you will run directly the script in the same folder (no need to manually delete files) + ${O2DPG_ROOT}/MC/bin/o2dpg_workflow_runner.py -f workflow.json -tt aod --rerun-from grpcreate + else + echo "Error: File does not exist or is empty: $INPUT" + exit 1 + fi +else + echo "Usage: $0 [--make] [-i|--input ]" + exit 1 +fi \ No newline at end of file diff --git a/UTILS/checkCorruptedAO2Ds.C b/UTILS/checkCorruptedAO2Ds.C new file mode 100644 index 000000000..85cbf0098 --- /dev/null +++ b/UTILS/checkCorruptedAO2Ds.C @@ -0,0 +1,52 @@ +#include +#include +#include +#include +#include +#include + +int checkCorruptedAO2Ds(TString infileName = "/alice/sim/2024/LHC24h2/535545/AOD/005/AO2D.root", bool fromAlien = true) { + + if (fromAlien) { + TGrid::Connect("alien://"); + if (!infileName.Contains("alien://")) { + infileName = "alien://" + infileName; + } + } + + auto inFile = TFile::Open(infileName.Data()); + if (!inFile || inFile->IsZombie()) { + return -1; + } + + // all VLA branches in the AO2Ds.root + std::map> branchesToCheck = { + {"O2mcparticle_001", std::vector{"fIndexArray_Mothers"}}, + {"O2ft0", std::vector{"fAmplitudeA", "fChannelA", "fAmplitudeC", "fChannelC"}}, + {"O2fv0a", std::vector{"fAmplitude", "fChannel"}}, + {"O2mccalolabel_001", std::vector{"fIndexArrayMcParticles", "fAmplitudeA"}}, + {"O2zdc_001", std::vector{"fEnergy", "fChannelE", "fAmplitude", "fTime", "fChannelT"}} + }; + + for (auto const& dirKey : *inFile->GetListOfKeys()) { + if (TString(dirKey->GetName()).Contains("DF")) { + auto df = static_cast(inFile->Get(dirKey->GetName())); + std::cout << dirKey->GetName() << std::endl; + for (auto const& pair : branchesToCheck) { + auto tree = static_cast(df->Get(pair.first.data())); + for (auto const& branchName : pair.second) { + auto leaf = static_cast(tree->GetLeaf(branchName.data())); + + for (int iEntry{0}; iEntryGetEntries(); ++iEntry) { + if (tree->GetEntry(iEntry) < 0) { + std::cout << "Found corrupted file! DF: " << dirKey->GetName() << " Tree:" << pair.first.data() << " Branch:" << branchName.data() << std::endl; + return -1; + } + } + } + } + } + } + + return 0; +} \ No newline at end of file diff --git a/UTILS/findCorruptedAO2Ds.sh b/UTILS/findCorruptedAO2Ds.sh new file mode 100755 index 000000000..fd4ebfccb --- /dev/null +++ b/UTILS/findCorruptedAO2Ds.sh @@ -0,0 +1,44 @@ +#!/bin/bash + +# Simple script to find corrupted AO2Ds using the checkCorruptedAO2Ds.C macro + +PRODUCTION=LHC24h2 +RUN=* # use * for all runs +NJOBS=90 + +OUTPUTFILE=corrupted_files_$PRODUCTION.txt +if [ -e "$OUTPUTFILE" ]; then + rm $OUTPUTFILE +fi + +# find all files in alien +if [ "$variable" == "*" ]; then + alien_find alien:///alice/sim/2024/${PRODUCTION} 5*/AOD/*/AO2D.root > files_to_check.txt +else + alien_find alien:///alice/sim/2024/${PRODUCTION} ${RUN}/AOD/*/AO2D.root > files_to_check.txt +fi +mapfile -t FILESTOCHECK < files_to_check.txt + +# process AO2Ds +process_file() { + IFS='/' read -a num <<< "$1" + INPUT=$1 + echo '.x checkCorruptedAO2Ds.C("'${INPUT}'", true)' | root -l -b > log_${num[5]}_${num[7]} + echo '.q' +} +export -f process_file + +parallel -j $NJOBS process_file ::: "${FILESTOCHECK[@]}" + +# create list of corrupted files +touch $OUTPUTFILE +ERRORSTR="Found corrupted file!" +for FILE in "${FILESTOCHECK[@]}"; do + IFS='/' read -a num <<< "$FILE" + if grep -q "$ERRORSTR" log_${num[5]}_${num[7]}; then + echo $FILE >> $OUTPUTFILE + fi +done + +rm files_to_check.txt +rm log* \ No newline at end of file diff --git a/UTILS/get_cherrypick_commit_list.sh b/UTILS/get_cherrypick_commit_list.sh new file mode 100755 index 000000000..1f452b268 --- /dev/null +++ b/UTILS/get_cherrypick_commit_list.sh @@ -0,0 +1,134 @@ +#!/bin/bash + +# given a commit cp_commit on branch source_branch as well +# as another branch target_branch, this script finds the list +# of all commits that are needed to succesfully cherry-pick cp_commit +# onto the target branch + +cp_commit=$1 # Commit you want to cherry-pick +source_branch=$2 # Branch containing commit cp_commit (e.g., master) +target_branch=$3 # Branch you want to cherry-pick onto (e.g., foo) + +# Check if all required arguments are provided +if [ -z "$cp_commit" ] || [ -z "$source_branch" ] || [ -z "$target_branch" ]; then + echo "Usage: $0 " + exit 1 +fi + +# Function to check if two git commits modify at least one common file +modifies_common_files() { + if [ "$#" -ne 2 ]; then + echo "Usage: check_common_files " + return 1 + fi + + local commit1="$1" + local commit2="$2" + + # Get the list of modified files for each commit + local files_commit1 + files_commit1=$(git diff --name-only "${commit1}^" "${commit1}") + + local files_commit2 + files_commit2=$(git diff --name-only "${commit2}^" "${commit2}") + + # Check for common files + local common_files + common_files=$(echo -e "${files_commit1}\n${files_commit2}" | sort | uniq -d) + + # Output result + if [ -n "$common_files" ]; then + return 1 + fi + return 0 +} + +# function to check if 2 commits can be swapped. This can determine if a commit needs +# to come stricly before another commit. +can_swap_commits() { + local commitA="$1" + shift + local commitB=("$@") # this is the list of commits that should swap (as a whole) with commitA + + reverseCommitList=() + + # Loop through the original array in reverse order + for ((i=${#commitB[@]}-1; i>=0; i--)); do + reverseCommitList+=("${commitB[i]}") + done + + # Create a temporary branch for testing + local temp_branch="temp_swap_test_branch" + + # record current state + GIT_CUR=$(git branch --show-current) + # Create a new temporary branch from the current HEAD + git checkout ${commitA}^ -b "$temp_branch" &> /dev/null + + RC=1 + for commit in "${reverseCommitList[@]}"; do + # Cherry-pick commit B onto a branch without commitA + if git cherry-pick "$commit" &> /dev/null; then + # RC=1 # Commits can be swapped + RC_local=1 + else + RC=0 # Cannot swap due to conflict when cherry-picking B + git cherry-pick --abort + fi + done + + # Cleanup: Reset to the original branch and delete the temp branch + git checkout ${GIT_CUR} &> /dev/null + git branch -D "$temp_branch" &>/dev/null + return ${RC} +} + +# Step 1: Identify branch-break off point +BRANCHPOINT=$(git merge-base "$source_branch" "$target_branch") + + +COMMITLIST=() +# Collect the initial set of commits to consider using a while loop +while IFS= read -r line; do + COMMITLIST+=("$line") +done < <(git log ${cp_commit}^...${BRANCHPOINT} --pretty=format:"%H") + +# filter out commits not touching the same files +FILTERED_COMMITS1=() +for commit_hash in "${COMMITLIST[@]}"; do + modifies_common_files ${commit_hash} ${cp_commit} + RC=$? + if [ ${RC} -eq 1 ]; then + FILTERED_COMMITS1+=(${commit_hash}) + fi +done + +# Next, filter out commits which are irrelevant for ${cp_commit} +CP_COMMIT_LIST=(${cp_commit}) # The list of CP=cherry_pick commits to keep/construct + +for commit_hash in "${FILTERED_COMMITS1[@]}"; do + if [ ! "${commit_hash}" == "${cp_commit}" ]; then + can_swap_commits "${commit_hash}" "${CP_COMMIT_LIST[@]}" + if [ $? -eq 0 ]; then + # echo "COMMIT ${commit_hash} is needed" + # in this case we need to record it to the list of relevant commits + # and also trace it's dependencies in turn + CP_COMMIT_LIST+=(${commit_hash}) + fi + fi +done + +# reverse the final list to have correct cherry-pick order + +CP_COMMITS_REVERSED=() +for ((i=${#CP_COMMIT_LIST[@]}-1; i>=0; i--)); do + CP_COMMITS_REVERSED+=("${CP_COMMIT_LIST[i]}") +done + +# List the commits +echo "To cherry-pick ${cp_commit} onto branch ${target_branch}, we need to apply:" +for ((i=0;i<${#CP_COMMITS_REVERSED[@]}; i++)); do + echo "${i}: ${CP_COMMITS_REVERSED[i]}" +done + +exit 0 diff --git a/UTILS/root_merger.py b/UTILS/root_merger.py new file mode 100755 index 000000000..0f673e2f2 --- /dev/null +++ b/UTILS/root_merger.py @@ -0,0 +1,41 @@ +#!/usr/bin/env python3 + +# Simple ROOT files merger + +from ROOT import TFile, TFileMerger +import sys +import os +import argparse + +output_file = '' +input_files = [] +# defining command line options + +parser = argparse.ArgumentParser(description='Simple ROOT files merger', + formatter_class=argparse.ArgumentDefaultsHelpFormatter) + +parser.add_argument('-o','--output', help='Output ROOT filename', required=True) +parser.add_argument('-i','--input', help='Input ROOT files to be merged, separated by a comma', required=True) + +args = parser.parse_args() + +output_file = args.output +input_files = args.input.split(',') + +merger = TFileMerger(False) +merger.OutputFile(output_file) + +for input_file in input_files: + if os.path.exists(input_file): + merger.AddFile(input_file) + else: + print(f"Fatal: {input_file} does not exist.") + sys.exit(1) + +if not merger.Merge(): + print("Error: Merging failed.") + sys.exit(2) +else: + print(f"Successfully merged files into {output_file}") + +sys.exit(0) \ No newline at end of file