Skip to content

Commit

Permalink
Truth matching updates for the Pandora LArRecoND CAFs
Browse files Browse the repository at this point in the history
  • Loading branch information
@jback08
jback08 committed Oct 9, 2024
1 parent 611c011 commit 27b5638
Show file tree
Hide file tree
Showing 3 changed files with 104 additions and 35 deletions.
25 changes: 24 additions & 1 deletion cfg/pandora.fcl
View file
Original file line number Diff line number Diff line change
@@ -1,6 +1,29 @@
#include "NDCAFMaker.fcl"
nd_cafmaker: @local::standard_nd_cafmaker
nd_cafmaker.CAFMakerSettings.GHEPFiles: ["/exp/dune/data/users/jback/ND_CAFs/MiniRun5/MiniRun5_1E19_RHC.genie.nu.0000001.GHEP.root"]

nd_cafmaker.CAFMakerSettings.GHEPFiles: [
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.nu.0000010.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.nu.0000011.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.nu.0000012.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.nu.0000013.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.nu.0000014.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.nu.0000015.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.nu.0000016.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.nu.0000017.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.nu.0000018.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.nu.0000019.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.rock.0000010.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.rock.0000011.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.rock.0000012.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.rock.0000013.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.rock.0000014.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.rock.0000015.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.rock.0000016.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.rock.0000017.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.rock.0000018.GHEP.root",
"/exp/dune/data/users/jback/ND_CAFs/MiniRun6/MiniRun5_1E19_RHC.genie.rock.0000019.GHEP.root"
]

nd_cafmaker.CAFMakerSettings.PandoraLArRecoNDFile: "/exp/dune/data/users/jback/ND_CAFs/MiniRun6/LArRecoND_MR6_0000001.root"
nd_cafmaker.CAFMakerSettings.OutputFile: "CAF_MR6_0000001.root"
nd_cafmaker.CAFMakerSettings.Verbosity: VERBOSE
Expand Down
108 changes: 77 additions & 31 deletions src/reco/PandoraLArRecoNDBranchFiller.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -49,7 +49,7 @@ namespace cafmaker
m_LArRecoNDTree->SetBranchAddress("dirZ", &m_dirZVect);
m_LArRecoNDTree->SetBranchAddress("energy", &m_energyVect);
m_LArRecoNDTree->SetBranchAddress("n3DHits", &m_n3DHitsVect);
m_LArRecoNDTree->SetBranchAddress("mcNuId", &m_mcNuIdVect);
m_LArRecoNDTree->SetBranchAddress("mcVertexId", &m_mcVertexIdVect);
m_LArRecoNDTree->SetBranchAddress("isPrimary", &m_isPrimaryVect);
m_LArRecoNDTree->SetBranchAddress("mcId", &m_mcIdVect);
m_LArRecoNDTree->SetBranchAddress("completeness", &m_completenessVect);
Expand Down Expand Up @@ -96,12 +96,13 @@ namespace cafmaker

// Fill track and shower info. Both use the same clusters (PFOs), and no
// distinction is made (yet) to identify which are tracks or showers
FillTracks(sr, nClusters);
FillShowers(sr, nClusters);
FillTracks(sr, nClusters, truthMatch);
FillShowers(sr, nClusters, truthMatch);
}

// ------------------------------------------------------------------------------
void PandoraLArRecoNDBranchFiller::FillTracks(caf::StandardRecord& sr, const int nClusters) const
void PandoraLArRecoNDBranchFiller::FillTracks(caf::StandardRecord &sr, const int nClusters,
const TruthMatcher *truthMatch) const
{
// Create tracks for each PFO (cluster) in the event
LOG.VERBOSE() << " Pandora LArRecoND FillTracks using " << nClusters <<" PFO clusters\n";
Expand Down Expand Up @@ -157,46 +158,69 @@ namespace cafmaker

// Use truth matching info from Pandora's LArContent hierarchy tools.
// For LArRecoND MC SpacePoints, we offset the MCId's to make them all unique:
// mcNuId = vertex_id + 10^8, so vertex_id = mcNuId - 10^8
// mcId = traj_id + nuIndex*10^6, where nuIndex = 0 to N-1 neutrinos
// nuIndex = int(mcId/10^6), so traj_id = mcId - nuIndex*10^6
const int mcNuId = (m_mcNuIdVect != nullptr) ? (*m_mcNuIdVect)[i] : 0;
const int vertex_id = (mcNuId > m_nuIdOffset) ? mcNuId - m_nuIdOffset : mcNuId;
// The vertex_id's are the same as those in the HDF5 files

const long vertex_id = (m_mcVertexIdVect != nullptr) ? (*m_mcVertexIdVect)[i] : 0;
const int isPrimary = (m_isPrimaryVect != nullptr) ? (*m_isPrimaryVect)[i] : -1;

const int mcId = (m_mcIdVect != nullptr) ? (*m_mcIdVect)[i] : 0;
const int nuIndex = int(mcId/m_maxMCId);
const int traj_id = mcId - nuIndex*m_maxMCId;
caf::TrueParticleID trueID;
trueID.ixn = vertex_id;

caf::TrueParticleID truePartID;

if (isPrimary == 1) {
trueID.type = caf::TrueParticleID::kPrimary;
truePartID.type = caf::TrueParticleID::kPrimary;
truePartID.part = traj_id;
} else if (isPrimary == -1) {
trueID.type = caf::TrueParticleID::kUnknown;
truePartID.type = caf::TrueParticleID::kUnknown;
} else {
trueID.type = caf::TrueParticleID::kSecondary;
truePartID.type = caf::TrueParticleID::kSecondary;
}

try
{
// Get the true interaction in the stack
caf::SRTrueInteraction &srTrueInt = truthMatch->GetTrueInteraction(sr, vertex_id);
const auto predicate = [&srTrueInt](const caf::SRTrueInteraction& ixn) { return ixn.id == srTrueInt.id; };
// Get the truth interaction index
const int srTrueIntIdx = std::distance(sr.mc.nu.begin(), std::find_if(sr.mc.nu.begin(), sr.mc.nu.end(), predicate));
truePartID.ixn = srTrueIntIdx;

// If the particle is not a primary, we might want to create a new particle if it wasn't created originally
if (isPrimary != 1) {
const auto pred = [&traj_id](const caf::SRTrueParticle& part) { return part.G4ID == traj_id; };
truePartID.part = std::distance(srTrueInt.sec.begin(),
std::find_if(srTrueInt.sec.begin(), srTrueInt.sec.end(), pred));
}
}
catch (...)
{
LOG.VERBOSE() << "Could not fill some of the truth information" << "\n";
}
trueID.part = traj_id;

// Just store the best MC match
std::vector<caf::TrueParticleID> trueIDVect;
trueIDVect.emplace_back(trueID);
track.truth = trueIDVect;
std::vector<caf::TrueParticleID> truePartIDVect;
truePartIDVect.emplace_back(truePartID);
track.truth = truePartIDVect;

// Fraction of true MC hits that are captured by the reconstructed cluster
const float completeness = (m_completenessVect != nullptr) ? (*m_completenessVect)[i] : 0.0;
std::vector<float> truthOverlap;
truthOverlap.emplace_back(completeness);
track.truthOverlap = truthOverlap;

// Add track to the record
sr.nd.lar.pandora[sliceId].tracks.emplace_back(std::move(track));
sr.nd.lar.pandora[sliceId].ntracks++;
}
}

// ------------------------------------------------------------------------------
void PandoraLArRecoNDBranchFiller::FillShowers(caf::StandardRecord& sr, const int nClusters) const
void PandoraLArRecoNDBranchFiller::FillShowers(caf::StandardRecord &sr, const int nClusters,
const TruthMatcher *truthMatch) const
{
// Create showers for each PFO (cluster) in the event
LOG.VERBOSE() << " Pandora LArRecoND FillShowers using " << nClusters <<" PFO clusters\n";
Expand Down Expand Up @@ -231,31 +255,53 @@ namespace cafmaker

// Use truth matching info from Pandora's LArContent hierarchy tools.
// For LArRecoND MC SpacePoints, we offset the MCId's to make them all unique:
// mcNuId = vertex_id + 10^8, so vertex_id = mcNuId - 10^8
// mcId = traj_id + nuIndex*10^6, where nuIndex = 0 to N-1 neutrinos
// nuIndex = int(mcId/10^6), so traj_id = mcId - nuIndex*10^6
const int mcNuId = (m_mcNuIdVect != nullptr) ? (*m_mcNuIdVect)[i] : 0;
const int vertex_id = (mcNuId > m_nuIdOffset) ? mcNuId - m_nuIdOffset : mcNuId;
// The vertex_id's are the same as those in the HDF5 files

const long vertex_id = (m_mcVertexIdVect != nullptr) ? (*m_mcVertexIdVect)[i] : 0;
const int isPrimary = (m_isPrimaryVect != nullptr) ? (*m_isPrimaryVect)[i] : -1;

const int mcId = (m_mcIdVect != nullptr) ? (*m_mcIdVect)[i] : 0;
const int nuIndex = int(mcId/m_maxMCId);
const int traj_id = mcId - nuIndex*m_maxMCId;

caf::TrueParticleID trueID;
trueID.ixn = vertex_id;
caf::TrueParticleID truePartID;

if (isPrimary == 1) {
trueID.type = caf::TrueParticleID::kPrimary;
truePartID.type = caf::TrueParticleID::kPrimary;
truePartID.part = traj_id;
} else if (isPrimary == -1) {
trueID.type = caf::TrueParticleID::kUnknown;
truePartID.type = caf::TrueParticleID::kUnknown;
} else {
trueID.type = caf::TrueParticleID::kSecondary;
truePartID.type = caf::TrueParticleID::kSecondary;
}

try
{
// Get the true interaction in the stack
caf::SRTrueInteraction &srTrueInt = truthMatch->GetTrueInteraction(sr, vertex_id);
const auto predicate = [&srTrueInt](const caf::SRTrueInteraction& ixn) { return ixn.id == srTrueInt.id; };
// Get the truth interaction index
const int srTrueIntIdx = std::distance(sr.mc.nu.begin(), std::find_if(sr.mc.nu.begin(), sr.mc.nu.end(), predicate));
truePartID.ixn = srTrueIntIdx;

// If the particle is not a primary, we might want to create a new particle if it wasn't created originally
if (isPrimary != 1) {
const auto pred = [&traj_id](const caf::SRTrueParticle& part) { return part.G4ID == traj_id; };
truePartID.part = std::distance(srTrueInt.sec.begin(),
std::find_if(srTrueInt.sec.begin(), srTrueInt.sec.end(), pred));
}
}
catch (...)
{
LOG.VERBOSE() << "Could not fill some of the truth information" << "\n";
}
trueID.part = traj_id;

// Just store the best MC match
std::vector<caf::TrueParticleID> trueIDVect;
trueIDVect.emplace_back(trueID);
shower.truth = trueIDVect;
std::vector<caf::TrueParticleID> truePartIDVect;
truePartIDVect.emplace_back(truePartID);
shower.truth = truePartIDVect;

// Fraction of true MC hits that are captured by the reconstructed cluster
const float completeness = (m_completenessVect != nullptr) ? (*m_completenessVect)[i] : 0.0;
Expand Down
6 changes: 3 additions & 3 deletions src/reco/PandoraLArRecoNDBranchFiller.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -40,8 +40,8 @@ namespace cafmaker
const cafmaker::Params &par,
const TruthMatcher *truthMatch= nullptr) const override;

void FillTracks(caf::StandardRecord &sr, const int nClusters) const;
void FillShowers(caf::StandardRecord &sr, const int nClusters) const;
void FillTracks(caf::StandardRecord &sr, const int nClusters, const TruthMatcher *truthMatch) const;
void FillShowers(caf::StandardRecord &sr, const int nClusters, const TruthMatcher *truthMatch) const;

TFile *m_LArRecoNDFile;
TTree *m_LArRecoNDTree;
Expand All @@ -64,7 +64,7 @@ namespace cafmaker
std::vector<float> *m_dirZVect = nullptr;
std::vector<float> *m_energyVect = nullptr;
std::vector<int> *m_n3DHitsVect = nullptr;
std::vector<int> *m_mcNuIdVect = nullptr;
std::vector<long> *m_mcVertexIdVect = nullptr;
std::vector<int> *m_isPrimaryVect = nullptr;
std::vector<int> *m_mcIdVect = nullptr;
std::vector<float> *m_completenessVect = nullptr;
Expand Down

0 comments on commit 27b5638

Please sign in to comment.