diff --git a/PhysicsAnalysis/NtupleDumper/CMakeLists.txt b/PhysicsAnalysis/NtupleDumper/CMakeLists.txt
index fab78641d7ae2a3842904122968a80170bca9a28..607c58d892d201731098f64829a82d556f1a302a 100644
--- a/PhysicsAnalysis/NtupleDumper/CMakeLists.txt
+++ b/PhysicsAnalysis/NtupleDumper/CMakeLists.txt
@@ -1,11 +1,14 @@
atlas_subdir( NtupleDumper )
+find_package( nlohmann_json )
+
atlas_add_component(
NtupleDumper
src/NtupleDumperAlg.h
src/NtupleDumperAlg.cxx
src/component/NtupleDumper_entries.cxx
LINK_LIBRARIES AthenaBaseComps StoreGateLib xAODFaserWaveform xAODFaserCalorimeter xAODFaserTrigger xAODFaserLHC ScintIdentifier FaserCaloIdentifier GeneratorObjects FaserActsGeometryLib TrackerSimEvent TrackerSimData TrackerIdentifier TrackerReadoutGeometry TrkTrack GeoPrimitives TrackerRIO_OnTrack TrackerSpacePoint FaserActsKalmanFilterLib
+ PRIVATE_LINK_LIBRARIES nlohmann_json::nlohmann_json
)
atlas_install_python_modules( python/*.py )
diff --git a/PhysicsAnalysis/NtupleDumper/scripts/faser_ntuple_maker.py b/PhysicsAnalysis/NtupleDumper/scripts/faser_ntuple_maker.py
index 8ea26ebb33fef5490110206cab8c940eecfa716b..1f4b990c2fb7bbdef7e109898e91bf5e0171d49d 100755
--- a/PhysicsAnalysis/NtupleDumper/scripts/faser_ntuple_maker.py
+++ b/PhysicsAnalysis/NtupleDumper/scripts/faser_ntuple_maker.py
@@ -55,6 +55,9 @@ parser.add_argument("--unblind", action='store_true',
help="Don't apply signal blinding (default: do)")
parser.add_argument("--onlyblind", action='store_true',
help="Only store events that were blinded (will override unblind arg)")
+
+parser.add_argument("--grl", default="",
+ help="Specify GRL to apply")
parser.add_argument("--fluka", action='store_true',
help="Add FLUKA weights to ntuple")
@@ -148,11 +151,6 @@ if filepath.is_dir():
print(f"First = {firstseg}")
print(f"Last = {lastseg}")
print(f"Tag = {args.tag}")
- print(f"Trigger Filter = {args.trigFilt}")
- print(f"Scintillator Filter = {args.scintFilt}")
- print(f"Track Filter = {not args.NoTrackFilt}")
- print(f"Blind = {not args.unblind}")
- print(f"OnlyBlinded = {args.onlyblind}")
# Find any tags
tagstr = firststem.replace(f"{firstfaser}-{firstshort}-{runstr}-{firstseg}", "")
@@ -212,6 +210,13 @@ if len(args.outfile) > 0:
print("Output file:")
print(outfile)
+print(f"Trigger Filter = {args.trigFilt}")
+print(f"Scintillator Filter = {args.scintFilt}")
+print(f"Track Filter = {not args.NoTrackFilt}")
+print(f"Blind = {not args.unblind}")
+print(f"OnlyBlinded = {args.onlyblind}")
+print(f"GRL = {args.grl}")
+
# OK, lets run the job here
from AthenaCommon.Logging import log, logging
from AthenaCommon.Constants import DEBUG, VERBOSE, INFO
@@ -245,18 +250,25 @@ ConfigFlags.lock()
acc = MainServicesCfg(ConfigFlags)
acc.merge(PoolReadCfg(ConfigFlags))
+# Create kwargs for NtupleDumper
+grl_kwargs = {}
+if args.grl:
+ grl_kwargs['GoodRunsList'] = args.grl
+ grl_kwargs['ApplyGoodRunsList'] = True
+
+mc_kwargs = {}
+if args.genie:
+ mc_kwargs['UseGenieWeights'] = True
+if args.fluka:
+ mc_kwargs['UseFlukaWeights'] = True
+
# algorithm
from NtupleDumper.NtupleDumperConfig import NtupleDumperAlgCfg
if args.isMC:
- if args.genie:
- acc.merge(NtupleDumperAlgCfg(ConfigFlags, outfile, UseGenieWeights=True))
- elif args.fluka:
- acc.merge(NtupleDumperAlgCfg(ConfigFlags, outfile, UseFlukaWeights=True))
- else:
- acc.merge(NtupleDumperAlgCfg(ConfigFlags, outfile))
+ acc.merge(NtupleDumperAlgCfg(ConfigFlags, outfile, **mc_kwargs))
else:
- acc.merge(NtupleDumperAlgCfg(ConfigFlags, outfile, DoBlinding=(not args.unblind), OnlyBlinded=args.onlyblind, DoScintFilter = args.scintFilt, DoTrackFilter = (not args.NoTrackFilt), DoTrigFilter = args.trigFilt))
+ acc.merge(NtupleDumperAlgCfg(ConfigFlags, outfile, DoBlinding=(not args.unblind), OnlyBlinded=args.onlyblind, DoScintFilter = args.scintFilt, DoTrackFilter = (not args.NoTrackFilt), DoTrigFilter = args.trigFilt, **grl_kwargs))
if not args.verbose:
from AthenaConfiguration.ComponentFactory import CompFactory
@@ -264,6 +276,10 @@ if not args.verbose:
AthenaEventLoopMgr.EventPrintoutInterval=1000
acc.addService(AthenaEventLoopMgr)
+else:
+ nd = acc.getEventAlgo("NtupleDumperAlg")
+ nd.OutputLevel = VERBOSE
+
# Hack to avoid problem with our use of MC databases when isMC = False
if not args.isMC:
replicaSvc = acc.getService("DBReplicaSvc")
diff --git a/PhysicsAnalysis/NtupleDumper/src/NtupleDumperAlg.cxx b/PhysicsAnalysis/NtupleDumper/src/NtupleDumperAlg.cxx
index 8b6aba7341e4c4c01aa752d8e6acbda13eb053cd..59a2be962c55923976f578c7e403cd2aae64a899 100644
--- a/PhysicsAnalysis/NtupleDumper/src/NtupleDumperAlg.cxx
+++ b/PhysicsAnalysis/NtupleDumper/src/NtupleDumperAlg.cxx
@@ -120,6 +120,19 @@ StatusCode NtupleDumperAlg::initialize()
ATH_CHECK(m_spacePointContainerKey.initialize());
+ // Read GRL if requested
+ if (!m_goodRunsList.value().empty()) {
+ ATH_MSG_INFO("Opening GRL file " << m_goodRunsList);
+ ATH_MSG_INFO(m_applyGoodRunsList);
+ std::ifstream f(m_goodRunsList);
+ m_grl = nlohmann::json::parse(f);
+ f.close();
+ ATH_MSG_INFO("Read GRL with " << m_grl.size() << " records");
+ } else {
+ // Make sure this is empty
+ m_grl.clear();
+ }
+
if (m_useFlukaWeights)
{
m_baseEventCrossSection = (m_flukaCrossSection * kfemtoBarnsPerMilliBarn)/m_flukaCollisions;
@@ -140,6 +153,7 @@ StatusCode NtupleDumperAlg::initialize()
m_tree->Branch("eventID", &m_event_number, "eventID/I");
m_tree->Branch("eventTime", &m_event_time, "eventTime/I");
m_tree->Branch("BCID", &m_bcid, "BCID/I");
+ m_tree->Branch("inGRL", &m_in_grl, "inGRL/I");
m_tree->Branch("fillNumber", &m_fillNumber, "fillNumber/I");
m_tree->Branch("betaStar", &m_betaStar, "betaStar/F");
@@ -417,6 +431,71 @@ StatusCode NtupleDumperAlg::execute(const EventContext &ctx) const
isMC = true;
}
+ // EventInfo data
+ m_run_number = ctx.eventID().run_number();
+ m_event_number = ctx.eventID().event_number();
+ m_event_time = ctx.eventID().time_stamp();
+ m_bcid = ctx.eventID().bunch_crossing_id();
+
+ // For real data, find if data is in GoodRunsList
+ if (!isMC && !m_grl.empty()) {
+ ATH_MSG_DEBUG("Checking GRL for run " << m_run_number << " event " << m_event_number);
+
+ // JSON keys are always strings
+ std::string run_string = std::to_string(m_run_number);
+ std::string reason("unknown");
+
+ if (m_grl.contains(run_string)) {
+ // m_in_grl set to 0 by clearTree
+ auto jrun = m_grl[run_string];
+
+ // Make sure this is in the stable beams range
+ if (jrun.contains("stable_list")) {
+ for ( auto& slist : jrun["stable_list"]) {
+ // Check if event falls in this range
+ if (slist["start_utime"] > m_event_time) continue;
+ if (slist["stop_utime"] <= m_event_time) continue;
+ m_in_grl = 1;
+ break;
+ }
+
+ // Check exclude list (not all runs will have this)
+ if (jrun.contains("excluded_list")) {
+ ATH_MSG_DEBUG("Checking excluded list");
+ for ( auto& blist : jrun["excluded_list"]) {
+ ATH_MSG_DEBUG("Checking excluded list start " << blist["start_utime"] << " stop " << blist["stop_utime"] << " event " << m_event_time);
+
+ // Check if event falls in this range
+ if (blist["start_utime"] > m_event_time) continue;
+ if (blist["stop_utime"] <= m_event_time) continue;
+ if (blist.contains("reason")) {
+ reason = blist["reason"];
+ }
+ m_in_grl = 0;
+ break;
+ }
+ }
+
+ if (m_in_grl == 0) {
+ ATH_MSG_DEBUG("Run " << run_string << " event " << m_event_number << " time " << m_event_time << " excluded for " << reason);
+ }
+
+ } else {
+ ATH_MSG_WARNING("Run " << run_string << " has no stable_list!");
+ }
+ } else {
+ ATH_MSG_WARNING("Run " << run_string << " not found in GRL!");
+ }
+ }
+
+ // Skip this event?
+ if (m_applyGoodRunsList && (m_in_grl == 0)) {
+ ATH_MSG_DEBUG("Skipping run " << m_run_number
+ << " event " << m_event_number << " - outside GRL" );
+ m_eventsFailedGRL += 1;
+ return StatusCode::SUCCESS;
+ }
+
// if real data, correct for diigitzer timing jitter with clock phase
if (!isMC) {
// correct waveform time with clock phase
@@ -577,11 +656,6 @@ StatusCode NtupleDumperAlg::execute(const EventContext &ctx) const
m_scintHit = m_scintHit | 128;
}
- m_run_number = ctx.eventID().run_number();
- m_event_number = ctx.eventID().event_number();
- m_event_time = ctx.eventID().time_stamp();
- m_bcid = ctx.eventID().bunch_crossing_id();
-
if (isMC) { // if simulation find MC cross section and primary lepton
// Work out effective cross section for MC
if (m_useFlukaWeights) {
@@ -1161,7 +1235,8 @@ StatusCode NtupleDumperAlg::execute(const EventContext &ctx) const
StatusCode NtupleDumperAlg::finalize()
{
- ATH_MSG_INFO("Number of events passed Ntuple selectioon = " << m_eventsPassed);
+ ATH_MSG_INFO("Number of events passed Ntuple selection = " << m_eventsPassed);
+ ATH_MSG_INFO("Number of events failing GRL selection = " << m_eventsFailedGRL);
return StatusCode::SUCCESS;
}
@@ -1173,6 +1248,7 @@ NtupleDumperAlg::clearTree() const
m_event_number = 0;
m_event_time = 0;
m_bcid = 0;
+ m_in_grl = 0;
m_fillNumber = 0;
m_betaStar = 0;
diff --git a/PhysicsAnalysis/NtupleDumper/src/NtupleDumperAlg.h b/PhysicsAnalysis/NtupleDumper/src/NtupleDumperAlg.h
index 322e38c90a68e917951a763b092c507efdb54c0a..339c0838b22047dd146d71163acebdc472f7d1f3 100644
--- a/PhysicsAnalysis/NtupleDumper/src/NtupleDumperAlg.h
+++ b/PhysicsAnalysis/NtupleDumper/src/NtupleDumperAlg.h
@@ -25,6 +25,7 @@
#include "StoreGate/ReadDecorHandle.h"
#include <vector>
+#include <nlohmann/json.hpp>
class TTree;
class TH1;
@@ -120,6 +121,12 @@ private:
DoubleProperty m_blindingCaloThreshold {this, "BlindingCaloThreshold", 100000.0, "Blind events with Ecal energy above threshold (in MeV)"};
DoubleProperty m_caloMC_FudgeFactor {this, "CaloFudgeFactorMC", 1.088, "Correct the MC energy calibration by this fudge factor"};
+ BooleanProperty m_applyGoodRunsList {this, "ApplyGoodRunsList", false, "Only write out events passing GRL (data only)"};
+ StringProperty m_goodRunsList {this, "GoodRunsList", "", "GoodRunsList in json format"};
+
+ // json object to hold data read from GRL file (or empty if not)
+ nlohmann::json m_grl;
+
double m_baseEventCrossSection {1.0};
const double kfemtoBarnsPerMilliBarn {1.0e12};
@@ -131,6 +138,7 @@ private:
mutable unsigned int m_event_number;
mutable unsigned int m_event_time;
mutable unsigned int m_bcid;
+ mutable unsigned int m_in_grl;
mutable unsigned int m_fillNumber;
mutable float m_betaStar;
@@ -344,7 +352,7 @@ private:
mutable double m_crossSection;
mutable int m_eventsPassed = 0;
-
+ mutable int m_eventsFailedGRL = 0;
};
inline const ServiceHandle <ITHistSvc> &NtupleDumperAlg::histSvc() const {