diff --git a/.gitignore b/.gitignore index e43b0f9..74d54ac 100644 --- a/.gitignore +++ b/.gitignore @@ -1 +1,2 @@ .DS_Store +tmpdir diff --git a/environment.yml b/environment.yml index 161c0e2..8cf2fae 100644 --- a/environment.yml +++ b/environment.yml @@ -2,6 +2,9 @@ name: eventdisplay-release-tests channels: - conda-forge dependencies: + - python>=3.8 + - numpy + - matplotlib - pandoc - pre-commit - towncrier @@ -11,3 +14,4 @@ dependencies: # activate: conda activate eventdisplay-release-tests # update (conda/mamba): conda env update -f environment.yml --prune # update (micromamba): micromamba update -f environment.yml +# update (micromamba): micromamba update -n eventdisplay-release-tests -f environment.yml diff --git a/release_tests/montecarlo/irf_plotting/irf_plotting.sh b/release_tests/montecarlo/irf_plotting/irf_plotting.sh index 02ddf03..7a41c35 100755 --- a/release_tests/montecarlo/irf_plotting/irf_plotting.sh +++ b/release_tests/montecarlo/irf_plotting/irf_plotting.sh @@ -18,6 +18,7 @@ CUT="NTel3-PointSource-Hard-TMVA-BDT" CUT="NTel2-PointSource-Soft-TMVA-BDT" CUT="NTel2-PointSource-Moderate-TMVA-BDT" CUT="NTel2-PointSource-Moderate" +echo "WARNING: CUT is hardwired to $CUT (for now)" # Comparison plots - version and simtype hardwired COMPAREVERSION="v492" COMPARESIMTYPE="CARE_202404" diff --git a/release_tests/montecarlo/mc_data_comparison/.gitignore b/release_tests/montecarlo/mc_data_comparison/.gitignore deleted file mode 100644 index 15c6bf7..0000000 --- a/release_tests/montecarlo/mc_data_comparison/.gitignore +++ /dev/null @@ -1 +0,0 @@ -tmpdir diff --git a/release_tests/montecarlo/mc_data_comparison/README.md b/release_tests/montecarlo/mc_data_comparison/README.md index a89720b..5596bb7 100644 --- a/release_tests/montecarlo/mc_data_comparison/README.md +++ b/release_tests/montecarlo/mc_data_comparison/README.md @@ -1,20 +1,97 @@ -# MC / data comparison +# MC / Data Comparison -Scripts and tools to compare MC distributions with results from the Crab Nebula. +**No apptainer usage!** -Fill and plot distributions with +Scripts and tools to compare Monte Carlo (MC) distributions with observational results from the Crab Nebula. + +## Quick Start ```console -./compareDatawithMC.sh +./compareDatawithMC.sh ``` -Requires as input: +Where `` specifies the observation conditions: + +- `SZE`: Small zenith angles (0-25°, 0.5° wobble) +- `MZE`: Medium zenith angles (40-50°, 0.5° wobble) +- `LZE`: Large zenith angles (50-70°, 0.5° wobble) +- `WOBBLE`: Large wobble offsets (MC at 1° wobble, summer atmosphere only) + +NOTE! The zenith angle bins are different compared to those used in the Crab Nebula analysis. + +## Prerequisites + +### Input Data Requirements + +- MC simulation files for each minor epoch (organized by atmosphere type) +- Crab Nebula mscw results for each minor epoch +- Run lists for Crab observations (generated separately) + +## Preparation Steps + +Before running the comparison, you must generate the Crab run lists: + +```console +# From release_tests/sources/Crab/ +./runlist_generator_from_anasum_log.sh +``` + +This generates run lists for minor epochs, zenith angle ranges, and atmosphere types. +The run lists are saved in `EventDisplay_Release_/Crab/runlists/`. + +## Running the Comparison + +### Step 1: Prepare Run Parameter File + +The `` must contain the following fields (one per line): -- MC files for each minor epoch -- Crab mscw results for each minor epoch +```text +* VERSION # e.g., V6 +* SIMTYPE # e.g., CARE, CARE_RedHV +* ATMOSPHERE # e.g., 61, 62 (repeat line for multiple values) +* EPOCH [nsb] # e.g., V6_2019_2020s 0 +* CRAB_NSB # e.g., 0, 1, or NOTSET +``` + +The leading `*` is part of the parameter-file syntax. To provide multiple atmospheres or epochs, repeat the corresponding line. If `CRAB_NSB` is `NOTSET`, the script reads the per-epoch NSB value from the optional 4th field of each `* EPOCH ...` line. + +### Step 2: Execute Comparison + +```console +cd release_tests/montecarlo/mc_data_comparison/ +./compareDatawithMC.sh SZE +``` + +This will: + +1. Parse the run parameter file +2. Locate MC simulation files based on epoch and atmosphere +3. Find corresponding Crab observational data +4. Submit Condor jobs for each combination +5. Generate comparison plots as PDF files + +### Step 3: Monitor Progress + +- Condor jobs are submitted with resource requests: 4000MB RAM, 10GB disk +- Temporary symbolic links are created in: `{output_dir}/tmp/{uuid}/` +- **Note**: Temporary directories must be cleaned up manually after completion + +## Output Structure + +Results are written to: + +```text +EventDisplay_Release_/mc_data_comparison/// +├── _ATM___/ +│ ├── mcdatacomparison.runparameter +│ ├── mcdatacomparison.root +│ ├── mcdatacomparison.log +│ └── *.pdf # Comparison plots +└── tmp/ # Temporary symlinks (manual cleanup required) +``` -**important:** -first run linking of epochs from `../../sources/Crab/`: `./runlist_generator.sh V6` -(or for any other epoch) +## Special Cases -Output and plots are written as PDFs into the `../../../../EventDisplay_Release_/mc_data_comparison/` directories +- **CARE_RedHV simulations**: Handled with modified BDT cut (BDT=0 instead of 1) +- **Summer vs Winter Atmosphere**: Automatically selected based on epoch suffix ('s' or 'w') +- **Wobble Mode**: Forces atmosphere 61 and 1.0° wobble offset diff --git a/release_tests/montecarlo/mc_data_comparison/compareDatawithMC.sh b/release_tests/montecarlo/mc_data_comparison/compareDatawithMC.sh index 440978b..7335918 100755 --- a/release_tests/montecarlo/mc_data_comparison/compareDatawithMC.sh +++ b/release_tests/montecarlo/mc_data_comparison/compareDatawithMC.sh @@ -4,13 +4,13 @@ # # requires: # - MC files for each minor epoch -# - Crab results for each minor epoch (read from Crab runlists; used pre-processed data) +# - Crab results for each minor epoch (read from Crab run lists; used pre-processed data) set -e if [[ $# < 2 ]]; then echo " - ./compareDatawithMC.sh + ./compareDatawithMC.sh --> choose zenith angle / wobble range SZE: small zenith angles (0.5 deg wobble) MZE: medium large zenith angles (0.5 deg wobble) @@ -23,23 +23,24 @@ exit fi ########################### -# read runparameter file +# read run parameter file if [[ ! -e ${1} ]]; then - echo "Error, runparameter file not found: ${1}" + echo "Error, run parameter file not found: ${1}" exit fi +RUNPAR="${1}" # Eventdisplay version -VERSION=$(grep VERSION ${1} | awk '{print $3}') +VERSION=$(awk '$1 == "*" && $2 == "VERSION" {print $3; exit}' "${RUNPAR}") # Simulation type -SIMTYPE=$(grep SIMTYPE ${1} | awk '{print $3}') +SIMTYPE=$(awk '$1 == "*" && $2 == "SIMTYPE" {print $3; exit}' "${RUNPAR}") # Atmosphere -ATMOS=($(grep ATMOSPHERE ${1} | grep "*" | awk '{print $3}')) +ATMOS=($(awk '$1 == "*" && $2 == "ATMOSPHERE" {print $3}' "${RUNPAR}")) # Epochs -EPOCH=($(grep EPOCH ${1} | grep "*" | grep -v MAJOR | awk '{print $3}')) +EPOCH=($(awk '$1 == "*" && $2 == "EPOCH" {print $3}' "${RUNPAR}")) # Wobble (default observing; WOBBLE mode fixed to 1 deg) MCWOFF="0.5" # Crab NSB level -CRABNSB=$(grep CRAB_NSB ${1} | awk '{print $3}') +CRABNSB=$(awk '$1 == "*" && $2 == "CRAB_NSB" {print $3; exit}' "${RUNPAR}") # Analysis type ANALYSISTYPE="AP" DIRRECOTYPE="_DISP" @@ -51,17 +52,30 @@ if [[ ! -z $VERITAS_ANALYSIS_TYPE ]]; then DIRRECOTYPE="" fi fi +# Stereo reconstruction method +# (0: dispBDT, 2: XGB) +RECOMETHOD=0 +STEREOMETHOD="" +if [[ $RECOMETHOD == 2 ]]; then + STEREOMETHOD="_XGB" +fi ########################### # elevation range [[ "$2" ]] && ELE=$2 || ELE="SZE" +case "$ELE" in + SZE|MZE|LZE|WOBBLE) + ;; + *) + echo "Error: invalid elevation/wobble selection '${ELE}' (expected SZE, MZE, LZE, or WOBBLE)" + exit 1 + ;; +esac # Directory for simulations SIMDIR=${VERITAS_IRFPRODUCTION_DIR}/${VERSION}/${ANALYSISTYPE}/$SIMTYPE/ if [[ ! -e ${SIMDIR} ]]; then - if [[ ! -e ${SIMDIR} ]]; then - echo "Error: simulation directory not found: $SIMDIR" - exit - fi + echo "Error: simulation directory not found: $SIMDIR" + exit fi # Directory for mscw data files DDIR="$VERITAS_PREPROCESSED_DATA_DIR/${ANALYSISTYPE}/mscw/" @@ -78,8 +92,9 @@ if [[ ! -e ${CDIR} ]]; then fi # output directory for MC/Data comparison -BDIR=$(readlink -f "../../../../EventDisplay_Release_${VERSION}/mc_data_comparison/${ANALYSISTYPE}${DIRRECOTYPE}/${SIMTYPE}/") +BDIR="../../../../EventDisplay_Release_${VERSION}/mc_data_comparison/${ANALYSISTYPE}${DIRRECOTYPE}${STEREOMETHOD}/${SIMTYPE}/" mkdir -p ${BDIR} +BDIR=$(readlink -f "../../../../EventDisplay_Release_${VERSION}/mc_data_comparison/${ANALYSISTYPE}${DIRRECOTYPE}${STEREOMETHOD}/${SIMTYPE}/") echo "Results will be written to $BDIR" PWDIR=$(pwd) @@ -88,14 +103,27 @@ get_mscw_file() { data_dir="${1}" runn="${2}" - if [ ! -e ${data_dir}/$runn.mscw.root ]; then - if [[ ${runn} -lt 100000 ]]; then - EDIR="${data_dir}/${runn:0:1}/" - else - EDIR="${data_dir}/${runn:0:2}/" - fi + local direct_file="${data_dir}/${runn}.mscw.root" + local edir="" + local nested_file="" + + if [[ -e "${direct_file}" ]]; then + echo "${direct_file}" + return 0 + fi + + if [[ ${runn} -lt 100000 ]]; then + edir="${data_dir}/${runn:0:1}" + else + edir="${data_dir}/${runn:0:2}" + fi + nested_file="${edir}/${runn}.mscw.root" + if [[ -e "${nested_file}" ]]; then + echo "${nested_file}" + return 0 fi - echo "$EDIR/$runn.mscw.root" + + return 1 } for I in "${EPOCH[@]}" @@ -106,11 +134,15 @@ do fi # Crab NSB level (depends on epoch) if [[ $CRABNSB == "NOTSET" ]]; then - NSB=$(grep ${I} ${1} | grep -v MAJOR | awk '{print $4}') + NSB=$(awk -v epoch="${I}" '$1 == "*" && $2 == "EPOCH" && $3 == epoch {print $4; exit}' "${RUNPAR}") else NSB=${CRABNSB} fi - # Simulation file (elevation depedent) + if [[ -z "${NSB}" ]]; then + echo "Error: no NSB value found for epoch ${I} in ${RUNPAR}" + exit 1 + fi + # Simulation file (elevation dependent) if [[ $ELE = "SZE" ]] then simfile="20deg_${MCWOFF}wob_NOISE${NSB}.mscw.root" @@ -123,7 +155,7 @@ do ZEMAX="50." elif [[ $ELE = "LZE" ]] then - ZEMIN="50" + ZEMIN="50." ZEMAX="70." simfile="55deg_${MCWOFF}wob_NOISE${NSB}.mscw.root" # wobble set (everything not 0.5 deg) @@ -151,6 +183,7 @@ do if [[ $SIMTYPE == "CARE_RedHV"* ]]; then REDHV="_redHV" fi + SIMATM="${atm}" echo "Processing $I $A ${atm} $REDHV" # Crab run list @@ -178,7 +211,19 @@ do mkdir -p "$DMSCWDIR" DFILES=$(cat $RUNLIST) for D in $DFILES; do - ln -s -f $(get_mscw_file $DDIR $D) "$DMSCWDIR"/"$D".mscw.root + if ! MSCWFILE=$(get_mscw_file "$DDIR" "$D"); then + echo "Error: mscw file not found for run ${D} in ${DDIR}" + exit 1 + fi + ln -s -f "${MSCWFILE}" "$DMSCWDIR"/"$D".mscw.root + if [[ ${RECOMETHOD} -eq 2 ]]; then + XGBFILE="${MSCWFILE%.mscw.root}.mscw.xgb_stereo.root" + if [[ ! -e "${XGBFILE}" ]]; then + echo "Error: XGB stereo file not found for run ${D}: ${XGBFILE}" + exit 1 + fi + ln -s -f "${XGBFILE}" "$DMSCWDIR"/"$D".mscw.xgb_stereo.root + fi done echo "TEMP DIRECTORY (to be deleted by hand): $DMSCWDIR" @@ -191,14 +236,11 @@ do SIMMSCW="MSCW_RECID0${DIRRECOTYPE}" # write run parameter file - if [[ $SIMTYPE == "CARE_RedHV" ]]; then - echo "* SIMS $SIMDIR/${I}_ATM61_gamma/${SIMMSCW}/${simfile} 4 ${MCWOFF} 0. 0. 0. ${ZEMIN} ${ZEMAX}" > $ODIR/mcdatacomparison.runparameter + if [[ $SIMTYPE == "CARE_RedHV"* ]]; then + SIMATM="61" + echo "* SIMS $SIMDIR/${I}_ATM${SIMATM}_gamma/${SIMMSCW}/${simfile} 4 ${MCWOFF} 0. 0. 0. ${ZEMIN} ${ZEMAX}" > $ODIR/mcdatacomparison.runparameter else - if [[ $ELE == "WOBBLE" ]] && [[ ${I: -1} == "s" ]]; then - echo "* SIMS $SIMDIR/${I}_ATM62_gamma/${SIMMSCW}/${simfile} 4 ${MCWOFF} 0. 0. 0. ${ZEMIN} ${ZEMAX}" > $ODIR/mcdatacomparison.runparameter - else - echo "* SIMS $SIMDIR/${I}_ATM${atm}_gamma/${SIMMSCW}/${simfile} 4 ${MCWOFF} 0. 0. 0. ${ZEMIN} ${ZEMAX}" > $ODIR/mcdatacomparison.runparameter - fi + echo "* SIMS $SIMDIR/${I}_ATM${SIMATM}_gamma/${SIMMSCW}/${simfile} 4 ${MCWOFF} 0. 0. 0. ${ZEMIN} ${ZEMAX}" > $ODIR/mcdatacomparison.runparameter fi echo "* ON ${DMSCWDIR}/[0-9]*.mscw.root 4 -99. -99. 0. 360. ${ZEMIN} ${ZEMAX}" >> $ODIR/mcdatacomparison.runparameter echo "* OFF ${DMSCWDIR}/[0-9]*.mscw.root 4 -99. -99. 0. 360. ${ZEMIN} ${ZEMAX}" >> $ODIR/mcdatacomparison.runparameter @@ -207,8 +249,9 @@ do FSCRIPT="$DMSCWDIR/compareDatawithMC_qsub_${SIMTYPE}_${I}${A}_${ELE}_${MCWOFF}_${NSB}" rm -f ${FSCRIPT}.sh sed -e "s|OUTDIR|$ODIR|" \ - -e "s|EEPOCHTM|${I}_ATM${atm}|" \ - -e "s|CURRENTDIR|$PWDIR|" compareDatawithMC_qsub.sh > ${FSCRIPT}.sh + -e "s|EEPOCHTM|${I}_ATM${SIMATM}|" \ + -e "s|CURRENTDIR|$PWDIR|" \ + -e "s|METHODRECO|$RECOMETHOD|" compareDatawithMC_qsub.sh > ${FSCRIPT}.sh echo "Run script: $FSCRIPT" chmod u+x $FSCRIPT.sh @@ -216,8 +259,5 @@ do $EVNDISPSCRIPTS/helper_scripts/UTILITY.condorSubmission.sh ${FSCRIPT}.sh 4000M 10G condor_submit ${FSCRIPT}.sh.condor - continue - - done done diff --git a/release_tests/montecarlo/mc_data_comparison/compareDatawithMC_qsub.sh b/release_tests/montecarlo/mc_data_comparison/compareDatawithMC_qsub.sh index 080398f..ad1b0f0 100755 --- a/release_tests/montecarlo/mc_data_comparison/compareDatawithMC_qsub.sh +++ b/release_tests/montecarlo/mc_data_comparison/compareDatawithMC_qsub.sh @@ -7,8 +7,9 @@ source $EVNDISPSYS/setObservatory.sh VTS ODIR=OUTDIR CDIR=CURRENTDIR EPOCHATM=EEPOCHTM +RECOMETHOD=METHODRECO -BDT="1" +BDT="0" if [[ ${ODIR} == *"CARE_RedHV"* ]]; then BDT="0" fi @@ -20,7 +21,7 @@ $EVNDISPSYS/bin/compareDatawithMC \ $ODIR/mcdatacomparison.runparameter \ -3 \ $ODIR/mcdatacomparison.root \ - ${BDT} $EPOCHATM \ + ${BDT} $EPOCHATM $RECOMETHOD \ > $ODIR/mcdatacomparison.log # prepare all plots diff --git a/release_tests/sources/Crab/README.md b/release_tests/sources/Crab/README.md index 34f56af..b1c9330 100644 --- a/release_tests/sources/Crab/README.md +++ b/release_tests/sources/Crab/README.md @@ -43,7 +43,7 @@ for Combine files using pre-processed anasum files and run list generated in step before: ```bash -./anasum_yearly.sh +./anasum_from_runlists.sh ``` ## Plotting diff --git a/release_tests/sources/Crab/plot_all.sh b/release_tests/sources/Crab/plot_all.sh index 5486d36..8564882 100755 --- a/release_tests/sources/Crab/plot_all.sh +++ b/release_tests/sources/Crab/plot_all.sh @@ -19,8 +19,6 @@ ODIR="${2}/$CUT" echo "Output directory: $ODIR" mkdir -p "$ODIR" -PDIR=$(pwd) - # copy anasum log files to release test directory for F in $LFIL; do TF=$(basename $F .combined.root) diff --git a/release_tests/sources/Crab/plot_energy_spectra.C b/release_tests/sources/Crab/plot_energy_spectra.C index 3fc3e6d..882d484 100644 --- a/release_tests/sources/Crab/plot_energy_spectra.C +++ b/release_tests/sources/Crab/plot_energy_spectra.C @@ -14,6 +14,8 @@ #include "TF1.h" +R__LOAD_LIBRARY($EVNDISPSYS/lib/libVAnaSum.so) + #include "../../utilities/parameters.C" #include "../../utilities/printutilities.C" diff --git a/release_tests/sources/Crab/plot_lightcurves.C b/release_tests/sources/Crab/plot_lightcurves.C index 458f068..5b341a7 100644 --- a/release_tests/sources/Crab/plot_lightcurves.C +++ b/release_tests/sources/Crab/plot_lightcurves.C @@ -13,6 +13,8 @@ #include #include +R__LOAD_LIBRARY($EVNDISPSYS/lib/libVAnaSum.so) + #include "../../utilities/parameters.C" #include "../../utilities/printutilities.C" diff --git a/release_tests/sources/Crab/plot_sensitivity.C b/release_tests/sources/Crab/plot_sensitivity.C index b8c4fd8..84de9dd 100644 --- a/release_tests/sources/Crab/plot_sensitivity.C +++ b/release_tests/sources/Crab/plot_sensitivity.C @@ -4,6 +4,7 @@ */ #include +R__LOAD_LIBRARY($EVNDISPSYS/lib/libVAnaSum.so) #include "../../utilities/parameters.C" #include "../../utilities/printutilities.C" diff --git a/release_tests/sources/Crab/plot_sensitivity_from_anasum_log.py b/release_tests/sources/Crab/plot_sensitivity_from_anasum_log.py new file mode 100644 index 0000000..1409d6a --- /dev/null +++ b/release_tests/sources/Crab/plot_sensitivity_from_anasum_log.py @@ -0,0 +1,187 @@ +"""Plot sensitivity diagnostics from anasum log files. + +Reads anasum log files from an input directory and writes a PNG summary figure. + +Usage: python plot_sensitivity_from_anasum_log.py [input_dir] [output_png] +""" + +import argparse +import glob +import os +import re + +import matplotlib.pyplot as plt +import numpy as np + +# Parse input directory argument +parser = argparse.ArgumentParser(description="Plot sensitivity metrics from anasum log files.") +parser.add_argument( + "input_dir", + nargs="?", + default=".", + help="Directory containing anasum log files (default: current directory)", +) +parser.add_argument( + "output_png", + nargs="?", + default="anasum_analysis_plots.png", + help="Output PNG filename (default: anasum_analysis_plots.png)", +) +args = parser.parse_args() + +if not os.path.isdir(args.input_dir): + raise SystemExit(f"Input directory does not exist: {args.input_dir}") + +log_files = sorted( + glob.glob(os.path.join(args.input_dir, "*anasum*V6_20*ATM*SZE_0.5deg.combined.log")) +) + +if not log_files: + raise SystemExit(f"No anasum log files found in: {args.input_dir}") + +epochs = [] +on_off_summary_rates = [] +off_summary_rates = [] +crab_001_times = [] +all_elevations = [] +all_on_rates = [] +all_off_rates = [] + +for filepath in log_files: + print(f"Processing file: {filepath}") + with open(filepath, "r") as f: + content = f.read() + + # Extract epoch from filename (e.g., V6_2023_2023s_ATM62 or V6_2014_2015_ATM61) + basename = os.path.basename(filepath) + epoch_match = re.search(r"(V6_[^_]+_\d+[a-z]*)_ATM\d+", basename) + if not epoch_match: + # Try alternative pattern for ATM61 files + epoch_match = re.search( + r"(V6_[^_]+_\d+[a-z]*)_ATM\d+", basename.replace("_SZE_0.5deg.combined", "") + ) + if not epoch_match: + continue + epoch = epoch_match.group(1) + + # Extract summary on/off rates + rate_match = re.search(r"Rates \(on/Off\): (\d+\.\d+)\s+(\d+\.\d+)", content) + + # Extract 0.01 Crab time in hours from the sensitivity table + lines = content.split("\n") + crab_time = None + in_table = False + for line in lines: + if "[h]" in line: + in_table = True + continue + if in_table and line.strip().startswith("0.01"): + parts = line.split() + if len(parts) >= 3: + try: + crab_time = float(parts[-1]) + except ValueError: + pass + break + + # Extract all elevation and rate values from RUN lines + el_values = [] + on_rates = [] + off_rates = [] + for line in lines: + el_match = re.search(r"RUN\s+\d+.*at\s+(\d+),\s+(-?\d+)\s+deg\s+El\.,\s+Az", line) + rate_match_line = re.search(r"Rates:\s+([\d.]+).*background:\s+([\d.]+)", line) + if el_match: + el_values.append(int(el_match.group(1))) + if rate_match_line: + on_rates.append(float(rate_match_line.group(1))) + off_rates.append(float(rate_match_line.group(2))) + + # Only add data if we have all values and non-empty arrays + if rate_match and crab_time is not None and el_values and on_rates and off_rates: + epochs.append(epoch) + on_off_summary_rates.append(float(rate_match.group(1))) + off_summary_rates.append(float(rate_match.group(2))) + crab_001_times.append(crab_time) + all_elevations.append(np.array(el_values)) + all_on_rates.append(np.array(on_rates)) + all_off_rates.append(np.array(off_rates)) + +# Calculate average elevation per epoch for sensitivity plot +avg_elevations = [np.mean(el) for el in all_elevations] + +# Create plots - 4 rows: 2x2 + 2 extra +fig, axs = plt.subplots(4, 2, figsize=(14, 20)) + +# Plot 1: On Rate vs Epoch +axs[0, 0].plot(epochs, on_off_summary_rates, "o-", color="blue") +axs[0, 0].set_xlabel("Epoch") +axs[0, 0].set_ylabel("On Rate") +axs[0, 0].set_title("On Rate vs Epoch") +axs[0, 0].tick_params(axis="x", rotation=45) + +# Plot 2: Off Rate vs Epoch +axs[0, 1].plot(epochs, off_summary_rates, "s-", color="orange") +axs[0, 1].set_xlabel("Epoch") +axs[0, 1].set_ylabel("Off Rate") +axs[0, 1].set_title("Off Rate vs Epoch") +axs[0, 1].tick_params(axis="x", rotation=45) + +# Plot 3: Violin plot of Elevations per Epoch +axs[1, 0].violinplot(all_elevations, positions=range(len(epochs)), widths=0.8, showmeans=True) +axs[1, 0].set_xticks(range(len(epochs))) +axs[1, 0].set_xticklabels(epochs, rotation=45, ha="right") +axs[1, 0].set_xlabel("Epoch") +axs[1, 0].set_ylabel("Elevation [deg]") +axs[1, 0].set_title("Elevation Distribution per Epoch") + +# Plot 4: Violin plot of On Rates per Epoch +axs[1, 1].violinplot(all_on_rates, positions=range(len(epochs)), widths=0.8, showmeans=True) +axs[1, 1].set_xticks(range(len(epochs))) +axs[1, 1].set_xticklabels(epochs, rotation=45, ha="right") +axs[1, 1].set_xlabel("Epoch") +axs[1, 1].set_ylabel("On Rate") +axs[1, 1].set_title("On Rate Distribution per Epoch") + +# Plot 5: Violin plot of Off Rates per Epoch +axs[2, 0].violinplot(all_off_rates, positions=range(len(epochs)), widths=0.8, showmeans=True) +axs[2, 0].set_xticks(range(len(epochs))) +axs[2, 0].set_xticklabels(epochs, rotation=45, ha="right") +axs[2, 0].set_xlabel("Epoch") +axs[2, 0].set_ylabel("Off Rate") +axs[2, 0].set_title("Off Rate Distribution per Epoch") + +# Plot 6: 0.01 Crab Sensitivity Time vs Epoch +axs[2, 1].plot(epochs, crab_001_times, "s-", color="purple") +axs[2, 1].set_xlabel("Epoch") +axs[2, 1].set_ylabel("0.01 Crab Time [h]") +axs[2, 1].set_title("0.01 Crab Sensitivity Time vs Epoch") +axs[2, 1].tick_params(axis="x", rotation=45) + +# Plot 7: Sensitivity vs Average Elevation +# Color points red if epoch year >= 2020 +colors = [ + ( + "red" + if "2020" in epoch + or "2021" in epoch + or "2022" in epoch + or "2023" in epoch + or "2024" in epoch + or "2025" in epoch + or "2026" in epoch + else "blue" + ) + for epoch in epochs +] +axs[3, 0].scatter(avg_elevations, crab_001_times, c=colors, alpha=0.7) +axs[3, 0].set_xlabel("Average Elevation [deg]") +axs[3, 0].set_ylabel("0.01 Crab Time [h]") +axs[3, 0].set_title("Sensitivity vs Elevation") + +# Hide the empty subplot +axs[3, 1].axis("off") + +plt.tight_layout() +plt.savefig(args.output_png, dpi=300, bbox_inches="tight") +print(f"Processed {len(epochs)} files. Plots saved as {args.output_png}") diff --git a/release_tests/sources/Crab/plot_skymaps.C b/release_tests/sources/Crab/plot_skymaps.C index 98d905a..5c8fda3 100644 --- a/release_tests/sources/Crab/plot_skymaps.C +++ b/release_tests/sources/Crab/plot_skymaps.C @@ -11,6 +11,8 @@ #include #include +R__LOAD_LIBRARY($EVNDISPSYS/lib/libVAnaSum.so) + #include "../../utilities/parameters.C" #include "../../utilities/printutilities.C" diff --git a/release_tests/sources/Crab/runlist_generator_from_anasum_log.sh b/release_tests/sources/Crab/runlist_generator_from_anasum_log.sh index 913b3b1..95470a8 100755 --- a/release_tests/sources/Crab/runlist_generator_from_anasum_log.sh +++ b/release_tests/sources/Crab/runlist_generator_from_anasum_log.sh @@ -99,6 +99,11 @@ do echo "Run $R - log file not found: ${ANASUMLOG}" exit fi + # Skip runs without 4-telescope cuts + if ! grep -q "VGammaHadronCuts::printCutSummary() (ntel=4" "${ANASUMLOG}"; then + echo "Run $R is not a 4-telescope run; skipping" + continue + fi echo "DATADIR ${ANASUMLOG}" # read and extract run info from files INSTRUMENT_EPOCH=$(grep "Instrument epoch selected" "${ANASUMLOG}" | head -n 1 | awk '{print $NF}') @@ -117,7 +122,7 @@ do OBSL="stdHV" fi ELEV=$(grep "mean elevation" "${ANASUMLOG}" | head -n 1 | awk '{print $3}') - EL=$(echo $ELEV | awk -v e=$ELEV '{if (e > 50 ) {print "SZE"} else if (e > 40 ) {print "MZE"} else if (e > 30 ) {print "LZE"} else {print "BZE"}}') + EL=$(echo $ELEV | awk -v e=$ELEV '{if (e > 60 ) {print "SZE"} else if (e > 45 ) {print "MZE"} else if (e > 35 ) {print "LZE"} else {print "BZE"}}') WOBBLESTRING=$(grep "Wobble offsets (currE)" "${ANASUMLOG}") n_offset=$(echo "$WOBBLESTRING" | head -n 1 | awk '{print $5}') w_offset=$(echo "$WOBBLESTRING" | head -n 1 | awk '{print $7}') diff --git a/release_tests/sources/Crab/runlist_releaseTestingV6.dat b/release_tests/sources/Crab/runlist_releaseTestingV6.dat index da6d0e3..e5d6a73 100644 --- a/release_tests/sources/Crab/runlist_releaseTestingV6.dat +++ b/release_tests/sources/Crab/runlist_releaseTestingV6.dat @@ -814,7 +814,6 @@ 104174 104199 104201 -104202 104251 104252 104253 @@ -992,7 +991,6 @@ 113299 113303 113352 -113375 113515 113637 113689 diff --git a/release_tests/utilities/parameters.C b/release_tests/utilities/parameters.C index a6880a8..5b2f35d 100644 --- a/release_tests/utilities/parameters.C +++ b/release_tests/utilities/parameters.C @@ -17,8 +17,9 @@ * Helper for ROOT macros: load libVAnaSum from environment. * Search order: * 1) $VERITAS_VANASUM_LIBRARY - * 2) $EVNDISP/lib/libVAnaSum.so - * 3) ROOT library path via "libVAnaSum.so" + * 2) $EVNDISPSYS/lib/libVAnaSum.so + * 3) $EVNDISP/lib/libVAnaSum.so + * 4) ROOT library path via "libVAnaSum.so" */ bool loadVAnaSumLibrary() { @@ -35,6 +36,14 @@ bool loadVAnaSumLibrary() iLibPath = iEnvLib; } if( iLibPath.size() == 0 ) + { + const char* iEvndispSys = gSystem->Getenv( "EVNDISPSYS" ); + if( iEvndispSys ) + { + iLibPath = string( iEvndispSys ) + "/lib/libVAnaSum.so"; + } + } + if( iLibPath.size() == 0 ) { const char* iEvndisp = gSystem->Getenv( "EVNDISP" ); if( iEvndisp ) @@ -45,21 +54,26 @@ bool loadVAnaSumLibrary() if( iLibPath.size() > 0 && !gSystem->AccessPathName( iLibPath.c_str() ) ) { - if( gSystem->Load( iLibPath.c_str() ) >= 0 ) + // Use gROOT->ProcessLine to properly load the library and its dictionaries + string loadCmd = "R__LOAD_LIBRARY(" + iLibPath + ");"; + Longptr_t result = gROOT->ProcessLine( loadCmd.c_str() ); + if( result >= 0 ) { iLibraryLoaded = true; return true; } } - if( gSystem->Load( "libVAnaSum.so" ) >= 0 ) + // Try with just the library name + Longptr_t result = gROOT->ProcessLine( "R__LOAD_LIBRARY(libVAnaSum.so);" ); + if( result >= 0 ) { iLibraryLoaded = true; return true; } cout << "Error: unable to load libVAnaSum.so. " - << "Set VERITAS_VANASUM_LIBRARY or EVNDISP." << endl; + << "Set VERITAS_VANASUM_LIBRARY, EVNDISPSYS, or EVNDISP." << endl; return false; }