diff --git a/Calorimeter/CaloRecAlgs/src/CaloRecAlg.cxx b/Calorimeter/CaloRecAlgs/src/CaloRecAlg.cxx
index 12ddca3a44a5ba4691ff3feffa2adae37d19e2d5..f4f78d488d6c542637003e38ea0d8959f89b6b1c 100644
--- a/Calorimeter/CaloRecAlgs/src/CaloRecAlg.cxx
+++ b/Calorimeter/CaloRecAlgs/src/CaloRecAlg.cxx
@@ -23,6 +23,12 @@ StatusCode CaloRecAlg::initialize() {
// Initalize tools
ATH_CHECK( m_recoCalibTool.retrieve() );
+ // Store calibrattiion factos in a vector for ease of access
+ m_EM_mu_Map[0] = m_calo_ch0_EM_mu;
+ m_EM_mu_Map[1] = m_calo_ch1_EM_mu;
+ m_EM_mu_Map[2] = m_calo_ch2_EM_mu;
+ m_EM_mu_Map[3] = m_calo_ch3_EM_mu;
+
return StatusCode::SUCCESS;
}
//----------------------------------------------------------------------
@@ -81,31 +87,34 @@ StatusCode CaloRecAlg::execute(const EventContext& ctx) const {
float MIPcharge_ref = m_recoCalibTool->getMIPcharge_ref(hit->channel()); // get reference MIP charge from database
- float hit_Edep = 0;
+ float charge = hit->integral()/50.0; // divide by 50 ohms to get charge
+ ATH_MSG_DEBUG("calo_hit filled has charge of " << charge << " pC");
+
+ float gainRatio = 1.0;
if (m_isMC) {
- hit_Edep = hit->integral() / MIPcharge_ref; // MC MIPcharge_ref from database is really digitization scaling factors that take simulated Edep in MeV to mV*ns
- ATH_MSG_DEBUG("isMC == True: calo hit Edep = " << hit_Edep << " MeV");
+ gainRatio = 1.0; // put dummy value for now, this will end up being ratio of digi scale factors
} else {
- float hit_charge = hit->integral()/50.0; // divide by 50 ohms to get charge
- ATH_MSG_DEBUG("calo_hit filled has charge of " << hit_charge << " pC");
-
- float HVgain_extrap = extrapolateHVgain(hit->channel());
- ATH_MSG_DEBUG("HV gain = " << HVgain_extrap );
+ gainRatio = extrapolateHVgain(hit->channel());
+ ATH_MSG_DEBUG("HV gain ratio = " << gainRatio );
+ }
- float hit_MIP = (hit_charge * HVgain_extrap) / MIPcharge_ref;
- ATH_MSG_DEBUG("Edep in MIPs = " << hit_MIP );
+ float Nmip = (charge * gainRatio) / MIPcharge_ref;
+ ATH_MSG_DEBUG("Nmip = " << Nmip );
+ calo_hit->set_Nmip(Nmip); // set Nmip value
- hit_Edep = hit_MIP * m_MIP_sim_Edep_calo;
- ATH_MSG_DEBUG("Edep in MeV = " << hit_Edep );
- }
+ float E_dep = Nmip * m_MIP_sim_Edep_calo;
+ ATH_MSG_DEBUG("E_dep in MeV = " << E_dep );
+ calo_hit->set_E_dep(E_dep); // set calibrated E_dep value
- calo_hit->set_Edep(hit_Edep); // set calibrated Edep value
+ float E_EM = Nmip * m_EM_mu_Map[hit->channel()];
+ ATH_MSG_DEBUG("Em E in MeV = " << E_EM );
+ calo_hit->set_E_EM(E_EM); // set calibrated E_EM value
- float hit_raw_Edep = 0.0;
+ float fit_to_raw_ratio = 1.0;
if (hit->integral() != 0.0) { // avoid possibility of division by zero error
- hit_raw_Edep = hit_Edep * (hit->raw_integral() / hit->integral());
+ fit_to_raw_ratio = hit->raw_integral() / hit->integral();
}
- calo_hit->set_raw_Edep(hit_raw_Edep); // set calibrated raw_Edep value
+ calo_hit->set_fit_to_raw_ratio(fit_to_raw_ratio); // set fit-to-raw-ratio that can be used to take any of the calibrated values to what they would be if we used the raw integral instead of the fit integral
calo_hit->set_channel(hit->channel()); // set channel number
@@ -126,28 +135,22 @@ StatusCode CaloRecAlg::execute(const EventContext& ctx) const {
float MIPcharge_ref = m_recoCalibTool->getMIPcharge_ref(hit->channel()); // get reference MIP charge from database
- float hit_Edep = 0;
- if (m_isMC) {
- hit_Edep = hit->integral() / MIPcharge_ref; // MC MIPcharge_ref from database is really digitization scaling factors that take simulated Edep in MeV to mV*ns
- ATH_MSG_DEBUG("isMC == True: preshower hit Edep = " << hit_Edep << " MeV");
- } else {
- float hit_charge = hit->integral()/50.0; // divide by 50 ohms to get charge
- ATH_MSG_DEBUG("preshower_hit filled has charge of " << hit_charge << " pC");
+ float charge = hit->integral()/50.0; // divide by 50 ohms to get charge
+ ATH_MSG_DEBUG("preshower_hit filled has charge of " << charge << " pC");
- float hit_MIP = hit_charge / MIPcharge_ref;
- ATH_MSG_DEBUG("Edep in MIPs = " << hit_MIP );
-
- hit_Edep = hit_MIP * m_MIP_sim_Edep_preshower;
- ATH_MSG_DEBUG("Edep in GeV = " << hit_Edep );
- }
+ float Nmip = charge / MIPcharge_ref;
+ ATH_MSG_DEBUG("Nmip = " << Nmip );
+ preshower_hit->set_Nmip(Nmip); // set Nmip value
- preshower_hit->set_Edep(hit_Edep); // set calibrated Edep value
+ float E_dep = Nmip * m_MIP_sim_Edep_preshower;
+ ATH_MSG_DEBUG("E_dep in GeV = " << E_dep );
+ preshower_hit->set_E_dep(E_dep); // set calibrated E_dep value
- float hit_raw_Edep = 0.0;
+ float fit_to_raw_ratio = 1.0;
if (hit->integral() != 0.0) { // avoid possibility of division by zero error
- hit_raw_Edep = hit_Edep * (hit->raw_integral() / hit->integral());
+ fit_to_raw_ratio = hit->raw_integral() / hit->integral();
}
- preshower_hit->set_raw_Edep(hit_raw_Edep); // set calibrated rsw_Edep value
+ preshower_hit->set_fit_to_raw_ratio(fit_to_raw_ratio); // set fit-to-raw-ratio that can be used to take any of the calibrated values to what they would be if we used the raw integral instead of the fit integral
preshower_hit->set_channel(hit->channel()); // set channel number
diff --git a/Calorimeter/CaloRecAlgs/src/CaloRecAlg.h b/Calorimeter/CaloRecAlgs/src/CaloRecAlg.h
index c9ff4c5146d9e59c94a279345f6bb5977519c997..721f521621ae519778a8ddadd8b5711f1beeac9b 100644
--- a/Calorimeter/CaloRecAlgs/src/CaloRecAlg.h
+++ b/Calorimeter/CaloRecAlgs/src/CaloRecAlg.h
@@ -88,6 +88,13 @@ class CaloRecAlg : public AthReentrantAlgorithm {
FloatProperty m_MIP_sim_Edep_calo {this, "MIP_sim_Edep_calo", 58.5}; // MIP deposits 5.85 MeV of energy in calo
FloatProperty m_MIP_sim_Edep_preshower {this, "MIP_sim_Edep_preshower", 4.894}; // MIP deposits 4.894 MeV of energy in a preshower layer
+ FloatProperty m_calo_ch0_EM_mu {this, "m_calo_ch0_EM_mu", 330.0}; // factor used to do rough calibration of calo ch0 to EM energy: 0.33 GeV or 330 MeV
+ FloatProperty m_calo_ch1_EM_mu {this, "m_calo_ch1_EM_mu", 330.0}; // factor used to do rough calibration of calo ch1 to EM energy: 0.33 GeV or 330 MeV
+ FloatProperty m_calo_ch2_EM_mu {this, "m_calo_ch2_EM_mu", 330.0}; // factor used to do rough calibration of calo ch2 to EM energy: 0.33 GeV or 330 MeV
+ FloatProperty m_calo_ch3_EM_mu {this, "m_calo_ch3_EM_mu", 330.0}; // factor used to do rough calibration of calo ch3 to EM energy: 0.33 GeV or 330 MeV
+
+ float m_EM_mu_Map[4]; // vector that holds EM_mu calibration factors for calo channels
+
Gaudi::Property<bool> m_isMC {this, "isMC", false};
};
diff --git a/Calorimeter/CaloRecTools/python/CaloRecToolConfig.py b/Calorimeter/CaloRecTools/python/CaloRecToolConfig.py
index e87ac3658ec8f92caf00bc19e061d3d9d95c8e5d..65be875da12b6deee3cbf9f35217cee6f40784ff 100644
--- a/Calorimeter/CaloRecTools/python/CaloRecToolConfig.py
+++ b/Calorimeter/CaloRecTools/python/CaloRecToolConfig.py
@@ -20,8 +20,7 @@ def CaloRecToolCfg(flags, name="CaloRecTool", **kwargs):
dbInstance = kwargs.get("dbInstance", "TRIGGER_OFL")
acc.merge(addFolders(flags, "/WAVE/Calibration/HV", dbInstance, className="CondAttrListCollection"))
- acc.merge(addFolders(flags, "/WAVE/Calibration/HV_ref", dbInstance, className="CondAttrListCollection"))
- acc.merge(addFolders(flags, "/WAVE/Calibration/MIP_charge_ref", dbInstance, className="CondAttrListCollection"))
+ acc.merge(addFolders(flags, "/WAVE/Calibration/MIP_ref", dbInstance, className="CondAttrListCollection"))
return acc
diff --git a/Calorimeter/CaloRecTools/src/CaloRecTool.cxx b/Calorimeter/CaloRecTools/src/CaloRecTool.cxx
index a6bb42c5140d10c1d2efdc486201c2bab81c0365..916011019805ab9087c57427859086dec456eb80 100644
--- a/Calorimeter/CaloRecTools/src/CaloRecTool.cxx
+++ b/Calorimeter/CaloRecTools/src/CaloRecTool.cxx
@@ -23,8 +23,7 @@ CaloRecTool::initialize() {
// Set keys to read calibratiion variables from data base
ATH_CHECK( m_PMT_HV_ReadKey.initialize() );
- ATH_CHECK( m_PMT_HV_ref_ReadKey.initialize() );
- ATH_CHECK( m_MIPcharge_ref_ReadKey.initialize() );
+ ATH_CHECK( m_MIP_ref_ReadKey.initialize() );
// access and store calo PMT HV gain curves
HVgaincurves_rootFile = TFile::Open(m_PMT_HV_Gain_Curve_file.value().c_str(),"read"); // open root file that has TF1 gain curves stored
@@ -112,7 +111,7 @@ float CaloRecTool::getHV_ref(const EventContext& ctx, int channel) const {
float HV_ref=0.;
// Read Cond Handle
- SG::ReadCondHandle<CondAttrListCollection> readHandle{m_PMT_HV_ref_ReadKey, ctx};
+ SG::ReadCondHandle<CondAttrListCollection> readHandle{m_MIP_ref_ReadKey, ctx};
const CondAttrListCollection* readCdo{*readHandle};
if (readCdo==nullptr) {
ATH_MSG_FATAL("Null pointer to the read conditions object");
@@ -154,7 +153,7 @@ float CaloRecTool::getMIPcharge_ref(const EventContext& ctx, int channel) const
float MIP_charge_ref =0.;
// Read Cond Handle
- SG::ReadCondHandle<CondAttrListCollection> readHandle{m_MIPcharge_ref_ReadKey, ctx};
+ SG::ReadCondHandle<CondAttrListCollection> readHandle{m_MIP_ref_ReadKey, ctx};
const CondAttrListCollection* readCdo{*readHandle};
if (readCdo==nullptr) {
ATH_MSG_FATAL("Null pointer to the read conditions object");
@@ -172,8 +171,8 @@ float CaloRecTool::getMIPcharge_ref(const EventContext& ctx, int channel) const
// Read offset for specific channel
const CondAttrListCollection::AttributeList& payload{readCdo->attributeList(channel)};
- if (payload.exists("MIP_charge_ref") and not payload["MIP_charge_ref"].isNull()) {
- MIP_charge_ref = payload["MIP_charge_ref"].data<float>();
+ if (payload.exists("charge_ref") and not payload["charge_ref"].isNull()) {
+ MIP_charge_ref = payload["charge_ref"].data<float>();
ATH_MSG_DEBUG("Found digitizer channel " << channel << ", MIP_charge_ref as " << MIP_charge_ref);
} else {
ATH_MSG_WARNING("No valid MIP_charge_ref found for channel "<<channel<<"!");
diff --git a/Calorimeter/CaloRecTools/src/CaloRecTool.h b/Calorimeter/CaloRecTools/src/CaloRecTool.h
index c6dc951818cb59472fe6a1ff077e97e6753ec509..6ab891ad42784e3388d933f42de9267545fab3d9 100644
--- a/Calorimeter/CaloRecTools/src/CaloRecTool.h
+++ b/Calorimeter/CaloRecTools/src/CaloRecTool.h
@@ -78,8 +78,7 @@ class CaloRecTool: public extends<AthAlgTool, ICaloRecTool> {
// Read Cond Handle
SG::ReadCondHandleKey<CondAttrListCollection> m_PMT_HV_ReadKey{this, "PMT_HV_ReadKey", "/WAVE/Calibration/HV", "Key of folder for PMT HV reading"};
- SG::ReadCondHandleKey<CondAttrListCollection> m_PMT_HV_ref_ReadKey{this, "PMT_HV_ref_ReadKey", "/WAVE/Calibration/HV_ref", "Key of folder for PMT HV at time of MIP calibration measurement"};
- SG::ReadCondHandleKey<CondAttrListCollection> m_MIPcharge_ref_ReadKey{this, "MIPcharge_ref_ReadKey", "/WAVE/Calibration/MIP_charge_ref", "Key of folder for MIPcharge calibration measurment"};
+ SG::ReadCondHandleKey<CondAttrListCollection> m_MIP_ref_ReadKey{this, "MIP_ref_ReadKey", "/WAVE/Calibration/MIP_ref", "Key of folder for MIP charge calibration measurment, also stores PMT HV used to measure the reference MIP charge"};
};
diff --git a/PhysicsAnalysis/NeutrinoSearch/CMakeLists.txt b/PhysicsAnalysis/NeutrinoSearch/CMakeLists.txt
index cff2ab8575553261c18eb0736c211528e44f22bb..29ea0c53b37e81e5cf4b4ec63e958ac0a344f299 100644
--- a/PhysicsAnalysis/NeutrinoSearch/CMakeLists.txt
+++ b/PhysicsAnalysis/NeutrinoSearch/CMakeLists.txt
@@ -5,7 +5,7 @@ atlas_add_component(
src/NeutrinoSearchAlg.h
src/NeutrinoSearchAlg.cxx
src/component/NeutrinoSearch_entries.cxx
- LINK_LIBRARIES AthenaBaseComps StoreGateLib xAODFaserWaveform ScintIdentifier FaserCaloIdentifier GeneratorObjects FaserActsGeometryLib TrackerSimEvent TrackerSimData TrackerIdentifier TrackerReadoutGeometry TrkTrack GeoPrimitives TrackerRIO_OnTrack
+ LINK_LIBRARIES AthenaBaseComps StoreGateLib xAODFaserWaveform ScintIdentifier ScintSimEvent FaserCaloIdentifier GeneratorObjects FaserActsGeometryLib TrackerSimEvent TrackerSimData TrackerIdentifier TrackerReadoutGeometry TrkTrack GeoPrimitives TrackerRIO_OnTrack
)
atlas_install_python_modules(python/*.py)
diff --git a/PhysicsAnalysis/NeutrinoSearch/python/FlukaSearchConfig.py b/PhysicsAnalysis/NeutrinoSearch/python/FlukaSearchConfig.py
index b26b82df91d117c889fead8e5630843620159488..4e81f6bcf7ba871622430d5ae246045f9fc043f6 100644
--- a/PhysicsAnalysis/NeutrinoSearch/python/FlukaSearchConfig.py
+++ b/PhysicsAnalysis/NeutrinoSearch/python/FlukaSearchConfig.py
@@ -47,50 +47,50 @@ if __name__ == "__main__":
# Configure
ConfigFlags.Input.Files = [
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00000-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00001-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00002-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00003-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00004-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00005-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00006-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00007-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00008-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00009-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00010-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00011-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00012-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00013-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00014-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00015-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00016-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00017-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00018-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00019-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00020-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00021-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00022-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00023-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00024-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00025-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00026-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210001/rec/r0008/FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00027-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00000-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00001-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00002-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00003-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00004-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00005-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00006-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00007-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00008-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00009-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00010-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00011-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00012-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00013-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00014-s0005-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/fluka/210002/rec/r0008/FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00015-s0005-r0008-xAOD.root'
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00011-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00003-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00001-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00015-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00010-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00014-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00009-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00000-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00008-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00012-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00002-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00004-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00013-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00007-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00005-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210002/rec/r0009/./FaserMC-MDC_Fluka_unit30_Pm_71m_m3750_v3-210002-00006-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00017-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00020-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00010-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00012-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00001-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00018-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00004-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00019-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00007-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00003-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00008-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00005-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00022-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00011-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00006-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00013-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00021-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00015-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00026-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00027-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00014-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00024-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00000-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00009-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00025-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00016-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00023-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/fluka/210001/rec/r0009/./FaserMC-MDC_Fluka_unit30_Nm_71m_m3750_v3-210001-00002-s0007-r0009-xAOD.root'
]
# Update this for samples with new field map
ConfigFlags.IOVDb.GlobalTag = "OFLCOND-FASER-02" # Always needed; must match FaserVersionS
diff --git a/PhysicsAnalysis/NeutrinoSearch/python/GenieSearchConfig.py b/PhysicsAnalysis/NeutrinoSearch/python/GenieSearchConfig.py
index a80de32bbc0643dbaca1ac12e28688535951b709..94d860f75196b82b1c5ca3cdd6b3c84e1e06f171 100644
--- a/PhysicsAnalysis/NeutrinoSearch/python/GenieSearchConfig.py
+++ b/PhysicsAnalysis/NeutrinoSearch/python/GenieSearchConfig.py
@@ -47,14 +47,27 @@ if __name__ == "__main__":
# Configure
ConfigFlags.Input.Files = [
- '/run/media/dcasper/Data/faser/genie/rec/FaserMC-MDC_Genie_all_150invfb_v1-200001-00000-00007-s0006-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/genie/rec/FaserMC-MDC_Genie_all_150invfb_v1-200001-00008-00015-s0006-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/genie/rec/FaserMC-MDC_Genie_all_150invfb_v1-200001-00016-00023-s0006-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/genie/rec/FaserMC-MDC_Genie_all_150invfb_v1-200001-00024-00031-s0006-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/genie/rec/FaserMC-MDC_Genie_all_150invfb_v1-200001-00032-00039-s0006-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/genie/rec/FaserMC-MDC_Genie_all_150invfb_v1-200001-00040-00047-s0006-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/genie/rec/FaserMC-MDC_Genie_all_150invfb_v1-200001-00048-00055-s0006-r0008-xAOD.root',
- '/run/media/dcasper/Data/faser/genie/rec/FaserMC-MDC_Genie_all_150invfb_v1-200001-00056-00063-s0006-r0008-xAOD.root'
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00040-00047-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00056-00063-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00120-00127-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00112-00119-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00000-00007-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00208-00210-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00200-00207-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00064-00071-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00080-00087-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00024-00031-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00072-00079-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00008-00015-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00176-00183-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00096-00103-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00104-00111-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00032-00039-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00048-00055-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00192-00199-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00088-00095-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00016-00023-s0007-r0009-xAOD.root',
+ '/run/media/dcasper/Data/faser/genie/r0009/rec/./FaserMC-MDC_Genie_all_600fbInv_v1-200003-00168-00175-s0007-r0009-xAOD.root'
]
ConfigFlags.IOVDb.GlobalTag = "OFLCOND-FASER-02" # Always needed; must match FaserVersionS
ConfigFlags.IOVDb.DatabaseInstance = "OFLP200" # Use MC conditions for now
diff --git a/PhysicsAnalysis/NeutrinoSearch/src/NeutrinoSearchAlg.cxx b/PhysicsAnalysis/NeutrinoSearch/src/NeutrinoSearchAlg.cxx
index e782dda9e3eb7d97da42edb1b5a79528268dac9f..e3d33d683e0b47d064216515cf83fb280cb9edd4 100644
--- a/PhysicsAnalysis/NeutrinoSearch/src/NeutrinoSearchAlg.cxx
+++ b/PhysicsAnalysis/NeutrinoSearch/src/NeutrinoSearchAlg.cxx
@@ -28,6 +28,8 @@ StatusCode NeutrinoSearchAlg::initialize()
{
ATH_CHECK(m_truthEventContainer.initialize());
ATH_CHECK(m_truthParticleContainer.initialize());
+ ATH_CHECK( m_vetoNuHitKey.initialize() );
+
ATH_CHECK(m_trackCollection.initialize());
ATH_CHECK(m_vetoNuContainer.initialize());
ATH_CHECK(m_vetoContainer.initialize());
@@ -46,6 +48,7 @@ StatusCode NeutrinoSearchAlg::initialize()
ATH_CHECK(detStore()->retrieve(m_detMgr, "SCT"));
ATH_CHECK(m_extrapolationTool.retrieve());
+ ATH_CHECK(m_trackingGeometryTool.retrieve());
if (m_useFlukaWeights)
{
@@ -64,6 +67,9 @@ StatusCode NeutrinoSearchAlg::initialize()
m_tree->Branch("run_number", &m_run_number, "run_number/I");
m_tree->Branch("event_number", &m_event_number, "event_number/I");
+// m_tree->Branch("VetoNuTruthX", &m_vetoNuHitsMeanX, "vetoNuTruthX/D");
+// m_tree->Branch("VetoNuTruthY", &m_vetoNuHitsMeanY, "vetoNuTruthY/D");
+
m_tree->Branch("VetoNuPmt0", &m_vetoNu0, "vetoNu0/D");
m_tree->Branch("VetoNuPmt1", &m_vetoNu1, "vetoNu1/D");
@@ -106,9 +112,17 @@ StatusCode NeutrinoSearchAlg::initialize()
m_tree->Branch("ndof", &m_ndof, "ndof/I");
m_tree->Branch("longTracks", &m_longTracks, "longTracks/I");
m_tree->Branch("pTruthLepton", &m_truthLeptonMomentum, "pTruthLepton/D");
+ m_tree->Branch("xTruthLepton", &m_truthLeptonX, "xTruthLepton/D");
+ m_tree->Branch("yTruthLepton", &m_truthLeptonY, "yTruthLepton/D");
m_tree->Branch("truthBarcode", &m_truthBarcode, "truthBarcode/I");
m_tree->Branch("truthPdg", &m_truthPdg, "truthPdg/I");
m_tree->Branch("CrossSection", &m_crossSection, "crossSection/D");
+ m_tree->Branch("xVetoNu", &m_xVetoNu, "xVetoNu/D");
+ m_tree->Branch("yVetoNu", &m_yVetoNu, "yVetoNu/D");
+ m_tree->Branch("sxVetoNu", &m_sxVetoNu, "sxVetoNu/D");
+ m_tree->Branch("syVetoNu", &m_syVetoNu, "syVetoNu/D");
+ m_tree->Branch("thetaxVetoNu", &m_thetaxVetoNu, "thetaxVetoNu/D");
+ m_tree->Branch("thetayVetoNu", &m_thetayVetoNu, "thetayVetoNu/D");
ATH_CHECK(histSvc()->regTree("/HIST2/tree", m_tree));
@@ -145,6 +159,23 @@ StatusCode NeutrinoSearchAlg::execute(const EventContext &ctx) const
double flukaWeight = truthEventContainer->at(0)->weights()[0];
ATH_MSG_ALWAYS("Found fluka weight = " << flukaWeight);
m_crossSection = m_baseEventCrossSection * flukaWeight;
+
+ // // Find crossing position of VetoNu so we can check extrapolation
+ // SG::ReadHandle<ScintHitCollection> h_vetoNuHits {m_vetoNuHitKey, ctx};
+ // if (h_vetoNuHits.isValid() && h_vetoNuHits->size()!=0){
+ // double totalEnergyLoss {0.0};
+ // for (const ScintHit& hit : *h_vetoNuHits)
+ // {
+ // m_vetoNuHitsMeanX += hit.energyLoss() * (hit.localStartPosition().x() + hit.localEndPosition().x()) / 2;
+ // m_vetoNuHitsMeanY += hit.energyLoss() * (hit.localStartPosition().y() + hit.localEndPosition().y()) / 2;
+ // totalEnergyLoss += hit.energyLoss();
+ // }
+ // if (totalEnergyLoss > 0)
+ // {
+ // m_vetoNuHitsMeanX /= totalEnergyLoss;
+ // m_vetoNuHitsMeanY /= totalEnergyLoss;
+ // }
+ // }
}
else if (m_useGenieWeights)
{
@@ -166,8 +197,14 @@ StatusCode NeutrinoSearchAlg::execute(const EventContext &ctx) const
if ( particle->absPdgId() == 11 || particle->absPdgId() == 13 || particle->absPdgId() == 15 )
{
if (particle->status() == 1 && (particle->nParents() == 0 || particle->nParents() == 2) )
+ {
m_truthLeptonMomentum = particle->p4().P();
-
+ double deltaZ = m_zExtrapolate - particle->prodVtx()->z();
+ double thetaX = particle->px()/particle->pz();
+ double thetaY = particle->py()/particle->pz();
+ m_truthLeptonX = particle->prodVtx()->x() + thetaX * deltaZ;
+ m_truthLeptonY = particle->prodVtx()->y() + thetaY * deltaZ;
+ }
break;
}
}
@@ -435,6 +472,7 @@ StatusCode NeutrinoSearchAlg::execute(const EventContext &ctx) const
const Trk::TrackParameters* upstreamParameters {nullptr};
for (auto params : *(track->trackParameters()))
{
+ if (params->position().z() < 0) continue; // Ignore IFT hits
if (upstreamParameters == nullptr || params->position().z() < upstreamParameters->position().z()) upstreamParameters = params;
}
if (candidateParameters == nullptr || upstreamParameters->momentum().mag() > candidateParameters->momentum().mag())
@@ -457,98 +495,98 @@ StatusCode NeutrinoSearchAlg::execute(const EventContext &ctx) const
if (cluster != nullptr)
{
// ATH_MSG_INFO("ClusterOnTrack is OK");
- cluster->dump(msg());
+ // cluster->dump(msg());
// Hack to work around issue with cluster navigation
- auto idRDO = cluster->identify();
+ // auto idRDO = cluster->identify();
+
+ // if (simDataCollection->count(idRDO) > 0)
+ // {
+ // // ATH_MSG_INFO("rdo entry found");
+ // const auto& simdata = simDataCollection->find(idRDO)->second;
+ // const auto& deposits = simdata.getdeposits();
+ // //loop through deposits and record contributions
+ // HepMcParticleLink primary{};
+ // for( const auto& depositPair : deposits)
+ // {
+ // // ATH_MSG_INFO("Deposit found");
+ // float eDep = depositPair.second;
+ // int barcode = depositPair.first->barcode();
+ // // if( depositPair.second > highestDep)
+ // // {
+ // // highestDep = depositPair.second;
+ // // barcode = depositPair.first->barcode();
+ // // primary = depositPair.first;
+ // // depositPair.first->print(std::cout);
+ // // ATH_MSG_INFO("pdg id "<<depositPair.first->pdg_id());
+ // // }
+ // if (truthMap.count(barcode) > 0)
+ // {
+ // truthMap[barcode] += eDep;
+ // }
+ // else
+ // {
+ // truthMap[barcode] = eDep;
+ // }
+ // }
+ // }
+
+
+
- if (simDataCollection->count(idRDO) > 0)
- {
- // ATH_MSG_INFO("rdo entry found");
- const auto& simdata = simDataCollection->find(idRDO)->second;
- const auto& deposits = simdata.getdeposits();
- //loop through deposits and record contributions
- HepMcParticleLink primary{};
- for( const auto& depositPair : deposits)
+ // const Tracker::FaserSCT_Cluster* origCluster = dynamic_cast<const Tracker::FaserSCT_Cluster*>(cluster->prepRawData());
+ auto origCluster = cluster->prepRawData();
+ if (origCluster != nullptr)
{
- // ATH_MSG_INFO("Deposit found");
- float eDep = depositPair.second;
- int barcode = depositPair.first->barcode();
- // if( depositPair.second > highestDep)
- // {
- // highestDep = depositPair.second;
- // barcode = depositPair.first->barcode();
- // primary = depositPair.first;
- // depositPair.first->print(std::cout);
- // ATH_MSG_INFO("pdg id "<<depositPair.first->pdg_id());
- // }
- if (truthMap.count(barcode) > 0)
+ // ATH_MSG_INFO("Orig Cluster is OK");
+ auto rdoList = origCluster->rdoList();
+ for (auto idRDO : rdoList)
{
- truthMap[barcode] += eDep;
- }
- else
- {
- truthMap[barcode] = eDep;
+ // ATH_MSG_INFO("rdoList not empty");
+ if (simDataCollection->count(idRDO) > 0)
+ {
+ // ATH_MSG_INFO("rdo entry found");
+ const auto& simdata = simDataCollection->find(idRDO)->second;
+ const auto& deposits = simdata.getdeposits();
+ //loop through deposits and record contributions
+ HepMcParticleLink primary{};
+ for( const auto& depositPair : deposits)
+ {
+ // ATH_MSG_INFO("Deposit found");
+ float eDep = depositPair.second;
+ int barcode = depositPair.first->barcode();
+ // if( depositPair.second > highestDep)
+ // {
+ // highestDep = depositPair.second;
+ // barcode = depositPair.first->barcode();
+ // primary = depositPair.first;
+ // depositPair.first->print(std::cout);
+ // ATH_MSG_INFO("pdg id "<<depositPair.first->pdg_id());
+ // }
+ if (truthMap.count(barcode) > 0)
+ {
+ truthMap[barcode] += eDep;
+ }
+ else
+ {
+ truthMap[barcode] = eDep;
+ }
+ }
+ }
}
}
- }
-
-
-
-
- // // const Tracker::FaserSCT_Cluster* origCluster = dynamic_cast<const Tracker::FaserSCT_Cluster*>(cluster->prepRawData());
- // auto origCluster = cluster->prepRawData();
- // if (origCluster != nullptr)
- // {
- // ATH_MSG_INFO("Orig Cluster is OK");
- // auto rdoList = origCluster->rdoList();
- // for (auto idRDO : rdoList)
- // {
- // ATH_MSG_INFO("rdoList not empty");
- // if (simDataCollection->count(idRDO) > 0)
- // {
- // ATH_MSG_INFO("rdo entry found");
- // const auto& simdata = simDataCollection->find(idRDO)->second;
- // const auto& deposits = simdata.getdeposits();
- // //loop through deposits and record contributions
- // HepMcParticleLink primary{};
- // for( const auto& depositPair : deposits)
- // {
- // ATH_MSG_INFO("Deposit found");
- // float eDep = depositPair.second;
- // int barcode = depositPair.first->barcode();
- // // if( depositPair.second > highestDep)
- // // {
- // // highestDep = depositPair.second;
- // // barcode = depositPair.first->barcode();
- // // primary = depositPair.first;
- // // depositPair.first->print(std::cout);
- // // ATH_MSG_INFO("pdg id "<<depositPair.first->pdg_id());
- // // }
- // if (truthMap.count(barcode) > 0)
- // {
- // truthMap[barcode] += eDep;
- // }
- // else
- // {
- // truthMap[barcode] = eDep;
- // }
- // }
- // }
- // }
- // }
}
}
std::vector<std::pair<int, float>> truth(truthMap.begin(), truthMap.end());
std::sort(truth.begin(), truth.end(), [](auto v1, auto v2) { return v1.second > v2.second; });
- if (truth.size()>0) ATH_MSG_ALWAYS("Selected track truth info:");
+ // if (truth.size()>0) ATH_MSG_ALWAYS("Selected track truth info:");
for (auto v : truth)
{
auto truthParticle = (*(std::find_if(truthParticleContainer->cbegin(), truthParticleContainer->cend(), [v](const xAOD::TruthParticle* p){ return p->barcode() == v.first; })));
if (m_truthPdg == 0) m_truthPdg = truthParticle->pdgId();
if (m_truthBarcode == 0) m_truthBarcode = v.first;
- ATH_MSG_ALWAYS("truth info: barcode = " << v.first << " ( " << truthParticle->p4().P()/1000 << " GeV/c, Id code = " << truthParticle->pdgId() << ") -> deposited energy: " << v.second/1000);
+ // ATH_MSG_ALWAYS("truth info: barcode = " << v.first << " ( " << truthParticle->p4().P()/1000 << " GeV/c, Id code = " << truthParticle->pdgId() << ") -> deposited energy: " << v.second/1000);
}
}
@@ -562,6 +600,63 @@ StatusCode NeutrinoSearchAlg::execute(const EventContext &ctx) const
m_pz = candidateParameters->momentum().z();
m_p = sqrt(m_px * m_px + m_py * m_py + m_pz * m_pz);
m_charge = (int) candidateParameters->charge();
+
+ FaserActsGeometryContext faserGeometryContext = m_trackingGeometryTool->getNominalGeometryContext();
+ auto gctx = faserGeometryContext.context();
+
+
+ Amg::Vector3D position = candidateParameters->position();
+ Amg::Vector3D momentum = candidateParameters->momentum();
+ auto covariance = *candidateParameters->covariance();
+ Acts::BoundVector params = Acts::BoundVector::Zero();
+ params[Acts::eBoundLoc0] = -position.y();
+ params[Acts::eBoundLoc1] = position.x();
+ params[Acts::eBoundPhi] = momentum.phi();
+ params[Acts::eBoundTheta] = momentum.theta();
+ params[Acts::eBoundQOverP] = candidateParameters->charge() / momentum.mag();
+ params[Acts::eBoundTime] = 0;
+
+ using namespace Acts::UnitLiterals;
+ std::optional<AmgSymMatrix(6)> cov = std::nullopt;
+ Acts::BoundMatrix newCov = Acts::BoundMatrix::Zero();
+ for (size_t i = 0; i < 5; i++)
+ for (size_t j = 0; j < 5; j++)
+ newCov(i, j) = covariance(i, j);
+ // ATH_MSG_ALWAYS("Covariance: " << covariance(1,1) << " newCov: " << newCov(1,1));
+ // Convert the covariance matrix from GeV
+ for(int i=0; i < newCov.rows(); i++){
+ newCov(i, 4) = newCov(i, 4)/1_MeV;
+ }
+ for(int i=0; i < newCov.cols(); i++){
+ newCov(4, i) = newCov(4, i)/1_MeV;
+ }
+ cov = std::optional<AmgSymMatrix(6)>(newCov);
+ // ATH_MSG_ALWAYS("cov: " << cov.value()(1,1));
+ auto startSurface = Acts::Surface::makeShared<Acts::PlaneSurface>(Acts::Vector3(0, 0, position.z()), Acts::Vector3(0, 0, 1));
+ Acts::BoundTrackParameters startParameters(std::move(startSurface), params, candidateParameters->charge(), cov);
+ auto targetSurface_VetoNu = Acts::Surface::makeShared<Acts::PlaneSurface>(Acts::Vector3(0, 0, m_zExtrapolate), Acts::Vector3(0, 0, 1)); // -3112 mm is z position of VetoNu planes touching
+ std::unique_ptr<const Acts::BoundTrackParameters> targetParameters_VetoNu =m_extrapolationTool->propagate(ctx, startParameters, *targetSurface_VetoNu, Acts::backward);
+ if (targetParameters_VetoNu != nullptr)
+ {
+ auto targetPosition_VetoNu = targetParameters_VetoNu->position(gctx);
+ auto targetMomentum_VetoNu = targetParameters_VetoNu->momentum();
+ // Acts::BoundSymMatrix targetCovariance_VetoNu { targetParameters_VetoNu->covariance() };
+ m_xVetoNu = targetPosition_VetoNu.x();
+ m_yVetoNu = targetPosition_VetoNu.y();
+ if (targetParameters_VetoNu->covariance().has_value())
+ {
+ auto targetCovariance_VetoNu = targetParameters_VetoNu->covariance().value() ;
+ m_sxVetoNu = sqrt(targetCovariance_VetoNu(Acts::eBoundLoc0, Acts::eBoundLoc0));
+ m_syVetoNu = sqrt(targetCovariance_VetoNu(Acts::eBoundLoc1, Acts::eBoundLoc1));
+ }
+ m_thetaxVetoNu = atan(targetMomentum_VetoNu[0]/targetMomentum_VetoNu[2]);
+ m_thetayVetoNu = atan(targetMomentum_VetoNu[1]/targetMomentum_VetoNu[2]);
+ ATH_MSG_INFO("vetoNu good targetParameters xV,yV = (" << m_xVetoNu << ", " << m_yVetoNu << ")");
+ }
+ else
+ {
+ ATH_MSG_INFO("vetoNu null targetParameters with p = " << momentum.mag() << ", theta = " << momentum.theta() << ", x,y = (" << position.x() << ", " << position.y() << ")");
+ }
}
// Here we apply the signal selection
@@ -571,7 +666,18 @@ StatusCode NeutrinoSearchAlg::execute(const EventContext &ctx) const
m_preshower0 == 0 || m_preshower1 == 0 ||
// m_ecalTotal == 0 ||
candidateParameters == nullptr)
- return StatusCode::SUCCESS;
+ {
+ // if (m_longTracks > 0)
+ // {
+ // ATH_MSG_ALWAYS(m_vetoUpstream);
+ // ATH_MSG_ALWAYS(m_vetoDownstream);
+ // ATH_MSG_ALWAYS(m_triggerTotal);
+ // ATH_MSG_ALWAYS(m_preshower0);
+ // ATH_MSG_ALWAYS(m_preshower1);
+ // ATH_MSG_ALWAYS((candidateParameters == nullptr ? "Null" : "Not Null"));
+ // }
+ return StatusCode::SUCCESS;
+ }
m_tree->Fill();
@@ -595,6 +701,8 @@ NeutrinoSearchAlg::clearTree() const
{
m_run_number = 0;
m_event_number = 0;
+// m_vetoNuHitsMeanX = 0;
+// m_vetoNuHitsMeanY = 0;
m_vetoNu0 = 0;
m_vetoNu1 = 0;
m_veto00 = 0;
@@ -631,7 +739,15 @@ NeutrinoSearchAlg::clearTree() const
m_y = 0;
m_z = 0;
m_longTracks = 0;
+ m_truthLeptonX = 0;
+ m_truthLeptonY = 0;
m_truthLeptonMomentum = 0;
m_truthBarcode = 0;
m_truthPdg = 0;
+ m_xVetoNu = 0;
+ m_yVetoNu = 0;
+ m_sxVetoNu = 0;
+ m_syVetoNu = 0;
+ m_thetaxVetoNu = 0;
+ m_thetayVetoNu = 0;
}
\ No newline at end of file
diff --git a/PhysicsAnalysis/NeutrinoSearch/src/NeutrinoSearchAlg.h b/PhysicsAnalysis/NeutrinoSearch/src/NeutrinoSearchAlg.h
index d1732eb1238a15bb7b2fb397869b3c4486aa50b3..4c8ca653e2b7832dac95645135a10c83c57a846a 100644
--- a/PhysicsAnalysis/NeutrinoSearch/src/NeutrinoSearchAlg.h
+++ b/PhysicsAnalysis/NeutrinoSearch/src/NeutrinoSearchAlg.h
@@ -8,10 +8,11 @@
#include "xAODFaserWaveform/WaveformHit.h"
#include "xAODTruth/TruthEventContainer.h"
#include "xAODTruth/TruthParticleContainer.h"
+#include "ScintSimEvent/ScintHitCollection.h"
#include "TrackerPrepRawData/FaserSCT_ClusterContainer.h"
#include "TrackerSimData/TrackerSimDataCollection.h"
#include "FaserActsGeometryInterfaces/IFaserActsExtrapolationTool.h"
-
+#include "FaserActsGeometryInterfaces/IFaserActsTrackingGeometryTool.h"
class TTree;
class FaserSCT_ID;
@@ -44,6 +45,7 @@ private:
SG::ReadHandleKey<xAOD::TruthEventContainer> m_truthEventContainer { this, "EventContainer", "TruthEvents", "Truth event container name." };
SG::ReadHandleKey<xAOD::TruthParticleContainer> m_truthParticleContainer { this, "ParticleContainer", "TruthParticles", "Truth particle container name." };
SG::ReadHandleKey<TrackerSimDataCollection> m_simDataCollection {this, "TrackerSimDataCollection", "SCT_SDO_Map"};
+ SG::ReadHandleKey<ScintHitCollection> m_vetoNuHitKey { this, "VetoNuHitCollection", "VetoNuHits" };
SG::ReadHandleKey<TrackCollection> m_trackCollection { this, "TrackCollection", "CKFTrackCollection", "Input track collection name" };
SG::ReadHandleKey<xAOD::WaveformHitContainer> m_vetoNuContainer { this, "VetoNuContainer", "VetoNuWaveformHits", "VetoNu hit container name" };
@@ -54,6 +56,7 @@ private:
SG::ReadHandleKey<Tracker::FaserSCT_ClusterContainer> m_clusterContainer { this, "ClusterContainer", "SCT_ClusterContainer", "Tracker cluster container name" };
ToolHandle<IFaserActsExtrapolationTool> m_extrapolationTool { this, "ExtrapolationTool", "FaserActsExtrapolationTool" };
+ ToolHandle<IFaserActsTrackingGeometryTool> m_trackingGeometryTool {this, "TrackingGeometryTool", "FaserActsTrackingGeometryTool"};
const TrackerDD::SCT_DetectorManager* m_detMgr {nullptr};
const FaserSCT_ID* m_sctHelper;
@@ -75,7 +78,8 @@ private:
BooleanProperty m_useGenieWeights { this, "UseGenieWeights", false, "Flag to weight events according to Genie luminosity" };
IntegerProperty m_flukaCollisions { this, "FlukaCollisions", 137130000, "Number of proton-proton collisions in FLUKA sample." };
DoubleProperty m_flukaCrossSection { this, "FlukaCrossSection", 80.0, "Fluka p-p inelastic cross-section in millibarns." };
- DoubleProperty m_genieLuminosity { this, "GenieLuminosity", 150.0, "Genie luminosity in inverse fb." };
+ DoubleProperty m_genieLuminosity { this, "GenieLuminosity", 600.0, "Genie luminosity in inverse fb." };
+ DoubleProperty m_zExtrapolate { this, "ZExtrapolationPlane", -3112.0, "Plane to extrapolate tracks for VetoNu" };
// BooleanProperty m_enforceBlinding { this, "EnforceBlinding", true, "Ignore data events with no VetoNu signals." };
const bool m_enforceBlinding {true};
@@ -111,6 +115,9 @@ private:
mutable int m_station2Clusters;
mutable int m_station3Clusters;
+ // mutable double m_vetoNuHitsMeanX;
+ // mutable double m_vetoNuHitsMeanY;
+
mutable double m_x;
mutable double m_y;
mutable double m_z;
@@ -123,10 +130,17 @@ private:
mutable int m_ndof;
mutable int m_longTracks;
mutable double m_truthLeptonMomentum;
+ mutable double m_truthLeptonX;
+ mutable double m_truthLeptonY;
mutable int m_truthBarcode;
mutable int m_truthPdg;
mutable double m_crossSection;
-
+ mutable double m_xVetoNu;
+ mutable double m_yVetoNu;
+ mutable double m_sxVetoNu;
+ mutable double m_syVetoNu;
+ mutable double m_thetaxVetoNu;
+ mutable double m_thetayVetoNu;
};
inline const ServiceHandle <ITHistSvc> &NeutrinoSearchAlg::histSvc() const {
diff --git a/xAOD/xAODFaserCalorimeter/Root/CalorimeterHitAuxContainer_v1.cxx b/xAOD/xAODFaserCalorimeter/Root/CalorimeterHitAuxContainer_v1.cxx
index f73a83dd127cc0d24d999c0e8f140e913a329e76..277ecd25c0d1e2056a6c2a1eeef101f5570f80c2 100644
--- a/xAOD/xAODFaserCalorimeter/Root/CalorimeterHitAuxContainer_v1.cxx
+++ b/xAOD/xAODFaserCalorimeter/Root/CalorimeterHitAuxContainer_v1.cxx
@@ -11,8 +11,10 @@ namespace xAOD {
: AuxContainerBase() {
AUX_VARIABLE(channel);
- AUX_VARIABLE(Edep);
- AUX_VARIABLE(raw_Edep);
+ AUX_VARIABLE(Nmip);
+ AUX_VARIABLE(E_dep);
+ AUX_VARIABLE(E_EM);
+ AUX_VARIABLE(fit_to_raw_ratio);
AUX_VARIABLE(WaveformLink);
}
diff --git a/xAOD/xAODFaserCalorimeter/Root/CalorimeterHit_v1.cxx b/xAOD/xAODFaserCalorimeter/Root/CalorimeterHit_v1.cxx
index f085cb0c92acca759173260ba56c26813fd72935..ff592dff67ca5c37bf0ea954dcc9eb63f3a26a28 100644
--- a/xAOD/xAODFaserCalorimeter/Root/CalorimeterHit_v1.cxx
+++ b/xAOD/xAODFaserCalorimeter/Root/CalorimeterHit_v1.cxx
@@ -15,9 +15,13 @@ namespace xAOD {
AUXSTORE_PRIMITIVE_SETTER_AND_GETTER( CalorimeterHit_v1, unsigned int, channel, set_channel )
- AUXSTORE_PRIMITIVE_SETTER_AND_GETTER( CalorimeterHit_v1, float, Edep, set_Edep )
+ AUXSTORE_PRIMITIVE_SETTER_AND_GETTER( CalorimeterHit_v1, float, Nmip, set_Nmip )
- AUXSTORE_PRIMITIVE_SETTER_AND_GETTER( CalorimeterHit_v1, float, raw_Edep, set_raw_Edep )
+ AUXSTORE_PRIMITIVE_SETTER_AND_GETTER( CalorimeterHit_v1, float, E_dep, set_E_dep )
+
+ AUXSTORE_PRIMITIVE_SETTER_AND_GETTER( CalorimeterHit_v1, float, E_EM, set_E_EM )
+
+ AUXSTORE_PRIMITIVE_SETTER_AND_GETTER( CalorimeterHit_v1, float, fit_to_raw_ratio, set_fit_to_raw_ratio )
// setters and getters for the Calo WaveformHit links
@@ -58,7 +62,7 @@ namespace xAOD {
std::ostream& operator<<(std::ostream& s, const xAOD::CalorimeterHit_v1& hit) {
s << "xAODCalorimeterHit:"
<< " channel = " << hit.channel()
- << ", Edep = " << hit.Edep()
+ << ", E_dep = " << hit.E_dep()
<< std::endl;
return s;
diff --git a/xAOD/xAODFaserCalorimeter/xAODFaserCalorimeter/versions/CalorimeterHitAuxContainer_v1.h b/xAOD/xAODFaserCalorimeter/xAODFaserCalorimeter/versions/CalorimeterHitAuxContainer_v1.h
index 41a89a0325e9dd5ce78dab40ad6ae08458eb37c3..230da2716246909fcf85653c288dc043aa9929a5 100644
--- a/xAOD/xAODFaserCalorimeter/xAODFaserCalorimeter/versions/CalorimeterHitAuxContainer_v1.h
+++ b/xAOD/xAODFaserCalorimeter/xAODFaserCalorimeter/versions/CalorimeterHitAuxContainer_v1.h
@@ -30,8 +30,10 @@ namespace xAOD {
/// @name Basic variables
///@ {
std::vector<unsigned int> channel;
- std::vector<float> Edep;
- std::vector<float> raw_Edep;
+ std::vector<float> Nmip;
+ std::vector<float> E_dep;
+ std::vector<float> E_EM;
+ std::vector<float> fit_to_raw_ratio;
typedef std::vector< ElementLink< WaveformHitContainer > > WaveformHitLink_t;
std::vector< WaveformHitLink_t > WaveformLink;
diff --git a/xAOD/xAODFaserCalorimeter/xAODFaserCalorimeter/versions/CalorimeterHit_v1.h b/xAOD/xAODFaserCalorimeter/xAODFaserCalorimeter/versions/CalorimeterHit_v1.h
index fefc0d5ae3dac013ee734e34ed82c600f8fbde5d..c27d9f033efc6eebaeb2120440a8855362b77700 100644
--- a/xAOD/xAODFaserCalorimeter/xAODFaserCalorimeter/versions/CalorimeterHit_v1.h
+++ b/xAOD/xAODFaserCalorimeter/xAODFaserCalorimeter/versions/CalorimeterHit_v1.h
@@ -38,11 +38,17 @@ namespace xAOD {
unsigned int channel() const;
void set_channel(unsigned int value);
- float Edep() const;
- void set_Edep(float value);
+ float Nmip() const;
+ void set_Nmip(float value);
- float raw_Edep() const;
- void set_raw_Edep(float value);
+ float E_dep() const;
+ void set_E_dep(float value);
+
+ float E_EM() const;
+ void set_E_EM(float value);
+
+ float fit_to_raw_ratio() const;
+ void set_fit_to_raw_ratio(float value);
// Waveform Hits
typedef std::vector< ElementLink< xAOD::WaveformHitContainer > > WaveformHitLinks_t;