diff --git a/.github/workflows/mirror.yaml b/.github/workflows/mirror.yaml new file mode 100644 index 00000000..eed20c12 --- /dev/null +++ b/.github/workflows/mirror.yaml @@ -0,0 +1,41 @@ +name: Mirror and Trigger EICweb + +on: + push: + workflow_dispatch: + +concurrency: + group: mirror + cancel-in-progress: false + +jobs: + build: + name: Mirror and Trigger EICweb + runs-on: ubuntu-latest + permissions: + actions: write + contents: read + steps: + - name: Checkout + uses: actions/checkout@v4 + with: + fetch-depth: 0 + - name: Push to EICweb + uses: eic/gitlab-sync@master + with: + url: https://eicweb.phy.anl.gov/EIC/benchmarks/detector_benchmarks.git/ + token: ${{ secrets.GITLAB_TOKEN }} + username: ${{ secrets.GITLAB_USERNAME }} + ciskip: true + - name: Trigger EICweb + uses: eic/trigger-gitlab-ci@v2 + with: + url: https://eicweb.phy.anl.gov + project_id: 399 + token: ${{ secrets.EICWEB_DETECTOR_BENCHMARKS_TRIGGER }} + ref_name: ${{ github.event.pull_request.head.ref || github.ref }} + variables: + GITHUB_REPOSITORY=${{ github.repository }} + GITHUB_SHA=${{ github.event.pull_request.head.sha || github.sha }} + GITHUB_PR=${{ github.event.pull_request.number }} + PIPELINE_NAME=${{ github.event.pull_request.title }} diff --git a/.gitignore b/.gitignore index fc58c11b..bf34adef 100644 --- a/.gitignore +++ b/.gitignore @@ -40,3 +40,7 @@ __pycache__/ calorimeters/test/ *.d *.pcm + +# org2py output +benchmarks/backgrounds/ecal_backwards.py +benchmarks/ecal_gaps/ecal_gaps.py diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index 1347d6ea..a6c80b27 100644 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -1,10 +1,19 @@ image: ${BENCHMARKS_REGISTRY}/${BENCHMARKS_CONTAINER}:${BENCHMARKS_TAG} +variables: + DETECTOR: epic + DETECTOR_CONFIG: epic_craterlake + DETECTOR_REPOSITORYURL: 'https://github.com/eic/epic.git' + GITHUB_SHA: '' + GITHUB_REPOSITORY: '' + workflow: + name: '$PIPELINE_NAME' rules: - if: '$CI_PIPELINE_SOURCE == "merge_request_event"' - if: '$CI_PIPELINE_SOURCE == "web"' - if: '$CI_PIPELINE_SOURCE == "webide"' + - if: '$CI_PIPELINE_SOURCE == "trigger"' - if: '$CI_COMMIT_BRANCH == "master"' - if: '$CI_COMMIT_TAG' @@ -21,10 +30,10 @@ default: - .local/share - results - config - - juggler.env + - .env - summary.txt reports: - dotenv: juggler.env + dotenv: .env stages: - status-pending @@ -40,6 +49,7 @@ stages: - status-report .status: + image: curlimages/curl:latest script: - | if [ -n "${GITHUB_SHA}" ] ; then @@ -67,49 +77,17 @@ benchmarks:detector:pending: common:setup: stage: config before_script: + script: - | if [[ "${COMMON_BENCH_VERSION}" == "" ]] ; then export COMMON_BENCH_VERSION="master" fi echo "COMMON_BENCH_VERSION = ${COMMON_BENCH_VERSION}" + echo "COMMON_BENCH_VERSION=${COMMON_BENCH_VERSION}" >> .env git clone -b "${COMMON_BENCH_VERSION}" https://eicweb.phy.anl.gov/EIC/benchmarks/common_bench.git setup - script: - - | - if [[ "x${CI_PIPELINE_SOURCE}" == "xmerge_request_event" || "$CI_COMMIT_BRANCH" == "master" ]]; then - echo "DETECTOR = ${DETECTOR}" - echo "DETECTOR_CONFIG = ${DETECTOR_CONFIG}" - echo "DETECTOR_VERSION = ${DETECTOR_VERSION}" - echo "DETECTOR_REPOSITORYURL = ${DETECTOR_REPOSITORYURL}" - echo "DETECTOR=$DETECTOR" >> juggler.env - echo "DETECTOR_CONFIG=$DETECTOR_CONFIG" >> juggler.env - echo "DETECTOR_VERSION=$DETECTOR_VERSION" >> juggler.env - echo "DETECTOR_REPOSITORYURL=$DETECTOR_REPOSITORYURL" >> juggler.env - echo "DETECTOR_DEPLOY_TOKEN_USERNAME=${DEPLOY_TOKEN_USERNAME}" >> juggler.env - echo "DETECTOR_DEPLOY_TOKEN_PASSWORD=${DEPLOY_TOKEN_PASSWORD}" >> juggler.env - echo "COMMON_BENCH_VERSION=$COMMON_BENCH_VERSION" >> juggler.env - else - if [[ "${DETECTOR}" == "" ]] ; then - export DETECTOR="epic" - fi - if [[ "${DETECTOR_CONFIG}" == "" ]] ; then - export DETECTOR_CONFIG="${DETECTOR}" - fi - if [[ "${DETECTOR_VERSION}" == "" ]] ; then - export DETECTOR_VERSION="main" - fi - if [[ "${DETECTOR_REPOSITORYURL}" == "" ]] ; then - export DETECTOR_REPOSITORYURL="https://github.com/eic/${DETECTOR}" - fi - echo "DETECTOR = ${DETECTOR}" - echo "DETECTOR_CONFIG = ${DETECTOR_CONFIG}" - echo "DETECTOR_VERSION = ${DETECTOR_VERSION}" - echo "DETECTOR_REPOSITORYURL = ${DETECTOR_REPOSITORYURL}" - echo "DETECTOR=$DETECTOR" >> juggler.env - echo "DETECTOR_CONFIG=$DETECTOR_CONFIG" >> juggler.env - echo "DETECTOR_VERSION=$DETECTOR_VERSION" >> juggler.env - echo "DETECTOR_REPOSITORYURL=$DETECTOR_REPOSITORYURL" >> juggler.env - echo "COMMON_BENCH_VERSION=$COMMON_BENCH_VERSION" >> juggler.env - fi + echo "BENCHMARKS_TAG: ${BENCHMARKS_TAG}" + echo "BENCHMARKS_CONTAINER: ${BENCHMARKS_CONTAINER}" + echo "BENCHMARKS_REGISTRY: ${BENCHMARKS_REGISTRY}" - source setup/bin/env.sh && ./setup/bin/install_common.sh @@ -140,6 +118,10 @@ get_data: - ls -lrtha - ln -s "${LOCAL_DATA_PATH}/sim_output" sim_output - ln -s "${LOCAL_DATA_PATH}/datasets/data" data + # snakemake support + - mkdir "${DETECTOR_CONFIG}" + - ln -s "${LOCAL_DATA_PATH}/sim_output" "${DETECTOR_CONFIG}/sim_output" + - ln -s "../results" "${DETECTOR_CONFIG}/results" - ls -lrtha retry: max: 2 @@ -147,27 +129,30 @@ get_data: - runner_system_failure include: +# - local: 'benchmarks/backgrounds/config.yml' + - local: 'benchmarks/ecal_gaps/config.yml' - local: 'benchmarks/tracking_detectors/config.yml' - local: 'benchmarks/barrel_ecal/config.yml' - local: 'benchmarks/barrel_hcal/config.yml' - local: 'benchmarks/zdc/config.yml' - local: 'benchmarks/material_maps/config.yml' + - local: 'benchmarks/material_scan/config.yml' - local: 'benchmarks/pid/config.yml' - local: 'benchmarks/timing/config.yml' - local: 'benchmarks/b0_tracker/config.yml' - local: 'benchmarks/LOWQ2/config.yml' - - local: 'benchmarks/others/config.yml' deploy_results: stage: deploy needs: - - ["collect_results:zdc","collect_results:barrel_ecal","collect_results:barrel_hcal","collect_results:materialscan"] + - ["collect_results:zdc","collect_results:barrel_ecal","collect_results:barrel_hcal","collect_results:material_scan"] script: - echo "deploy results!" - find results -print | sort | tee summary.txt benchmarks:detector:success: stage: status-report + dependencies: [] extends: .status variables: STATE: "success" @@ -176,6 +161,7 @@ benchmarks:detector:success: benchmarks:detector:failure: stage: status-report + dependencies: [] extends: .status variables: STATE: "failure" @@ -194,6 +180,9 @@ benchmarks:reconstruction: DETECTOR_DEPLOY_TOKEN_USERNAME: "${DETECTOR_DEPLOY_TOKEN_USERNAME}" DETECTOR_DEPLOY_TOKEN_PASSWORD: "${DETECTOR_DEPLOY_TOKEN_PASSWORD}" COMMON_BENCH_VERSION: "$COMMON_BENCH_VERSION" + BENCHMARKS_TAG: "${BENCHMARKS_TAG}" + BENCHMARKS_CONTAINER: "${BENCHMARKS_CONTAINER}" + BENCHMARKS_REGISTRY: "${BENCHMARKS_REGISTRY}" trigger: project: EIC/benchmarks/reconstruction_benchmarks strategy: depend @@ -211,6 +200,9 @@ benchmarks:physics: DETECTOR_DEPLOY_TOKEN_USERNAME: "${DETECTOR_DEPLOY_TOKEN_USERNAME}" DETECTOR_DEPLOY_TOKEN_PASSWORD: "${DETECTOR_DEPLOY_TOKEN_PASSWORD}" COMMON_BENCH_VERSION: "$COMMON_BENCH_VERSION" + BENCHMARKS_TAG: "${BENCHMARKS_TAG}" + BENCHMARKS_CONTAINER: "${BENCHMARKS_CONTAINER}" + BENCHMARKS_REGISTRY: "${BENCHMARKS_REGISTRY}" trigger: project: EIC/benchmarks/physics_benchmarks strategy: depend diff --git a/.rootlogon.C b/.rootlogon.C index 5ce62268..b63b6b05 100644 --- a/.rootlogon.C +++ b/.rootlogon.C @@ -3,6 +3,7 @@ gROOT->ProcessLine(".include include"); R__LOAD_LIBRARY(fmt) + R__LOAD_LIBRARY(HepMC3) // Setting for Graphs gROOT->SetStyle("Plain"); diff --git a/Snakefile b/Snakefile new file mode 100644 index 00000000..a543b210 --- /dev/null +++ b/Snakefile @@ -0,0 +1,59 @@ +configfile: "config.yaml" + +if config["remote"] == "S3": + from snakemake.remote.S3 import RemoteProvider as S3RemoteProvider + provider = S3RemoteProvider( + endpoint_url="https://eics3.sdcc.bnl.gov:9000", + access_key_id=os.environ["S3_ACCESS_KEY"], + secret_access_key=os.environ["S3_SECRET_KEY"], + ) + remote_path = lambda path: f"eictest/{path}" +elif config["remote"] == "XRootD": + from snakemake.remote.XRootD import RemoteProvider as XRootDRemoteProvider + provider = XRootDRemoteProvider( + stay_on_remote=False, + ) + remote_path = lambda path: f"root://dtn-eic.jlab.org//work/eic2/{path}" +else: + raise ValueError(f"Unexpected config[\"remote\"] = {config['remote']}") + +include: "benchmarks/backgrounds/Snakefile" +include: "benchmarks/barrel_ecal/Snakefile" +include: "benchmarks/ecal_gaps/Snakefile" +include: "benchmarks/material_scan/Snakefile" + + +rule warmup_run: + output: + "warmup/{DETECTOR_CONFIG}.edm4hep.root", + message: "Ensuring that calibrations/fieldmaps are available for {wildcards.DETECTOR_CONFIG}" + shell: """ +ddsim \ + --runType batch \ + --numberOfEvents 1 \ + --compactFile "$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml" \ + --outputFile "{output}" \ + --enableGun +""" + + +rule matplotlibrc: + output: + ".matplotlibrc", + run: + with open(output[0], "wt") as fp: + fp.write("backend: Agg\n") + # interactive mode prevents plt.show() from blocking + fp.write("interactive : True\n") + + +rule org2py: + input: + notebook=workflow.basedir + "/{NOTEBOOK}.org", + converter=workflow.source_path("benchmarks/common/org2py.awk"), + output: + "{NOTEBOOK}.py" + shell: + """ +awk -f {input.converter} {input.notebook} > {output} +""" diff --git a/benchmarks/backgrounds/Snakefile b/benchmarks/backgrounds/Snakefile new file mode 100644 index 00000000..300e38b0 --- /dev/null +++ b/benchmarks/backgrounds/Snakefile @@ -0,0 +1,91 @@ +import os +import shutil + + +rule backgrounds_get_beam_gas_electron: + input: + provider.remote(remote_path("EPIC/EVGEN/BACKGROUNDS/BEAMGAS/electron/beam_gas_ep_10GeV_foam_emin10keV_10Mevt_vtx.hepmc")), + output: + "input/backgrounds/beam_gas_electron.hepmc", + run: + shutil.move(input[0], output[0]) + + +rule backgrounds_get_beam_gas_proton: + input: + provider.remote(remote_path("EPIC/EVGEN/BACKGROUNDS/BEAMGAS/proton/ProtonBeamGasEvents/100GeV/100GeV_1.hepmc")), + output: + "input/backgrounds/beam_gas_proton.hepmc", + run: + shutil.move(input[0], output[0]) + + +rule backgrounds_get_DIS: + input: + provider.remote(remote_path("EPIC/EVGEN/DIS/NC/{BEAM}/minQ2={MINQ2}/pythia8NCDIS_{BEAM}_minQ2={MINQ2}_{SUFFIX}.hepmc")), + wildcard_constraints: + BEAM="\d+x\d+", + MINQ2="\d+", + output: + "input/backgrounds/pythia8NCDIS_{BEAM}_minQ2={MINQ2}_{SUFFIX}.hepmc", + run: + shutil.move(input[0], output[0]) + + +rule backgrounds_sim: + input: + hepmc="input/backgrounds/{NAME}.hepmc", + warmup="warmup/{DETECTOR_CONFIG}.edm4hep.root", + output: + "sim_output/{DETECTOR_CONFIG}/backgrounds/{NAME}.edm4hep.root", + log: + "sim_output/{DETECTOR_CONFIG}/backgrounds/{NAME}.edm4hep.root.log", + params: + N_EVENTS=100 + shell: + """ +ddsim \ + --runType batch \ + --part.minimalKineticEnergy 100*GeV \ + --filter.tracker edep0 \ + -v WARNING \ + --numberOfEvents {params.N_EVENTS} \ + --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ + --inputFiles {input.hepmc} \ + --outputFile {output} +""" + + +DETECTOR_CONFIG=os.environ["DETECTOR_CONFIG"] + +rule backgrounds_ecal_backwards: + input: + matplotlibrc=".matplotlibrc", + script="benchmarks/backgrounds/ecal_backwards.py", + electron_beam_gas_gen="input/backgrounds/beam_gas_electron.hepmc", + electron_beam_gas_sim="sim_output/" + DETECTOR_CONFIG + "/backgrounds/beam_gas_electron.edm4hep.root", + physics_signal_sim="sim_output/" + DETECTOR_CONFIG + "/backgrounds/pythia8NCDIS_10x100_minQ2=1_beamEffects_xAngle=-0.025_hiDiv_1.edm4hep.root", + proton_beam_gas_gen="input/backgrounds/beam_gas_proton.hepmc", + proton_beam_gas_sim="sim_output/" + DETECTOR_CONFIG + "/backgrounds/beam_gas_proton.edm4hep.root", + output: + directory("results/backgrounds/backwards_ecal") + threads: workflow.cores + shell: + """ +PORT=$RANDOM +dask scheduler --port $PORT & +export DASK_SCHEDULER=localhost:$PORT +SCHEDULER_PID=$! +dask worker tcp://$DASK_SCHEDULER --nworkers {threads} --nthreads 1 & +WORKER_PID=$! +env \ +MATPLOTLIBRC={input.matplotlibrc} \ +ELECTRON_BEAM_GAS_GEN=$(realpath {input.electron_beam_gas_gen}) \ +ELECTRON_BEAM_GAS_SIM=$(realpath {input.electron_beam_gas_sim}) \ +PHYSICS_PROCESS_SIM=$(realpath {input.physics_signal_sim}) \ +PROTON_BEAM_GAS_GEN=$(realpath {input.proton_beam_gas_gen}) \ +PROTON_BEAM_GAS_SIM=$(realpath {input.proton_beam_gas_sim}) \ +OUTPUT_DIR={output} \ +python {input.script} +kill $WORKER_PID $SCHEDULER_PID +""" diff --git a/benchmarks/backgrounds/config.yml b/benchmarks/backgrounds/config.yml new file mode 100644 index 00000000..acf45c56 --- /dev/null +++ b/benchmarks/backgrounds/config.yml @@ -0,0 +1,26 @@ +sim:backgrounds: + extends: .det_benchmark + stage: simulate + script: + - mkdir -p $LOCAL_DATA_PATH/input + - ln -s $LOCAL_DATA_PATH/input input + - snakemake --cores 2 sim_output/$DETECTOR_CONFIG/beam_gas_{electron,proton}.edm4hep.root sim_output/$DETECTOR_CONFIG/pythia8NCDIS_10x100_minQ2=1_beamEffects_xAngle=-0.025_hiDiv_vtxfix_1.edm4hep.root + +bench:backgrounds_emcal_backwards: + extends: .det_benchmark + stage: benchmarks + needs: + - ["sim:backgrounds"] + script: + - ln -s $LOCAL_DATA_PATH/input input + - export PYTHONUSERBASE=$LOCAL_DATA_PATH/deps + - pip install -r benchmarks/backgrounds/requirements.txt + - snakemake --cores 8 backgrounds_ecal_backwards + +collect_results:backgrounds: + extends: .det_benchmark + stage: collect + needs: + - "bench:backgrounds_emcal_backwards" + script: + - ls -lrht diff --git a/benchmarks/backgrounds/ecal_backwards.org b/benchmarks/backgrounds/ecal_backwards.org new file mode 100644 index 00000000..a420358e --- /dev/null +++ b/benchmarks/backgrounds/ecal_backwards.org @@ -0,0 +1,447 @@ +#+PROPERTY: header-args:jupyter-python :session /jpy:localhost#8888:eeemcal :async yes :results drawer :exports both + +#+TITLE: ePIC EEEMCal background rates study +#+AUTHOR: Dmitry Kalinkin +#+OPTIONS: d:t + +#+begin_src jupyter-python :results silent +import os +from pathlib import Path + +import awkward as ak +import boost_histogram as bh +import dask +import dask_awkward as dak +import dask_histogram as dh +import numpy as np +import uproot +import pyhepmc +#+end_src + +#+begin_src jupyter-python :results silent +from dask.distributed import Client + +client = Client(os.environ.get("DASK_SCHEDULER", "localhost:8786")) +#+end_src + +* Plotting setup + +#+begin_src jupyter-python :results silent +import matplotlib as mpl +import matplotlib.pyplot as plt + +def setup_presentation_style(): + mpl.rcParams.update(mpl.rcParamsDefault) + plt.style.use('ggplot') + plt.rcParams.update({ + 'axes.labelsize': 8, + 'axes.titlesize': 9, + 'figure.titlesize': 9, + 'figure.figsize': (4, 3), + 'legend.fontsize': 7, + 'legend.loc': 'upper right', + 'xtick.labelsize': 8, + 'ytick.labelsize': 8, + 'pgf.rcfonts': False, + }) + +setup_presentation_style() +#+end_src + +* Analysis + +** Input files + +#+begin_src jupyter-python :results silent +ELECTRON_BEAM_GAS_GEN=os.environ.get("ELECTRON_BEAM_GAS_GEN", "../beam_gas_ep_10GeV_foam_emin10keV_10Mevt_vtx.hepmc") +ELECTRON_BEAM_GAS_SIM=os.environ.get("ELECTRON_BEAM_GAS_SIM", "electron_beam_gas.edm4hep.root") +PHYSICS_PROCESS_SIM=os.environ.get("PHYSICS_PROCESS_SIM", "pythia8NCDIS_10x100_minQ2=1_beamEffects_xAngle=-0.025_hiDiv_1.edm4hep.root") +PROTON_BEAM_GAS_GEN=os.environ.get("PROTON_BEAM_GAS_GEN", "100GeV.hepmc") +PROTON_BEAM_GAS_SIM=os.environ.get("PROTON_BEAM_GAS_SIM", "proton+beam_gas_ep.edm4hep.root") + +output_dir=Path(os.environ.get("OUTPUT_DIR", "./")) +output_dir.mkdir(parents=True, exist_ok=True) +#+end_src + +#+begin_src jupyter-python :results silent +builder = ak.ArrayBuilder() +n = 0 +with pyhepmc.open(PROTON_BEAM_GAS_GEN) as f: + it = iter(f) + for event in f: + builder.begin_list() + for vertex in event.vertices: + assert event.length_unit.name == "MM" + with builder.record("Vector4D"): + builder.field("x") + builder.real(vertex.position.x) + builder.field("y") + builder.real(vertex.position.y) + builder.field("z") + builder.real(vertex.position.z) + builder.field("t") + builder.real(vertex.position.t) + builder.end_list() + n += 1 + if n > 10000: break + +vertices_proton_beam_gas = builder.snapshot() + +builder = ak.ArrayBuilder() +n = 0 +with pyhepmc.open(ELECTRON_BEAM_GAS_GEN) as f: + it = iter(f) + for event in f: + builder.begin_list() + for vertex in event.vertices: + assert event.length_unit.name == "MM" + with builder.record("Vector3D"): + builder.field("x") + builder.real(vertex.position.x) + builder.field("y") + builder.real(vertex.position.y) + builder.field("z") + builder.real(vertex.position.z) + builder.field("t") + builder.real(vertex.position.t) + builder.end_list() + n += 1 + if n > 10000: break + +vertices_electron_beam_gas = builder.snapshot() +#+end_src + +#+begin_src jupyter-python :results silent +def filter_name(name): + return "Hits." in name or "MCParticles." in name + +datasets = { + # https://wiki.bnl.gov/EPIC/index.php?title=Electron_Beam_Gas + "electron beam gas 10 GeV": { + "vertices": vertices_electron_beam_gas, + "events": uproot.dask({ELECTRON_BEAM_GAS_SIM: "events"}, filter_name=filter_name, open_files=False, steps_per_file=32), + "cmap": "cool", + "rate": 3.2e6, # Hz + }, + "DIS 10x100, $Q^2 > 1$ GeV${}^2$": { + "events": uproot.dask({PHYSICS_PROCESS_SIM: "events"}, filter_name=filter_name, open_files=False, steps_per_file=32), + "rate": 184e3, # Hz + }, + # https://wiki.bnl.gov/EPIC/index.php?title=Hadron_Beam_Gas + "proton beam gas 100 GeV": { + "vertices": vertices_proton_beam_gas, + "events": uproot.dask({PROTON_BEAM_GAS_SIM: "events"}, filter_name=filter_name, open_files=False, steps_per_file=32), + "cmap": "summer", + "rate": 31.9e3, # Hz + }, +} + +for ds in datasets.values(): + ds["events"].eager_compute_divisions() +#+end_src + +** Vertex distributions + +#+begin_src jupyter-python +for label, ds in datasets.items(): + if "vertices" not in ds: continue + vs = ds["vertices"] + weight = ds["rate"] / ak.num(ds["vertices"], axis=0) + plt.hist(ak.ravel(vs.t[:,0]), bins=100, histtype="step", label=label) +plt.minorticks_on() +plt.xlabel("vertex[0].t, mm") +plt.legend() +plt.savefig(output_dir / "vertex_time_distribution.png", bbox_inches="tight") +plt.show() + +for label, ds in datasets.items(): + if "vertices" not in ds: continue + vs = ds["vertices"] + weight = ds["rate"] / ak.num(ds["vertices"], axis=0) + plt.hist(ak.ravel(vs.z[:,0]), bins=100, histtype="step", label=label) +plt.minorticks_on() +plt.xlabel("vertex[0].z, mm") +plt.legend() +plt.savefig(output_dir / "vertex_z_distribution.png", bbox_inches="tight") +plt.show() + +for label, ds in datasets.items(): + if "vertices" not in ds: continue + vs = ds["vertices"] + cmap = ds["cmap"] + weight = ds["rate"] / ak.num(ds["vertices"], axis=0) + plt.hist2d(vs.z[:,0].to_numpy(), vs.x[:,0].to_numpy(), bins=(100, np.linspace(-130, 130, 160)), cmin=1e-30, label=label, cmap=cmap) + plt.plot([], color=mpl.colormaps[cmap](0.5), label=label) +plt.minorticks_on() +plt.xlabel("vertex[0].z, mm") +plt.ylabel("vertex[0].x, mm") +plt.legend() +plt.savefig(output_dir / "vertex_xz_distribution.png", bbox_inches="tight") +plt.show() + +for ix, (label, ds) in enumerate(datasets.items()): + if "vertices" not in ds: continue + vs = ds["vertices"] + cmap = ds["cmap"] + weight = ds["rate"] / ak.num(ds["vertices"], axis=0) + plt.hist2d(vs.z[:,0].to_numpy(), vs.y[:,0].to_numpy(), bins=(100, 100), cmin=1e-30, cmap=cmap) + plt.colorbar() + plt.minorticks_on() + plt.xlabel("vertex[0].z, mm") + plt.ylabel("vertex[0].y, mm") + plt.title(label) + plt.savefig(output_dir / f"vertex_yz_distribution_{ix}.png", bbox_inches="tight") + plt.show() +#+end_src + +** Simulation results + +#+begin_src jupyter-python +for collection_name in ["EcalEndcapNHits", "EcalEndcapPHits"]: + for dataset_ix, (label, ds) in enumerate(datasets.items()): + events = ds["events"] + + energy_sums = ak.sum(events[f"{collection_name}.energy"].head(10000), axis=1) + event_id = ak.argmax(energy_sums) + xs = events[f"{collection_name}.position.x"].head(event_id + 1)[event_id].to_numpy() + ys = events[f"{collection_name}.position.y"].head(event_id + 1)[event_id].to_numpy() + + bin_widths = [None, None] + for ix, vals in enumerate([xs, ys]): + centers = np.unique(vals) + diffs = centers[1:] - centers[:-1] + bin_widths[ix] = np.min(diffs[diffs > 0]) if np.sum(diffs > 0) > 0 else 1. + print(f"bin_widths[{ix}]", bin_widths[ix]) + + bins = { + "EcalEndcapNHits": [np.arange(-750., 750., bin_width) for bin_width in bin_widths], + "EcalEndcapPHits": [np.arange(-1800., 1800., bin_width) for bin_width in bin_widths], + }[collection_name] + + plt.hist2d( + xs, + ys, + weights=events[f"{collection_name}.energy"].head(event_id + 1)[event_id].to_numpy(), + bins=bins, + cmin=1e-10, + ) + plt.colorbar().set_label("energy, GeV", loc="top") + plt.title(f"{label}, event_id={event_id}\n{collection_name}") + plt.xlabel("hit x, mm", loc="right") + plt.ylabel("hit y, mm", loc="top") + plt.savefig(output_dir / f"{collection_name}_event_display_{dataset_ix}.png", bbox_inches="tight") + plt.show() +#+end_src + +** Discovering number of cells + +Using HyperLogLog algorithm would be faster here, or actually load +DD4hep geometry and count sensitive volumes. + +#+begin_src jupyter-python +def unique(array): + if ak.backend(array) == "typetracer": + ak.typetracer.touch_data(array) + return array + return ak.from_numpy(np.unique(ak.to_numpy(ak.ravel(array)))) +unique_delayed = dask.delayed(unique) +len_delayed = dask.delayed(len) + +cellID_for_r = dict() + +for collection_name in ["EcalEndcapNHits", "EcalEndcapPHits"]: + r_axis = { + "EcalEndcapNHits": bh.axis.Regular(75, 0., 750.), + "EcalEndcapPHits": bh.axis.Regular(90, 0., 1800.), + }[collection_name] + ds = datasets["DIS 10x100, $Q^2 > 1$ GeV${}^2$"] + events = ds["events"] + + r = np.hypot( + ak.ravel(events[f"{collection_name}.position.x"]), + ak.ravel(events[f"{collection_name}.position.y"]), + ) + cellID = ak.ravel(events[f"{collection_name}.cellID"]) + + cellID_for_r[collection_name] = np.array(client.gather(client.compute([ + len_delayed(unique_delayed( + cellID[(r >= r_min) & (r < r_max)].map_partitions(unique) + )) + for r_min, r_max in zip(r_axis.edges[:-1], r_axis.edges[1:]) + ]))) + + print(cellID_for_r[collection_name]) + print(sum(cellID_for_r[collection_name])) + + plt.stairs( + cellID_for_r[collection_name], + r_axis.edges, + ) + + plt.title(f"{collection_name}") + plt.legend() + plt.xlabel("r, mm", loc="right") + dr = (r_axis.edges[1] - r_axis.edges[0]) + plt.ylabel(f"Number of towers per {dr} mm slice in $r$", loc="top") + plt.savefig(output_dir / f"{collection_name}_num_towers.png", bbox_inches="tight") + plt.show() +#+end_src + +** Plotting the rates + +#+begin_src jupyter-python +for collection_name in ["EcalEndcapNHits", "EcalEndcapPHits"]: + r_axis = { + "EcalEndcapNHits": bh.axis.Regular(75, 0., 750.), + "EcalEndcapPHits": bh.axis.Regular(90, 0., 1800.), + }[collection_name] + for edep_min in [0.005, 0.015, 0.050]: # GeV + for label, ds in datasets.items(): + events = ds["events"] + weight = ds["rate"] / len(events) + + r = np.hypot( + ak.ravel(events[f"{collection_name}.position.x"]), + ak.ravel(events[f"{collection_name}.position.y"]), + ) + edep = ak.ravel(events[f"{collection_name}.energy"]) + r = r[edep > edep_min] + + hist = dh.factory( + r, + axes=(r_axis,), + ).compute() + plt.stairs( + hist.values() * weight / cellID_for_r[collection_name], + hist.axes[0].edges, + label=f"{label}", + ) + + plt.title(f"for $E_{{dep.}} >$ {edep_min * 1000} MeV\n{collection_name}") + plt.legend() + plt.xlabel("r, mm", loc="right") + plt.ylabel("rate per tower, Hz", loc="top") + plt.yscale("log") + plt.savefig(output_dir / f"{collection_name}_hit_rate_vs_r_edep_min_{edep_min:.3f}.png", bbox_inches="tight") + plt.show() +#+end_src + +#+begin_src jupyter-python +for collection_name in ["EcalEndcapNHits", "EcalEndcapPHits"]: + for totedep_min in [-1, 0, 0.1, 0.5, 1.0, 5.0, 10.]: # GeV + for label, ds in datasets.items(): + events = ds["events"] + weight = ds["rate"] / len(events) + + z = ds["events"]["MCParticles.vertex.z"][:,1] + totedep = ak.sum(events[f"{collection_name}.energy"], axis=1) + z = z[totedep > totedep_min] + + hist = dh.factory( + z, + axes=(bh.axis.Regular(250, -7500., 17500.),), + ).compute() + plt.stairs( + hist.values() * weight, + hist.axes[0].edges, + label=f"{label}", + ) + + plt.title(rf"for events with $E_{{\mathrm{{dep. tot.}}}}$ $>$ {totedep_min} GeV" + f"\n{collection_name}") + plt.legend() + plt.xlabel("$z$ of the first interaction vertex, mm", loc="right") + plt.ylabel("rate, Hz", loc="top") + plt.yscale("log") + plt.savefig(output_dir / f"{collection_name}_hit_rate_vs_z_totedep_min_{totedep_min:.1f}.png", bbox_inches="tight") + plt.show() +#+end_src + +#+begin_src jupyter-python +num_towers_cache = { + "LumiSpecCAL": 200, + "LumiDirectPCAL": 1, + "ZDCHcal": 1470, + "LFHCAL": 578338, + "ZDC_WSi_": 187043, + "EcalBarrelScFi": 124205483, + "EcalEndcapP": 15037, + "ZDCEcal": 400, + "EcalEndcapPInsert": 536, + "HcalEndcapPInsert": 20251, + "B0ECal": 131, + "HcalEndcapN": 13800, + "HcalBarrel": 7680, + "EcalBarrelImaging": 5765469, + "EcalEndcapN": 2988, + "ZDC_PbSi_": 44344, +} + +fig_cmb = plt.figure() +ax_cmb = fig_cmb.gca() + +for edep_min in [0]: # GeV + for dataset_ix, (x_offset, (ds_label, ds)) in enumerate(zip(np.linspace(-0.3, 0.3, len(datasets)), datasets.items())): + events = ds["events"] + weight = ds["rate"] / len(events) + + labels = [] + values = [] + norms = [] + + for branch_name in events.fields: + if ".energy" not in branch_name: continue + if "ZDC_SiliconPix_Hits" in branch_name: continue + + edep = ak.ravel(events[branch_name]) + + #cellID = ak.ravel(events[branch_name.replace(".energy", ".cellID")]) + #num_towers = len(unique_delayed( + # cellID.map_partitions(unique) + #).compute()) + + key = branch_name.replace("Hits.energy", "") + if key not in num_towers_cache: + print(f"The \"{key}\" not in num_towers_cache. Skipping.") + continue + num_towers = num_towers_cache[key] + + labels.append(branch_name.replace("Hits.energy", "")) + values.append(ak.count(edep[edep > edep_min])) + norms.append(num_towers if num_towers != 0 else np.nan) + + fig_cur = plt.figure() + ax_cur = fig_cur.gca() + + values, = dask.compute(values) + for ax, per_tower, offset, width in [ + (ax_cmb, True, x_offset, 2 * 0.3 / (len(datasets) - 1)), + (ax_cur, False, 0, 2 * 0.3), + ]: + ax.bar( + np.arange(len(labels)) + offset, + weight * np.array(values) / (np.array(norms) if per_tower else np.ones_like(norms)), + width=width, + tick_label=labels, + label=ds_label, + color=f"C{dataset_ix}", + ) + + plt.sca(ax_cur) + plt.legend() + plt.title(f"for $E_{{dep.}} >$ {edep_min * 1000} MeV") + plt.ylabel("rate, Hz", loc="top") + plt.yscale("log") + plt.xticks(rotation=90, ha='right') + fig_cur.savefig(f"rates_edep_min_{edep_min}_{dataset_ix}.png", bbox_inches="tight") + plt.show() + plt.close(fig_cur) + + plt.sca(ax_cmb) + plt.legend() + plt.title(f"for $E_{{dep.}} >$ {edep_min * 1000} MeV") + plt.ylabel("rate per tower, Hz", loc="top") + plt.yscale("log") + plt.xticks(rotation=90, ha='right') + fig_cmb.savefig(f"rates_edep_min_{edep_min}.png", bbox_inches="tight") + plt.show() +#+end_src diff --git a/benchmarks/backgrounds/requirements.txt b/benchmarks/backgrounds/requirements.txt new file mode 100644 index 00000000..738152dd --- /dev/null +++ b/benchmarks/backgrounds/requirements.txt @@ -0,0 +1,7 @@ +awkward >= 2.4.0 +dask >= 2023 +dask_awkward +dask_histogram +distributed >= 2023 +pyhepmc +uproot >= 5.2.0 diff --git a/benchmarks/barrel_ecal/Snakefile b/benchmarks/barrel_ecal/Snakefile new file mode 100644 index 00000000..707c2a89 --- /dev/null +++ b/benchmarks/barrel_ecal/Snakefile @@ -0,0 +1,197 @@ +DETECTOR_PATH = os.environ["DETECTOR_PATH"] + + +rule emcal_barrel_particles_gen: + input: + workflow.source_path("scripts/emcal_barrel_common_functions.h"), + script=workflow.source_path("scripts/emcal_barrel_particles_gen.cxx"), + params: + JUGGLER_N_EVENTS = 100, + output: + "data/emcal_barrel_{PARTICLE}_energies{E_MIN}_{E_MAX}.hepmc", + shell: + """ +env \ + JUGGLER_GEN_FILE=../{output} \ + root -l -b -q '{input.script}+({params.JUGGLER_N_EVENTS}, {wildcards.E_MIN}, {wildcards.E_MAX}, \"{wildcards.PARTICLE}\")' +""" + + +rule emcal_barrel_particles: + input: + "data/emcal_barrel_{PARTICLE}_energies{E_MIN}_{E_MAX}.hepmc" + params: + JUGGLER_N_EVENTS = 100, + output: + "{DETECTOR_CONFIG}/sim_output/sim_emcal_barrel_{PARTICLE}_energies{E_MIN}_{E_MAX}.edm4hep.root" + shell: + """ +ddsim \ + --runType batch \ + -v WARNING \ + --part.minimalKineticEnergy 0.5*GeV \ + --filter.tracker edep0 \ + --numberOfEvents {params.JUGGLER_N_EVENTS} \ + --compactFile """ + DETECTOR_PATH + """/{wildcards.DETECTOR_CONFIG}.xml \ + --inputFiles {input} \ + --outputFile {output} +""" + + +# This is needed to bridge snakemake rule files to "normal" benchmarks +rule emcal_barrel_particles_compat_normal: + input: + "{DETECTOR_CONFIG}/sim_output/sim_emcal_barrel_{PARTICLE}_energies5.0_5.0.edm4hep.root", + wildcard_constraints: + PARTICLE="[^_]+", + output: + "{DETECTOR_CONFIG}/sim_output/sim_emcal_barrel_{PARTICLE}.edm4hep.root", + shell: + """ +ln {input} {output} +""" + + +# This is needed to bridge snakemake rule files to "energy_scan" benchmarks +rule emcal_barrel_particles_compat_energy_scan: + input: + "{DETECTOR_CONFIG}/sim_output/sim_emcal_barrel_{PARTICLE}_energies{E}_{E}.edm4hep.root", + wildcard_constraints: + PARTICLE="[^_]+", + output: + "{DETECTOR_CONFIG}/sim_output/energy_scan/{E}/sim_emcal_barrel_{PARTICLE}.edm4hep.root", + shell: + """ +ln {input} {output} +""" + + +rule emcal_barrel_particles_analysis: + input: + workflow.source_path("scripts/emcal_barrel_common_functions.h"), + script=workflow.source_path("scripts/emcal_barrel_particles_analysis.cxx"), + sim="{DETECTOR_CONFIG}/sim_output/sim_emcal_barrel_{PARTICLE}.edm4hep.root", + wildcard_constraints: + PARTICLE="(electron|photon|piplus|piminus)", # avoid clash with "pions" + output: + "{DETECTOR_CONFIG}/results/emcal_barrel_{PARTICLE}_calibration.json", + "{DETECTOR_CONFIG}/results/emcal_barrel_{PARTICLE}_Ethr.png", + "{DETECTOR_CONFIG}/results/emcal_barrel_{PARTICLE}_Ethr.pdf", + "{DETECTOR_CONFIG}/results/emcal_barrel_{PARTICLE}_nhits.png", + "{DETECTOR_CONFIG}/results/emcal_barrel_{PARTICLE}_nhits.pdf", + "{DETECTOR_CONFIG}/results/emcal_barrel_{PARTICLE}_Esim.png", + "{DETECTOR_CONFIG}/results/emcal_barrel_{PARTICLE}_Esim.pdf", + "{DETECTOR_CONFIG}/results/emcal_barrel_{PARTICLE}_fsam.png", + "{DETECTOR_CONFIG}/results/emcal_barrel_{PARTICLE}_fsam.pdf", + "{DETECTOR_CONFIG}/results/emcal_barrel_{PARTICLE}_fsamImg.png", + "{DETECTOR_CONFIG}/results/emcal_barrel_{PARTICLE}_fsamImg.pdf", + "{DETECTOR_CONFIG}/results/emcal_barrel_{PARTICLE}_fsamScFi.png", + "{DETECTOR_CONFIG}/results/emcal_barrel_{PARTICLE}_fsamScFi.pdf", + shell: + """ +cd {wildcards.DETECTOR_CONFIG} +root -l -b -q '{input.script}+("{wildcards.PARTICLE}", true)' +""" + + +rule emcal_barrel_pions_analysis: + input: + script=workflow.source_path("scripts/emcal_barrel_pions_analysis.cxx"), + sim="{DETECTOR_CONFIG}/sim_output/sim_emcal_barrel_piplus_energies5.0_5.0.edm4hep.root", + output: + expand( + "{{DETECTOR_CONFIG}}/results/emcal_barrel_pions_{var_name}.{extension}", + var_name=["Ethr", "nhits", "Esim", "fsam", "pid"], + extension=["pdf", "png"], + ), + shell: + """ +cd {wildcards.DETECTOR_CONFIG} +root -l -b -q '{input.script}+("../{input.sim}")' +""" + + +rule emcal_barrel_pi0_analysis: + input: + script=workflow.source_path("scripts/emcal_barrel_pi0_analysis.cxx"), + sim="{DETECTOR_CONFIG}/sim_output/sim_emcal_barrel_pi0_energies5.0_5.0.edm4hep.root", + fsam="{DETECTOR_CONFIG}/results/emcal_barrel_electron_calibration.json", + output: + expand( + "{{DETECTOR_CONFIG}}/results/emcal_barrel_pi0_{var_name}.{extension}", + var_name=["Ethr", "nhits", "Esim", "dE_rel"], + extension=["pdf", "png"], + ), + "{DETECTOR_CONFIG}/results/Barrel_emcal_pi0.json" + shell: + """ +cd {wildcards.DETECTOR_CONFIG} +root -l -b -q '{input.script}+("../{input.sim}")' +""" + + +ENERGY_SCAN_ENERGIES = [0.5, 1., 2., 5., 10.] +rule emcal_energy_scan: + input: + # Require all simulations produced for this rule + expand("{{DETECTOR_CONFIG}}/sim_output/energy_scan/{energy}/sim_emcal_barrel_{{PARTICLE}}.edm4hep.root", energy=ENERGY_SCAN_ENERGIES), + output: + "{DETECTOR_CONFIG}/sim_output/emcal_barrel_energy_scan_points_{PARTICLE}.txt", + run: + with open(output[0], "wt") as fp: + for energy in ENERGY_SCAN_ENERGIES: + fp.write(f"{energy}\n") + + +rule emcal_barrel_particles_energy_scan_analysis: + input: + script=workflow.source_path("scripts/emcal_barrel_energy_scan_analysis.cxx"), + scan_points="{DETECTOR_CONFIG}/sim_output/emcal_barrel_energy_scan_points_{PARTICLE}.txt", + output: + "{DETECTOR_CONFIG}/results/energy_scan/emcal_barrel_{PARTICLE}_fsam_scan.png", + "{DETECTOR_CONFIG}/results/energy_scan/emcal_barrel_{PARTICLE}_fsam_scan_res.png", + expand( + "{{DETECTOR_CONFIG}}/results/energy_scan/{energy}/emcal_barrel_{{PARTICLE}}_{plot}.png", + energy=ENERGY_SCAN_ENERGIES, + plot=["Esim_layer", "Layer_nodep", "Layer_Esim_mean", "Ethr", "nhits", "Esim", "fsam"], + ), + shell: + """ +cd {wildcards.DETECTOR_CONFIG} +root -l -b -q '{input.script}+("{wildcards.PARTICLE}")' +""" + + +rule emcal_barrel_pion_rejection_analysis: + input: + workflow.source_path("scripts/emcal_barrel_common_functions.h"), + script=workflow.source_path("scripts/emcal_barrel_pion_rejection_analysis.cxx"), + electron="{DETECTOR_CONFIG}/sim_output/sim_emcal_barrel_electron_energies1.0_18.0.edm4hep.root", + piminus="{DETECTOR_CONFIG}/sim_output/sim_emcal_barrel_piminus_energies1.0_18.0.edm4hep.root", + output: + "{DETECTOR_CONFIG}/results/emcal_barrel_pion_rej_RatioRej.png", + "{DETECTOR_CONFIG}/results/Barrel_emcal_pion_rej.json", + expand( + "{{DETECTOR_CONFIG}}/results/emcal_barrel_pion_rej_uncut_comb_{var_save}.png", + var_save=["Esim", "EsimTot", "EDep6", "EDep6OverP", "pT", "eta", "EsimScFi", "EsimScFiOverP"], + ), + expand( + "{{DETECTOR_CONFIG}}/results/emcal_barrel_pion_rej_uncut_comb_E{i}Eta{j}.png", + i=range(6), + j=range(2), + ), + expand( + "{{DETECTOR_CONFIG}}/results/emcal_barrel_pion_rej_{tag}_E{energy}_eta{eta_bin}.{extension}", + energy=[5, 10, 18], + eta_bin=range(2, 4), + tag=( + ["cut_mom_ele", "cut_mom_pim", "cut_ratio_pim"] + + sum([[f"cut_{var}_ele", f"cut_{var}_pim", f"cut_{var}_comb"] for var in ["pT", "EDep6OverP"]], []) + ), + extension=["pdf", "png"], + ), + shell: + """ +cd {wildcards.DETECTOR_CONFIG} +root -l -b -q '{input.script}+("../{input.electron}", "../{input.piminus}")' +""" diff --git a/benchmarks/barrel_ecal/config.yml b/benchmarks/barrel_ecal/config.yml index 100037cc..a236f68b 100644 --- a/benchmarks/barrel_ecal/config.yml +++ b/benchmarks/barrel_ecal/config.yml @@ -2,36 +2,34 @@ sim:emcal_barrel_pions: extends: .det_benchmark stage: simulate script: - - bash benchmarks/barrel_ecal/run_emcal_barrel_particles.sh piplus - - bash benchmarks/barrel_ecal/run_emcal_barrel_particles.sh piminus + - snakemake --cores 1 $DETECTOR_CONFIG/sim_output/sim_emcal_barrel_{piplus,piminus}_energies5.0_5.0.edm4hep.root sim:emcal_barrel_pi0: extends: .det_benchmark stage: simulate script: - - bash benchmarks/barrel_ecal/run_emcal_barrel_particles.sh pi0 + - snakemake --cores 1 $DETECTOR_CONFIG/sim_output/sim_emcal_barrel_pi0_energies5.0_5.0.edm4hep.root sim:emcal_barrel_electrons: extends: .det_benchmark stage: simulate script: - - if [[ "$RUN_EXTENDED_BENCHMARK" == "true" ]] ; then bash benchmarks/barrel_ecal/run_emcal_barrel_energy_scan.sh electron ; fi + - if [[ "$RUN_EXTENDED_BENCHMARK" == "true" ]] ; then snakemake --cores $DETECTOR_CONFIG/results/energy_scan/emcal_barrel_electron_fsam_scan.png; fi - export JUGGLER_N_EVENTS=400 - - bash benchmarks/barrel_ecal/run_emcal_barrel_particles.sh electron + - snakemake --cores 1 $DETECTOR_CONFIG/sim_output/sim_emcal_barrel_electron_energies5.0_5.0.edm4hep.root sim:emcal_barrel_photons: extends: .det_benchmark stage: simulate script: - - if [[ "$RUN_EXTENDED_BENCHMARK" == "true" ]] ; then bash benchmarks/barrel_ecal/run_emcal_barrel_energy_scan.sh photon ; fi - - bash benchmarks/barrel_ecal/run_emcal_barrel_particles.sh photon + - if [[ "$RUN_EXTENDED_BENCHMARK" == "true" ]] ; then snakemake --cores $DETECTOR_CONFIG/results/energy_scan/emcal_barrel_proton_fsam_scan.png; fi + - snakemake --cores 1 $DETECTOR_CONFIG/sim_output/sim_emcal_barrel_photon_energies5.0_5.0.edm4hep.root sim:emcal_barrel_pion_rejection: extends: .det_benchmark stage: simulate script: - - bash benchmarks/barrel_ecal/run_emcal_barrel_pion_rejection.sh electron - - bash benchmarks/barrel_ecal/run_emcal_barrel_pion_rejection.sh piminus + - snakemake --cores 1 $DETECTOR_CONFIG/sim_output/sim_emcal_barrel_{piminus,electron}_energies1.0_18.0.edm4hep.root calib:emcal_barrel_electrons: extends: .det_benchmark @@ -40,8 +38,8 @@ calib:emcal_barrel_electrons: - ["sim:emcal_barrel_electrons"] script: - ls -lhtR sim_output/ - - rootls -t sim_output/sim_emcal_barrel_electron.edm4hep.root - - root -b -q 'benchmarks/barrel_ecal/scripts/emcal_barrel_particles_analysis.cxx+("electron", true)' + - rootls -t sim_output/sim_emcal_barrel_electron_energies5.0_5.0.edm4hep.root + - snakemake --cores 1 $DETECTOR_CONFIG/results/emcal_barrel_electron_calibration.json - mv sim_output/sim_emcal_barrel_electron.edm4hep.root results/. - echo "JSON file(s) from analysis:" ; cat results/*.json @@ -51,7 +49,7 @@ bench:emcal_barrel_pions: needs: - ["sim:emcal_barrel_pions"] script: - - root -b -q benchmarks/barrel_ecal/scripts/emcal_barrel_pions_analysis.cxx+ + - snakemake --cores 1 $DETECTOR_CONFIG/results/emcal_barrel_pions_Ethr.png bench:emcal_barrel_electrons_scan: extends: .det_benchmark @@ -59,7 +57,7 @@ bench:emcal_barrel_electrons_scan: needs: - ["sim:emcal_barrel_electrons"] script: - - if [[ "$RUN_EXTENDED_BENCHMARK" == "true" ]] ; then root -b -q 'benchmarks/barrel_ecal/scripts/emcal_barrel_energy_scan_analysis.cxx+("electron")' ; fi + - if [[ "$RUN_EXTENDED_BENCHMARK" == "true" ]] ; then snakemake --cores 1 $DETECTOR_CONFIG/results/energy_scan/emcal_barrel_electron_fsam_scan.png; fi bench:emcal_barrel_pi0: extends: .det_benchmark @@ -68,7 +66,7 @@ bench:emcal_barrel_pi0: - ["sim:emcal_barrel_pi0", "calib:emcal_barrel_electrons"] script: - echo "JSON file(s) from analysis:" ; cat results/*.json - - root -b -q benchmarks/barrel_ecal/scripts/emcal_barrel_pi0_analysis.cxx+ + - snakemake --cores 1 epic_craterlake/results/Barrel_emcal_pi0.json bench:emcal_barrel_photons: extends: .det_benchmark @@ -77,10 +75,10 @@ bench:emcal_barrel_photons: - ["sim:emcal_barrel_photons"] script: - ls -lhtR sim_output/ - - rootls -t sim_output/sim_emcal_barrel_photon.edm4hep.root - - root -b -q 'benchmarks/barrel_ecal/scripts/emcal_barrel_particles_analysis.cxx+("photon", false)' + - rootls -t sim_output/sim_emcal_barrel_photon_energies5.0_5.0.edm4hep.root + - snakemake --cores 1 $DETECTOR_CONFIG/results/emcal_barrel_photon_calibration.json - mv sim_output/sim_emcal_barrel_photon.edm4hep.root results/. - - if [[ "$RUN_EXTENDED_BENCHMARK" == "true" ]] ; then root -b -q 'benchmarks/barrel_ecal/scripts/emcal_barrel_energy_scan_analysis.cxx+("photon")' ; fi + - if [[ "$RUN_EXTENDED_BENCHMARK" == "true" ]] ; then snakemake --cores $DETECTOR_CONFIG/results/energy_scan/emcal_barrel_proton_fsam_scan.png; fi bench:emcal_barrel_pion_rejection: extends: .det_benchmark @@ -89,9 +87,9 @@ bench:emcal_barrel_pion_rejection: - ["sim:emcal_barrel_pion_rejection"] script: - ls -lhtR sim_output/ - - rootls -t sim_output/sim_emcal_barrel_piRej_piminus.edm4hep.root - - rootls -t sim_output/sim_emcal_barrel_piRej_electron.edm4hep.root - - root -b -q benchmarks/barrel_ecal/scripts/emcal_barrel_pion_rejection_analysis.cxx+ + - rootls -t $DETECTOR_CONFIG/sim_output/sim_emcal_barrel_piminus_energies1.0_18.0.edm4hep.root + - rootls -t $DETECTOR_CONFIG/sim_output/sim_emcal_barrel_electron_energies1.0_18.0.edm4hep.root + - snakemake --cores 1 $DETECTOR_CONFIG/results/Barrel_emcal_pion_rej.json collect_results:barrel_ecal: extends: .det_benchmark diff --git a/benchmarks/barrel_ecal/scripts/emcal_barrel_particles_analysis.cxx b/benchmarks/barrel_ecal/scripts/emcal_barrel_particles_analysis.cxx index ad6e057a..683ae5be 100644 --- a/benchmarks/barrel_ecal/scripts/emcal_barrel_particles_analysis.cxx +++ b/benchmarks/barrel_ecal/scripts/emcal_barrel_particles_analysis.cxx @@ -257,7 +257,7 @@ void emcal_barrel_particles_analysis(std::string particle_name = "electron", boo {"sampling_fraction_error_scfi", fSam_scfi_mean_err} }; if (save_calib) { - std::string calib_output_path = "results/emcal_barrel_calibration.json"; + std::string calib_output_path = fmt::format("results/emcal_barrel_{}_calibration.json", particle_name); std::cout << "Saving calibration results to " << calib_output_path << std::endl; std::ofstream o(calib_output_path); o << std::setw(4) << j << std::endl; diff --git a/benchmarks/barrel_ecal/scripts/emcal_barrel_pi0_analysis.cxx b/benchmarks/barrel_ecal/scripts/emcal_barrel_pi0_analysis.cxx index 6d4a51ef..55354fe1 100644 --- a/benchmarks/barrel_ecal/scripts/emcal_barrel_pi0_analysis.cxx +++ b/benchmarks/barrel_ecal/scripts/emcal_barrel_pi0_analysis.cxx @@ -49,9 +49,8 @@ void emcal_barrel_pi0_analysis( ROOT::RDataFrame d0("events", input_fname); // Sampling Fraction grabbed from json file - // Note that this value is derived from electrons json j; - std::ifstream prev_steps_ifstream("results/emcal_barrel_calibration.json"); + std::ifstream prev_steps_ifstream("results/emcal_barrel_electron_calibration.json"); prev_steps_ifstream >> j; double samp_frac = j["electron"]["sampling_fraction"]; diff --git a/benchmarks/barrel_ecal/scripts/emcal_barrel_pion_rejection_analysis.cxx b/benchmarks/barrel_ecal/scripts/emcal_barrel_pion_rejection_analysis.cxx index 27950554..836810f3 100644 --- a/benchmarks/barrel_ecal/scripts/emcal_barrel_pion_rejection_analysis.cxx +++ b/benchmarks/barrel_ecal/scripts/emcal_barrel_pion_rejection_analysis.cxx @@ -87,9 +87,8 @@ void emcal_barrel_pion_rejection_analysis( /* // Sampling Fraction grabbed from json file - // Note that this value is derived from electrons json j; - std::ifstream prev_steps_ifstream("results/emcal_barrel_calibration.json"); + std::ifstream prev_steps_ifstream("results/emcal_barrel_electron_calibration.json"); prev_steps_ifstream >> j; // Sampling Fraction diff --git a/benchmarks/barrel_hcal/config.yml b/benchmarks/barrel_hcal/config.yml index 1674c8e8..fa40a22f 100644 --- a/benchmarks/barrel_hcal/config.yml +++ b/benchmarks/barrel_hcal/config.yml @@ -1,106 +1,53 @@ -sim:hcal_barrel_pions: +sim:hcal_barrel: extends: .det_benchmark stage: simulate + parallel: + matrix: + - PARTICLE: ["piplus", "piminus", "kplus", "kminus", "kshort", "klong", "muon", "antimuon", "proton"] script: - - bash benchmarks/barrel_hcal/run_hcal_barrel_particles.sh piplus - - bash benchmarks/barrel_hcal/run_hcal_barrel_particles.sh piminus + - bash benchmarks/barrel_hcal/run_hcal_barrel_particles.sh ${PARTICLE} -sim:hcal_barrel_kaons: +sim:hcal_barrel:scan: extends: .det_benchmark stage: simulate + parallel: + matrix: + - PARTICLE: ["muon", "antimuon", "proton"] + E: ["0.25", "0.5", "1", "2", "3", "4", "7", "15", "20"] script: - - bash benchmarks/barrel_hcal/run_hcal_barrel_particles.sh kplus - - bash benchmarks/barrel_hcal/run_hcal_barrel_particles.sh kminus - - bash benchmarks/barrel_hcal/run_hcal_barrel_particles.sh kshort - - bash benchmarks/barrel_hcal/run_hcal_barrel_particles.sh klong + - bash benchmarks/barrel_hcal/run_hcal_barrel_energy_scan.sh ${PARTICLE} ${E} -sim:hcal_barrel_muons: - extends: .det_benchmark - stage: simulate - script: - - bash benchmarks/barrel_hcal/run_hcal_barrel_energy_scan.sh muon - - bash benchmarks/barrel_hcal/run_hcal_barrel_particles.sh muon - -sim:hcal_barrel_antimuons: - extends: .det_benchmark - stage: simulate - script: - - bash benchmarks/barrel_hcal/run_hcal_barrel_energy_scan.sh antimuon - - bash benchmarks/barrel_hcal/run_hcal_barrel_particles.sh antimuon - -sim:hcal_barrel_protons: - extends: .det_benchmark - stage: simulate - script: - - bash benchmarks/barrel_hcal/run_hcal_barrel_energy_scan.sh proton - - bash benchmarks/barrel_hcal/run_hcal_barrel_particles.sh proton - -bench:hcal_barrel_protons: - extends: .det_benchmark - stage: benchmarks - needs: - - ["sim:hcal_barrel_protons"] - script: - - ls -lhtR sim_output/ - - rootls -t sim_output/sim_hcal_barrel_proton.edm4hep.root - - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_particles_analysis.cxx+("proton")' - - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_energy_scan_analysis.cxx+("proton")' - -bench:hcal_barrel_muons: - extends: .det_benchmark - stage: benchmarks - needs: - - ["sim:hcal_barrel_muons"] - script: - - ls -lhtR sim_output/ - - rootls -t sim_output/sim_hcal_barrel_muon.edm4hep.root - - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_particles_analysis.cxx+("muon")' - - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_energy_scan_analysis.cxx+("muon")' - -bench:hcal_barrel_antimuons: - extends: .det_benchmark - stage: benchmarks - needs: - - ["sim:hcal_barrel_antimuons"] - script: - - ls -lhtR sim_output/ - - rootls -t sim_output/sim_hcal_barrel_antimuon.edm4hep.root - - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_particles_analysis.cxx+("antimuon")' - - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_energy_scan_analysis.cxx+("antimuon")' - -bench:hcal_barrel_kaons: +bench:hcal_barrel: extends: .det_benchmark stage: benchmarks needs: - - ["sim:hcal_barrel_kaons"] + - ["sim:hcal_barrel"] + parallel: + matrix: + - PARTICLE: ["piplus", "piminus", "kplus", "kminus", "kshort", "klong", "muon", "antimuon", "proton"] script: - ls -lhtR sim_output/ - - rootls -t sim_output/sim_hcal_barrel_kplus.edm4hep.root - - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_particles_analysis.cxx+("kplus")' - - rootls -t sim_output/sim_hcal_barrel_kminus.edm4hep.root - - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_particles_analysis.cxx+("kminus")' - - rootls -t sim_output/sim_hcal_barrel_kshort.edm4hep.root - - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_particles_analysis.cxx+("kshort")' - - rootls -t sim_output/sim_hcal_barrel_klong.edm4hep.root - - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_particles_analysis.cxx+("klong")' + - rootls -t sim_output/sim_hcal_barrel_${PARTICLE}.edm4hep.root + - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_particles_analysis.cxx+("'${PARTICLE}'")' -bench:hcal_barrel_pions: +bench:hcal_barrel:scan: extends: .det_benchmark stage: benchmarks needs: - - ["sim:hcal_barrel_pions"] + - ["sim:hcal_barrel:scan"] + parallel: + matrix: + - PARTICLE: ["muon", "antimuon", "proton"] script: - ls -lhtR sim_output/ - - rootls -t sim_output/sim_hcal_barrel_piplus.edm4hep.root - - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_particles_analysis.cxx+("piplus")' - - rootls -t sim_output/sim_hcal_barrel_piminus.edm4hep.root - - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_particles_analysis.cxx+("piminus")' + - sort -n sim_output/hcal_barrel_energy_scan_points_${PARTICLE}_*.txt > sim_output/hcal_barrel_energy_scan_points_${PARTICLE}.txt + - root -b -q 'benchmarks/barrel_hcal/scripts/hcal_barrel_energy_scan_analysis.cxx+("'${PARTICLE}'")' collect_results:barrel_hcal: extends: .det_benchmark stage: collect needs: - - ["bench:hcal_barrel_muons", "bench:hcal_barrel_protons", "bench:hcal_barrel_kaons", "bench:hcal_barrel_pions"] + - ["bench:hcal_barrel", "bench:hcal_barrel:scan"] script: - ls -lrht - echo " FIX ME" diff --git a/benchmarks/barrel_hcal/run_hcal_barrel_energy_scan.sh b/benchmarks/barrel_hcal/run_hcal_barrel_energy_scan.sh index 92388286..74282724 100755 --- a/benchmarks/barrel_hcal/run_hcal_barrel_energy_scan.sh +++ b/benchmarks/barrel_hcal/run_hcal_barrel_energy_scan.sh @@ -1,21 +1,28 @@ #!/bin/bash export PARTICLE=$1 -E_file="sim_output/hcal_barrel_energy_scan_points_${PARTICLE}.txt" +shift +E_file="sim_output/hcal_barrel_energy_scan_points_${PARTICLE}_${CI_JOB_ID}.txt" -#for E in 1 2 6 10 -for E in 0.25 0.5 1 2 3 4 7 15 20 +if [ $# -gt 0 ] ; then + E_VALUES=("$@") +else + E_VALUES=(0.25 0.5 1 2 3 4 7 15 20) +fi + +for E in ${E_VALUES[@]} do export E_START="$E" export E_END="$E" + export JUGGLER_FILE_NAME_TAG=hcal_barrel_${PARTICLE}_${E} bash benchmarks/barrel_hcal/run_hcal_barrel_particles.sh "${PARTICLE}" && echo "$E" >> "$E_file" || exit 1 path_rootfiles="sim_output/energy_scan/${E}/" path_plots="results/energy_scan/${E}/" mkdir -p "$path_rootfiles" mkdir -p "$path_plots" ls -lthaR sim_output/ - mv "sim_output/sim_hcal_barrel_${PARTICLE}.edm4hep.root" "$path_rootfiles" + mv "sim_output/sim_${JUGGLER_FILE_NAME_TAG}.edm4hep.root" "$path_rootfiles/sim_hcal_barrel_${PARTICLE}.edm4hep.root" done ls -lthaR sim_output diff --git a/benchmarks/barrel_hcal/run_hcal_barrel_particles.sh b/benchmarks/barrel_hcal/run_hcal_barrel_particles.sh index ddae1318..d9f459f8 100755 --- a/benchmarks/barrel_hcal/run_hcal_barrel_particles.sh +++ b/benchmarks/barrel_hcal/run_hcal_barrel_particles.sh @@ -21,10 +21,13 @@ if [[ ! -n "${PARTICLE}" ]] ; then export PARTICLE="electron" fi -export JUGGLER_FILE_NAME_TAG="hcal_barrel_${PARTICLE}" -export JUGGLER_GEN_FILE="${JUGGLER_FILE_NAME_TAG}.hepmc" +if [[ ! -n "${JUGGLER_FILE_NAME_TAG}" ]] ; then + export JUGGLER_FILE_NAME_TAG="hcal_barrel_${PARTICLE}" +fi + +export JUGGLER_GEN_FILE="data/${JUGGLER_FILE_NAME_TAG}.hepmc" -export JUGGLER_SIM_FILE="sim_${JUGGLER_FILE_NAME_TAG}.edm4hep.root" +export JUGGLER_SIM_FILE="sim_output/sim_${JUGGLER_FILE_NAME_TAG}.edm4hep.root" export JUGGLER_REC_FILE="rec_${JUGGLER_FILE_NAME_TAG}.root" echo "JUGGLER_N_EVENTS = ${JUGGLER_N_EVENTS}" @@ -51,8 +54,8 @@ ddsim --runType batch \ --filter.tracker edep0 \ --numberOfEvents ${JUGGLER_N_EVENTS} \ --compactFile ${DETECTOR_PATH}/${DETECTOR_CONFIG}.xml \ - --inputFiles data/${JUGGLER_FILE_NAME_TAG}.hepmc \ - --outputFile sim_output/${JUGGLER_SIM_FILE} + --inputFiles ${JUGGLER_GEN_FILE} \ + --outputFile ${JUGGLER_SIM_FILE} if [[ "$?" -ne "0" ]] ; then echo "ERROR running npdet" @@ -61,7 +64,3 @@ fi # Directory for plots mkdir -p results - -# Move ROOT output file -#mv ${JUGGLER_REC_FILE} sim_output/ - diff --git a/benchmarks/barrel_hcal/scripts/hcal_barrel_particles_gen.cxx b/benchmarks/barrel_hcal/scripts/hcal_barrel_particles_gen.cxx index 297a74f5..a35829db 100644 --- a/benchmarks/barrel_hcal/scripts/hcal_barrel_particles_gen.cxx +++ b/benchmarks/barrel_hcal/scripts/hcal_barrel_particles_gen.cxx @@ -18,6 +18,7 @@ #include #include #include +#include #include #include "hcal_barrel_common_functions.h" @@ -25,7 +26,12 @@ using namespace HepMC3; void hcal_barrel_particles_gen(int n_events = 1e6, double e_start = 0.0, double e_end = 20.0, std::string particle_name = "muon") { - std::string out_fname = fmt::format("./data/hcal_barrel_{}.hepmc", particle_name); + std::string out_fname; + auto env_fname = getenv("JUGGLER_GEN_FILE"); + if (env_fname != nullptr) + out_fname = env_fname; + else + out_fname = fmt::format("./data/hcal_barrel_{}.hepmc", particle_name); WriterAscii hepmc_output(out_fname); int events_parsed = 0; GenEvent evt(Units::GEV, Units::MM); diff --git a/benchmarks/barrel_hcal/scripts/hcal_barrel_particles_reader.cxx b/benchmarks/barrel_hcal/scripts/hcal_barrel_particles_reader.cxx index c21db831..f40f7c29 100644 --- a/benchmarks/barrel_hcal/scripts/hcal_barrel_particles_reader.cxx +++ b/benchmarks/barrel_hcal/scripts/hcal_barrel_particles_reader.cxx @@ -19,6 +19,7 @@ #include #include +#include #include #include "hcal_barrel_common_functions.h" @@ -50,7 +51,12 @@ void hcal_barrel_particles_reader(std::string particle_name = "electron") gStyle->SetPadLeftMargin(0.14); gStyle->SetPadRightMargin(0.17); - std::string in_fname = fmt::format("./data/hcal_barrel_{}.hepmc",particle_name); + std::string in_fname; + auto env_fname = getenv("JUGGLER_GEN_FILE"); + if (env_fname != nullptr) + in_fname = env_fname; + else + in_fname = fmt::format("./data/hcal_barrel_{}.hepmc", particle_name); ReaderAscii hepmc_input(in_fname); int events_parsed = 0; GenEvent evt(Units::GEV, Units::MM); diff --git a/benchmarks/common/org2py.awk b/benchmarks/common/org2py.awk new file mode 100644 index 00000000..0d23cf25 --- /dev/null +++ b/benchmarks/common/org2py.awk @@ -0,0 +1,22 @@ +BEGIN { + in_src = 0 + IGNORECASE=1 +} + +/^#\+begin_src\s+[^\s]*python/ { + in_src = 1 + match($0, /^ */) + spaces = RLENGTH + next +} + +/^#\+end_src/ { + in_src = 0 + next +} + +in_src { + re = "^ {" spaces "}" + gsub(re,"") + print +} diff --git a/benchmarks/ecal_gaps/Snakefile b/benchmarks/ecal_gaps/Snakefile new file mode 100644 index 00000000..a312ac1e --- /dev/null +++ b/benchmarks/ecal_gaps/Snakefile @@ -0,0 +1,85 @@ +import os + + +rule ecal_gaps_sim: + input: + steering_file=provider.remote(remote_path("EPIC/EVGEN/SINGLE/{PARTICLE}/{ENERGY}/{PHASE_SPACE}/{PARTICLE}_{ENERGY}_{PHASE_SPACE}.steer")), + warmup="warmup/{DETECTOR_CONFIG}.edm4hep.root", + output: + "sim_output/{DETECTOR_CONFIG}/{PARTICLE}/{ENERGY}/{PHASE_SPACE}/{PARTICLE}_{ENERGY}_{PHASE_SPACE}.{INDEX}.edm4hep.root", + log: + "sim_output/{DETECTOR_CONFIG}/{PARTICLE}/{ENERGY}/{PHASE_SPACE}/{PARTICLE}_{ENERGY}_{PHASE_SPACE}.{INDEX}.edm4hep.root.log", + wildcard_constraints: + PARTICLE="e-", + ENERGY="(500MeV|5GeV|20GeV)", + PHASE_SPACE="(3to50|45to135|130to177)deg", + INDEX="\d{4}", + params: + N_EVENTS=1000 + shell: + """ +ddsim \ + --runType batch \ + --enableGun \ + --steeringFile "{input.steering_file}" \ + --random.seed 1{wildcards.INDEX} \ + --filter.tracker edep0 \ + -v WARNING \ + --numberOfEvents {params.N_EVENTS} \ + --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ + --outputFile {output} +""" + + +rule ecal_gaps_recon: + input: + "sim_output/{DETECTOR_CONFIG}/{PARTICLE}/{ENERGY}/{PHASE_SPACE}/{PARTICLE}_{ENERGY}_{PHASE_SPACE}.{INDEX}.edm4hep.root", + output: + "sim_output/{DETECTOR_CONFIG}/{PARTICLE}/{ENERGY}/{PHASE_SPACE}/{PARTICLE}_{ENERGY}_{PHASE_SPACE}.{INDEX}.eicrecon.tree.edm4eic.root", + log: + "sim_output/{DETECTOR_CONFIG}/{PARTICLE}/{ENERGY}/{PHASE_SPACE}/{PARTICLE}_{ENERGY}_{PHASE_SPACE}.{INDEX}.eicrecon.tree.edm4eic.root.log", + wildcard_constraints: + INDEX="\d{4}", + shell: """ +eicrecon {input} -Ppodio:output_file={output} \ + -Ppodio:output_include_collections=EcalEndcapNRecHits,EcalBarrelScFiRecHits,EcalBarrelImagingRecHits,EcalEndcapPRecHits,MCParticles +""" + + +DETECTOR_CONFIG=os.environ["DETECTOR_CONFIG"] + +rule ecal_gaps: + input: + matplotlibrc=".matplotlibrc", + script="benchmarks/ecal_gaps/ecal_gaps.py", + # TODO pass as a file list? + _=expand( + "sim_output/{DETECTOR_CONFIG}/{PARTICLE}/{ENERGY}/{PHASE_SPACE}/{PARTICLE}_{ENERGY}_{PHASE_SPACE}.{INDEX:04d}.eicrecon.tree.edm4eic.root", + DETECTOR_CONFIG=DETECTOR_CONFIG, + PARTICLE=["e-"], + ENERGY=["500MeV", "5GeV", "20GeV"], + PHASE_SPACE=["3to50deg", "45to135deg", "130to177deg"], + INDEX=range(1), + ), + output: + directory("results/ecal_gaps"), + threads: workflow.cores + shell: + """ +cleanup() {{ + echo Cleaning up + kill $WORKER_PID $SCHEDULER_PID +}} +trap cleanup EXIT + +PORT=$RANDOM +dask scheduler --port $PORT & +export DASK_SCHEDULER=localhost:$PORT +SCHEDULER_PID=$! +dask worker tcp://$DASK_SCHEDULER --nworkers {threads} --nthreads 1 & +WORKER_PID=$! +env \ +MATPLOTLIBRC={input.matplotlibrc} \ +OUTPUT_DIR={output} \ +python {input.script} +""" diff --git a/benchmarks/ecal_gaps/config.yml b/benchmarks/ecal_gaps/config.yml new file mode 100644 index 00000000..e854964e --- /dev/null +++ b/benchmarks/ecal_gaps/config.yml @@ -0,0 +1,26 @@ +sim:ecal_gaps: + extends: .det_benchmark + stage: simulate + script: + - mkdir -p $LOCAL_DATA_PATH/input + - ln -s $LOCAL_DATA_PATH/input input + - snakemake --cores 10 ecal_gaps --omit-from ecal_gaps + +bench:ecal_gaps: + extends: .det_benchmark + stage: benchmarks + needs: + - ["sim:ecal_gaps"] + script: + - ln -s $LOCAL_DATA_PATH/input input + - export PYTHONUSERBASE=$LOCAL_DATA_PATH/deps + - pip install -r benchmarks/ecal_gaps/requirements.txt + - snakemake --cores 8 ecal_gaps + +collect_results:ecal_gaps: + extends: .det_benchmark + stage: collect + needs: + - "bench:ecal_gaps" + script: + - ls -lrht diff --git a/benchmarks/ecal_gaps/ecal_gaps.org b/benchmarks/ecal_gaps/ecal_gaps.org new file mode 100644 index 00000000..36b64e03 --- /dev/null +++ b/benchmarks/ecal_gaps/ecal_gaps.org @@ -0,0 +1,221 @@ +#+PROPERTY: header-args:jupyter-python :session /jpy:localhost#8888:gap :async yes :results drawer :exports both + +#+TITLE: ePIC ECal gap study +#+AUTHOR: Dmitry Kalinkin +#+OPTIONS: d:t + +#+begin_src jupyter-python :results silent +import os +from pathlib import Path + +import awkward as ak +import boost_histogram as bh +import dask +import dask_awkward as dak +import dask_histogram as dh +import numpy as np +import uproot +#+end_src + +#+begin_src jupyter-python :results slient +from dask.distributed import Client +client = Client(os.environ.get("DASK_SCHEDULER", "localhost:8786")) +#+end_src + +* Plotting setup + +#+begin_src jupyter-python :results silent +import matplotlib as mpl +import matplotlib.pyplot as plt + +def setup_presentation_style(): + mpl.rcParams.update(mpl.rcParamsDefault) + plt.style.use('ggplot') + plt.rcParams.update({ + 'axes.labelsize': 8, + 'axes.titlesize': 9, + 'figure.titlesize': 9, + 'figure.figsize': (4, 3), + 'legend.fontsize': 7, + 'xtick.labelsize': 8, + 'ytick.labelsize': 8, + 'pgf.rcfonts': False, + }) + +setup_presentation_style() +#+end_src + +** Settings + +#+begin_src jupyter-python :results silent +DETECTOR_CONFIG=os.environ.get("DETECTOR_CONFIG") + +output_dir=Path(os.environ.get("OUTPUT_DIR", "./")) +output_dir.mkdir(parents=True, exist_ok=True) +#+end_src + +* Analysis + +#+begin_src jupyter-python :results silent +def filter_name(name): + return ( + "PARAMETERS" not in name + and all(coll not in name for coll in ["_intMap", "_floatMap", "_stringMap", "_doubleMap"]) + ) + +import functools + +axis_eta = bh.axis.Regular(200, -4, 4) +axis_eta_coarse = bh.axis.Regular(100, -4, 4) + +@functools.cache +def get_events(particle="e-", energy="20GeV", num_files=1): + events = uproot.dask( + {} + | {f"sim_output/{DETECTOR_CONFIG}/{particle}/{energy}/3to50deg/{particle}_{energy}_3to50deg.{INDEX:04d}.eicrecon.tree.edm4eic.root": "events" for INDEX in range(num_files)} + | {f"sim_output/{DETECTOR_CONFIG}/{particle}/{energy}/45to135deg/{particle}_{energy}_45to135deg.{INDEX:04d}.eicrecon.tree.edm4eic.root": "events" for INDEX in range(num_files)} + | {f"sim_output/{DETECTOR_CONFIG}/{particle}/{energy}/130to177deg/{particle}_{energy}_130to177deg.{INDEX:04d}.eicrecon.tree.edm4eic.root": "events" for INDEX in range(num_files)} + , + filter_name=filter_name, open_files=False, steps_per_file=1, + ) + + pt = np.hypot(events["MCParticles.momentum.x"][:,0], events["MCParticles.momentum.y"][:,0]) + theta = np.arctan2(pt, events["MCParticles.momentum.z"][:,0]) + eta = -np.log(np.tan(theta / 2)) + p = np.hypot(pt, events["MCParticles.momentum.z"][:,0]) + + hist_norm = dh.factory( + eta, + axes=( + axis_eta, + ), + ).compute() + weight_lut = 1 / hist_norm.values(flow=True) + + def get_weight(array): + if ak.backend(array) == "typetracer": + ak.typetracer.touch_data(array) + return array + return ak.from_numpy(weight_lut[axis_eta.index(ak.to_numpy(array))]) + + weights = eta.map_partitions(get_weight) + + events["p_thrown"] = p + events["eta_thrown"] = eta + events["weights"] = weights + + return events +#+end_src + +#+begin_src jupyter-python +particle = "e-" + +for energy in ["500MeV", "5GeV", "20GeV"]: + events = get_events(particle=particle, energy=energy) + + for calo_name in ["EcalEndcapN", "EcalBarrelScFi", "EcalBarrelImaging", "EcalEndcapP"]: + cond = ak.sum(events[f"{calo_name}RecHits.energy"], axis=1) > 10e-3 # GeV + + hist, profile = client.gather(client.compute([ + dh.factory( + events["eta_thrown"][cond], + (ak.sum(events[f"{calo_name}RecHits.energy"], axis=1) / events["p_thrown"])[cond], + axes=( + axis_eta, + bh.axis.Regular(100, 0.1, 1.1), + ), + weights=events["weights"][cond], + ), + dh.factory( + events["eta_thrown"][cond], + sample=(ak.sum(events[f"{calo_name}RecHits.energy"], axis=1) / events["p_thrown"])[cond], + storage=bh.storage.WeightedMean(), + axes=( + axis_eta_coarse, + ), + weights=events["weights"][cond], + ), + ])) + + cmap = plt.get_cmap(name="rainbow", lut=None).copy() + cmap.set_under("none") + + plt.pcolormesh( + hist.axes[0].edges, + hist.axes[1].edges, + hist.values().T, + cmap=cmap, + norm=mpl.colors.LogNorm( + vmin=np.min(hist.values()[hist.values() > 0]), + ), + ) + plt.colorbar() + std = np.sqrt(profile.variances()) + cond = profile.values() > std + plt.errorbar(profile.axes[0].centers[cond], profile.values()[cond], yerr=std[cond], marker=".", markersize=2, color="black", ls="none", lw=0.6, capsize=1.) + plt.xlabel(r"$\eta_{thrown}$", loc="right") + plt.ylabel(r"$\sum E_{\mathrm{dep.}} / p_{\mathrm{thrown}}$", loc="top") + plt.title(f"{energy} {particle} in {calo_name}") + plt.minorticks_on() + plt.savefig(output_dir / f"ecal_gap_{particle}_{energy}_{calo_name}.png", bbox_inches="tight") + plt.show() + plt.clf() +#+end_src + +#+begin_src jupyter-python +particle = "e-" + +for energy in ["500MeV", "5GeV", "20GeV"]: + events = get_events(particle=particle, energy=energy) + + calos = ["EcalEndcapN", "EcalBarrelScFi", "EcalEndcapP"] + total_energy = sum([ + ak.sum(events[f"{calo_name}RecHits.energy"], axis=1) + for calo_name in calos + ]) + + hist, profile = client.gather(client.compute([ + dh.factory( + events["eta_thrown"], + total_energy / events["p_thrown"], + axes=( + axis_eta, + bh.axis.Regular(100, 0.1, 1.1), + ), + weights=events["weights"], + ), + dh.factory( + events["eta_thrown"], + sample=total_energy / events["p_thrown"], + storage=bh.storage.WeightedMean(), + axes=( + axis_eta_coarse, + ), + weights=events["weights"], + ), + ])) + + cmap = plt.get_cmap(name="rainbow", lut=None).copy() + cmap.set_under("none") + + plt.pcolormesh( + hist.axes[0].edges, + hist.axes[1].edges, + hist.values().T, + cmap=cmap, + norm=mpl.colors.LogNorm( + vmin=np.min(hist.values()[hist.values() > 0]), + ), + ) + plt.colorbar() + std = np.sqrt(profile.variances()) + cond = profile.values() > std + plt.errorbar(profile.axes[0].centers[cond], profile.values()[cond], yerr=std[cond], marker=".", markersize=2, color="black", ls="none", lw=0.6, capsize=1.) + plt.xlabel(r"$\eta_{thrown}$", loc="right") + plt.ylabel(r"$\sum E_{\mathrm{dep.}} / p_{\mathrm{thrown}}$", loc="top") + plt.title(f"{energy} {particle}\n" + "+".join(calos)) + plt.minorticks_on() + plt.savefig(output_dir / f"ecal_gap_{particle}_{energy}_sum_all.png", bbox_inches="tight") + plt.show() + plt.clf() +#+end_src diff --git a/benchmarks/ecal_gaps/requirements.txt b/benchmarks/ecal_gaps/requirements.txt new file mode 100644 index 00000000..dc376fe2 --- /dev/null +++ b/benchmarks/ecal_gaps/requirements.txt @@ -0,0 +1,7 @@ +awkward >= 2.4.0 +dask >= 2023 +dask_awkward +dask_histogram +distributed >= 2023 +pyhepmc +uproot ~= 5.2.0 diff --git a/benchmarks/material_maps/config.yml b/benchmarks/material_maps/config.yml index 8db84363..dcfb04fe 100644 --- a/benchmarks/material_maps/config.yml +++ b/benchmarks/material_maps/config.yml @@ -20,3 +20,4 @@ material_maps: artifacts: paths: - results/material_maps + allow_failure: true diff --git a/benchmarks/material_scan/Snakefile b/benchmarks/material_scan/Snakefile new file mode 100644 index 00000000..b295715a --- /dev/null +++ b/benchmarks/material_scan/Snakefile @@ -0,0 +1,21 @@ +rule material_scan_fetch_script: + output: + "material_scan.py", + shell: """ +curl -L --output {output} https://github.com/eic/epic/raw/main/scripts/subdetector_tests/material_scan.py +""" + +rule material_scan: + input: + script="material_scan.py", + output: + "{DETECTOR_CONFIG}/results/material_scan.png", + "{DETECTOR_CONFIG}/results/material_scan_agg.csv", + "{DETECTOR_CONFIG}/results/material_scan_details.pdf", + log: + "{DETECTOR_CONFIG}/results/material_scan.log", + shadow: "full" # avoid putting calibrations/fieldmaps to results + shell: """ +cd {wildcards.DETECTOR_CONFIG}/results +python ../../{input.script} $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml +""" diff --git a/benchmarks/material_scan/config.yml b/benchmarks/material_scan/config.yml new file mode 100644 index 00000000..f8bc0822 --- /dev/null +++ b/benchmarks/material_scan/config.yml @@ -0,0 +1,13 @@ +bench:material_scan: + extends: .det_benchmark + stage: benchmarks + script: + - snakemake --cores 1 epic_craterlake/results/material_scan_details.pdf + +collect_results:material_scan: + extends: .det_benchmark + stage: collect + needs: + - ["bench:material_scan"] + script: + - ls -lrht diff --git a/benchmarks/others/.gitignore b/benchmarks/others/.gitignore deleted file mode 100644 index e69de29b..00000000 diff --git a/benchmarks/others/config.yml b/benchmarks/others/config.yml deleted file mode 100644 index e3c6f147..00000000 --- a/benchmarks/others/config.yml +++ /dev/null @@ -1,81 +0,0 @@ -variables: - ETAMIN: "-4.5" - ETAMAX: "+4.5" - ETASTEP: "0.01" - PHIMIN: "0.0" - PHIMAX: "6.28318530718" - PHISTEP: "0.01" - TRACKING_RHOMAX: "103." - TRACKING_ZNMAX: "191." - TRACKING_ZPMAX: "350." - ECAL_RHOMAX: "230." - ECAL_ZNMAX: "355." - ECAL_ZPMAX: "380." - -bench:materialScanEta: - stage: benchmarks - extends: .det_benchmark - script: - - echo ".x benchmarks/others/materialScanEta.cxx($ETAMIN, $ETAMAX, $ETASTEP, $PHIMIN)" | materialScan ${DETECTOR_PATH}/${DETECTOR_CONFIG}.xml -interactive - - mkdir results/images/ - - mv materialScanEta.png results/images/materialScanEta.png - - mv materialScanEta.pdf results/images/materialScanEta.pdf - allow_failure: true - -bench:materialScanEta:tracking: - stage: benchmarks - extends: .det_benchmark - script: - - echo ".x benchmarks/others/materialScanEta.cxx($ETAMIN, $ETAMAX, $ETASTEP, $PHIMIN, $TRACKING_RHOMAX, $TRACKING_ZNMAX, $TRACKING_ZPMAX)" | materialScan ${DETECTOR_PATH}/${DETECTOR_CONFIG}.xml -interactive - - mkdir results/images/ - - mv materialScanEta.png results/images/materialScanEtaTracking.png - - mv materialScanEta.pdf results/images/materialScanEtaTracking.pdf - allow_failure: true - -bench:materialScanEta:ecal: - stage: benchmarks - extends: .det_benchmark - script: - - echo ".x benchmarks/others/materialScanEta.cxx($ETAMIN, $ETAMAX, $ETASTEP, $PHIMIN, $ECAL_RHOMAX, $ECAL_ZNMAX, $ECAL_ZPMAX)" | materialScan ${DETECTOR_PATH}/${DETECTOR_CONFIG}.xml -interactive - - mkdir results/images/ - - mv materialScanEta.png results/images/materialScanEtaEcal.png - - mv materialScanEta.pdf results/images/materialScanEtaEcal.pdf - allow_failure: true - -bench:materialScanEtaPhi: - stage: benchmarks - extends: .det_benchmark - script: - - echo ".x benchmarks/others/materialScanEtaPhi.cxx($ETAMIN, $ETAMAX, $ETASTEP, $PHIMIN, $PHIMAX, $PHISTEP)" | materialScan ${DETECTOR_PATH}/${DETECTOR_CONFIG}.xml -interactive - - mkdir results/images/ - - mv materialScanEtaPhi.png results/images/materialScanEtaPhi.png - - mv materialScanEtaPhi.pdf results/images/materialScanEtaPhi.pdf - allow_failure: true - -bench:materialScanEtaPhi:tracking: - stage: benchmarks - extends: .det_benchmark - script: - - echo ".x benchmarks/others/materialScanEtaPhi.cxx($ETAMIN, $ETAMAX, $ETASTEP, $PHIMIN, $PHIMAX, $PHISTEP, $TRACKING_RHOMAX, $TRACKING_ZNMAX, $TRACKING_ZPMAX)" | materialScan ${DETECTOR_PATH}/${DETECTOR_CONFIG}.xml -interactive - - mkdir results/images/ - - mv materialScanEtaPhi.png results/images/materialScanEtaPhiTracking.png - - mv materialScanEtaPhi.pdf results/images/materialScanEtaPhiTracking.pdf - allow_failure: true - -bench:materialScanEtaPhi:ecal: - stage: benchmarks - extends: .det_benchmark - script: - - echo ".x benchmarks/others/materialScanEtaPhi.cxx($ETAMIN, $ETAMAX, $ETASTEP, $PHIMIN, $PHIMAX, $PHISTEP, $ECAL_RHOMAX, $ECAL_ZNMAX, $ECAL_ZPMAX)" | materialScan ${DETECTOR_PATH}/${DETECTOR_CONFIG}.xml -interactive - - mkdir results/images/ - - mv materialScanEtaPhi.png results/images/materialScanEtaPhiEcal.png - - mv materialScanEtaPhi.pdf results/images/materialScanEtaPhiEcal.pdf - allow_failure: true - -collect_results:materialscan: - extends: .det_benchmark - stage: collect - needs: - - ["bench:materialScanEta", "bench:materialScanEta:tracking", "bench:materialScanEta:ecal", "bench:materialScanEtaPhi", "bench:materialScanEtaPhi:tracking", "bench:materialScanEtaPhi:ecal"] - script: - - ls -lrht diff --git a/benchmarks/others/materialScanEta.cxx b/benchmarks/others/materialScanEta.cxx deleted file mode 100644 index 1dc43c3f..00000000 --- a/benchmarks/others/materialScanEta.cxx +++ /dev/null @@ -1,132 +0,0 @@ -#import -#import - -#include -#include - -#include -#include - -Color_t color(const Material& m) { - if (m.name() == std::string("Silicon")) return kGray; - else if (m.name() == std::string("Aluminum")) return kAzure; - else if (m.name() == std::string("CarbonFiber")) return kGray; - else if (m.name() == std::string("Beryllium")) return kGreen; - else if (m.name() == std::string("Gold")) return kYellow; - else if (m.name() == std::string("Mylar")) return kGreen; - else if (m.name() == std::string("Kapton")) return kGreen; - else if (m.name() == std::string("Copper")) return kGreen; - else if (m.name() == std::string("C2F6_DRICH")) return kOrange; - else if (m.name() == std::string("Ar10CO2")) return kOrange; - else if (m.name() == std::string("Aerogel")) return kPink; - else if (m.name() == std::string("AerogelOptical")) return kPink; - else if (m.name() == std::string("Aerogel_DRICH")) return kPink; - else if (m.name() == std::string("Lead")) return kBlack; - else if (m.name() == std::string("Steel235")) return kGray+2; - else if (m.name() == std::string("TungstenDens24")) return kBlack; - else if (m.name() == std::string("Polystyrene")) return kGray; - else if (m.name() == std::string("PolystyreneFoam")) return kGray; - else if (m.name() == std::string("Epoxy")) return kGray; - else if (m.name() == std::string("PlasticScint")) return kGray; - else if (m.name() == std::string("AcrylicOptical")) return kGray; - else if (m.name() == std::string("Acrylic_DRICH")) return kGray; - else if (m.name() == std::string("Quartz")) return kViolet; - else if (m.name() == std::string("Air")) return kBlue; - else if (m.name() == std::string("AirOptical")) return kBlue; - else if (m.name() == std::string("Vacuum")) return kWhite; - else { - std::cout << "Unknown material: " << m.name() << std::endl; - return kRed; - } -} - -void materialScanEta( - double etamin = -2.0, // minimum eta - double etamax = +1.0, // maximum eta - double etastep = 0.2, // steps in eta - double phi = 0.0, // phi angle for material scan - double rhomax = 10000.0, // maximum distance from z axis - double znmax = 10000.0, // maximum negative endcap z plane (positive number) - double zpmax = 10000.0 // maximum positive endcap z plane (positive number) -) { - // check inputs - if (etamin > etamax) { - std::cout << "Error: ordered eta range required" << std::endl; - return; - } - if (rhomax <= 0.0) { - std::cout << "Error: positive rhomax required" << std::endl; - return; - } - if (znmax <= 0.0) { - std::cout << "Error: positive znmax required" << std::endl; - return; - } - if (zpmax <= 0.0) { - std::cout << "Error: positive zpmax required" << std::endl; - return; - } - - // get material scans - size_t total{0}; - std::vector scan; - double x0{0}, y0{0}, z0{0}; - for (double eta = etamin; eta <= etamax + 0.5*etastep; eta += etastep) { - double theta = 2.0 * (atan(1) - atan(exp(-eta))); // |theta| < 90 deg, cos(theta) > 0, sin(theta) can be 0 - double r = min((theta < 0? -znmax: zpmax) / sin(theta), rhomax / cos(theta)); // theta == 0 results in min(+inf, rhomax) - double x = r * cos(theta) * cos(phi); - double y = r * cos(theta) * sin(phi); - double z = r * sin(theta); - scan.emplace_back(gMaterialScan->scan(x0,y0,z0,x,y,z)); - total += scan.back().size(); - } - - // start creating histograms for stacking: - // - start with first material layer at central eta bin - // - pop front layers first positive - size_t layer = 0; - std::vector histograms; - while (total > 0) { - // find next layer, starting from center bin outwards - size_t eta0 = static_cast(std::rint(abs((etamax-etamin)/etastep/2))); - for (size_t i = 0, j = eta0; - i < scan.size(); - ++i, j += (2*eta0 < scan.size()? -1: +1) * (i <= 2*min(eta0,scan.size()-eta0-1)? (i%2==0? -i: +i): -1)) { - if (scan.at(j).size() == 0) continue; - // define layer - auto layer_anchor = scan.at(j).at(0); - histograms.emplace_back(Form("h%zu",layer++),layer_anchor.first.name(),scan.size(),etamin,etamax); - //histograms.back().SetLineColor(color(layer_anchor.first)); - histograms.back().SetFillColor(color(layer_anchor.first)); - // add all bins to this layer - for (size_t bin = 0; bin < scan.size(); bin++) { - double X0{0}; - for (auto& mat_len: scan.at(bin)) { - if (mat_len.first.name() == layer_anchor.first.name()) { - X0 += mat_len.second / mat_len.first.radLength(); - scan.at(bin).erase(scan.at(bin).begin()); total--; - } else { - break; - } - } - histograms.back().SetBinContent(bin+1, X0); // bins start at 1 - } - } - } - std::cout << histograms.size() << " histograms created" << std::endl; - - // plot histograms as stack - THStack hs("hs",Form("Material Scan (rho < %.0f cm, -%.0f cm < z < %.0f cm)", rhomax, znmax, zpmax)); - for (auto& h: histograms) { - hs.Add(&h); - } - TCanvas cs("cs","Material Scan",1920,1080); - auto pad = cs.cd(); - pad->SetLogy(); - hs.Draw(); - hs.GetXaxis()->SetTitle("eta"); - hs.GetYaxis()->SetTitle("Fraction X0"); - hs.SetMinimum(2.5e-3); - cs.SaveAs("materialScanEta.png"); - cs.SaveAs("materialScanEta.pdf"); -} diff --git a/benchmarks/others/materialScanEtaPhi.cxx b/benchmarks/others/materialScanEtaPhi.cxx deleted file mode 100644 index 5051ca1b..00000000 --- a/benchmarks/others/materialScanEtaPhi.cxx +++ /dev/null @@ -1,74 +0,0 @@ -#import -#import - -#include -#include - -#include -#include - -void materialScanEtaPhi( - double etamin = -2.0, // minimum eta - double etamax = +1.0, // maximum eta - double etastep = 0.2, // steps in eta - double phimin = 0.0, // minimum phi - double phimax = 2.0 * M_PI, // maximum phi - double phistep = M_PI / 2, // steps in phi - double rhomax = 10000.0, // maximum distance from z axis - double znmax = 10000.0, // maximum negative endcap z plane (positive number) - double zpmax = 10000.0 // maximum positive endcap z plane (positive number) -) { - // check inputs - if (etamin > etamax) { - std::cout << "Error: ordered eta range required" << std::endl; - return; - } - if (phimin > phimax) { - std::cout << "Error: ordered phi range required" << std::endl; - return; - } - if (rhomax <= 0.0) { - std::cout << "Error: positive rhomax required" << std::endl; - return; - } - if (znmax <= 0.0) { - std::cout << "Error: positive znmax required" << std::endl; - return; - } - if (zpmax <= 0.0) { - std::cout << "Error: positive zpmax required" << std::endl; - return; - } - - // get material scans - size_t total{0}; - double x0{0}, y0{0}, z0{0}; - TH2F h2("h2","Material Scan Eta vs Phi", - (etamax-etamin)/etastep+1,etamin,etamax, - (phimax-phimin)/phistep+1,phimin,phimax - ); - for (double eta = etamin; eta <= etamax + 0.5*etastep; eta += etastep) { - for (double phi = phimin; phi <= phimax + 0.5*phistep; phi += phistep) { - double theta = 2.0 * (atan(1) - atan(exp(-eta))); // |theta| < 90 deg, cos(theta) > 0, sin(theta) can be 0 - double r = min((theta < 0? -znmax: zpmax) / sin(theta), rhomax / cos(theta)); // theta == 0 results in min(+inf, rhomax) - double x = r * cos(theta) * cos(phi); - double y = r * cos(theta) * sin(phi); - double z = r * sin(theta); - auto scan = gMaterialScan->scan(x0,y0,z0,x,y,z); - double X0{0}; - for (auto& mat_len: scan) - X0 += mat_len.second / mat_len.first.radLength(); - h2.Fill(eta,phi,X0); - } - } - - // plot histograms as stack - TCanvas cs("cs","Material Scan Eta vs Phi",1920,1080); - auto pad = cs.cd(); - cs.SetLogz(); - h2.Draw("colz"); - h2.GetXaxis()->SetTitle("eta"); - h2.GetYaxis()->SetTitle("phi"); - cs.SaveAs("materialScanEtaPhi.png"); - cs.SaveAs("materialScanEtaPhi.pdf"); -} diff --git a/benchmarks/zdc/config.yml b/benchmarks/zdc/config.yml index 2c299366..bee6a909 100644 --- a/benchmarks/zdc/config.yml +++ b/benchmarks/zdc/config.yml @@ -6,7 +6,6 @@ sim:zdc: parallel: matrix: - PARTICLE: ["neutron", "photon"] - allow_failure: true bench:zdc_benchmark: extends: .det_benchmark @@ -14,11 +13,10 @@ bench:zdc_benchmark: needs: - ["sim:zdc"] script: - - root -b -q benchmarks/zdc/scripts/zdc_analysis.cxx+ + - root -b -q -x 'benchmarks/zdc/scripts/analysis_zdc_particles.cxx+("sim_output/sim_zdc_'${PARTICLE}'.edm4hep.root", "results/far_forward/zdc/'${PARTICLE}'")' parallel: matrix: - PARTICLE: ["neutron", "photon"] - allow_failure: true collect_results:zdc: extends: .det_benchmark @@ -27,4 +25,3 @@ collect_results:zdc: - ["bench:zdc_benchmark"] script: - ls -lrht - allow_failure: true diff --git a/benchmarks/zdc/scripts/analysis_zdc_particles.cxx b/benchmarks/zdc/scripts/analysis_zdc_particles.cxx index c3339ce7..50026918 100644 --- a/benchmarks/zdc/scripts/analysis_zdc_particles.cxx +++ b/benchmarks/zdc/scripts/analysis_zdc_particles.cxx @@ -89,14 +89,14 @@ void analysis_zdc_particles( // Define variables auto d1 = d0 .Define("Ethr", Ethr, {"MCParticles"}) - .Define("nhits_Ecal", nhits, {"ZDCEcalHits"}) - .Define("Esim_Ecal", Esim, {"ZDCEcalHits"}) + .Define("nhits_Ecal", nhits, {"EcalFarForwardZDCHits"}) + .Define("Esim_Ecal", Esim, {"EcalFarForwardZDCHits"}) .Define("fsam_Ecal", fsam, {"Esim_Ecal", "Ethr"}) - .Define("nhits_Hcal", nhits, {"ZDCHcalHits"}) - .Define("Esim_Hcal", Esim, {"ZDCHcalHits"}) + .Define("nhits_Hcal", nhits, {"HcalFarForwardZDCHits"}) + .Define("Esim_Hcal", Esim, {"HcalFarForwardZDCHits"}) .Define("fsam_Hcal", fsam, {"Esim_Hcal", "Ethr"}) - .Define("hit_x_position", hit_x_position, {"ZDCHcalHits"}) - .Define("hit_y_position", hit_y_position, {"ZDCHcalHits"}) + .Define("hit_x_position", hit_x_position, {"HcalFarForwardZDCHits"}) + .Define("hit_y_position", hit_y_position, {"HcalFarForwardZDCHits"}) ; // Define Histograms diff --git a/config.yaml b/config.yaml new file mode 100644 index 00000000..75d4ff68 --- /dev/null +++ b/config.yaml @@ -0,0 +1 @@ +remote: S3