2121#include " FWCore/Framework/interface/MakerMacros.h"
2222#include " FWCore/ParameterSet/interface/ParameterSet.h"
2323#include " FWCore/ServiceRegistry/interface/Service.h"
24+ #include " FWCore/Utilities/interface/transform.h"
2425
2526#include " DataFormats/Provenance/interface/EventID.h"
2627#include " DataFormats/CaloRecHit/interface/CaloCluster.h"
3334#include " DataFormats/Math/interface/Point3D.h"
3435#include " DataFormats/GeometrySurface/interface/BoundDisk.h"
3536#include " DataFormats/HGCalReco/interface/Common.h"
37+ #include " DataFormats/HGCRecHit/interface/HGCRecHitCollections.h"
3638#include " SimDataFormats/CaloAnalysis/interface/CaloParticle.h"
3739#include " SimDataFormats/CaloAnalysis/interface/SimCluster.h"
40+ #include " SimDataFormats/CaloHit/interface/PCaloHit.h"
3841#include " DataFormats/EgammaReco/interface/SuperClusterFwd.h"
3942#include " DataFormats/EgammaReco/interface/SuperCluster.h"
4043
@@ -631,6 +634,11 @@ class TICLDumper : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::Sh
631634
632635 const edm::EDGetTokenT<std::vector<SimCluster>> simclusters_token_;
633636 const edm::EDGetTokenT<std::vector<CaloParticle>> caloparticles_token_;
637+ const std::vector<edm::InputTag> label_rechits;
638+ const std::vector<edm::EDGetTokenT<HGCRecHitCollection>> rechits_tokens_;
639+ const std::vector<edm::InputTag> label_simhits;
640+ const std::vector<edm::EDGetTokenT<std::vector<PCaloHit>>> simhits_tokens_;
641+ const edm::EDGetTokenT<std::unordered_map<DetId, const unsigned int >> hitMapToken_;
634642
635643 const edm::ESGetToken<CaloGeometry, CaloGeometryRecord> geometry_token_;
636644 const std::string detector_;
@@ -646,6 +654,7 @@ class TICLDumper : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::Sh
646654 bool saveTICLCandidate_;
647655 bool saveSimTICLCandidate_;
648656 bool saveTracks_;
657+ bool saveHits_;
649658
650659 // Output tree
651660 TTree* tree_;
@@ -723,6 +732,7 @@ class TICLDumper : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::Sh
723732 std::vector<float > cluster_time;
724733 std::vector<float > cluster_timeErr;
725734 std::vector<uint32_t > cluster_number_of_hits;
735+ std::vector<std::vector<uint32_t >> rechits_inLC;
726736
727737 // Tracks
728738 std::vector<unsigned int > track_id;
@@ -752,11 +762,34 @@ class TICLDumper : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::Sh
752762 std::vector<int > track_isMuon;
753763 std::vector<int > track_isTrackerMuon;
754764
765+ // rechits
766+ std::vector<uint32_t > rechit_ID;
767+ std::vector<float > rechit_energy;
768+ std::vector<float > rechit_x;
769+ std::vector<float > rechit_y;
770+ std::vector<float > rechit_z;
771+ std::vector<float > rechit_time;
772+ std::vector<float > rechit_radius;
773+ std::vector<float > rechit_simEnergy;
774+ std::vector<float > rechit_simEnergyEM;
775+ std::vector<float > rechit_simEnergyHad;
776+
777+ std::vector<uint32_t > simhit_ID;
778+ std::vector<float > simhit_energy;
779+ std::vector<float > simhit_energyEM;
780+ std::vector<float > simhit_energyHad;
781+ std::vector<float > simhit_x;
782+ std::vector<float > simhit_y;
783+ std::vector<float > simhit_z;
784+ std::vector<float > simhit_time;
785+
755786 TTree* cluster_tree_;
756787 TTree* candidate_tree_;
757788 TTree* superclustering_tree_;
758789 TTree* tracks_tree_;
759790 TTree* simTICLCandidate_tree;
791+ TTree* rechits_tree_;
792+ TTree* simhits_tree_;
760793};
761794
762795void TICLDumper::clearVariables () {
@@ -833,6 +866,7 @@ void TICLDumper::clearVariables() {
833866 cluster_time.clear ();
834867 cluster_timeErr.clear ();
835868 cluster_number_of_hits.clear ();
869+ rechits_inLC.clear ();
836870
837871 track_id.clear ();
838872 track_hgcal_x.clear ();
@@ -860,6 +894,26 @@ void TICLDumper::clearVariables() {
860894 track_nhits.clear ();
861895 track_isMuon.clear ();
862896 track_isTrackerMuon.clear ();
897+
898+ rechit_ID.clear ();
899+ rechit_energy.clear ();
900+ rechit_x.clear ();
901+ rechit_y.clear ();
902+ rechit_z.clear ();
903+ rechit_time.clear ();
904+ rechit_radius.clear ();
905+ rechit_simEnergy.clear ();
906+ rechit_simEnergyEM.clear ();
907+ rechit_simEnergyHad.clear ();
908+
909+ simhit_ID.clear ();
910+ simhit_energy.clear ();
911+ simhit_energyEM.clear ();
912+ simhit_energyHad.clear ();
913+ simhit_x.clear ();
914+ simhit_y.clear ();
915+ simhit_z.clear ();
916+ simhit_time.clear ();
863917};
864918
865919TICLDumper::TICLDumper (const edm::ParameterSet& ps)
@@ -900,6 +954,14 @@ TICLDumper::TICLDumper(const edm::ParameterSet& ps)
900954 associations_dumperHelpers_(associations_parameterSets_.size()),
901955 simclusters_token_(consumes(ps.getParameter<edm::InputTag>(" simclusters" ))),
902956 caloparticles_token_(consumes(ps.getParameter<edm::InputTag>(" caloparticles" ))),
957+ label_rechits(ps.getParameter<std::vector<edm::InputTag>>(" label_rechits" )),
958+ rechits_tokens_{edm::vector_transform (
959+ label_rechits, [this ](const edm::InputTag& lab) { return consumes<HGCRecHitCollection>(lab); })},
960+ label_simhits (ps.getParameter<std::vector<edm::InputTag>>(" label_simhits" )),
961+ simhits_tokens_{edm::vector_transform (
962+ label_simhits, [this ](const edm::InputTag& lab) { return consumes<std::vector<PCaloHit>>(lab); })},
963+ hitMapToken_ (
964+ consumes<std::unordered_map<DetId, const unsigned int >>(ps.getParameter<edm::InputTag>(" hitMapTag" ))),
903965 geometry_token_(esConsumes<CaloGeometry, CaloGeometryRecord, edm::Transition::BeginRun>()),
904966 detector_(ps.getParameter<std::string>(" detector" )),
905967 propName_(ps.getParameter<std::string>(" propagator" )),
@@ -912,7 +974,8 @@ TICLDumper::TICLDumper(const edm::ParameterSet& ps)
912974 saveRecoSuperclusters_(ps.getParameter<bool >(" saveRecoSuperclusters" )),
913975 saveTICLCandidate_(ps.getParameter<bool >(" saveSimTICLCandidate" )),
914976 saveSimTICLCandidate_(ps.getParameter<bool >(" saveSimTICLCandidate" )),
915- saveTracks_(ps.getParameter<bool >(" saveTracks" )) {
977+ saveTracks_(ps.getParameter<bool >(" saveTracks" )),
978+ saveHits_(ps.getParameter<bool >(" saveHits" )) {
916979 if (saveSuperclustering_) {
917980 superclustering_linkedResultTracksters_token =
918981 consumes<std::vector<std::vector<unsigned int >>>(ps.getParameter <edm::InputTag>(" superclustering" ));
@@ -966,6 +1029,29 @@ void TICLDumper::beginJob() {
9661029 tracksters_trees.push_back (tree);
9671030 tracksters_dumperHelpers_[i].initTree (tree, &eventId_);
9681031 }
1032+ if (saveHits_) {
1033+ rechits_tree_ = fs->make <TTree>(" rechits" , " HGCAL rechits" );
1034+ rechits_tree_->Branch (" ID" , &rechit_ID);
1035+ rechits_tree_->Branch (" energy" , &rechit_energy);
1036+ rechits_tree_->Branch (" position_x" , &rechit_x);
1037+ rechits_tree_->Branch (" position_y" , &rechit_y);
1038+ rechits_tree_->Branch (" position_z" , &rechit_z);
1039+ rechits_tree_->Branch (" time" , &rechit_time);
1040+ rechits_tree_->Branch (" radiusToSide" , &rechit_radius);
1041+ rechits_tree_->Branch (" simEnergy" , &rechit_simEnergy);
1042+ rechits_tree_->Branch (" simEnergyEM" , &rechit_simEnergyEM);
1043+ rechits_tree_->Branch (" simEnergyHad" , &rechit_simEnergyHad);
1044+
1045+ simhits_tree_ = fs->make <TTree>(" simhits" , " HGCAL simhits" );
1046+ simhits_tree_->Branch (" ID" , &simhit_ID);
1047+ simhits_tree_->Branch (" energy" , &simhit_energy);
1048+ simhits_tree_->Branch (" energyEM" , &simhit_energyEM);
1049+ simhits_tree_->Branch (" energyHad" , &simhit_energyHad);
1050+ simhits_tree_->Branch (" position_x" , &simhit_x);
1051+ simhits_tree_->Branch (" position_y" , &simhit_y);
1052+ simhits_tree_->Branch (" position_z" , &simhit_z);
1053+ simhits_tree_->Branch (" time" , &simhit_time);
1054+ }
9691055 if (saveLCs_) {
9701056 cluster_tree_ = fs->make <TTree>(" clusters" , " TICL tracksters" );
9711057 cluster_tree_->Branch (" event" , &eventId_);
@@ -983,6 +1069,7 @@ void TICLDumper::beginJob() {
9831069 cluster_tree_->Branch (" cluster_time" , &cluster_time);
9841070 cluster_tree_->Branch (" cluster_timeErr" , &cluster_timeErr);
9851071 cluster_tree_->Branch (" cluster_number_of_hits" , &cluster_number_of_hits);
1072+ cluster_tree_->Branch (" rechits" , &rechits_inLC);
9861073 }
9871074 if (saveTICLCandidate_) {
9881075 candidate_tree_ = fs->make <TTree>(" candidates" , " TICL candidates" );
@@ -1202,6 +1289,73 @@ void TICLDumper::analyze(const edm::Event& event, const edm::EventSetup& setup)
12021289
12031290 nclusters_ = clusters.size ();
12041291
1292+ if (saveHits_) {
1293+ edm::Handle<std::unordered_map<DetId, const unsigned int >> hitMap;
1294+ event.getByToken (hitMapToken_, hitMap);
1295+
1296+ struct ThreeFloat {
1297+ ThreeFloat () : energy(0 .f), energyEM(0 .f), energyHad(0 .f) {};
1298+ ThreeFloat (double e, double eEM, double eHad) : energy((float )e), energyEM((float )eEM), energyHad((float )eHad) {};
1299+ ThreeFloat (float e, float eEM, float eHad) : energy(e), energyEM(eEM), energyHad(eHad) {};
1300+
1301+ ThreeFloat& operator +=(const ThreeFloat& other) {
1302+ energy += other.energy ;
1303+ energyEM += other.energyEM ;
1304+ energyHad += other.energyHad ;
1305+ return *this ;
1306+ }
1307+
1308+ float energy;
1309+ float energyEM;
1310+ float energyHad;
1311+ };
1312+
1313+ std::vector<ThreeFloat> hitIdToEnergies (hitMap->size ());
1314+
1315+ for (auto const & sh_token : simhits_tokens_) {
1316+ edm::Handle<std::vector<PCaloHit>> simhit_handle;
1317+ event.getByToken (sh_token, simhit_handle);
1318+ const auto & shColl = *simhit_handle;
1319+ for (auto const & sh : shColl) {
1320+ simhit_energy.push_back (sh.energy ());
1321+ simhit_energyEM.push_back (sh.energyEM ());
1322+ simhit_energyHad.push_back (sh.energyHad ());
1323+ auto const shPosition = detectorTools_->rhtools .getPosition (sh.id ());
1324+ simhit_x.push_back (shPosition.x ());
1325+ simhit_y.push_back (shPosition.y ());
1326+ simhit_z.push_back (shPosition.z ());
1327+ simhit_ID.push_back (sh.id ());
1328+ simhit_time.push_back (sh.time ());
1329+ const auto hitId = hitMap->find (DetId (sh.id ()));
1330+ if (hitId != hitMap->end ()) {
1331+ hitIdToEnergies[hitId->second ] += {sh.energy (), sh.energyEM (), sh.energyHad ()};
1332+ }
1333+ }
1334+ }
1335+
1336+ for (auto const & rh_token : rechits_tokens_) {
1337+ edm::Handle<HGCRecHitCollection> rechit_handle;
1338+ event.getByToken (rh_token, rechit_handle);
1339+ const auto & rhColl = *rechit_handle;
1340+ for (auto const & rh : rhColl) {
1341+ rechit_energy.push_back (rh.energy ());
1342+ auto const rhPosition = detectorTools_->rhtools .getPosition (rh.detid ());
1343+ rechit_x.push_back (rhPosition.x ());
1344+ rechit_y.push_back (rhPosition.y ());
1345+ rechit_z.push_back (rhPosition.z ());
1346+ rechit_ID.push_back (rh.detid ());
1347+ rechit_time.push_back (rh.time ());
1348+ rechit_radius.push_back (detectorTools_->rhtools .getRadiusToSide (rh.detid ()));
1349+ const auto hitId = hitMap->find (DetId (rh.detid ()));
1350+ if (hitId != hitMap->end ()) {
1351+ rechit_simEnergy.push_back (hitIdToEnergies[hitId->second ].energy );
1352+ rechit_simEnergyEM.push_back (hitIdToEnergies[hitId->second ].energyEM );
1353+ rechit_simEnergyHad.push_back (hitIdToEnergies[hitId->second ].energyHad );
1354+ }
1355+ }
1356+ }
1357+ }
1358+
12051359 // Save all the trackster collections
12061360 for (unsigned int i = 0 ; i < tracksters_dumperHelpers_.size (); i++) {
12071361 edm::Handle<std::vector<ticl::Trackster>> tracksters_handle;
@@ -1304,6 +1458,11 @@ void TICLDumper::analyze(const edm::Event& event, const edm::EventSetup& setup)
13041458 cluster_timeErr.push_back (layerClustersTimes.get (c_id).second );
13051459 cluster_time.push_back (layerClustersTimes.get (c_id).first );
13061460 c_id += 1 ;
1461+ std::vector<uint32_t > hits_detid;
1462+ for (auto const & handf : cluster_iterator->hitsAndFractions ()) {
1463+ hits_detid.push_back (handf.first );
1464+ }
1465+ rechits_inLC.push_back (hits_detid);
13071466 }
13081467
13091468 tracksters_in_candidate.resize (ticlcandidates.size ());
@@ -1407,6 +1566,10 @@ void TICLDumper::analyze(const edm::Event& event, const edm::EventSetup& setup)
14071566 tracks_tree_->Fill ();
14081567 if (saveSimTICLCandidate_)
14091568 simTICLCandidate_tree->Fill ();
1569+ if (saveHits_) {
1570+ rechits_tree_->Fill ();
1571+ simhits_tree_->Fill ();
1572+ }
14101573}
14111574
14121575void TICLDumper::endJob () {}
@@ -1456,6 +1619,15 @@ void TICLDumper::fillDescriptions(edm::ConfigurationDescriptions& descriptions)
14561619 desc.add <edm::InputTag>(" simtrackstersSC" , edm::InputTag (" ticlSimTracksters" ))
14571620 ->setComment (" SimTrackster from CaloParticle collection to use for simTICLcandidates" );
14581621 desc.add <edm::InputTag>(" simTICLCandidates" , edm::InputTag (" ticlSimTracksters" ));
1622+ desc.add <std::vector<edm::InputTag>>(" label_rechits" ,
1623+ {edm::InputTag (" HGCalRecHit" , " HGCEERecHits" ),
1624+ edm::InputTag (" HGCalRecHit" , " HGCHEFRecHits" ),
1625+ edm::InputTag (" HGCalRecHit" , " HGCHEBRecHits" )});
1626+ desc.add <std::vector<edm::InputTag>>(" label_simhits" ,
1627+ {edm::InputTag (" g4SimHits" , " HGCHitsEE" ),
1628+ edm::InputTag (" g4SimHits" , " HGCHitsHEfront" ),
1629+ edm::InputTag (" g4SimHits" , " HGCHitsHEback" )});
1630+ desc.add <edm::InputTag>(" hitMapTag" , edm::InputTag (" recHitMapProducer" , " hgcalRecHitMap" ));
14591631
14601632 // Settings for dumping trackster associators (recoToSim & simToReco)
14611633 edm::ParameterSetDescription associatorDescValidator;
@@ -1479,6 +1651,7 @@ void TICLDumper::fillDescriptions(edm::ConfigurationDescriptions& descriptions)
14791651 desc.add <bool >(" saveSuperclustering" , true );
14801652 desc.add <bool >(" saveRecoSuperclusters" , true )
14811653 ->setComment (" Save superclustering Egamma collections (as reco::SuperCluster)" );
1654+ desc.add <bool >(" saveHits" , false );
14821655 descriptions.add (" ticlDumper" , desc);
14831656}
14841657
0 commit comments