Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add the detector table to our hdf5 files #13

Merged
merged 11 commits into from
May 15, 2024
40 changes: 40 additions & 0 deletions numl/NeutrinoML/HDF5Maker_module.cc
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,10 @@
#include "canvas/Persistency/Common/FindManyP.h"
#include "messagefacility/MessageLogger/MessageLogger.h"

#include "larcore/Geometry/Geometry.h"
#include "larcorealg/Geometry/GeometryCore.h"
#include "larcorealg/Geometry/TPCGeo.h"

#include "nusimdata/SimulationBase/MCTruth.h"
#include "lardataobj/RecoBase/Hit.h"
#include "lardataobj/RecoBase/SpacePoint.h"
Expand Down Expand Up @@ -103,6 +107,13 @@ class HDF5Maker : public art::EDAnalyzer {
hep_hpc::hdf5::Column<int, 1>, // g4 id
hep_hpc::hdf5::Column<float, 1> // deposited energy [ MeV ]
>* fEnergyDepNtuple; ///< energy deposition ntuple

hep_hpc::hdf5::Ntuple<hep_hpc::hdf5::Column<float, 1>, // TPC module center x
hep_hpc::hdf5::Column<float, 1>, // TPC module center y
hep_hpc::hdf5::Column<float, 1>, // TPC module center z
hep_hpc::hdf5::Column<int, 1>, // drift direction
hep_hpc::hdf5::Column<int, 1> // tpc id
>* fDetectorNtuple;
};


Expand All @@ -120,6 +131,7 @@ HDF5Maker::HDF5Maker(fhicl::ParameterSet const& p)
throw art::Exception(art::errors::Configuration)
<< "EventInfo must be \"none\" or \"nu\", not " << fEventInfo;
}
std::vector<double> tpc_ids_checked = {}; // add the TPC ids here so we only check them once

void HDF5Maker::analyze(art::Event const& e)
{
Expand Down Expand Up @@ -257,6 +269,24 @@ void HDF5Maker::analyze(art::Event const& e)
<< ", local wire " << wireid.Wire
<< ", local time " << hit->PeakTime();

// if we haven't checked a tpc id add the info to the table
if(std::find(tpc_ids_checked.begin(), tpc_ids_checked.end(), wireid.TPC ) == tpc_ids_checked.end()){
auto const& tpcgeom = art::ServiceHandle<geo::Geometry>()->TPC(geo::TPCID{0,wireid.TPC});
geo::Point_t center = tpcgeom.GetCenter();
fDetectorNtuple->insert((float)geo::vect::Xcoord(center),(float)geo::vect::Ycoord(center),(float)geo::vect::Zcoord(center),tpcgeom.DetectDriftDirection(),wireid.TPC
);
mf::LogInfo("HDF5Maker") << "Filling detector table"
<< "\nrun " << evtID[0] << ", subrun " << evtID[1]
<< ", event " << evtID[2]
<< "\n tpc center x" << (float)geo::vect::Xcoord(center)
<< "\n tpc center y" << (float)geo::vect::Ycoord(center)
<< "\n tpc center z" << (float)geo::vect::Zcoord(center)
<< "\n drift direction " << tpcgeom.DetectDriftDirection()
<< "\n tpc " << wireid.TPC;

tpc_ids_checked.push_back(wireid.TPC); // once we check a tpc id add it to the list of ids we checked
vhewes marked this conversation as resolved.
Show resolved Hide resolved
}

// Fill energy deposit table
if (fUseMap) {
std::vector<art::Ptr<simb::MCParticle>> particle_vec = hittruth->at(hit.key());
Expand Down Expand Up @@ -406,6 +436,15 @@ void HDF5Maker::beginSubRun(art::SubRun const& sr) {
hep_hpc::hdf5::make_scalar_column<int>("g4_id"),
hep_hpc::hdf5::make_scalar_column<float>("energy")
));

fDetectorNtuple = new hep_hpc::hdf5::Ntuple(
hep_hpc::hdf5::make_ntuple({fFile, "detector_table", 1000},
vhewes marked this conversation as resolved.
Show resolved Hide resolved
hep_hpc::hdf5::make_scalar_column<float>("tpc_center_x"),
hep_hpc::hdf5::make_scalar_column<float>("tpc_center_y"),
hep_hpc::hdf5::make_scalar_column<float>("tpc_center_z"),
hep_hpc::hdf5::make_scalar_column<int>("drift_direction"),
hep_hpc::hdf5::make_scalar_column<int>("tpc")
));
}

void HDF5Maker::endSubRun(art::SubRun const& sr) {
Expand All @@ -415,6 +454,7 @@ void HDF5Maker::endSubRun(art::SubRun const& sr) {
delete fHitNtuple;
delete fParticleNtuple;
delete fEnergyDepNtuple;
delete fDetectorNtuple;
fFile.close();
}

Expand Down