Skip to content
Snippets Groups Projects
Commit 622c9f34 authored by Whitney Armstrong's avatar Whitney Armstrong
Browse files

Added tracking scripts and mod options file .

-modified:   tracking/options/tracker_reconstruction.py
-modified:   tracking/scripts/rec_central_electrons.cxx
parent 0625c883
Branches
No related tags found
1 merge request!31Added tracking scripts and mod options file .
...@@ -29,6 +29,8 @@ from Configurables import Jug__Reco__TrackerHitReconstruction as TrackerHitRecon ...@@ -29,6 +29,8 @@ from Configurables import Jug__Reco__TrackerHitReconstruction as TrackerHitRecon
from Configurables import Jug__Reco__TrackerSourceLinker as TrackerSourceLinker from Configurables import Jug__Reco__TrackerSourceLinker as TrackerSourceLinker
from Configurables import Jug__Reco__TrackParamTruthInit as TrackParamTruthInit from Configurables import Jug__Reco__TrackParamTruthInit as TrackParamTruthInit
from Configurables import Jug__Reco__TrackParamClusterInit as TrackParamClusterInit from Configurables import Jug__Reco__TrackParamClusterInit as TrackParamClusterInit
from Configurables import Jug__Reco__TrackParamVertexClusterInit as TrackParamVertexClusterInit
from Configurables import Jug__Reco__TrackFindingAlgorithm as TrackFindingAlgorithm from Configurables import Jug__Reco__TrackFindingAlgorithm as TrackFindingAlgorithm
from Configurables import Jug__Reco__ParticlesFromTrackFit as ParticlesFromTrackFit from Configurables import Jug__Reco__ParticlesFromTrackFit as ParticlesFromTrackFit
from Configurables import Jug__Reco__EMCalReconstruction as EMCalReconstruction from Configurables import Jug__Reco__EMCalReconstruction as EMCalReconstruction
...@@ -38,77 +40,106 @@ from Configurables import Jug__Reco__SimpleClustering as SimpleClustering ...@@ -38,77 +40,106 @@ from Configurables import Jug__Reco__SimpleClustering as SimpleClustering
podioinput = PodioInput("PodioReader", podioinput = PodioInput("PodioReader",
collections=["mcparticles","SiTrackerBarrelHits","EcalBarrelHits"])#, OutputLevel=DEBUG) collections=["mcparticles","SiTrackerBarrelHits","SiVertexBarrelHits","EcalBarrelHits"])#, OutputLevel=DEBUG)
## copiers to get around input --> output copy bug. Note the "2" appended to the output collection. ## copiers to get around input --> output copy bug. Note the "2" appended to the output collection.
copier = MCCopier("MCCopier", inputCollection="mcparticles", outputCollection="mcparticles2",OutputLevel=DEBUG) copier = MCCopier("MCCopier",
trkcopier = TrkCopier("TrkCopier", inputCollection="SiTrackerBarrelHits", outputCollection="SiTrackerBarrelHits2",OutputLevel=DEBUG) inputCollection="mcparticles",
#caldigi = ExampleCaloDigi(inputHitCollection="FAEC_ShHits",outputHitCollection="RawFAECShowerHits") outputCollection="mcparticles2")
#ufsd_digi = UFSDTrackerDigi("ufsd_digi", inputHitCollection="SiVertexBarrelHits",outputHitCollection="VertexRawHits",timeResolution=8) trkcopier = TrkCopier("TrkCopier",
inputCollection="SiTrackerBarrelHits",
outputCollection="SiTrackerBarrelHits2")
ecal_digi = EMCalorimeterDigi("ecal_digi", ecal_digi = EMCalorimeterDigi("ecal_digi",
inputHitCollection="EcalBarrelHits", inputHitCollection="EcalBarrelHits",
outputHitCollection="RawEcalBarrelHits") outputHitCollection="RawEcalBarrelHits")
ufsd_digi = UFSDTrackerDigi("ufsd_digi",
inputHitCollection="SiTrackerBarrelHits",
outputHitCollection="SiTrackerBarrelRawHits",
timeResolution=8)
vtx_digi = UFSDTrackerDigi("vtx_digi",
inputHitCollection="SiVertexBarrelHits",
outputHitCollection="SiVertexBarrelRawHits",
timeResolution=8)
ecal_reco = EMCalReconstruction("ecal_reco", ecal_reco = EMCalReconstruction("ecal_reco",
inputHitCollection="RawEcalBarrelHits", inputHitCollection="RawEcalBarrelHits",
outputHitCollection="RecEcalBarrelHits", outputHitCollection="RecEcalBarrelHits",
minModuleEdep=1.0*units.MeV,OutputLevel=DEBUG) minModuleEdep=1.0*units.MeV)
simple_cluster = SimpleClustering("simple_cluster", simple_cluster = SimpleClustering("simple_cluster",
inputHitCollection="RecEcalBarrelHits", inputHitCollection="RecEcalBarrelHits",
outputClusters="SimpleClusters", outputClusters="SimpleClusters",
minModuleEdep=1.0*units.MeV,OutputLevel=DEBUG) minModuleEdep=1.0*units.MeV,OutputLevel=DEBUG)
cluster_init = TrackParamClusterInit("trk_clust_init", trk_barrel_reco = TrackerHitReconstruction("trk_barrel_reco",
inputClusters="SimpleClusters", inputHitCollection="SiTrackerBarrelRawHits",
outputInitialTrackParameters="InitTrackParamsFromClusters", outputHitCollection="TrackerBarrelRecHits")
OutputLevel=DEBUG)
ufsd_digi = UFSDTrackerDigi("ufsd_digi", vtx_barrel_reco = TrackerHitReconstruction("vtx_barrel_reco",
inputHitCollection="SiTrackerBarrelHits", inputHitCollection = vtx_digi.outputHitCollection,
outputHitCollection="SiTrackerBarrelRawHits", outputHitCollection="VertexBarrelRecHits")
timeResolution=8)
trackpartruth = TrackParamTruthInit("trk_par_init", # Source linker
sourcelinker = TrackerSourceLinker("trk_srclinker",
inputHitCollection="TrackerBarrelRecHits",
outputSourceLinks="BarrelTrackSourceLinks")
## Track param init
truth_trk_init = TrackParamTruthInit("truth_trk_init",
inputMCParticles="mcparticles", inputMCParticles="mcparticles",
outputInitialTrackParameters="InitTrackParams", outputInitialTrackParameters="InitTrackParams",
OutputLevel=DEBUG) OutputLevel=DEBUG)
trackerhit = TrackerHitReconstruction("trk_hit_reco", clust_trk_init = TrackParamClusterInit("clust_trk_init",
inputHitCollection="SiTrackerBarrelRawHits", inputClusters="SimpleClusters",
outputHitCollection="TrackerBarrelRecHits", outputInitialTrackParameters="InitTrackParamsFromClusters",
OutputLevel=DEBUG) OutputLevel=DEBUG)
sourcelinker = TrackerSourceLinker("trk_srclinker", vtxcluster_trk_init = TrackParamVertexClusterInit("vtxcluster_trk_init",
inputHitCollection="TrackerBarrelRecHits", inputVertexHits="VertexBarrelRecHits",
outputSourceLinks="BarrelTrackSourceLinks", inputClusters="SimpleClusters",
outputInitialTrackParameters="InitTrackParamsFromVtxClusters",
OutputLevel=DEBUG) OutputLevel=DEBUG)
# Tracking algorithms
trk_find_alg = TrackFindingAlgorithm("trk_find_alg", trk_find_alg = TrackFindingAlgorithm("trk_find_alg",
inputSourceLinks="BarrelTrackSourceLinks", inputSourceLinks="BarrelTrackSourceLinks",
inputInitialTrackParameters= "InitTrackParams",#"InitTrackParamsFromClusters", inputInitialTrackParameters= "InitTrackParams",#"InitTrackParamsFromClusters",
outputTrajectories="trajectories", outputTrajectories="trajectories",
OutputLevel=DEBUG) OutputLevel=DEBUG)
trk_find_alg1 = TrackFindingAlgorithm("trk_find_alg1",
inputSourceLinks="BarrelTrackSourceLinks",
inputInitialTrackParameters= "InitTrackParamsFromClusters",
outputTrajectories="trajectories1",
OutputLevel=DEBUG)
parts_from_fit = ParticlesFromTrackFit("parts_from_fit", parts_from_fit = ParticlesFromTrackFit("parts_from_fit",
inputTrajectories="trajectories", inputTrajectories="trajectories",
outputParticles="ReconstructedParticles", outputParticles="ReconstructedParticles",
outputTrackParameters="outputTrackParameters", outputTrackParameters="outputTrackParameters",
OutputLevel=DEBUG) OutputLevel=DEBUG)
trk_find_alg1 = TrackFindingAlgorithm("trk_find_alg1",
inputSourceLinks="BarrelTrackSourceLinks",
inputInitialTrackParameters= "InitTrackParamsFromClusters",
outputTrajectories="trajectories1",
OutputLevel=DEBUG)
parts_from_fit1 = ParticlesFromTrackFit("parts_from_fit1", parts_from_fit1 = ParticlesFromTrackFit("parts_from_fit1",
inputTrajectories="trajectories1", inputTrajectories="trajectories1",
outputParticles="ReconstructedParticles1", outputParticles="ReconstructedParticles1",
outputTrackParameters="outputTrackParameters1", outputTrackParameters="outputTrackParameters1",
OutputLevel=DEBUG) OutputLevel=DEBUG)
trk_find_alg2 = TrackFindingAlgorithm("trk_find_alg2",
inputSourceLinks="BarrelTrackSourceLinks",
inputInitialTrackParameters= "InitTrackParamsFromVtxClusters",
outputTrajectories="trajectories2",
OutputLevel=DEBUG)
parts_from_fit2 = ParticlesFromTrackFit("parts_from_fit2",
inputTrajectories="trajectories2",
outputParticles="ReconstructedParticles2",
outputTrackParameters="outputTrackParameters2",
OutputLevel=DEBUG)
#types = [] #types = []
## this printout is useful to check that the type information is passed to python correctly ## this printout is useful to check that the type information is passed to python correctly
#print("---------------------------------------\n") #print("---------------------------------------\n")
...@@ -140,15 +171,18 @@ out.outputCommands = ["keep *", ...@@ -140,15 +171,18 @@ out.outputCommands = ["keep *",
ApplicationMgr( ApplicationMgr(
TopAlg = [podioinput, TopAlg = [podioinput,
copier, trkcopier, copier, trkcopier,
ecal_digi,ecal_reco, ecal_digi, ufsd_digi, vtx_digi,
ecal_reco,
simple_cluster, simple_cluster,
ufsd_digi, trackerhit, trk_barrel_reco,
vtx_barrel_reco,
sourcelinker, sourcelinker,
cluster_init, trackpartruth, clust_trk_init,
trk_find_alg, truth_trk_init,
trk_find_alg1, vtxcluster_trk_init,
parts_from_fit, trk_find_alg, parts_from_fit,
parts_from_fit1, trk_find_alg1, parts_from_fit1,
trk_find_alg2, parts_from_fit2,
out out
], ],
EvtSel = 'NONE', EvtSel = 'NONE',
......
...@@ -81,6 +81,16 @@ auto fourvec = [](ROOT::VecOps::RVec<dd4pod::Geant4ParticleData> const& in) { ...@@ -81,6 +81,16 @@ auto fourvec = [](ROOT::VecOps::RVec<dd4pod::Geant4ParticleData> const& in) {
return result; return result;
}; };
auto delta_p = [](const std::vector<double>& tracks, const std::vector<double>& thrown) {
std::vector<double> res;
for (const auto& p1 : thrown) {
for (const auto& p2 : tracks) {
res.push_back(p1 - p2);
}
}
return res;
};
int rec_central_electrons(const char* fname = "topside/rec_central_electrons.root") int rec_central_electrons(const char* fname = "topside/rec_central_electrons.root")
{ {
...@@ -94,30 +104,10 @@ int rec_central_electrons(const char* fname = "topside/rec_central_electrons.roo ...@@ -94,30 +104,10 @@ int rec_central_electrons(const char* fname = "topside/rec_central_electrons.roo
.Define("nTracks", "outputTrackParameters.size()") .Define("nTracks", "outputTrackParameters.size()")
.Define("p_track", p_track, {"outputTrackParameters"}) .Define("p_track", p_track, {"outputTrackParameters"})
.Define("p_track1", p_track, {"outputTrackParameters1"}) .Define("p_track1", p_track, {"outputTrackParameters1"})
.Define("delta_p", .Define("p_track2", p_track, {"outputTrackParameters2"})
[](const std::vector<double>& tracks, .Define("delta_p",delta_p, {"p_track", "p_thrown"})
const std::vector<double>& thrown) { .Define("delta_p1",delta_p, {"p_track1", "p_thrown"})
std::vector<double> res; .Define("delta_p2",delta_p, {"p_track2", "p_thrown"});
for (const auto &p1 : thrown) {
for (const auto &p2 : tracks) {
res.push_back(p1 - p2);
}
}
return res;
},
{"p_track", "p_thrown"})
.Define("delta_p1",
[](const std::vector<double>& tracks,
const std::vector<double>& thrown) {
std::vector<double> res;
for (const auto &p1 : thrown) {
for (const auto &p2 : tracks) {
res.push_back(p1 - p2);
}
}
return res;
},
{"p_track1", "p_thrown"});
auto h_nTracks = df0.Histo1D({"h_nTracks", "; N tracks ", 10, 0, 10}, "nTracks"); auto h_nTracks = df0.Histo1D({"h_nTracks", "; N tracks ", 10, 0, 10}, "nTracks");
auto h_pTracks = df0.Histo1D({"h_pTracks", "; GeV/c ", 100, 0, 10}, "p_track"); auto h_pTracks = df0.Histo1D({"h_pTracks", "; GeV/c ", 100, 0, 10}, "p_track");
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment