1 #include "art/Framework/Core/ModuleMacros.h"
2 #include "art/Framework/Core/EDFilter.h"
3 #include "art/Framework/Core/EDAnalyzer.h"
4 #include "art/Framework/Principal/Event.h"
5 #include "art/Framework/Principal/SubRun.h"
6 #include "art/Framework/Principal/Handle.h"
7 #include "art/Framework/Services/Registry/ServiceHandle.h"
9 #include "art_root_io/TFileService.h"
10 #include "art_root_io/TFileDirectory.h"
42 #include "nusimdata/SimulationBase/MCParticle.h"
43 #include "nusimdata/SimulationBase/MCTruth.h"
44 #include "nusimdata/SimulationBase/simb.h"
45 #include "nusimdata/SimulationBase/MCFlux.h"
46 #include "nusimdata/SimulationBase/GTruth.h"
48 #include "canvas/Utilities/ensurePointer.h"
49 #include "canvas/Persistency/Common/FindManyP.h"
50 #include "canvas/Persistency/Common/FindMany.h"
51 #include "canvas/Persistency/Common/FindOneP.h"
52 #include "canvas/Persistency/Common/FindOne.h"
54 #include "fhiclcpp/ParameterSet.h"
55 #include "messagefacility/MessageLogger/MessageLogger.h"
56 #include "cetlib_except/exception.h"
70 #include "TGraphDelaunay.h"
72 #include "TGeoPolygon.h"
100 namespace single_photon
146 bool endSubRun(art::SubRun& sr)
override;
153 std::cout<<
"SinglePhoton::"<<__FUNCTION__<<
" kicks off ------------------------------"<<std::endl;
159 paras.
s_SCE = lar::providerFrom<spacecharge::SpaceChargeService>();
161 std::cout<<
"SinglePhoton::"<<__FUNCTION__<<
" now can start jobs --------------------"<<std::endl;
169 std::cout<<
"SinglePhoton::"<<__FUNCTION__<<
" starts ------------------------------"<<std::endl;
193 std::cout<<
"SinglePhoton::reconfigure || whats configured? "<<std::endl;
217 paras.
s_pidLabel = pset.get<std::string>(
"ParticleIDLabel",
"pandoracalipidSCE");
307 std::cout<<
"SinglePhoton::"<<__FUNCTION__<<
" finishes ------------------------------"<<std::endl;
320 auto detClocks = art::ServiceHandle<detinfo::DetectorClocksService const>()->DataFor(evt);
321 auto theDetector = art::ServiceHandle<detinfo::DetectorPropertiesService const>()->DataFor(evt, detClocks);
322 std::cout<<
"SinglePhoton::"<<__FUNCTION__<<
" a new event ------------------------------"<<std::endl;
353 paras.
_time2cm =
sampling_rate(detClocks) / 1000.0 * theDetector.DriftVelocity( theDetector.Efield(), theDetector.Temperature() );
365 art::ValidHandle<std::vector<recob::PFParticle>>
const & pfParticleHandle = evt.getValidHandle<std::vector<recob::PFParticle>>(
paras.
s_pandoraLabel);
366 std::vector<art::Ptr<recob::PFParticle>> pfParticleVector;
367 art::fill_ptr_vector(pfParticleVector,pfParticleHandle);
369 if (!pfParticleHandle.isValid())
370 { mf::LogDebug(
"SinglePhoton") <<
" Failed to find the PFParticles.\n";
375 art::ValidHandle<std::vector<recob::Cluster>>
const & clusterHandle = evt.getValidHandle<std::vector<recob::Cluster>>(
paras.
s_pandoraLabel);
376 std::vector< art::Ptr<recob::Cluster> > clusterVector;
377 art::fill_ptr_vector(clusterVector,clusterHandle);
387 art::ValidHandle<std::vector<recob::Vertex>>
const & vertexHandle = evt.getValidHandle<std::vector<recob::Vertex>>(
paras.
s_pandoraLabel);
388 std::vector<art::Ptr<recob::Vertex>> vertexVector;
389 art::fill_ptr_vector(vertexVector,vertexHandle);
393 art::FindManyP<recob::Vertex> vertices_per_pfparticle(pfParticleHandle, evt,
paras.
s_pandoraLabel);
394 std::map< art::Ptr<recob::PFParticle>, std::vector<art::Ptr<recob::Vertex>> > pfParticlesToVerticesMap;
395 for(
size_t i=0; i< pfParticleVector.size(); ++i){
396 auto pfp = pfParticleVector[i];
397 pfParticlesToVerticesMap[pfp] =vertices_per_pfparticle.at(pfp.key());
401 art::FindOneP<recob::Shower> showerreco3D_per_pfparticle(pfParticleHandle, evt,
paras.
s_shower3dLabel);
405 art::FindManyP< larpandoraobj::PFParticleMetadata > pfPartToMetadataAssoc(pfParticleHandle, evt,
paras.
s_pandoraLabel);
408 art::FindManyP< recob::Track > pfPartToTrackAssoc(pfParticleHandle, evt,
paras.
s_trackLabel);
409 art::FindManyP< recob::Shower > pfPartToShowerAssoc(pfParticleHandle, evt,
paras.
s_showerLabel);
411 art::FindManyP<recob::Cluster> clusters_per_pfparticle(pfParticleHandle, evt,
paras.
s_pandoraLabel);
414 std::map<art::Ptr<recob::PFParticle>, std::vector<art::Ptr<larpandoraobj::PFParticleMetadata>> > pfParticleToMetadataMap;
415 for(
size_t i=0; i< pfParticleVector.size(); ++i){
416 const art::Ptr<recob::PFParticle> pfp = pfParticleVector[i];
417 pfParticleToMetadataMap[pfp] = pfPartToMetadataAssoc.at(pfp.key());
422 std::vector<PandoraPFParticle> allPPFParticles;
423 for(
size_t index=0; index< pfParticleVector.size(); ++index){
424 auto pfp = pfParticleVector[index];
425 int temp_key = pfp.key();
429 pfPartToMetadataAssoc.at(temp_key),
430 vertices_per_pfparticle.at(temp_key),
431 clusters_per_pfparticle.at(temp_key),
432 pfPartToShowerAssoc.at(temp_key),
433 pfPartToTrackAssoc.at(temp_key),
437 allPPFParticles.push_back(temp_pf);
445 art::ValidHandle<std::vector<recob::Slice>>
const & sliceHandle = evt.getValidHandle<std::vector<recob::Slice>>(
paras.
s_pandoraLabel);
446 std::vector<art::Ptr<recob::Slice>> sliceVector;
447 art::fill_ptr_vector(sliceVector,sliceHandle);
450 art::FindManyP<recob::PFParticle> pfparticles_per_slice(sliceHandle, evt,
paras.
s_pandoraLabel);
454 std::map< art::Ptr<recob::Slice>, std::vector<art::Ptr<recob::PFParticle>> > sliceToPFParticlesMap;
455 std::map<int, std::vector<art::Ptr<recob::PFParticle>> > sliceIDToPFParticlesMap;
456 for(
size_t i=0; i< sliceVector.size(); ++i){
457 auto slice = sliceVector[i];
458 PPFP_FindSliceIDandHits(allPPFParticles, slice, pfparticles_per_slice.at(slice.key()), hits_per_slice.at(slice.key()) );
459 sliceToPFParticlesMap[slice] =pfparticles_per_slice.at(slice.key());
460 sliceIDToPFParticlesMap[slice->ID()] = pfparticles_per_slice.at(slice.key());
470 std::map<int, std::vector<art::Ptr<recob::Hit>> > sliceIDToHitsMap;
471 for(
size_t i=0; i< sliceVector.size(); ++i){
473 auto slice = sliceVector[i];
474 sliceIDToHitsMap[slice->ID()] = hits_per_slice.at(slice.key());
483 std::vector< art::Ptr<recob::PFParticle> > nuParticles;
484 std::vector< art::Ptr<recob::PFParticle> > crParticles;
486 for(
size_t jndex=0; jndex< allPPFParticles.size(); ++jndex){
502 art::FindManyP<recob::SpacePoint> spacePoints_per_pfparticle(pfParticleHandle, evt,
paras.
s_pandoraLabel);
503 std::map<art::Ptr<recob::PFParticle>, std::vector<art::Ptr<recob::SpacePoint>> > pfParticleToSpacePointsMap;
504 for(
size_t i=0; i< nuParticles.size(); ++i){
505 const art::Ptr<recob::PFParticle> pfp = nuParticles[i];
506 pfParticleToSpacePointsMap[pfp] = spacePoints_per_pfparticle.at(pfp.key());
512 std::map<art::Ptr<recob::PFParticle>, std::vector<art::Ptr<recob::Cluster>> > pfParticleToClustersMap;
513 std::map<art::Ptr<recob::Cluster>, std::vector<art::Ptr<recob::Hit>> > clusterToHitsMap;
515 for(
size_t i=0; i< nuParticles.size(); ++i){
516 auto pfp = nuParticles[i];
517 pfParticleToClustersMap[pfp] = clusters_per_pfparticle.at(pfp.key());
520 for(
size_t i=0; i< clusterVector.size(); ++i){
521 auto cluster = clusterVector[i];
529 std::map<art::Ptr<recob::PFParticle>, std::vector<art::Ptr<recob::Hit>> > pfParticleToHitsMap;
534 for(
size_t i=0; i<nuParticles.size(); ++i){
535 auto pfp = nuParticles[i];
538 std::vector<art::Ptr<recob::Cluster>> clusters_vec = pfParticleToClustersMap[pfp] ;
541 std::vector<art::Ptr<recob::Hit>> hits_for_pfp = {};
545 for (art::Ptr<recob::Cluster>
cluster: clusters_vec){
546 std::vector<art::Ptr<recob::Hit>> hits_vec = clusterToHitsMap[
cluster];
549 hits_for_pfp.insert( hits_for_pfp.end(), hits_vec.begin(), hits_vec.end() );
553 pfParticleToHitsMap[pfp] = hits_for_pfp;
563 std::map<art::Ptr<recob::PFParticle>, std::vector<art::Ptr<recob::Cluster>> > cr_pfParticleToClustersMap;
564 std::map<art::Ptr<recob::PFParticle>, std::vector<art::Ptr<recob::Hit>> > cr_pfParticleToHitsMap;
580 for(
size_t i=0; i< crParticles.size(); ++i){
581 auto pfp = crParticles[i];
582 cr_pfParticleToClustersMap[pfp] = clusters_per_pfparticle.at(pfp.key());
585 for(
size_t i=0; i< crParticles.size(); ++i){
586 auto pfp = crParticles[i];
590 std::vector<art::Ptr<recob::Cluster>> clusters_vec = cr_pfParticleToClustersMap[pfp] ;
593 std::vector<art::Ptr<recob::Hit>> hits_for_pfp = {};
598 for (art::Ptr<recob::Cluster>
cluster: clusters_vec){
599 std::vector<art::Ptr<recob::Hit>> hits_vec = clusterToHitsMap[
cluster];
603 hits_for_pfp.insert( hits_for_pfp.end(), hits_vec.begin(), hits_vec.end() );
607 cr_pfParticleToHitsMap[pfp] = hits_for_pfp;
624 std::vector< art::Ptr<recob::Track> >
tracks;
625 std::vector< art::Ptr<recob::Shower> > showers;
626 std::map< art::Ptr<recob::Track> , art::Ptr<recob::PFParticle >> trackToNuPFParticleMap;
627 std::map< art::Ptr<recob::Shower> , art::Ptr<recob::PFParticle>> showerToNuPFParticleMap;
631 for(
size_t jndex=0; jndex< allPPFParticles.size(); ++jndex){
635 tracks.push_back(temp_pf->
pTrack);
639 showers.push_back(temp_pf->
pShower);
650 art::ValidHandle<std::vector<recob::Track>>
const & trackHandle = evt.getValidHandle<std::vector<recob::Track>>(
paras.
s_trackLabel);
652 art::FindManyP<anab::Calorimetry> calo_per_track(trackHandle, evt,
paras.
s_caloLabel);
655 if (!calo_per_track.isValid())
657 mf::LogDebug(
"SinglePhoton") <<
" Failed to get Assns between recob::Track and anab::Calorimetry.\n";
661 art::FindOneP<anab::ParticleID> pid_per_track(trackHandle, evt,
paras.
s_pidLabel);
662 std::map<art::Ptr<recob::Track>, art::Ptr<anab::ParticleID> > trackToPIDMap;
664 for(
size_t i=0; i< tracks.size(); ++i){
665 if(calo_per_track.at(tracks[i].key()).
size() ==0){
666 std::cerr<<
"Track Calorimetry Breaking! the vector of calo_per_track is of length 0 at this track."<<std::endl;
684 for(
size_t i=0; i< tracks.size(); ++i){
685 trackToPIDMap[tracks[i]] = pid_per_track.at(tracks[i].key());
698 std::vector<art::Ptr<simb::MCTruth>> mcTruthVector;
699 std::vector<art::Ptr<simb::MCParticle>> mcParticleVector;
702 std::map< art::Ptr<simb::MCParticle>, std::vector<art::Ptr<recob::Hit> > > mcParticleToHitsMap;
703 std::map< art::Ptr<recob::Hit>, art::Ptr<simb::MCParticle> > hitToMCParticleMap;
708 std::map< art::Ptr<simb::MCTruth>, std::vector<art::Ptr<simb::MCParticle>>> MCTruthToMCParticlesMap;
709 std::map< art::Ptr<simb::MCParticle>, art::Ptr<simb::MCTruth>> MCParticleToMCTruthMap;
710 std::map<int, art::Ptr<simb::MCParticle> > MCParticleToTrackIdMap;
712 std::vector<art::Ptr<sim::MCTrack>> mcTrackVector;
713 std::vector<art::Ptr<sim::MCShower>> mcShowerVector;
715 std::vector<art::Ptr<simb::MCParticle>> matchedMCParticleVector;
716 std::map<art::Ptr<recob::Track>, art::Ptr<simb::MCParticle> > trackToMCParticleMap;
717 std::map<art::Ptr<recob::Shower>, art::Ptr<simb::MCParticle> > showerToMCParticleMap;
720 std::map< art::Ptr<simb::MCParticle>, art::Ptr<sim::MCTrack> > MCParticleToMCTrackMap;
721 std::map< art::Ptr<simb::MCParticle>, art::Ptr<sim::MCShower> > MCParticleToMCShowerMap;
724 std::cout <<
"SinglePhoton::analyze()\t||\t Consolidated event summary:" <<
"\n";
725 std::cout <<
"SinglePhoton::analyze()\t||\t - Number of primary cosmic-ray PFParticles : " << crParticles.size() <<
"\n";
726 std::cout <<
"SinglePhoton::analyze()\t||\t - Number of neutrino final-state PFParticles : " << nuParticles.size() <<
"\n";
727 std::cout <<
"SinglePhoton::analyze()\t||\t ... of which are track-like : " << tracks.size() <<
"\n";
728 std::cout <<
"SinglePhoton::analyze()\t||\t ... of which are showers-like : " << showers.size() <<
"\n";
743 std::vector<std::pair<art::Ptr<recob::PFParticle>,
int>> primaryPFPSliceIdVec;
744 std::map<int, double> sliceIdToNuScoreMap;
745 std::map<art::Ptr<recob::PFParticle>,
bool> PFPToClearCosmicMap;
747 std::map<art::Ptr<recob::PFParticle>,
bool> PFPToNuSliceMap;
748 std::map<art::Ptr<recob::PFParticle>,
double> PFPToTrackScoreMap;
749 std::map<int, int> sliceIdToNumPFPsMap;
757 art::ValidHandle<std::vector<recob::OpFlash>>
const & flashHandle = evt.getValidHandle<std::vector<recob::OpFlash>>(
paras.
s_flashLabel);
758 std::vector<art::Ptr<recob::OpFlash>> flashVector;
759 art::fill_ptr_vector(flashVector,flashHandle);
761 std::map<art::Ptr<recob::OpFlash>, std::vector< art::Ptr<sbn::crt::CRTHit>>> crtvetoToFlashMap;
764 art::FindManyP<sbn::crt::CRTHit> crtveto_per_flash(flashHandle, evt,
paras.
s_CRTVetoLabel);
765 for(
size_t i=0; i< flashVector.size(); ++i){
766 crtvetoToFlashMap[flashVector[i]] = crtveto_per_flash.at(flashVector[i].key());
770 art::Handle<std::vector<sbn::crt::CRTHit>> crthit_h;
772 double evt_timeGPS_nsec = -999 ;
779 art::Timestamp evtTimeGPS = evt.time();
780 evt_timeGPS_nsec = evtTimeGPS.timeLow();
783 std::cout<<
"SinglePhoton::analyze \t||\t Got CRT hits"<<std::endl;
811 std::cout<<
"\nSinglePhoton::analyze \t||\t Start on Track Analysis "<<std::endl;
817 pfParticleToSpacePointsMap, MCParticleToTrackIdMap, sliceIdToNuScoreMap,
827 std::cout<<
"\nSinglePhoton::analyze \t||\t Start on Shower Analysis "<<std::endl;
840 std::cout<<
"\nSinglePhoton::analyze \t||\t Finish Shower Analysis "<<std::endl;
843 std::cout<<
"SinglePhoton::Reconstruction Summary ---------------------------------------------- "<<std::endl;
844 std::vector<int> spacers =
Printer_header({
" pfpID",
" Parent ",
"Vertex(x, ",
" y, ",
", z )",
"slice ",
" nu_score",
" trk_score ",
" prim? ",
" nu? ",
" nuSlice?",
" cos? ",
" S# ",
" T# "});
845 for(
size_t jndex=0; jndex< allPPFParticles.size(); ++jndex){
847 art::Ptr<recob::PFParticle> pfp = temp_pf.
pPFParticle;
858 (pfp->IsPrimary())?
"T":
" ",
874 for(
size_t i_shr = 0; i_shr<showers.size();i_shr++){
875 const art::Ptr<recob::Shower>
s = showers[i_shr];
876 const art::Ptr<recob::PFParticle> pfp = showerToNuPFParticleMap[
s];
877 const std::vector< art::Ptr<recob::SpacePoint> > shr_spacepoints = pfParticleToSpacePointsMap[pfp];
882 for(
auto &sp: shr_spacepoints){
883 std::vector<double> tmp_spt = {sp->XYZ()[0],sp->XYZ()[1] , sp->XYZ()[2]};
890 if(showers.size()==1){
901 art::ValidHandle<std::vector<recob::Hit>>
const & hitHandle = evt.getValidHandle<std::vector<recob::Hit>>(
paras.
s_hitfinderLabel);
902 std::vector<art::Ptr<recob::Hit>> hitVector;
903 art::fill_ptr_vector(hitVector,hitHandle);
928 art::ValidHandle<std::vector<simb::MCTruth>>
const & mcTruthHandle= evt.getValidHandle<std::vector<simb::MCTruth>>(
paras.
s_generatorLabel);
929 art::fill_ptr_vector(mcTruthVector,mcTruthHandle);
932 art::ValidHandle<std::vector<simb::MCParticle>>
const & mcParticleHandle= evt.getValidHandle<std::vector<simb::MCParticle>>(
paras.
s_geantModuleLabel);
933 art::fill_ptr_vector(mcParticleVector,mcParticleHandle);
938 if(
g_is_verbose)
std::cout<<
"SinglePhoton::AnalyzeGeant4()\t||\t Begininning recob::Geant4 analysis suite"<<std::endl;;
947 std::vector<art::Ptr<simb::MCParticle>> particle_vec;
948 std::vector<anab::BackTrackerHitMatchingData const *> match_vec;
951 for(
size_t j=0; j<hitVector.size();j++){
952 const art::Ptr<recob::Hit>
hit = hitVector[j];
954 particle_vec.clear(); match_vec.clear();
955 mcparticles_per_hit.get(hit.key(), particle_vec, match_vec);
956 if(particle_vec.size() > 0){
967 std::cout<<
"\nSinglePhoton\t||\t Starting backtracker on recob::track"<<std::endl;
969 trackToMCParticleMap,
970 trackToNuPFParticleMap,
973 matchedMCParticleVector,
976 std::cout<<
"\nSinglePhoton\t||\t Starting backtracker on recob::shower"<<std::endl;
981 showerToMCParticleMap,
983 matchedMCParticleVector,
984 MCParticleToTrackIdMap,
989 std::cout<<
"\nSinglePhoton\t||\tStarting RecoMCTracks "<<std::endl;
993 trackToMCParticleMap,
994 MCParticleToMCTruthMap,
996 MCParticleToTrackIdMap,
1000 std::cout<<
"\nSinglePhoton\t||\tStarting AnalyzeMCTruths"<<std::endl;
1009 std::cout<<
"\nSinglePhoton\t||\tStarting AnalyzeEventWeight"<<std::endl;
1013 std::cout<<
"\nSinglePhoton\t||\tStarting AnalyzeRecoMCSlices"<<std::endl;
1016 std::cout<<
"\nSinglePhoton\t||\tFinish AnalyzeRecoMCSlices"<<std::endl;
1022 art::Handle<std::vector<evwgh::MCEventWeight>> ev_evw_ph ;
1023 if( evt.getByLabel(
"eventweight",ev_evw_ph)){
1024 std::map<std::string, std::vector<double>>
const & weight_map = ev_evw_ph->front().fWeight;
1025 for (
auto const&
x : weight_map){
1028 <<
x.second.size() << std::endl ;
1029 if(
x.first ==
"photonuclear_photon_PhotoNuclear"){
1030 auto vec =
x.second;
1031 double ph_low = vec[1];
1032 double ph_high = vec[0];
1033 std::cout<<
"PhotoNuBit: "<<ph_low<<
" "<<ph_high<<std::endl;
1060 std::cout<<
"SinglePhoton::analyze\t||\t finnished loop for this event"<<std::endl;
1067 std::cout<<
"\nSinglePhoton::analyze \t||\t Start IsolationStudy "<<std::endl;
1069 std::cout<<
"\nSinglePhoton::analyze \t||\t Finish IsolationStudy "<<std::endl;
1078 std::cout <<
"----------------- Stub clustering --------------------------- " << std::endl;
1087 art::Ptr<recob::Shower> p_shr = showers.front();
1089 art::Ptr<recob::PFParticle> p_pfp = ppfp->
pPFParticle;
1090 std::vector<art::Ptr<recob::Hit>> p_hits = pfParticleToHitsMap[p_pfp];
1093 auto p_slice_hits = sliceIDToHitsMap[p_sliceid];
1124 for(
auto &trk: tracks){
1125 art::Ptr<recob::PFParticle> p_pfp_trk = trackToNuPFParticleMap[trk];
1126 std::vector<art::Ptr<recob::Hit>> p_hits_trk = pfParticleToHitsMap[p_pfp_trk];
1134 for(
auto &cr: crParticles){
1135 std::vector<art::Ptr<recob::Hit>> p_hits_cr = cr_pfParticleToHitsMap[cr];
1154 std::vector<seaview::cluster> vec_SEAclusters ;
1159 std::cout<<
"After SEAview we have "<<vec_SEAclusters.size()<<
" Stub clusters to chat about"<<std::endl;
1162 for(
size_t c=0; c< vec_SEAclusters.size(); c++){
1163 auto& clu = vec_SEAclusters.at(c);
1164 int pl = clu.getPlane();
1165 auto hitz = clu.getHits();
1167 int remerge = clu.getShowerRemerge();
1170 std::cout<<c<<
" "<<pl<<
" "<<Ep<<
" "<<clu.getMinHitImpactParam()<<
" "<<clu.getMinHitConvDist()<<
" "<<clu.getMinHitIOC()<<
" "<<clu.getMeanADC()<<
" "<<clu.getADCrms()<<
" "<< clu.getLinearChi() <<
" " << remerge<<std::endl;
1239 std::cout <<
"===============================================================" << std::endl;
1259 std::cout <<
"------------- Shower clustering --------------------" << std::endl;
1267 art::Ptr<recob::Shower> p_shr = showers.front();
1270 art::Ptr<recob::PFParticle> p_pfp = ppfp->
pPFParticle;
1271 std::vector<art::Ptr<recob::Hit>> p_hits = pfParticleToHitsMap[p_pfp];
1275 auto p_slice_hits = sliceIDToHitsMap[p_sliceid];
1304 for(
auto &trk: tracks){
1305 art::Ptr<recob::PFParticle> p_pfp_trk = trackToNuPFParticleMap[trk];
1306 std::vector<art::Ptr<recob::Hit>> p_hits_trk = pfParticleToHitsMap[p_pfp_trk];
1331 std::vector<seaview::cluster> vec_SEAclusters ;
1335 std::cout<<
"After SEAview we have "<<vec_SEAclusters.size()<<
" Shower clusters to chat about"<<std::endl;
1338 for(
size_t c=0; c< vec_SEAclusters.size(); c++){
1339 auto clu = vec_SEAclusters.at(c);
1340 int pl = clu.getPlane();
1341 auto hitz = clu.getHits();
1343 int remerge = clu.getShowerRemerge();
1346 std::cout<<c<<
" "<<pl<<
" "<<Ep<<
" "<<clu.getImpactParam()<<
" "<<clu.getFitSlope()<<
" "<<clu.getFitCons()<<
" "<<clu.getMeanADC() <<
" " << clu.getADCrms() <<
" "<<clu.getAngleWRTShower()<<
" "<<remerge<<std::endl;
1409 std::cout <<
"===============================================================" << std::endl;
1424 for(
int i =0; i<(int)showers.size(); i++){
1436 std::cout<<
"------------ Shower3D --------------"<<std::endl;
1483 mf::LogDebug(
"SinglePhoton") <<
" *** beginJob() *** " <<
"\n";
1485 art::ServiceHandle<art::TFileService>
tfs;
1489 vars.
pot_tree = tfs->make<TTree>(
"pot_tree",
"pot_tree");
1527 std::cout <<
"SinglePhoton \t||\t Running in selected-event only mode " << std::endl;
1538 while(std::getline(infile, line)){
1539 std::istringstream ss(line);
1541 std::vector<int> event_info;
1542 for(
int i; ss >> i; ) event_info.push_back(i);
1550 std::cout <<
"Selected Events: " << std::endl;
1551 std::cout <<
"Run \t SubRun \t Event" << std::endl;
1553 std::for_each(v.begin(), v.end(), [](
int n){
std::cout <<
n<<
" \t "; });
1567 double this_pot = 0;
1576 art::Handle<sumdata::POTSummary> gen_pot_hand;
1578 this_pot = gen_pot_hand->totgoodpot;
1584 art::Handle<sumdata::POTSummary> potSummaryHandlebnbETOR875;
1585 if (sr.getByLabel(
"beamdata",
"bnbETOR875",potSummaryHandlebnbETOR875)){
1586 this_pot =potSummaryHandlebnbETOR875->totpot;
1588 std::cout<<
"SinglePhoton::beginSubRun()\t||\t SubRun POT: "<<potSummaryHandlebnbETOR875->totpot<<
" . Current total POT this file: "<<
vars.
m_pot_count<<
" (label) "<<
paras.
s_potLabel<<std::endl;
std::vector< double > m_sss_candidate_PCA
std::vector< double > m_sss_candidate_overlay_fraction
void AnalyzeTrackCalo(const std::vector< art::Ptr< recob::Track >> &tracks, std::vector< PandoraPFParticle > all_PPFPs, var_all &vars, para_all ¶s)
TTree * ncdelta_slice_tree
void CreateShowerBranches(var_all &vars)
std::vector< double > m_trackstub_candidate_PCA
std::vector< double > SecondShowerMatching(std::vector< art::Ptr< recob::Hit >> &hitz, art::FindManyP< simb::MCParticle, anab::BackTrackerHitMatchingData > &mcparticles_per_hit, std::vector< art::Ptr< simb::MCParticle >> &mcParticleVector, std::map< int, art::Ptr< simb::MCParticle >> &MCParticleToTrackIdMap, var_all &vars)
std::string s_hitMCParticleAssnsLabel
int loadVertex(double m_vertex_pos_x, double m_vertex_pos_y, double m_vertex_pos_z)
std::vector< double > m_trackstub_candidate_mean_ADC_first_to_second_ratio
std::vector< double > m_trackstub_candidate_max_wire
const bool get_IsNuSlice() const
std::vector< double > m_sss_candidate_wire_tick_based_length
bool beginSubRun(art::SubRun &sr) override
: grab run, subrun number, and subrun POT, fill the TTree
bool endSubRun(art::SubRun &sr) override
double s_SEAviewPlotDistance
std::vector< double > m_reco_track_spacepoint_principal0
std::string s_showerLabel
The label for the shower producer from PFParticles.
std::vector< double > m_trackstub_candidate_linear_fit_chi2
std::vector< int > m_trackstub_candidate_in_nu_slice
void CreateEventWeightBranches(var_all &vars)
PandoraPFParticle * PPFP_GetPPFPFromPFID(std::vector< PandoraPFParticle > &PPFPs, int id)
void CollectMCParticles(const art::Event &evt, const std::string &label, std::map< art::Ptr< simb::MCTruth >, std::vector< art::Ptr< simb::MCParticle >>> &truthToParticles, std::map< art::Ptr< simb::MCParticle >, art::Ptr< simb::MCTruth >> &particlesToTruth, std::map< int, art::Ptr< simb::MCParticle > > &MCParticleToTrackIdMap, var_all &vars)
void CreateFlashBranches(var_all &vars)
Utilities related to art service access.
void ClearFlashes(var_all &vars)
std::vector< double > m_sss_candidate_ADC_RMS
const int get_HasTrack() const
std::set< std::vector< int > > m_selected_set
const int get_AncestorID() const
void CreateSecondShowerBranches3D(var_all &vars)
double s_beamgate_flash_start
void set_ParticleID(const art::Ptr< anab::ParticleID > input_ParticleID)
void ClearMCTruths(var_all &vars)
std::vector< double > m_trackstub_candidate_track_angle_wrt_shower_direction
ClusterModuleLabel join with tracks
std::vector< double > m_sss_candidate_fit_constant
int addPFParticleHits(std::vector< art::Ptr< recob::Hit >> &hits, std::string leg, double arg1=0.0, double arg2=0.0, double arg3=0.0)
void AnalyzeEventWeight(art::Event const &e, var_all &vars)
art::Ptr< recob::Track > pTrack
int s_SEAviewNumRecoTrack
std::vector< int > m_sss_candidate_remerge
std::vector< double > m_trackstub_candidate_mean_tick
std::vector< int > m_sss_candidate_plane
std::vector< double > m_trackstub_candidate_mean_ADC_first_half
process_name opflash particleana ie x
BEGIN_PROLOG could also be cerr
std::vector< double > m_sss_candidate_impact_parameter
std::vector< double > trackRecoMCmatching(std::vector< art::Ptr< recob::Track >> &objectVector, std::map< art::Ptr< recob::Track >, art::Ptr< simb::MCParticle >> &objectToMCParticleMap, std::map< art::Ptr< recob::Track >, art::Ptr< recob::PFParticle >> &objectToPFParticleMap, std::map< art::Ptr< recob::PFParticle >, std::vector< art::Ptr< recob::Hit >> > &pfParticleToHitsMap, art::FindManyP< simb::MCParticle, anab::BackTrackerHitMatchingData > &mcparticles_per_hit, std::vector< art::Ptr< simb::MCParticle >> &mcParticleVector, var_all &vars)
const double * get_Vertex_pos() const
double s_track_calo_min_dEdx
std::vector< int > m_sss_candidate_num_wires
std::vector< double > m_reco_shower_impact_parameter
std::vector< double > m_reco_shower_end_dist_to_SCB
std::vector< double > m_trackstub_candidate_group_timeoverlap_fraction
std::vector< int > m_trackstub_candidate_remerge
void ClearGeant4Branches(var_all &vars)
: fill event weight related variables
int m_number_of_events_in_subrun
std::vector< double > analyzeShowerLikeClusters(double eps, const std::map< art::Ptr< recob::Shower >, art::Ptr< recob::PFParticle >> &showerToPFParticleMap, const std::map< art::Ptr< recob::PFParticle >, std::vector< art::Ptr< recob::Hit >> > &pfParticleToHitsMap, std::vector< seaview::cluster > &vec_c)
std::vector< int > m_trackstub_candidate_plane
spacecharge::SpaceCharge const * s_SCE
void CreateTrackBranches(var_all &vars)
std::string s_caloLabel
The label for calorimetry associations producer.
int DefineNuSlice(std::vector< PandoraPFParticle > &PPFPs)
Declaration of signal hit object.
int m_trackstub_num_unassociated_hits
std::vector< double > m_sss_candidate_angle_to_shower
void BuildMCParticleHitMaps(const art::Event &evt, const std::string &label, const std::vector< art::Ptr< recob::Hit >> &hitVector, std::map< art::Ptr< simb::MCParticle >, std::vector< art::Ptr< recob::Hit > > > &particlesToHits, std::map< art::Ptr< recob::Hit >, art::Ptr< simb::MCParticle > > &hitsToParticles, const lar_pandora::LArPandoraHelper::DaughterMode daughterMode, std::map< int, art::Ptr< simb::MCParticle > > &MCParticleToTrackIdMap, var_all &vars)
void CreateSecondShowerBranches(var_all &vars)
double m_photonu_weight_low
void ClearSecondShowers(var_all &vars)
bool s_SEAviewStubMakePDF
art::Ptr< recob::Shower > pShower
std::vector< int > m_trackstub_candidate_pdg
int setHitThreshold(double)
double s_SEAviewStubDbscanEps
double s_track_calo_trunc_fraction
std::vector< int > m_sss_candidate_trackid
std::vector< int > m_sss_candidate_num_hits
void CreateGeant4Branches(var_all &vars)
std::vector< int > m_trackstub_candidate_parent_pdg
art::Ptr< recob::PFParticle > pPFParticle
std::string s_pandoraLabel
The label for the pandora producer.
std::pair< int, std::pair< std::vector< std::vector< double > >, std::vector< double > > > GroupClusterCandidate(int num_clusters, const std::vector< int > &cluster_planes, const std::vector< double > &cluster_max_ticks, const std::vector< double > &cluster_min_ticks)
std::size_t size(FixedBins< T, C > const &) noexcept
int addTrack(art::Ptr< recob::Track > &trk)
Contains data associated to particles from detector simulation.
int m_trackstub_num_candidates
int m_trackstub_num_candidate_groups
void AnalyzeGeant4(const std::vector< art::Ptr< simb::MCParticle >> &mcParticleVector, var_all &vars)
void RecoMCTracks(std::vector< PandoraPFParticle > all_PPFPs, const std::vector< art::Ptr< recob::Track >> &tracks, std::map< art::Ptr< recob::Track >, art::Ptr< simb::MCParticle > > &trackToMCParticleMap, std::map< art::Ptr< simb::MCParticle >, art::Ptr< simb::MCTruth >> &MCParticleToMCTruthMap, std::vector< art::Ptr< simb::MCParticle >> &mcParticleVector, std::map< int, art::Ptr< simb::MCParticle > > &MCParticleToTrackIdMap, std::vector< double > &vfrac, var_all &vars)
std::string s_hitfinderLabel
int m_trackstub_unassociated_hits_below_threshold
double s_beamgate_flash_end
std::vector< int > m_sss_candidate_pdg
void showerRecoMCmatching(std::vector< PandoraPFParticle > all_PPFPs, std::vector< art::Ptr< recob::Shower >> &showerVector, std::map< art::Ptr< recob::Shower >, art::Ptr< simb::MCParticle >> &showerToMCParticleMap, art::FindManyP< simb::MCParticle, anab::BackTrackerHitMatchingData > &mcparticles_per_hit, std::vector< art::Ptr< simb::MCParticle >> &mcParticleVector, std::map< int, art::Ptr< simb::MCParticle > > &MCParticleToTrackIdMap, var_all &vars)
std::string s_CRTHitProducer
int s_SEAviewStubNumRecoTrack
std::vector< double > m_trackstub_candidate_min_dist
void beginJob() override
Begin the job, setting up !
std::vector< double > m_reco_shower_reclustered_energy_plane0
void PPFP_FindAncestor(std::vector< PandoraPFParticle > &PPFPs)
void CreateStubBranches(var_all &vars)
void Save_EventMeta(art::Event &evt, var_all &vars)
std::vector< double > m_trackstub_candidate_ADC_RMS
std::vector< double > m_trackstub_candidate_min_tick
int s_SEAviewStubNumRecoShower
void ClearEventWeightBranches(var_all &vars)
double s_SEAviewDbscanEps
std::vector< double > s_gain_data
bool filter(art::Event &evt) override
Analyze an event!
std::vector< double > s_gain_mc
void setTPCGeom(para_all ¶s)
bool Pi0PreselectionFilter2g0p(var_all &vars, para_all ¶s)
std::vector< int > m_trackstub_candidate_trackid
double s_exiting_photon_energy_threshold
int m_sss_num_unassociated_hits_below_threshold
int s_SEAviewNumRecoShower
std::string s_trackLabel
The label for the track producer from PFParticles.
std::vector< double > m_reco_shower_reclustered_energy_plane2
std::vector< double > m_trackstub_candidate_mean_ADC_second_half
std::vector< double > m_sss_candidate_true_energy
double s_track_calo_min_dEdx_hits
SinglePhoton(fhicl::ParameterSet const &pset)
Constructor.
bool Pi0PreselectionFilter(var_all &vars, para_all ¶s)
std::vector< double > m_trackstub_candidate_matched_energy_fraction_best_plane
std::vector< int > Printer_header(std::vector< std::string > headings)
std::vector< double > m_reco_shower_spacepoint_x
std::string s_generatorLabel
std::vector< double > m_reco_shower_energy_max
void IsolationStudy(std::vector< PandoraPFParticle > all_PPFPs, const std::vector< art::Ptr< recob::Track >> &tracks, const std::vector< art::Ptr< recob::Shower >> &showers, detinfo::DetectorPropertiesData const &theDetector, var_all &vars, para_all ¶s)
std::vector< double > m_reco_shower_end_dist_to_active_TPC
void CreateMCTruthBranches(var_all &vars)
std::vector< double > m_sss_candidate_energy
std::vector< double > m_reco_shower_reclustered_energy_max
const double get_TrackScore() const
double s_recombination_factor
std::string s_CRTVetoLabel
std::string s_true_eventweight_label
bool s_runTrueEventweight
std::string s_selected_event_list
std::vector< double > m_trackstub_candidate_overlay_fraction
geo::GeometryCore const * s_geom
std::string s_shower3dLabel
void CreateSliceBranches(var_all &vars)
std::vector< double > m_sss_candidate_mean_tick
std::vector< double > m_trackstub_candidate_max_tick
std::vector< double > m_sss_candidate_max_wire
std::vector< double > m_sss_candidate_min_dist
std::vector< double > analyzeTrackLikeClusters(double eps, const std::map< art::Ptr< recob::Shower >, art::Ptr< recob::PFParticle >> &showerToPFParticleMap, const std::map< art::Ptr< recob::PFParticle >, std::vector< art::Ptr< recob::Hit >> > &pfParticleToHitsMap, std::vector< seaview::cluster > &vec_c)
int distToSCB(double &dist, std::vector< double > &vec, para_all ¶s)
int m_trackstub_associated_hits
std::vector< int > m_trackstub_candidate_num_ticks
std::vector< std::vector< double > > m_grouped_trackstub_candidate_indices
void ClearSecondShowers3D(var_all &vars)
Declaration of cluster object.
Class def header for mctrack data container.
int Print(double plot_distance)
PandoraPFParticle * PPFP_GetPPFPFromTrack(std::vector< PandoraPFParticle > &PPFPs, art::Ptr< recob::Track > pTrack)
void AnalyzeTracks(std::vector< PandoraPFParticle > all_PPFPs, const std::vector< art::Ptr< recob::Track >> &tracks, std::map< art::Ptr< recob::PFParticle >, std::vector< art::Ptr< recob::SpacePoint >>> &pfParticleToSpacePointsMap, std::map< int, art::Ptr< simb::MCParticle > > &MCParticleToTrackIdMap, std::map< int, double > &sliceIdToNuScoreMap, var_all &vars, para_all ¶s)
void set_Calorimetries(const std::vector< art::Ptr< anab::Calorimetry >> input_Calorimetries)
TTree * true_eventweight_tree
std::vector< double > m_sss_candidate_matched_energy_fraction_best_plane
std::vector< double > m_trackstub_candidate_wire_tick_based_length
Provides recob::Track data product.
int addShower(art::Ptr< recob::Shower > &shr)
double s_exiting_proton_energy_threshold
int addAllHits(std::vector< art::Ptr< recob::Hit >> &hits)
void PPFP_FindSliceIDandHits(std::vector< PandoraPFParticle > &PPFPs, art::Ptr< recob::Slice > slice, const std::vector< art::Ptr< recob::PFParticle > > PFP_in_slice, const std::vector< art::Ptr< recob::Hit > > Hit_inslice)
const int get_SliceID() const
std::string s_Spline_CV_label
void ResizeShowers(size_t size, var_all &vars)
std::vector< double > m_reco_shower_spacepoint_z
std::vector< double > m_sss_candidate_fit_slope
void reconfigure(fhicl::ParameterSet const &pset)
Configure memeber variables using FHiCL parameters.
double CalcEShowerPlane(const std::vector< art::Ptr< recob::Hit >> &hits, int this_plane, para_all ¶s)
void endJob() override
End the job, setting down !
void ClearMeta(var_all &vars)
: reset/clear data members
std::vector< int > m_sss_candidate_parent_pdg
std::vector< double > m_trackstub_candidate_ioc_based_length
void ClearStubs(var_all &vars)
double s_SEAviewStubDbscanMinPts
std::vector< double > m_trackstub_candidate_min_impact_parameter_to_shower
void ResizeMCTruths(size_t size, var_all &vars)
std::string s_CRTTzeroLabel
void Printer_content(std::vector< std::string > nums, std::vector< int > spacers)
std::vector< double > m_trackstub_candidate_true_energy
std::string to_string(WindowPattern const &pattern)
void ClearIsolation(var_all &vars)
then echo File list $list not found else cat $list while read file do echo $file sed s
bool s_run_pi0_filter_2g0p
void Save_PFParticleInfo(std::vector< PandoraPFParticle > PPFPs, var_all &vars, para_all ¶s)
double s_SEAviewStubPlotDistance
PandoraPFParticle * PPFP_GetPPFPFromShower(std::vector< PandoraPFParticle > &PPFPs, art::Ptr< recob::Shower > pShower)
void SetClusterLegend(int cluster, double energy, int is_matched, int matched_pdg, double overlay_fraction)
double s_SEAviewMaxPtsLinFit
double s_SEAviewStubHitThreshold
object containing MC truth information necessary for making RawDigits and doing back tracking ...
int runseaDBSCAN(double min_pts, double eps)
int m_sss_num_associated_hits
Class def header for MCShower data container.
const double get_NuScore() const
int trigger_offset(DetectorClocksData const &data)
std::vector< double > m_trackstub_candidate_min_ioc_to_shower_start
bool IsEventInList(int run, int subrun, int event, var_all &vars)
void SimpleSecondShowerCluster(var_all &vars, para_all ¶s)
double s_SEAviewDbscanMinPts
void set_HasPID(const bool input_bool)
double m_photonu_weight_high
std::vector< double > m_reco_shower_spacepoint_y
std::vector< int > m_sss_candidate_matched
const int get_HasShower() const
const bool get_IsClearCosmic() const
art::ServiceHandle< art::TFileService > tfs
std::vector< int > m_trackstub_candidate_num_wires
int addHitsToConsider(std::vector< art::Ptr< recob::Hit >> &hits)
void CollectPID(std::vector< art::Ptr< recob::Track >> &tracks, std::vector< PandoraPFParticle > all_PPFPs, var_all &vars)
void ResizeSlices(size_t size, var_all &vars)
std::vector< double > m_trackstub_candidate_min_wire
const bool get_IsNeutrino() const
double distToTPCActive(std::vector< double > &vec, para_all ¶s)
std::vector< int > m_sss_candidate_in_nu_slice
std::vector< double > m_reco_shower_reclustered_energy_plane1
std::string s_pidLabel
For PID stuff.
void AnalyzeShowers(std::vector< PandoraPFParticle > all_PPFPs, const std::vector< art::Ptr< recob::Shower >> &showers, std::map< art::Ptr< recob::Cluster >, std::vector< art::Ptr< recob::Hit >> > &clusterToHitMap, double triggeroffset, detinfo::DetectorPropertiesData const &theDetector, var_all &vars, para_all ¶s)
bool s_run_pi0_filter_2g1p
void ResizeFlashes(size_t size, var_all &vars)
std::string s_truthmatching_signaldef
bool s_use_PID_algorithms
double sampling_rate(DetectorClocksData const &data)
Returns the period of the TPC readout electronics clock.
std::vector< double > m_trackstub_candidate_energy
std::vector< int > calcUnassociatedHits()
std::vector< double > m_trackstub_candidate_mean_wire
std::vector< double > m_trackstub_candidate_min_conversion_dist_to_shower_start
void SecondShowerSearch3D(std::vector< art::Ptr< recob::Shower >> &showers, std::map< art::Ptr< recob::Shower >, art::Ptr< recob::PFParticle >> &NormalShowerToPFParticleMap, std::vector< art::Ptr< recob::Track >> &tracks, std::map< art::Ptr< recob::Track >, art::Ptr< recob::PFParticle >> &NormalTrackToPFParticleMap, art::Event const &evt, var_all &vars, para_all ¶s)
void CreateIsolationBranches(var_all &vars)
std::vector< int > m_trackstub_candidate_num_hits
std::vector< double > m_reco_shower_conversion_distance
void AnalyzeMCTruths(std::vector< art::Ptr< simb::MCTruth >> &mcTruthVector, std::vector< art::Ptr< simb::MCParticle >> &mcParticleVector, var_all &vars, para_all ¶s)
void AnalyzeRecoMCSlices(std::string signal_def, std::vector< PandoraPFParticle > all_PPFPs, std::map< int, art::Ptr< simb::MCParticle >> &MCParticleToTrackIDMap, std::map< art::Ptr< recob::Shower >, art::Ptr< simb::MCParticle > > &showerToMCParticleMap, std::map< art::Ptr< recob::Track >, art::Ptr< simb::MCParticle > > &trackToMCParticleMap, var_all &vars, para_all ¶s)
std::vector< double > m_sss_candidate_max_tick
std::vector< double > m_trackstub_candidate_mean_ADC
void AnalyzeFlashes(const std::vector< art::Ptr< recob::OpFlash >> &flashes, art::Handle< std::vector< sbn::crt::CRTHit >> crthit_h, double evt_timeGPS_nsec, std::map< art::Ptr< recob::OpFlash >, std::vector< art::Ptr< sbn::crt::CRTHit >>> crtvetoToFlashMap, var_all &vars, para_all ¶s)
art framework interface to geometry description
BEGIN_PROLOG could also be cout
std::vector< int > m_sss_candidate_num_ticks
int m_sss_num_unassociated_hits
std::vector< double > m_sss_candidate_mean_wire
std::vector< int > m_trackstub_candidate_matched
std::string s_geantModuleLabel
void ClearSlices(var_all &vars)
std::vector< double > m_reco_track_length
std::vector< double > m_sss_candidate_mean_ADC
double s_SEAviewHitThreshold
std::vector< double > m_sss_candidate_min_tick
void ClearShowers(var_all &vars)
int filterConsideredHits(double dist_to_vertex)
std::vector< double > m_sss_candidate_min_wire
void CreateMetaBranches(var_all &vars)
double s_track_calo_max_dEdx
void ClearTracks(var_all &vars)