Fix clustering benchmarks for missing ID.value

4ba38b23 · Wouter Deconinck · 16efa6e7 · 4ba38b23 · 4ba38b23 · 4ba38b23
Commit 4ba38b23 authored 2 years ago by Wouter Deconinck
--- a/benchmarks/clustering/scripts/cluster_plots.py
+++ b/benchmarks/clustering/scripts/cluster_plots.py
@@ -158,10 +158,14 @@ if __name__ == '__main__':
            print('{} do not have valid entries, skip it'.format(coll))
            continue
        df.rename(columns={c: c.replace(coll + '.', '') for c in df.columns}, inplace=True)
-        # print(df[['eta', 'polar.theta', 'position.x', 'position.y', 'position.z']])
+        df['r'] = np.sqrt(df['position.x'].values**2 + df['position.y'].values**2 + df['position.z'].values**2)
+        df['phi'] = np.arctan2(df['position.y'].values, df['position.x'].values)
+        df['theta'] = np.arccos(df['position.z'].values/df['r'].values)
+        df['eta'] = -np.log(np.tan(df['theta'].values/2.))
+        # print(df[['eta', 'theta', 'position.x', 'position.y', 'position.z']])
        fig, axs = plt.subplots(2, 2, figsize=(12, 8), dpi=160)
-        ncl = df.groupby('event')['ID.value'].nunique().values
-        axs[0][0].hist(ncl, weights=np.repeat(1./float(ncl.shape[0]), ncl.shape[0]),
+        ncl = df.groupby('event')['event'].count()
+        axs[0][0].hist(ncl, weights=np.repeat(1./float(ncl.count()), ncl.count()),
                bins=np.arange(0, 10), align='mid', ec='k')
        axs[0][0].set_xlabel('Number of Clusters', fontsize=16)


--- a/benchmarks/imaging_ecal/config.yml
+++ b/benchmarks/imaging_ecal/config.yml
@@ -4,7 +4,6 @@ imaging_ecal_electrons:
  stage: run
  script:
    - bash benchmarks/imaging_ecal/run_emcal_barrel.sh -t emcal_barrel_electrons -p "electron" -n 100
-  allow_failure: true

 imaging_ecal_photons:
  extends: .rec_benchmark
@@ -12,7 +11,6 @@ imaging_ecal_photons:
  stage: run
  script:
    - bash benchmarks/imaging_ecal/run_emcal_barrel.sh -t emcal_barrel_photons -p "photon" -n 100
-  allow_failure: true

 imaging_ecal_pions:
  extends: .rec_benchmark
@@ -20,7 +18,6 @@ imaging_ecal_pions:
  stage: run
  script:
    - bash benchmarks/imaging_ecal/run_emcal_barrel.sh -t emcal_barrel_pions -p "pion-" -n 100
-  allow_failure: true

 imaging_ecal_pion0:
  extends: .rec_benchmark
@@ -28,7 +25,6 @@ imaging_ecal_pion0:
  stage: run
  script:
    - bash benchmarks/imaging_ecal/run_imcal_pion0.sh -t imcal_barrel_pion0 -p "pion0" -n 100
-  allow_failure: true

 imaging_ecal_pion_rejection:
  extends: .rec_benchmark
@@ -38,7 +34,6 @@ imaging_ecal_pion_rejection:
    - compile_analyses.py imaging_ecal
    - bash benchmarks/imaging_ecal/run_emcal_barrel_pion_rej.sh -t emcal_barrel_pion_rej_electron -p "electron" -n 100
    - bash benchmarks/imaging_ecal/run_emcal_barrel_pion_rej.sh -t emcal_barrel_pion_rej_piminus -p "pion-" -n 100
-  allow_failure: true

 imaging_ecal_pion_rejection:bench:
  extends: .rec_benchmark
@@ -51,4 +46,3 @@ imaging_ecal_pion_rejection:bench:
    - compile_analyses.py imaging_ecal
    - root -b -q benchmarks/imaging_ecal/analysis/emcal_barrel_pion_rejection_analysis.cxx+
    #- bash run_pion_rej_analysis.sh
-  allow_failure: true
--- a/benchmarks/imaging_ecal/options/imaging_2dcluster.py
+++ b/benchmarks/imaging_ecal/options/imaging_2dcluster.py
@@ -80,7 +80,6 @@ imcal_barrel_clreco = RecoCoG('imcal_barrel_clreco',
        # OutputLevel=DEBUG,
        inputProtoClusterCollection = imcal_barrel_cl.outputProtoClusterCollection,
        outputClusterCollection='EcalBarrelImagingClusters',
-        mcHits="EcalBarrelHits",
        logWeightBase=6.2,
        samplingFraction=kwargs['img_sf'])


--- a/benchmarks/imaging_ecal/run_emcal_barrel.sh
+++ b/benchmarks/imaging_ecal/run_emcal_barrel.sh
@@ -142,13 +142,14 @@ for iev in "${ADDR[@]}"; do
    fi
 done

-python ${CB_EMCAL_SCRIPT_DIR}/energy_profile.py \
-    ${CB_EMCAL_REC_FILE} --type=EM --energy=${CB_EMCAL_ENERGY} -o results/${particle} \
-    --save=results/profile.csv --color=royalblue
-if [[ "$?" -ne "0" ]] ; then
-  echo "ERROR running analysis script: energy_profile"
-  exit 1
-fi
+# FIXME energy_profile disabled due to change in layer encoding, needs work
+#python ${CB_EMCAL_SCRIPT_DIR}/energy_profile.py \
+#    ${CB_EMCAL_REC_FILE} --type=EM --energy=${CB_EMCAL_ENERGY} -o results/${particle} \
+#    --save=results/profile.csv --color=royalblue
+#if [[ "$?" -ne "0" ]] ; then
+#  echo "ERROR running analysis script: energy_profile"
+#  exit 1
+#fi

 root_filesize=$(stat --format=%s "${CB_EMCAL_REC_FILE}")
 if [[ "${CB_EMCAL_NUMEV}" -lt "500" ]] ; then

--- a/benchmarks/imaging_ecal/run_emcal_barrel_pion_rej.sh
+++ b/benchmarks/imaging_ecal/run_emcal_barrel_pion_rej.sh
@@ -143,13 +143,14 @@ for iev in "${ADDR[@]}"; do
    fi
 done

-python ${CB_EMCAL_SCRIPT_DIR}/energy_profile.py \
-    ${CB_EMCAL_REC_FILE} --type=EM --energy=${CB_EMCAL_ENERGY} -o results/${particle} \
-    --save=results/profile.csv --color=royalblue
-if [[ "$?" -ne "0" ]] ; then
-  echo "ERROR running analysis script: energy_profile"
-  exit 1
-fi
+# FIXME energy_profile disabled due to change in layer encoding, needs work
+#python ${CB_EMCAL_SCRIPT_DIR}/energy_profile.py \
+#    ${CB_EMCAL_REC_FILE} --type=EM --energy=${CB_EMCAL_ENERGY} -o results/${particle} \
+#    --save=results/profile.csv --color=royalblue
+#if [[ "$?" -ne "0" ]] ; then
+#  echo "ERROR running analysis script: energy_profile"
+#  exit 1
+#fi

 echo "Reconstruction File"
 rootls -t "${CB_EMCAL_REC_FILE}"

--- a/benchmarks/imaging_ecal/scripts/draw_cluster.py
+++ b/benchmarks/imaging_ecal/scripts/draw_cluster.py
@@ -167,7 +167,7 @@ if __name__ == '__main__':
        # Calculate geometric variables of decaying particles
        dfdecaymcp['r'] = np.sqrt(dfdecaymcp['vex'].values**2 + dfdecaymcp['vey'].values**2 + dfdecaymcp['vez'].values**2)
        dfdecaymcp['phi'] = np.arctan2(dfdecaymcp['vey'].values, dfdecaymcp['vex'].values)*1000.
-        dfdecaymcp['theta'] = np.arccos(dfdecaymcp['vez'].values/dfdecaymcp['r'].values)*1000.
+        dfdecaymcp['theta'] = np.arctan2(dfdecaymcp['vez'].values, dfdecaymcp['r'].values)*1000.
        dfdecaymcp['eta'] = -np.log(np.tan(dfdecaymcp['theta'].values/1000./2.))

    # truth

--- a/benchmarks/imaging_ecal/scripts/draw_cluster_layers.py
+++ b/benchmarks/imaging_ecal/scripts/draw_cluster_layers.py
@@ -153,7 +153,7 @@ if __name__ == '__main__':
        # gif frames
        fig.savefig('ltmp.png')
        plt.close(fig)
-        frames.append(imageio.imread('ltmp.png'))
+        frames.append(imageio.v2.imread('ltmp.png'))
    pdf.close()
    os.remove('ltmp.png')

@@ -191,7 +191,7 @@ if __name__ == '__main__':
        # gif frames
        fig.savefig('ltmp.png')
        plt.close(fig)
-        frames.append(imageio.imread('ltmp.png'))
+        frames.append(imageio.v2.imread('ltmp.png'))
    pdf.close()
    os.remove('ltmp.png')


--- a/benchmarks/imaging_ecal/scripts/utils.py
+++ b/benchmarks/imaging_ecal/scripts/utils.py
@@ -126,8 +126,8 @@ def get_hits_data(path, evnums=None, branch='RecoEcalBarreImaginglHits'):
            continue

        events.GetEntry(iev)
-        for hit in getattr(events, branch):
-            dbuf[idb] = (iev, hit.ID.value, hit.layer, hit.position.x, hit.position.y,
+        for ihit, hit in enumerate(getattr(events, branch)):
+            dbuf[idb] = (iev, ihit, hit.layer, hit.position.x, hit.position.y,
                    hit.position.z, hit.energy*1000.)
            idb += 1

@@ -152,8 +152,8 @@ def get_layers_data(path, evnums=None, branch="EcalBarrelImagingClustersLayers")
            continue

        events.GetEntry(iev)
-        for layer in getattr(events, branch):
-            dbuf[idb] = (iev, layer.clusterID.value, layer.layer,
+        for icl, layer in enumerate(getattr(events, branch)):
+            dbuf[idb] = (iev, icl, layer.layer,
                         layer.position.x, layer.position.y, layer.position.z,
                         layer.energy*1000.)
            idb += 1