Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-29 11:58:07

0001 // Original Authors:  Philipp Zehetner, Wahid Redjeb
0002 
0003 #include "TTree.h"
0004 #include "TFile.h"
0005 
0006 #include <iostream>
0007 #include <fstream>
0008 #include <sstream>
0009 #include <variant>
0010 
0011 #include <memory>  // unique_ptr
0012 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0013 #include "FWCore/ParameterSet/interface/PluginDescription.h"
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015 #include "FWCore/Utilities/interface/ESGetToken.h"
0016 #include "FWCore/Framework/interface/ESHandle.h"
0017 #include "FWCore/Framework/interface/Frameworkfwd.h"
0018 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0019 #include "FWCore/Framework/interface/Event.h"
0020 #include "FWCore/Framework/interface/MakerMacros.h"
0021 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0022 #include "FWCore/ServiceRegistry/interface/Service.h"
0023 
0024 #include "DataFormats/CaloRecHit/interface/CaloCluster.h"
0025 #include "DataFormats/HGCalReco/interface/Trackster.h"
0026 #include "DataFormats/HGCalReco/interface/TICLCandidate.h"
0027 #include "DataFormats/TrackReco/interface/Track.h"
0028 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0029 #include "DataFormats/DetId/interface/DetId.h"
0030 #include "DataFormats/Math/interface/Point3D.h"
0031 #include "DataFormats/GeometrySurface/interface/BoundDisk.h"
0032 #include "DataFormats/HGCalReco/interface/Common.h"
0033 #include "SimDataFormats/CaloAnalysis/interface/CaloParticle.h"
0034 #include "SimDataFormats/CaloAnalysis/interface/SimCluster.h"
0035 
0036 #include "RecoLocalCalo/HGCalRecAlgos/interface/RecHitTools.h"
0037 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
0038 #include "TrackingTools/GeomPropagators/interface/Propagator.h"
0039 #include "TrackingTools/Records/interface/TrackingComponentsRecord.h"
0040 
0041 #include "MagneticField/Engine/interface/MagneticField.h"
0042 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0043 
0044 #include "Geometry/HGCalCommonData/interface/HGCalDDDConstants.h"
0045 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0046 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0047 
0048 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0049 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0050 #include "RecoLocalCalo/HGCalRecAlgos/interface/RecHitTools.h"
0051 
0052 #include "SimDataFormats/Associations/interface/TracksterToSimTracksterHitLCAssociator.h"
0053 #include "RecoHGCal/TICL/interface/commons.h"
0054 
0055 // TFileService
0056 #include "FWCore/ServiceRegistry/interface/Service.h"
0057 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0058 
0059 class TICLDumper : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::SharedResources> {
0060 public:
0061   explicit TICLDumper(const edm::ParameterSet&);
0062   ~TICLDumper() override;
0063   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0064   typedef math::XYZVector Vector;
0065   typedef std::vector<double> Vec;
0066 
0067 private:
0068   void beginJob() override;
0069   void beginRun(const edm::Run&, const edm::EventSetup&) override;
0070 
0071   void initialize(const HGCalDDDConstants* hgcons,
0072                   const hgcal::RecHitTools rhtools,
0073                   const edm::ESHandle<MagneticField> bfieldH,
0074                   const edm::ESHandle<Propagator> propH);
0075   void buildLayers();
0076 
0077   void analyze(const edm::Event&, const edm::EventSetup&) override;
0078   void endRun(edm::Run const& iEvent, edm::EventSetup const&) override{};
0079   void endJob() override;
0080 
0081   // Define Tokens
0082   const edm::EDGetTokenT<std::vector<ticl::Trackster>> tracksters_token_;
0083   const edm::EDGetTokenT<std::vector<reco::CaloCluster>> layer_clusters_token_;
0084   const edm::EDGetTokenT<std::vector<TICLCandidate>> ticl_candidates_token_;
0085   const edm::EDGetTokenT<std::vector<reco::Track>> tracks_token_;
0086   const edm::EDGetTokenT<std::vector<bool>> tracks_mask_token_;
0087   const edm::EDGetTokenT<edm::ValueMap<float>> tracks_time_token_;
0088   const edm::EDGetTokenT<edm::ValueMap<float>> tracks_time_quality_token_;
0089   const edm::EDGetTokenT<edm::ValueMap<float>> tracks_time_err_token_;
0090   const edm::EDGetTokenT<std::vector<double>> hgcaltracks_x_token_;
0091   const edm::EDGetTokenT<std::vector<double>> hgcaltracks_y_token_;
0092   const edm::EDGetTokenT<std::vector<double>> hgcaltracks_z_token_;
0093   const edm::EDGetTokenT<std::vector<double>> hgcaltracks_eta_token_;
0094   const edm::EDGetTokenT<std::vector<double>> hgcaltracks_phi_token_;
0095   const edm::EDGetTokenT<std::vector<double>> hgcaltracks_px_token_;
0096   const edm::EDGetTokenT<std::vector<double>> hgcaltracks_py_token_;
0097   const edm::EDGetTokenT<std::vector<double>> hgcaltracks_pz_token_;
0098   const edm::EDGetTokenT<std::vector<ticl::Trackster>> tracksters_merged_token_;
0099   const edm::EDGetTokenT<edm::ValueMap<std::pair<float, float>>> clustersTime_token_;
0100   const edm::EDGetTokenT<std::vector<int>> tracksterSeeds_token_;
0101   edm::ESGetToken<CaloGeometry, CaloGeometryRecord> caloGeometry_token_;
0102   const edm::EDGetTokenT<std::vector<ticl::Trackster>> simTracksters_SC_token_;
0103   const edm::EDGetTokenT<std::vector<ticl::Trackster>> simTracksters_CP_token_;
0104   const edm::EDGetTokenT<std::vector<ticl::Trackster>> simTracksters_PU_token_;
0105   const edm::EDGetTokenT<std::vector<TICLCandidate>> simTICLCandidate_token_;
0106   const edm::EDGetTokenT<hgcal::RecoToSimCollectionSimTracksters> tsRecoToSimSC_token_;
0107   const edm::EDGetTokenT<hgcal::SimToRecoCollectionSimTracksters> tsSimToRecoSC_token_;
0108   const edm::EDGetTokenT<hgcal::RecoToSimCollectionSimTracksters> tsRecoToSimCP_token_;
0109   const edm::EDGetTokenT<hgcal::SimToRecoCollectionSimTracksters> tsSimToRecoCP_token_;
0110   const edm::EDGetTokenT<hgcal::RecoToSimCollectionSimTracksters> MergeRecoToSimSC_token_;
0111   const edm::EDGetTokenT<hgcal::SimToRecoCollectionSimTracksters> MergeSimToRecoSC_token_;
0112   const edm::EDGetTokenT<hgcal::RecoToSimCollectionSimTracksters> MergeRecoToSimCP_token_;
0113   const edm::EDGetTokenT<hgcal::SimToRecoCollectionSimTracksters> MergeSimToRecoCP_token_;
0114   const edm::EDGetTokenT<hgcal::RecoToSimCollectionSimTracksters> MergeRecoToSimPU_token_;
0115   const edm::EDGetTokenT<hgcal::SimToRecoCollectionSimTracksters> MergeSimToRecoPU_token_;
0116   const edm::EDGetTokenT<std::vector<SimCluster>> simclusters_token_;
0117   const edm::EDGetTokenT<std::vector<CaloParticle>> caloparticles_token_;
0118 
0119   const edm::ESGetToken<CaloGeometry, CaloGeometryRecord> geometry_token_;
0120   const std::string detector_;
0121   const std::string propName_;
0122   const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> bfield_token_;
0123   const edm::ESGetToken<Propagator, TrackingComponentsRecord> propagator_token_;
0124   hgcal::RecHitTools rhtools_;
0125   edm::ESGetToken<HGCalDDDConstants, IdealGeometryRecord> hdc_token_;
0126   const HGCalDDDConstants* hgcons_;
0127   std::unique_ptr<GeomDet> firstDisk_[2];
0128   std::unique_ptr<GeomDet> interfaceDisk_[2];
0129   edm::ESHandle<MagneticField> bfield_;
0130   edm::ESHandle<Propagator> propagator_;
0131   bool saveLCs_;
0132   bool saveCLUE3DTracksters_;
0133   bool saveTrackstersMerged_;
0134   bool saveSimTrackstersSC_;
0135   bool saveSimTrackstersCP_;
0136   bool saveTICLCandidate_;
0137   bool saveSimTICLCandidate_;
0138   bool saveTracks_;
0139   bool saveAssociations_;
0140 
0141   // Output tree
0142   TTree* tree_;
0143 
0144   void clearVariables();
0145 
0146   // Variables for branches
0147   unsigned int ev_event_;
0148   unsigned int ntracksters_;
0149   unsigned int nclusters_;
0150   unsigned int stsSC_ntracksters_;
0151   unsigned int stsCP_ntracksters_;
0152   size_t nsimTrackstersSC;
0153   size_t nsimTrackstersCP;
0154 
0155   std::vector<float> trackster_time;
0156   std::vector<float> trackster_timeError;
0157   std::vector<float> trackster_regressed_energy;
0158   std::vector<float> trackster_raw_energy;
0159   std::vector<float> trackster_raw_em_energy;
0160   std::vector<float> trackster_raw_pt;
0161   std::vector<float> trackster_raw_em_pt;
0162   std::vector<float> trackster_barycenter_x;
0163   std::vector<float> trackster_barycenter_y;
0164   std::vector<float> trackster_barycenter_z;
0165   std::vector<float> trackster_EV1;
0166   std::vector<float> trackster_EV2;
0167   std::vector<float> trackster_EV3;
0168   std::vector<float> trackster_eVector0_x;
0169   std::vector<float> trackster_eVector0_y;
0170   std::vector<float> trackster_eVector0_z;
0171   std::vector<float> trackster_sigmaPCA1;
0172   std::vector<float> trackster_sigmaPCA2;
0173   std::vector<float> trackster_sigmaPCA3;
0174   std::vector<float> trackster_barycenter_eta;
0175   std::vector<float> trackster_barycenter_phi;
0176   std::vector<std::vector<float>> trackster_id_probabilities;
0177   std::vector<std::vector<uint32_t>> trackster_vertices_indexes;
0178   std::vector<std::vector<float>> trackster_vertices_x;
0179   std::vector<std::vector<float>> trackster_vertices_y;
0180   std::vector<std::vector<float>> trackster_vertices_z;
0181   std::vector<std::vector<float>> trackster_vertices_time;
0182   std::vector<std::vector<float>> trackster_vertices_timeErr;
0183   std::vector<std::vector<float>> trackster_vertices_energy;
0184   std::vector<std::vector<float>> trackster_vertices_correctedEnergy;
0185   std::vector<std::vector<float>> trackster_vertices_correctedEnergyUncertainty;
0186   std::vector<std::vector<float>> trackster_vertices_multiplicity;
0187 
0188   std::vector<float> stsSC_trackster_time;
0189   std::vector<float> stsSC_trackster_timeError;
0190   std::vector<float> stsSC_trackster_regressed_energy;
0191   std::vector<float> stsSC_trackster_regressed_pt;
0192   std::vector<float> stsSC_trackster_raw_energy;
0193   std::vector<float> stsSC_trackster_raw_em_energy;
0194   std::vector<float> stsSC_trackster_raw_pt;
0195   std::vector<float> stsSC_trackster_raw_em_pt;
0196   std::vector<float> stsSC_trackster_barycenter_x;
0197   std::vector<float> stsSC_trackster_barycenter_y;
0198   std::vector<float> stsSC_trackster_barycenter_z;
0199   std::vector<float> stsSC_trackster_barycenter_eta;
0200   std::vector<float> stsSC_trackster_barycenter_phi;
0201   std::vector<float> stsSC_trackster_EV1;
0202   std::vector<float> stsSC_trackster_EV2;
0203   std::vector<float> stsSC_trackster_EV3;
0204   std::vector<float> stsSC_trackster_eVector0_x;
0205   std::vector<float> stsSC_trackster_eVector0_y;
0206   std::vector<float> stsSC_trackster_eVector0_z;
0207   std::vector<float> stsSC_trackster_sigmaPCA1;
0208   std::vector<float> stsSC_trackster_sigmaPCA2;
0209   std::vector<float> stsSC_trackster_sigmaPCA3;
0210   std::vector<int> stsSC_pdgID;
0211   std::vector<int> stsSC_trackIdx;
0212   std::vector<float> stsSC_trackTime;
0213   std::vector<float> stsSC_boundaryX;
0214   std::vector<float> stsSC_boundaryY;
0215   std::vector<float> stsSC_boundaryZ;
0216   std::vector<float> stsSC_boundaryEta;
0217   std::vector<float> stsSC_boundaryPhi;
0218   std::vector<float> stsSC_boundaryPx;
0219   std::vector<float> stsSC_boundaryPy;
0220   std::vector<float> stsSC_boundaryPz;
0221   std::vector<float> stsSC_track_boundaryX;
0222   std::vector<float> stsSC_track_boundaryY;
0223   std::vector<float> stsSC_track_boundaryZ;
0224   std::vector<float> stsSC_track_boundaryEta;
0225   std::vector<float> stsSC_track_boundaryPhi;
0226   std::vector<float> stsSC_track_boundaryPx;
0227   std::vector<float> stsSC_track_boundaryPy;
0228   std::vector<float> stsSC_track_boundaryPz;
0229   std::vector<std::vector<float>> stsSC_trackster_id_probabilities;
0230   std::vector<std::vector<uint32_t>> stsSC_trackster_vertices_indexes;
0231   std::vector<std::vector<float>> stsSC_trackster_vertices_x;
0232   std::vector<std::vector<float>> stsSC_trackster_vertices_y;
0233   std::vector<std::vector<float>> stsSC_trackster_vertices_z;
0234   std::vector<std::vector<float>> stsSC_trackster_vertices_time;
0235   std::vector<std::vector<float>> stsSC_trackster_vertices_timeErr;
0236   std::vector<std::vector<float>> stsSC_trackster_vertices_energy;
0237   std::vector<std::vector<float>> stsSC_trackster_vertices_correctedEnergy;
0238   std::vector<std::vector<float>> stsSC_trackster_vertices_correctedEnergyUncertainty;
0239   std::vector<std::vector<float>> stsSC_trackster_vertices_multiplicity;
0240   std::vector<float> stsCP_trackster_time;
0241   std::vector<float> stsCP_trackster_timeError;
0242   std::vector<float> stsCP_trackster_regressed_energy;
0243   std::vector<float> stsCP_trackster_regressed_pt;
0244   std::vector<float> stsCP_trackster_raw_energy;
0245   std::vector<float> stsCP_trackster_raw_em_energy;
0246   std::vector<float> stsCP_trackster_raw_pt;
0247   std::vector<float> stsCP_trackster_raw_em_pt;
0248   std::vector<float> stsCP_trackster_barycenter_x;
0249   std::vector<float> stsCP_trackster_barycenter_y;
0250   std::vector<float> stsCP_trackster_barycenter_z;
0251   std::vector<float> stsCP_trackster_barycenter_eta;
0252   std::vector<float> stsCP_trackster_barycenter_phi;
0253   std::vector<float> stsCP_trackster_EV1;
0254   std::vector<float> stsCP_trackster_EV2;
0255   std::vector<float> stsCP_trackster_EV3;
0256   std::vector<float> stsCP_trackster_eVector0_x;
0257   std::vector<float> stsCP_trackster_eVector0_y;
0258   std::vector<float> stsCP_trackster_eVector0_z;
0259   std::vector<float> stsCP_trackster_sigmaPCA1;
0260   std::vector<float> stsCP_trackster_sigmaPCA2;
0261   std::vector<float> stsCP_trackster_sigmaPCA3;
0262   std::vector<int> stsCP_pdgID;
0263   std::vector<int> stsCP_trackIdx;
0264   std::vector<float> stsCP_trackTime;
0265   std::vector<float> stsCP_boundaryX;
0266   std::vector<float> stsCP_boundaryY;
0267   std::vector<float> stsCP_boundaryZ;
0268   std::vector<float> stsCP_boundaryEta;
0269   std::vector<float> stsCP_boundaryPhi;
0270   std::vector<float> stsCP_boundaryPx;
0271   std::vector<float> stsCP_boundaryPy;
0272   std::vector<float> stsCP_boundaryPz;
0273   std::vector<float> stsCP_track_boundaryX;
0274   std::vector<float> stsCP_track_boundaryY;
0275   std::vector<float> stsCP_track_boundaryZ;
0276   std::vector<float> stsCP_track_boundaryEta;
0277   std::vector<float> stsCP_track_boundaryPhi;
0278   std::vector<float> stsCP_track_boundaryPx;
0279   std::vector<float> stsCP_track_boundaryPy;
0280   std::vector<float> stsCP_track_boundaryPz;
0281   std::vector<std::vector<float>> stsCP_trackster_id_probabilities;
0282   std::vector<std::vector<uint32_t>> stsCP_trackster_vertices_indexes;
0283   std::vector<std::vector<float>> stsCP_trackster_vertices_x;
0284   std::vector<std::vector<float>> stsCP_trackster_vertices_y;
0285   std::vector<std::vector<float>> stsCP_trackster_vertices_z;
0286   std::vector<std::vector<float>> stsCP_trackster_vertices_time;
0287   std::vector<std::vector<float>> stsCP_trackster_vertices_timeErr;
0288   std::vector<std::vector<float>> stsCP_trackster_vertices_energy;
0289   std::vector<std::vector<float>> stsCP_trackster_vertices_correctedEnergy;
0290   std::vector<std::vector<float>> stsCP_trackster_vertices_correctedEnergyUncertainty;
0291   std::vector<std::vector<float>> stsCP_trackster_vertices_multiplicity;
0292 
0293   std::vector<float> simTICLCandidate_raw_energy;
0294   std::vector<float> simTICLCandidate_regressed_energy;
0295   std::vector<std::vector<int>> simTICLCandidate_simTracksterCPIndex;
0296   std::vector<float> simTICLCandidate_boundaryX;
0297   std::vector<float> simTICLCandidate_boundaryY;
0298   std::vector<float> simTICLCandidate_boundaryZ;
0299   std::vector<float> simTICLCandidate_boundaryPx;
0300   std::vector<float> simTICLCandidate_boundaryPy;
0301   std::vector<float> simTICLCandidate_boundaryPz;
0302   std::vector<float> simTICLCandidate_trackTime;
0303   std::vector<float> simTICLCandidate_trackBeta;
0304   std::vector<float> simTICLCandidate_caloParticleMass;
0305   std::vector<int> simTICLCandidate_pdgId;
0306   std::vector<int> simTICLCandidate_charge;
0307   std::vector<int> simTICLCandidate_track_in_candidate;
0308 
0309   // from TICLCandidate, product of linking
0310   size_t nCandidates;
0311   std::vector<int> candidate_charge;
0312   std::vector<int> candidate_pdgId;
0313   std::vector<float> candidate_energy;
0314   std::vector<double> candidate_px;
0315   std::vector<double> candidate_py;
0316   std::vector<double> candidate_pz;
0317   std::vector<float> candidate_time;
0318   std::vector<float> candidate_time_err;
0319   std::vector<std::vector<float>> candidate_id_probabilities;
0320   std::vector<std::vector<uint32_t>> tracksters_in_candidate;
0321   std::vector<int> track_in_candidate;
0322 
0323   // merged tracksters
0324   size_t nTrackstersMerged;
0325   std::vector<float> tracksters_merged_time;
0326   std::vector<float> tracksters_merged_timeError;
0327   std::vector<float> tracksters_merged_regressed_energy;
0328   std::vector<float> tracksters_merged_raw_energy;
0329   std::vector<float> tracksters_merged_raw_em_energy;
0330   std::vector<float> tracksters_merged_raw_pt;
0331   std::vector<float> tracksters_merged_raw_em_pt;
0332   std::vector<float> tracksters_merged_barycenter_x;
0333   std::vector<float> tracksters_merged_barycenter_y;
0334   std::vector<float> tracksters_merged_barycenter_z;
0335   std::vector<float> tracksters_merged_barycenter_eta;
0336   std::vector<float> tracksters_merged_barycenter_phi;
0337   std::vector<float> tracksters_merged_EV1;
0338   std::vector<float> tracksters_merged_EV2;
0339   std::vector<float> tracksters_merged_EV3;
0340   std::vector<float> tracksters_merged_eVector0_x;
0341   std::vector<float> tracksters_merged_eVector0_y;
0342   std::vector<float> tracksters_merged_eVector0_z;
0343   std::vector<float> tracksters_merged_sigmaPCA1;
0344   std::vector<float> tracksters_merged_sigmaPCA2;
0345   std::vector<float> tracksters_merged_sigmaPCA3;
0346   std::vector<std::vector<uint32_t>> tracksters_merged_vertices_indexes;
0347   std::vector<std::vector<float>> tracksters_merged_vertices_x;
0348   std::vector<std::vector<float>> tracksters_merged_vertices_y;
0349   std::vector<std::vector<float>> tracksters_merged_vertices_z;
0350   std::vector<std::vector<float>> tracksters_merged_vertices_time;
0351   std::vector<std::vector<float>> tracksters_merged_vertices_timeErr;
0352   std::vector<std::vector<float>> tracksters_merged_vertices_energy;
0353   std::vector<std::vector<float>> tracksters_merged_vertices_correctedEnergy;
0354   std::vector<std::vector<float>> tracksters_merged_vertices_correctedEnergyUncertainty;
0355   std::vector<std::vector<float>> tracksters_merged_vertices_multiplicity;
0356   std::vector<std::vector<float>> tracksters_merged_id_probabilities;
0357 
0358   // associations
0359   std::vector<std::vector<uint32_t>> trackstersCLUE3D_recoToSim_SC;
0360   std::vector<std::vector<float>> trackstersCLUE3D_recoToSim_SC_score;
0361   std::vector<std::vector<float>> trackstersCLUE3D_recoToSim_SC_sharedE;
0362   std::vector<std::vector<uint32_t>> trackstersCLUE3D_simToReco_SC;
0363   std::vector<std::vector<float>> trackstersCLUE3D_simToReco_SC_score;
0364   std::vector<std::vector<float>> trackstersCLUE3D_simToReco_SC_sharedE;
0365 
0366   std::vector<std::vector<uint32_t>> trackstersCLUE3D_recoToSim_CP;
0367   std::vector<std::vector<float>> trackstersCLUE3D_recoToSim_CP_score;
0368   std::vector<std::vector<float>> trackstersCLUE3D_recoToSim_CP_sharedE;
0369   std::vector<std::vector<uint32_t>> trackstersCLUE3D_simToReco_CP;
0370   std::vector<std::vector<float>> trackstersCLUE3D_simToReco_CP_score;
0371   std::vector<std::vector<float>> trackstersCLUE3D_simToReco_CP_sharedE;
0372 
0373   std::vector<std::vector<uint32_t>> MergeTracksters_recoToSim_SC;
0374   std::vector<std::vector<float>> MergeTracksters_recoToSim_SC_score;
0375   std::vector<std::vector<float>> MergeTracksters_recoToSim_SC_sharedE;
0376   std::vector<std::vector<uint32_t>> MergeTracksters_simToReco_SC;
0377   std::vector<std::vector<float>> MergeTracksters_simToReco_SC_score;
0378   std::vector<std::vector<float>> MergeTracksters_simToReco_SC_sharedE;
0379 
0380   std::vector<std::vector<uint32_t>> MergeTracksters_recoToSim_CP;
0381   std::vector<std::vector<float>> MergeTracksters_recoToSim_CP_score;
0382   std::vector<std::vector<float>> MergeTracksters_recoToSim_CP_sharedE;
0383   std::vector<std::vector<uint32_t>> MergeTracksters_simToReco_CP;
0384   std::vector<std::vector<float>> MergeTracksters_simToReco_CP_score;
0385   std::vector<std::vector<float>> MergeTracksters_simToReco_CP_sharedE;
0386 
0387   std::vector<std::vector<uint32_t>> MergeTracksters_recoToSim_PU;
0388   std::vector<std::vector<float>> MergeTracksters_recoToSim_PU_score;
0389   std::vector<std::vector<float>> MergeTracksters_recoToSim_PU_sharedE;
0390   std::vector<std::vector<uint32_t>> MergeTracksters_simToReco_PU;
0391   std::vector<std::vector<float>> MergeTracksters_simToReco_PU_score;
0392   std::vector<std::vector<float>> MergeTracksters_simToReco_PU_sharedE;
0393 
0394   std::vector<uint32_t> cluster_seedID;
0395   std::vector<float> cluster_energy;
0396   std::vector<float> cluster_correctedEnergy;
0397   std::vector<float> cluster_correctedEnergyUncertainty;
0398   std::vector<float> cluster_position_x;
0399   std::vector<float> cluster_position_y;
0400   std::vector<float> cluster_position_z;
0401   std::vector<float> cluster_position_eta;
0402   std::vector<float> cluster_position_phi;
0403   std::vector<unsigned int> cluster_layer_id;
0404   std::vector<int> cluster_type;
0405   std::vector<float> cluster_time;
0406   std::vector<float> cluster_timeErr;
0407   std::vector<uint32_t> cluster_number_of_hits;
0408 
0409   std::vector<unsigned int> track_id;
0410   std::vector<float> track_hgcal_x;
0411   std::vector<float> track_hgcal_y;
0412   std::vector<float> track_hgcal_z;
0413   std::vector<float> track_hgcal_px;
0414   std::vector<float> track_hgcal_py;
0415   std::vector<float> track_hgcal_pz;
0416   std::vector<float> track_hgcal_eta;
0417   std::vector<float> track_hgcal_phi;
0418   std::vector<float> track_hgcal_pt;
0419   std::vector<float> track_pt;
0420   std::vector<int> track_quality;
0421   std::vector<int> track_missing_outer_hits;
0422   std::vector<int> track_charge;
0423   std::vector<double> track_time;
0424   std::vector<float> track_time_quality;
0425   std::vector<float> track_time_err;
0426   std::vector<int> track_nhits;
0427 
0428   TTree* trackster_tree_;
0429   TTree* cluster_tree_;
0430   TTree* candidate_tree_;
0431   TTree* tracksters_merged_tree_;
0432   TTree* associations_tree_;
0433   TTree* simtrackstersSC_tree_;
0434   TTree* simtrackstersCP_tree_;
0435   TTree* tracks_tree_;
0436   TTree* simTICLCandidate_tree;
0437 };
0438 
0439 void TICLDumper::clearVariables() {
0440   // event info
0441   ntracksters_ = 0;
0442   nclusters_ = 0;
0443 
0444   trackster_time.clear();
0445   trackster_timeError.clear();
0446   trackster_regressed_energy.clear();
0447   trackster_raw_energy.clear();
0448   trackster_raw_em_energy.clear();
0449   trackster_raw_pt.clear();
0450   trackster_raw_em_pt.clear();
0451   trackster_barycenter_x.clear();
0452   trackster_barycenter_y.clear();
0453   trackster_barycenter_z.clear();
0454   trackster_EV1.clear();
0455   trackster_EV2.clear();
0456   trackster_EV3.clear();
0457   trackster_eVector0_x.clear();
0458   trackster_eVector0_y.clear();
0459   trackster_eVector0_z.clear();
0460   trackster_sigmaPCA1.clear();
0461   trackster_sigmaPCA2.clear();
0462   trackster_sigmaPCA3.clear();
0463   trackster_barycenter_eta.clear();
0464   trackster_barycenter_phi.clear();
0465   trackster_id_probabilities.clear();
0466   trackster_vertices_indexes.clear();
0467   trackster_vertices_x.clear();
0468   trackster_vertices_y.clear();
0469   trackster_vertices_z.clear();
0470   trackster_vertices_time.clear();
0471   trackster_vertices_timeErr.clear();
0472   trackster_vertices_energy.clear();
0473   trackster_vertices_correctedEnergy.clear();
0474   trackster_vertices_correctedEnergyUncertainty.clear();
0475   trackster_vertices_multiplicity.clear();
0476 
0477   stsSC_trackster_time.clear();
0478   stsSC_trackster_timeError.clear();
0479   stsSC_trackster_regressed_energy.clear();
0480   stsSC_trackster_regressed_pt.clear();
0481   stsSC_trackster_raw_energy.clear();
0482   stsSC_trackster_raw_em_energy.clear();
0483   stsSC_trackster_raw_pt.clear();
0484   stsSC_trackster_raw_em_pt.clear();
0485   stsSC_trackster_barycenter_x.clear();
0486   stsSC_trackster_barycenter_y.clear();
0487   stsSC_trackster_barycenter_z.clear();
0488   stsSC_trackster_EV1.clear();
0489   stsSC_trackster_EV2.clear();
0490   stsSC_trackster_EV3.clear();
0491   stsSC_trackster_eVector0_x.clear();
0492   stsSC_trackster_eVector0_y.clear();
0493   stsSC_trackster_eVector0_z.clear();
0494   stsSC_trackster_sigmaPCA1.clear();
0495   stsSC_trackster_sigmaPCA2.clear();
0496   stsSC_trackster_sigmaPCA3.clear();
0497   stsSC_trackster_barycenter_eta.clear();
0498   stsSC_trackster_barycenter_phi.clear();
0499   stsSC_pdgID.clear();
0500   stsSC_trackIdx.clear();
0501   stsSC_trackTime.clear();
0502   stsSC_boundaryX.clear();
0503   stsSC_boundaryY.clear();
0504   stsSC_boundaryZ.clear();
0505   stsSC_boundaryEta.clear();
0506   stsSC_boundaryPhi.clear();
0507   stsSC_boundaryPx.clear();
0508   stsSC_boundaryPy.clear();
0509   stsSC_boundaryPz.clear();
0510   stsSC_track_boundaryX.clear();
0511   stsSC_track_boundaryY.clear();
0512   stsSC_track_boundaryZ.clear();
0513   stsSC_track_boundaryEta.clear();
0514   stsSC_track_boundaryPhi.clear();
0515   stsSC_track_boundaryPx.clear();
0516   stsSC_track_boundaryPy.clear();
0517   stsSC_track_boundaryPz.clear();
0518   stsSC_trackster_id_probabilities.clear();
0519   stsSC_trackster_vertices_indexes.clear();
0520   stsSC_trackster_vertices_x.clear();
0521   stsSC_trackster_vertices_y.clear();
0522   stsSC_trackster_vertices_z.clear();
0523   stsSC_trackster_vertices_time.clear();
0524   stsSC_trackster_vertices_timeErr.clear();
0525   stsSC_trackster_vertices_energy.clear();
0526   stsSC_trackster_vertices_correctedEnergy.clear();
0527   stsSC_trackster_vertices_correctedEnergyUncertainty.clear();
0528   stsSC_trackster_vertices_multiplicity.clear();
0529 
0530   stsCP_trackster_time.clear();
0531   stsCP_trackster_timeError.clear();
0532   stsCP_trackster_regressed_energy.clear();
0533   stsCP_trackster_regressed_pt.clear();
0534   stsCP_trackster_raw_energy.clear();
0535   stsCP_trackster_raw_em_energy.clear();
0536   stsCP_trackster_raw_pt.clear();
0537   stsCP_trackster_raw_em_pt.clear();
0538   stsCP_trackster_barycenter_x.clear();
0539   stsCP_trackster_barycenter_y.clear();
0540   stsCP_trackster_barycenter_z.clear();
0541   stsCP_trackster_sigmaPCA1.clear();
0542   stsCP_trackster_sigmaPCA2.clear();
0543   stsCP_trackster_sigmaPCA3.clear();
0544   stsCP_trackster_barycenter_eta.clear();
0545   stsCP_trackster_barycenter_phi.clear();
0546   stsCP_pdgID.clear();
0547   stsCP_trackIdx.clear();
0548   stsCP_trackTime.clear();
0549   stsCP_boundaryX.clear();
0550   stsCP_boundaryY.clear();
0551   stsCP_boundaryZ.clear();
0552   stsCP_boundaryEta.clear();
0553   stsCP_boundaryPhi.clear();
0554   stsCP_boundaryPx.clear();
0555   stsCP_boundaryPy.clear();
0556   stsCP_boundaryPz.clear();
0557   stsCP_track_boundaryX.clear();
0558   stsCP_track_boundaryY.clear();
0559   stsCP_track_boundaryZ.clear();
0560   stsCP_track_boundaryEta.clear();
0561   stsCP_track_boundaryPhi.clear();
0562   stsCP_track_boundaryPx.clear();
0563   stsCP_track_boundaryPy.clear();
0564   stsCP_track_boundaryPz.clear();
0565   stsCP_trackster_id_probabilities.clear();
0566   stsCP_trackster_vertices_indexes.clear();
0567   stsCP_trackster_vertices_x.clear();
0568   stsCP_trackster_vertices_y.clear();
0569   stsCP_trackster_vertices_z.clear();
0570   stsCP_trackster_vertices_time.clear();
0571   stsCP_trackster_vertices_timeErr.clear();
0572   stsCP_trackster_vertices_energy.clear();
0573   stsCP_trackster_vertices_correctedEnergy.clear();
0574   stsCP_trackster_vertices_correctedEnergyUncertainty.clear();
0575   stsCP_trackster_vertices_multiplicity.clear();
0576 
0577   simTICLCandidate_raw_energy.clear();
0578   simTICLCandidate_regressed_energy.clear();
0579   simTICLCandidate_simTracksterCPIndex.clear();
0580   simTICLCandidate_boundaryX.clear();
0581   simTICLCandidate_boundaryY.clear();
0582   simTICLCandidate_boundaryZ.clear();
0583   simTICLCandidate_boundaryPx.clear();
0584   simTICLCandidate_boundaryPy.clear();
0585   simTICLCandidate_boundaryPz.clear();
0586   simTICLCandidate_trackTime.clear();
0587   simTICLCandidate_trackBeta.clear();
0588   simTICLCandidate_caloParticleMass.clear();
0589   simTICLCandidate_pdgId.clear();
0590   simTICLCandidate_charge.clear();
0591   simTICLCandidate_track_in_candidate.clear();
0592 
0593   nCandidates = 0;
0594   candidate_charge.clear();
0595   candidate_pdgId.clear();
0596   candidate_energy.clear();
0597   candidate_px.clear();
0598   candidate_py.clear();
0599   candidate_pz.clear();
0600   candidate_time.clear();
0601   candidate_time_err.clear();
0602   candidate_id_probabilities.clear();
0603   tracksters_in_candidate.clear();
0604   track_in_candidate.clear();
0605 
0606   nTrackstersMerged = 0;
0607   tracksters_merged_time.clear();
0608   tracksters_merged_timeError.clear();
0609   tracksters_merged_regressed_energy.clear();
0610   tracksters_merged_raw_energy.clear();
0611   tracksters_merged_raw_em_energy.clear();
0612   tracksters_merged_raw_pt.clear();
0613   tracksters_merged_raw_em_pt.clear();
0614   tracksters_merged_barycenter_x.clear();
0615   tracksters_merged_barycenter_y.clear();
0616   tracksters_merged_barycenter_z.clear();
0617   tracksters_merged_barycenter_eta.clear();
0618   tracksters_merged_barycenter_phi.clear();
0619   tracksters_merged_EV1.clear();
0620   tracksters_merged_EV2.clear();
0621   tracksters_merged_EV3.clear();
0622   tracksters_merged_eVector0_x.clear();
0623   tracksters_merged_eVector0_y.clear();
0624   tracksters_merged_eVector0_z.clear();
0625   tracksters_merged_sigmaPCA1.clear();
0626   tracksters_merged_sigmaPCA2.clear();
0627   tracksters_merged_sigmaPCA3.clear();
0628   tracksters_merged_id_probabilities.clear();
0629   tracksters_merged_time.clear();
0630   tracksters_merged_timeError.clear();
0631   tracksters_merged_regressed_energy.clear();
0632   tracksters_merged_raw_energy.clear();
0633   tracksters_merged_raw_em_energy.clear();
0634   tracksters_merged_raw_pt.clear();
0635   tracksters_merged_raw_em_pt.clear();
0636 
0637   tracksters_merged_vertices_indexes.clear();
0638   tracksters_merged_vertices_x.clear();
0639   tracksters_merged_vertices_y.clear();
0640   tracksters_merged_vertices_z.clear();
0641   tracksters_merged_vertices_time.clear();
0642   tracksters_merged_vertices_timeErr.clear();
0643   tracksters_merged_vertices_energy.clear();
0644   tracksters_merged_vertices_correctedEnergy.clear();
0645   tracksters_merged_vertices_correctedEnergyUncertainty.clear();
0646   tracksters_merged_vertices_multiplicity.clear();
0647 
0648   trackstersCLUE3D_recoToSim_SC.clear();
0649   trackstersCLUE3D_recoToSim_SC_score.clear();
0650   trackstersCLUE3D_recoToSim_SC_sharedE.clear();
0651   trackstersCLUE3D_simToReco_SC.clear();
0652   trackstersCLUE3D_simToReco_SC_score.clear();
0653   trackstersCLUE3D_simToReco_SC_sharedE.clear();
0654 
0655   trackstersCLUE3D_recoToSim_CP.clear();
0656   trackstersCLUE3D_recoToSim_CP_score.clear();
0657   trackstersCLUE3D_recoToSim_CP_sharedE.clear();
0658   trackstersCLUE3D_simToReco_CP.clear();
0659   trackstersCLUE3D_simToReco_CP_score.clear();
0660   trackstersCLUE3D_simToReco_CP_sharedE.clear();
0661 
0662   MergeTracksters_recoToSim_SC.clear();
0663   MergeTracksters_recoToSim_SC_score.clear();
0664   MergeTracksters_recoToSim_SC_sharedE.clear();
0665   MergeTracksters_simToReco_SC.clear();
0666   MergeTracksters_simToReco_SC_score.clear();
0667   MergeTracksters_simToReco_SC_sharedE.clear();
0668 
0669   MergeTracksters_recoToSim_CP.clear();
0670   MergeTracksters_recoToSim_CP_score.clear();
0671   MergeTracksters_recoToSim_CP_sharedE.clear();
0672   MergeTracksters_simToReco_CP.clear();
0673   MergeTracksters_simToReco_CP_score.clear();
0674   MergeTracksters_simToReco_CP_sharedE.clear();
0675 
0676   MergeTracksters_recoToSim_PU.clear();
0677   MergeTracksters_recoToSim_PU_score.clear();
0678   MergeTracksters_recoToSim_PU_sharedE.clear();
0679   MergeTracksters_simToReco_PU.clear();
0680   MergeTracksters_simToReco_PU_score.clear();
0681   MergeTracksters_simToReco_PU_sharedE.clear();
0682 
0683   nsimTrackstersSC = 0;
0684 
0685   cluster_seedID.clear();
0686   cluster_energy.clear();
0687   cluster_correctedEnergy.clear();
0688   cluster_correctedEnergyUncertainty.clear();
0689   cluster_position_x.clear();
0690   cluster_position_y.clear();
0691   cluster_position_z.clear();
0692   cluster_position_eta.clear();
0693   cluster_position_phi.clear();
0694   cluster_layer_id.clear();
0695   cluster_type.clear();
0696   cluster_time.clear();
0697   cluster_timeErr.clear();
0698   cluster_number_of_hits.clear();
0699 
0700   track_id.clear();
0701   track_hgcal_x.clear();
0702   track_hgcal_y.clear();
0703   track_hgcal_z.clear();
0704   track_hgcal_eta.clear();
0705   track_hgcal_phi.clear();
0706   track_hgcal_px.clear();
0707   track_hgcal_py.clear();
0708   track_hgcal_pz.clear();
0709   track_hgcal_pt.clear();
0710   track_quality.clear();
0711   track_pt.clear();
0712   track_missing_outer_hits.clear();
0713   track_charge.clear();
0714   track_time.clear();
0715   track_time_quality.clear();
0716   track_time_err.clear();
0717   track_nhits.clear();
0718 };
0719 
0720 TICLDumper::TICLDumper(const edm::ParameterSet& ps)
0721     : tracksters_token_(consumes<std::vector<ticl::Trackster>>(ps.getParameter<edm::InputTag>("trackstersclue3d"))),
0722       layer_clusters_token_(consumes<std::vector<reco::CaloCluster>>(ps.getParameter<edm::InputTag>("layerClusters"))),
0723       ticl_candidates_token_(consumes<std::vector<TICLCandidate>>(ps.getParameter<edm::InputTag>("ticlcandidates"))),
0724       tracks_token_(consumes<std::vector<reco::Track>>(ps.getParameter<edm::InputTag>("tracks"))),
0725       tracks_time_token_(consumes<edm::ValueMap<float>>(ps.getParameter<edm::InputTag>("tracksTime"))),
0726       tracks_time_quality_token_(consumes<edm::ValueMap<float>>(ps.getParameter<edm::InputTag>("tracksTimeQual"))),
0727       tracks_time_err_token_(consumes<edm::ValueMap<float>>(ps.getParameter<edm::InputTag>("tracksTimeErr"))),
0728       tracksters_merged_token_(
0729           consumes<std::vector<ticl::Trackster>>(ps.getParameter<edm::InputTag>("trackstersmerged"))),
0730       clustersTime_token_(
0731           consumes<edm::ValueMap<std::pair<float, float>>>(ps.getParameter<edm::InputTag>("layer_clustersTime"))),
0732       caloGeometry_token_(esConsumes<CaloGeometry, CaloGeometryRecord, edm::Transition::BeginRun>()),
0733       simTracksters_SC_token_(
0734           consumes<std::vector<ticl::Trackster>>(ps.getParameter<edm::InputTag>("simtrackstersSC"))),
0735       simTracksters_CP_token_(
0736           consumes<std::vector<ticl::Trackster>>(ps.getParameter<edm::InputTag>("simtrackstersCP"))),
0737       simTracksters_PU_token_(
0738           consumes<std::vector<ticl::Trackster>>(ps.getParameter<edm::InputTag>("simtrackstersPU"))),
0739       simTICLCandidate_token_(
0740           consumes<std::vector<TICLCandidate>>(ps.getParameter<edm::InputTag>("simTICLCandidates"))),
0741       tsRecoToSimSC_token_(
0742           consumes<hgcal::RecoToSimCollectionSimTracksters>(ps.getParameter<edm::InputTag>("recoToSimAssociatorSC"))),
0743       tsSimToRecoSC_token_(
0744           consumes<hgcal::SimToRecoCollectionSimTracksters>(ps.getParameter<edm::InputTag>("simToRecoAssociatorSC"))),
0745       tsRecoToSimCP_token_(
0746           consumes<hgcal::RecoToSimCollectionSimTracksters>(ps.getParameter<edm::InputTag>("recoToSimAssociatorCP"))),
0747       tsSimToRecoCP_token_(
0748           consumes<hgcal::SimToRecoCollectionSimTracksters>(ps.getParameter<edm::InputTag>("simToRecoAssociatorCP"))),
0749       MergeRecoToSimSC_token_(consumes<hgcal::RecoToSimCollectionSimTracksters>(
0750           ps.getParameter<edm::InputTag>("MergerecoToSimAssociatorSC"))),
0751       MergeSimToRecoSC_token_(consumes<hgcal::SimToRecoCollectionSimTracksters>(
0752           ps.getParameter<edm::InputTag>("MergesimToRecoAssociatorSC"))),
0753       MergeRecoToSimCP_token_(consumes<hgcal::RecoToSimCollectionSimTracksters>(
0754           ps.getParameter<edm::InputTag>("MergerecoToSimAssociatorCP"))),
0755       MergeSimToRecoCP_token_(consumes<hgcal::SimToRecoCollectionSimTracksters>(
0756           ps.getParameter<edm::InputTag>("MergesimToRecoAssociatorCP"))),
0757       MergeRecoToSimPU_token_(consumes<hgcal::RecoToSimCollectionSimTracksters>(
0758           ps.getParameter<edm::InputTag>("MergerecoToSimAssociatorPU"))),
0759       MergeSimToRecoPU_token_(consumes<hgcal::SimToRecoCollectionSimTracksters>(
0760           ps.getParameter<edm::InputTag>("MergesimToRecoAssociatorPU"))),
0761       simclusters_token_(consumes(ps.getParameter<edm::InputTag>("simclusters"))),
0762       caloparticles_token_(consumes(ps.getParameter<edm::InputTag>("caloparticles"))),
0763       geometry_token_(esConsumes<CaloGeometry, CaloGeometryRecord, edm::Transition::BeginRun>()),
0764       detector_(ps.getParameter<std::string>("detector")),
0765       propName_(ps.getParameter<std::string>("propagator")),
0766       bfield_token_(esConsumes<MagneticField, IdealMagneticFieldRecord, edm::Transition::BeginRun>()),
0767       propagator_token_(
0768           esConsumes<Propagator, TrackingComponentsRecord, edm::Transition::BeginRun>(edm::ESInputTag("", propName_))),
0769       saveLCs_(ps.getParameter<bool>("saveLCs")),
0770       saveCLUE3DTracksters_(ps.getParameter<bool>("saveCLUE3DTracksters")),
0771       saveTrackstersMerged_(ps.getParameter<bool>("saveTrackstersMerged")),
0772       saveSimTrackstersSC_(ps.getParameter<bool>("saveSimTrackstersSC")),
0773       saveSimTrackstersCP_(ps.getParameter<bool>("saveSimTrackstersCP")),
0774       saveTICLCandidate_(ps.getParameter<bool>("saveSimTICLCandidate")),
0775       saveSimTICLCandidate_(ps.getParameter<bool>("saveSimTICLCandidate")),
0776       saveTracks_(ps.getParameter<bool>("saveTracks")),
0777       saveAssociations_(ps.getParameter<bool>("saveAssociations")) {
0778   std::string detectorName_ = (detector_ == "HFNose") ? "HGCalHFNoseSensitive" : "HGCalEESensitive";
0779   hdc_token_ =
0780       esConsumes<HGCalDDDConstants, IdealGeometryRecord, edm::Transition::BeginRun>(edm::ESInputTag("", detectorName_));
0781 };
0782 
0783 TICLDumper::~TICLDumper() { clearVariables(); };
0784 
0785 void TICLDumper::beginRun(edm::Run const&, edm::EventSetup const& es) {
0786   const CaloGeometry& geom = es.getData(caloGeometry_token_);
0787   rhtools_.setGeometry(geom);
0788 
0789   edm::ESHandle<HGCalDDDConstants> hdc = es.getHandle(hdc_token_);
0790   hgcons_ = hdc.product();
0791   edm::ESHandle<MagneticField> bfield_ = es.getHandle(bfield_token_);
0792   edm::ESHandle<Propagator> propagator = es.getHandle(propagator_token_);
0793   initialize(hgcons_, rhtools_, bfield_, propagator);
0794 }
0795 
0796 // Define tree and branches
0797 void TICLDumper::beginJob() {
0798   edm::Service<TFileService> fs;
0799   if (saveCLUE3DTracksters_) {
0800     trackster_tree_ = fs->make<TTree>("tracksters", "TICL tracksters");
0801     trackster_tree_->Branch("event", &ev_event_);
0802     trackster_tree_->Branch("NClusters", &nclusters_);
0803     trackster_tree_->Branch("NTracksters", &ntracksters_);
0804     trackster_tree_->Branch("time", &trackster_time);
0805     trackster_tree_->Branch("timeError", &trackster_timeError);
0806     trackster_tree_->Branch("regressed_energy", &trackster_regressed_energy);
0807     trackster_tree_->Branch("raw_energy", &trackster_raw_energy);
0808     trackster_tree_->Branch("raw_em_energy", &trackster_raw_em_energy);
0809     trackster_tree_->Branch("raw_pt", &trackster_raw_pt);
0810     trackster_tree_->Branch("raw_em_pt", &trackster_raw_em_pt);
0811     trackster_tree_->Branch("barycenter_x", &trackster_barycenter_x);
0812     trackster_tree_->Branch("barycenter_y", &trackster_barycenter_y);
0813     trackster_tree_->Branch("barycenter_z", &trackster_barycenter_z);
0814     trackster_tree_->Branch("barycenter_eta", &trackster_barycenter_eta);
0815     trackster_tree_->Branch("barycenter_phi", &trackster_barycenter_phi);
0816     trackster_tree_->Branch("EV1", &trackster_EV1);
0817     trackster_tree_->Branch("EV2", &trackster_EV2);
0818     trackster_tree_->Branch("EV3", &trackster_EV3);
0819     trackster_tree_->Branch("eVector0_x", &trackster_eVector0_x);
0820     trackster_tree_->Branch("eVector0_y", &trackster_eVector0_y);
0821     trackster_tree_->Branch("eVector0_z", &trackster_eVector0_z);
0822     trackster_tree_->Branch("sigmaPCA1", &trackster_sigmaPCA1);
0823     trackster_tree_->Branch("sigmaPCA2", &trackster_sigmaPCA2);
0824     trackster_tree_->Branch("sigmaPCA3", &trackster_sigmaPCA3);
0825     trackster_tree_->Branch("id_probabilities", &trackster_id_probabilities);
0826     trackster_tree_->Branch("vertices_indexes", &trackster_vertices_indexes);
0827     trackster_tree_->Branch("vertices_x", &trackster_vertices_x);
0828     trackster_tree_->Branch("vertices_y", &trackster_vertices_y);
0829     trackster_tree_->Branch("vertices_z", &trackster_vertices_z);
0830     trackster_tree_->Branch("vertices_time", &trackster_vertices_time);
0831     trackster_tree_->Branch("vertices_timeErr", &trackster_vertices_timeErr);
0832     trackster_tree_->Branch("vertices_energy", &trackster_vertices_energy);
0833     trackster_tree_->Branch("vertices_correctedEnergy", &trackster_vertices_correctedEnergy);
0834     trackster_tree_->Branch("vertices_correctedEnergyUncertainty", &trackster_vertices_correctedEnergyUncertainty);
0835     trackster_tree_->Branch("vertices_multiplicity", &trackster_vertices_multiplicity);
0836   }
0837   if (saveLCs_) {
0838     cluster_tree_ = fs->make<TTree>("clusters", "TICL tracksters");
0839     cluster_tree_->Branch("seedID", &cluster_seedID);
0840     cluster_tree_->Branch("energy", &cluster_energy);
0841     cluster_tree_->Branch("correctedEnergy", &cluster_correctedEnergy);
0842     cluster_tree_->Branch("correctedEnergyUncertainty", &cluster_correctedEnergyUncertainty);
0843     cluster_tree_->Branch("position_x", &cluster_position_x);
0844     cluster_tree_->Branch("position_y", &cluster_position_y);
0845     cluster_tree_->Branch("position_z", &cluster_position_z);
0846     cluster_tree_->Branch("position_eta", &cluster_position_eta);
0847     cluster_tree_->Branch("position_phi", &cluster_position_phi);
0848     cluster_tree_->Branch("cluster_layer_id", &cluster_layer_id);
0849     cluster_tree_->Branch("cluster_type", &cluster_type);
0850     cluster_tree_->Branch("cluster_time", &cluster_time);
0851     cluster_tree_->Branch("cluster_timeErr", &cluster_timeErr);
0852     cluster_tree_->Branch("cluster_number_of_hits", &cluster_number_of_hits);
0853   }
0854   if (saveTICLCandidate_) {
0855     candidate_tree_ = fs->make<TTree>("candidates", "TICL candidates");
0856     candidate_tree_->Branch("NCandidates", &nCandidates);
0857     candidate_tree_->Branch("candidate_charge", &candidate_charge);
0858     candidate_tree_->Branch("candidate_pdgId", &candidate_pdgId);
0859     candidate_tree_->Branch("candidate_id_probabilities", &candidate_id_probabilities);
0860     candidate_tree_->Branch("candidate_time", &candidate_time);
0861     candidate_tree_->Branch("candidate_timeErr", &candidate_time_err);
0862     candidate_tree_->Branch("candidate_energy", &candidate_energy);
0863     candidate_tree_->Branch("candidate_px", &candidate_px);
0864     candidate_tree_->Branch("candidate_py", &candidate_py);
0865     candidate_tree_->Branch("candidate_pz", &candidate_pz);
0866     candidate_tree_->Branch("track_in_candidate", &track_in_candidate);
0867     candidate_tree_->Branch("tracksters_in_candidate", &tracksters_in_candidate);
0868   }
0869   if (saveTrackstersMerged_) {
0870     tracksters_merged_tree_ = fs->make<TTree>("trackstersMerged", "TICL tracksters merged");
0871     tracksters_merged_tree_->Branch("event", &ev_event_);
0872     tracksters_merged_tree_->Branch("time", &tracksters_merged_time);
0873     tracksters_merged_tree_->Branch("timeError", &tracksters_merged_timeError);
0874     tracksters_merged_tree_->Branch("regressed_energy", &tracksters_merged_regressed_energy);
0875     tracksters_merged_tree_->Branch("raw_energy", &tracksters_merged_raw_energy);
0876     tracksters_merged_tree_->Branch("raw_em_energy", &tracksters_merged_raw_em_energy);
0877     tracksters_merged_tree_->Branch("raw_pt", &tracksters_merged_raw_pt);
0878     tracksters_merged_tree_->Branch("raw_em_pt", &tracksters_merged_raw_em_pt);
0879     tracksters_merged_tree_->Branch("NTrackstersMerged", &nTrackstersMerged);
0880     tracksters_merged_tree_->Branch("barycenter_x", &tracksters_merged_barycenter_x);
0881     tracksters_merged_tree_->Branch("barycenter_y", &tracksters_merged_barycenter_y);
0882     tracksters_merged_tree_->Branch("barycenter_z", &tracksters_merged_barycenter_z);
0883     tracksters_merged_tree_->Branch("barycenter_eta", &tracksters_merged_barycenter_eta);
0884     tracksters_merged_tree_->Branch("barycenter_phi", &tracksters_merged_barycenter_phi);
0885     tracksters_merged_tree_->Branch("EV1", &tracksters_merged_EV1);
0886     tracksters_merged_tree_->Branch("EV2", &tracksters_merged_EV2);
0887     tracksters_merged_tree_->Branch("EV3", &tracksters_merged_EV3);
0888     tracksters_merged_tree_->Branch("eVector0_x", &tracksters_merged_eVector0_x);
0889     tracksters_merged_tree_->Branch("eVector0_y", &tracksters_merged_eVector0_y);
0890     tracksters_merged_tree_->Branch("eVector0_z", &tracksters_merged_eVector0_z);
0891     tracksters_merged_tree_->Branch("sigmaPCA1", &tracksters_merged_sigmaPCA1);
0892     tracksters_merged_tree_->Branch("sigmaPCA2", &tracksters_merged_sigmaPCA2);
0893     tracksters_merged_tree_->Branch("sigmaPCA3", &tracksters_merged_sigmaPCA3);
0894     tracksters_merged_tree_->Branch("id_probabilities", &tracksters_merged_id_probabilities);
0895     tracksters_merged_tree_->Branch("vertices_indexes", &tracksters_merged_vertices_indexes);
0896     tracksters_merged_tree_->Branch("vertices_x", &tracksters_merged_vertices_x);
0897     tracksters_merged_tree_->Branch("vertices_y", &tracksters_merged_vertices_y);
0898     tracksters_merged_tree_->Branch("vertices_z", &tracksters_merged_vertices_z);
0899     tracksters_merged_tree_->Branch("vertices_time", &tracksters_merged_vertices_time);
0900     tracksters_merged_tree_->Branch("vertices_timeErr", &tracksters_merged_vertices_timeErr);
0901     tracksters_merged_tree_->Branch("vertices_energy", &tracksters_merged_vertices_energy);
0902     tracksters_merged_tree_->Branch("vertices_correctedEnergy", &tracksters_merged_vertices_correctedEnergy);
0903     tracksters_merged_tree_->Branch("vertices_correctedEnergyUncertainty",
0904                                     &tracksters_merged_vertices_correctedEnergyUncertainty);
0905     tracksters_merged_tree_->Branch("vertices_multiplicity", &tracksters_merged_vertices_multiplicity);
0906   }
0907   if (saveAssociations_) {
0908     associations_tree_ = fs->make<TTree>("associations", "Associations");
0909     associations_tree_->Branch("tsCLUE3D_recoToSim_SC", &trackstersCLUE3D_recoToSim_SC);
0910     associations_tree_->Branch("tsCLUE3D_recoToSim_SC_score", &trackstersCLUE3D_recoToSim_SC_score);
0911     associations_tree_->Branch("tsCLUE3D_recoToSim_SC_sharedE", &trackstersCLUE3D_recoToSim_SC_sharedE);
0912     associations_tree_->Branch("tsCLUE3D_simToReco_SC", &trackstersCLUE3D_simToReco_SC);
0913     associations_tree_->Branch("tsCLUE3D_simToReco_SC_score", &trackstersCLUE3D_simToReco_SC_score);
0914     associations_tree_->Branch("tsCLUE3D_simToReco_SC_sharedE", &trackstersCLUE3D_simToReco_SC_sharedE);
0915 
0916     associations_tree_->Branch("tsCLUE3D_recoToSim_CP", &trackstersCLUE3D_recoToSim_CP);
0917     associations_tree_->Branch("tsCLUE3D_recoToSim_CP_score", &trackstersCLUE3D_recoToSim_CP_score);
0918     associations_tree_->Branch("tsCLUE3D_recoToSim_CP_sharedE", &trackstersCLUE3D_recoToSim_CP_sharedE);
0919     associations_tree_->Branch("tsCLUE3D_simToReco_CP", &trackstersCLUE3D_simToReco_CP);
0920     associations_tree_->Branch("tsCLUE3D_simToReco_CP_score", &trackstersCLUE3D_simToReco_CP_score);
0921     associations_tree_->Branch("tsCLUE3D_simToReco_CP_sharedE", &trackstersCLUE3D_simToReco_CP_sharedE);
0922 
0923     associations_tree_->Branch("Mergetstracksters_recoToSim_SC", &MergeTracksters_recoToSim_SC);
0924     associations_tree_->Branch("Mergetstracksters_recoToSim_SC_score", &MergeTracksters_recoToSim_SC_score);
0925     associations_tree_->Branch("Mergetstracksters_recoToSim_SC_sharedE", &MergeTracksters_recoToSim_SC_sharedE);
0926     associations_tree_->Branch("Mergetstracksters_simToReco_SC", &MergeTracksters_simToReco_SC);
0927     associations_tree_->Branch("Mergetstracksters_simToReco_SC_score", &MergeTracksters_simToReco_SC_score);
0928     associations_tree_->Branch("Mergetstracksters_simToReco_SC_sharedE", &MergeTracksters_simToReco_SC_sharedE);
0929 
0930     associations_tree_->Branch("Mergetracksters_recoToSim_CP", &MergeTracksters_recoToSim_CP);
0931     associations_tree_->Branch("Mergetracksters_recoToSim_CP_score", &MergeTracksters_recoToSim_CP_score);
0932     associations_tree_->Branch("Mergetracksters_recoToSim_CP_sharedE", &MergeTracksters_recoToSim_CP_sharedE);
0933     associations_tree_->Branch("Mergetracksters_simToReco_CP", &MergeTracksters_simToReco_CP);
0934     associations_tree_->Branch("Mergetracksters_simToReco_CP_score", &MergeTracksters_simToReco_CP_score);
0935     associations_tree_->Branch("Mergetracksters_simToReco_CP_sharedE", &MergeTracksters_simToReco_CP_sharedE);
0936 
0937     associations_tree_->Branch("Mergetracksters_recoToSim_PU", &MergeTracksters_recoToSim_PU);
0938     associations_tree_->Branch("Mergetracksters_recoToSim_PU_score", &MergeTracksters_recoToSim_PU_score);
0939     associations_tree_->Branch("Mergetracksters_recoToSim_PU_sharedE", &MergeTracksters_recoToSim_PU_sharedE);
0940     associations_tree_->Branch("Mergetracksters_simToReco_PU", &MergeTracksters_simToReco_PU);
0941     associations_tree_->Branch("Mergetracksters_simToReco_PU_score", &MergeTracksters_simToReco_PU_score);
0942     associations_tree_->Branch("Mergetracksters_simToReco_PU_sharedE", &MergeTracksters_simToReco_PU_sharedE);
0943   }
0944 
0945   if (saveSimTrackstersSC_) {
0946     simtrackstersSC_tree_ = fs->make<TTree>("simtrackstersSC", "TICL simTracksters SC");
0947     simtrackstersSC_tree_->Branch("event", &ev_event_);
0948     simtrackstersSC_tree_->Branch("NTracksters", &stsSC_ntracksters_);
0949     simtrackstersSC_tree_->Branch("time", &stsSC_trackster_time);
0950     simtrackstersSC_tree_->Branch("timeError", &stsSC_trackster_timeError);
0951     simtrackstersSC_tree_->Branch("regressed_energy", &stsSC_trackster_regressed_energy);
0952     simtrackstersSC_tree_->Branch("regressed_pt", &stsSC_trackster_regressed_pt);
0953     simtrackstersSC_tree_->Branch("raw_energy", &stsSC_trackster_raw_energy);
0954     simtrackstersSC_tree_->Branch("raw_em_energy", &stsSC_trackster_raw_em_energy);
0955     simtrackstersSC_tree_->Branch("raw_pt", &stsSC_trackster_raw_pt);
0956     simtrackstersSC_tree_->Branch("raw_em_pt", &stsSC_trackster_raw_em_pt);
0957     simtrackstersSC_tree_->Branch("barycenter_x", &stsSC_trackster_barycenter_x);
0958     simtrackstersSC_tree_->Branch("barycenter_y", &stsSC_trackster_barycenter_y);
0959     simtrackstersSC_tree_->Branch("barycenter_z", &stsSC_trackster_barycenter_z);
0960     simtrackstersSC_tree_->Branch("barycenter_eta", &stsSC_trackster_barycenter_eta);
0961     simtrackstersSC_tree_->Branch("barycenter_phi", &stsSC_trackster_barycenter_phi);
0962     simtrackstersSC_tree_->Branch("EV1", &stsSC_trackster_EV1);
0963     simtrackstersSC_tree_->Branch("EV2", &stsSC_trackster_EV2);
0964     simtrackstersSC_tree_->Branch("EV3", &stsSC_trackster_EV3);
0965     simtrackstersSC_tree_->Branch("eVector0_x", &stsSC_trackster_eVector0_x);
0966     simtrackstersSC_tree_->Branch("eVector0_y", &stsSC_trackster_eVector0_y);
0967     simtrackstersSC_tree_->Branch("eVector0_z", &stsSC_trackster_eVector0_z);
0968     simtrackstersSC_tree_->Branch("sigmaPCA1", &stsSC_trackster_sigmaPCA1);
0969     simtrackstersSC_tree_->Branch("sigmaPCA2", &stsSC_trackster_sigmaPCA2);
0970     simtrackstersSC_tree_->Branch("sigmaPCA3", &stsSC_trackster_sigmaPCA3);
0971     simtrackstersSC_tree_->Branch("pdgID", &stsSC_pdgID);
0972     simtrackstersSC_tree_->Branch("trackIdx", &stsSC_trackIdx);
0973     simtrackstersSC_tree_->Branch("trackTime", &stsSC_trackTime);
0974     simtrackstersSC_tree_->Branch("boundaryX", &stsSC_boundaryX);
0975     simtrackstersSC_tree_->Branch("boundaryY", &stsSC_boundaryY);
0976     simtrackstersSC_tree_->Branch("boundaryZ", &stsSC_boundaryZ);
0977     simtrackstersSC_tree_->Branch("boundaryEta", &stsSC_boundaryEta);
0978     simtrackstersSC_tree_->Branch("boundaryPhi", &stsSC_boundaryPhi);
0979     simtrackstersSC_tree_->Branch("boundaryPx", &stsSC_boundaryPx);
0980     simtrackstersSC_tree_->Branch("boundaryPy", &stsSC_boundaryPy);
0981     simtrackstersSC_tree_->Branch("boundaryPz", &stsSC_boundaryPz);
0982     simtrackstersSC_tree_->Branch("track_boundaryX", &stsSC_track_boundaryX);
0983     simtrackstersSC_tree_->Branch("track_boundaryY", &stsSC_track_boundaryY);
0984     simtrackstersSC_tree_->Branch("track_boundaryZ", &stsSC_track_boundaryZ);
0985     simtrackstersSC_tree_->Branch("track_boundaryEta", &stsSC_track_boundaryEta);
0986     simtrackstersSC_tree_->Branch("track_boundaryPhi", &stsSC_track_boundaryPhi);
0987     simtrackstersSC_tree_->Branch("track_boundaryPx", &stsSC_track_boundaryPx);
0988     simtrackstersSC_tree_->Branch("track_boundaryPy", &stsSC_track_boundaryPy);
0989     simtrackstersSC_tree_->Branch("track_boundaryPz", &stsSC_track_boundaryPz);
0990     simtrackstersSC_tree_->Branch("id_probabilities", &stsSC_trackster_id_probabilities);
0991     simtrackstersSC_tree_->Branch("vertices_indexes", &stsSC_trackster_vertices_indexes);
0992     simtrackstersSC_tree_->Branch("vertices_x", &stsSC_trackster_vertices_x);
0993     simtrackstersSC_tree_->Branch("vertices_y", &stsSC_trackster_vertices_y);
0994     simtrackstersSC_tree_->Branch("vertices_z", &stsSC_trackster_vertices_z);
0995     simtrackstersSC_tree_->Branch("vertices_time", &stsSC_trackster_vertices_time);
0996     simtrackstersSC_tree_->Branch("vertices_timeErr", &stsSC_trackster_vertices_timeErr);
0997     simtrackstersSC_tree_->Branch("vertices_energy", &stsSC_trackster_vertices_energy);
0998     simtrackstersSC_tree_->Branch("vertices_correctedEnergy", &stsSC_trackster_vertices_correctedEnergy);
0999     simtrackstersSC_tree_->Branch("vertices_correctedEnergyUncertainty",
1000                                   &stsSC_trackster_vertices_correctedEnergyUncertainty);
1001     simtrackstersSC_tree_->Branch("vertices_multiplicity", &stsSC_trackster_vertices_multiplicity);
1002     simtrackstersSC_tree_->Branch("NsimTrackstersSC", &nsimTrackstersSC);
1003   }
1004   if (saveSimTrackstersCP_) {
1005     simtrackstersCP_tree_ = fs->make<TTree>("simtrackstersCP", "TICL simTracksters CP");
1006     simtrackstersCP_tree_->Branch("event", &ev_event_);
1007     simtrackstersCP_tree_->Branch("NTracksters", &stsCP_ntracksters_);
1008     simtrackstersCP_tree_->Branch("time", &stsCP_trackster_time);
1009     simtrackstersCP_tree_->Branch("timeError", &stsCP_trackster_timeError);
1010     simtrackstersCP_tree_->Branch("regressed_energy", &stsCP_trackster_regressed_energy);
1011     simtrackstersCP_tree_->Branch("regressed_pt", &stsCP_trackster_regressed_pt);
1012     simtrackstersCP_tree_->Branch("raw_energy", &stsCP_trackster_raw_energy);
1013     simtrackstersCP_tree_->Branch("raw_em_energy", &stsCP_trackster_raw_em_energy);
1014     simtrackstersCP_tree_->Branch("raw_pt", &stsCP_trackster_raw_pt);
1015     simtrackstersCP_tree_->Branch("raw_em_pt", &stsCP_trackster_raw_em_pt);
1016     simtrackstersCP_tree_->Branch("barycenter_x", &stsCP_trackster_barycenter_x);
1017     simtrackstersCP_tree_->Branch("barycenter_y", &stsCP_trackster_barycenter_y);
1018     simtrackstersCP_tree_->Branch("barycenter_z", &stsCP_trackster_barycenter_z);
1019     simtrackstersCP_tree_->Branch("barycenter_eta", &stsCP_trackster_barycenter_eta);
1020     simtrackstersCP_tree_->Branch("barycenter_phi", &stsCP_trackster_barycenter_phi);
1021     simtrackstersCP_tree_->Branch("pdgID", &stsCP_pdgID);
1022     simtrackstersCP_tree_->Branch("trackIdx", &stsCP_trackIdx);
1023     simtrackstersCP_tree_->Branch("trackTime", &stsCP_trackTime);
1024     simtrackstersCP_tree_->Branch("boundaryX", &stsCP_boundaryX);
1025     simtrackstersCP_tree_->Branch("boundaryY", &stsCP_boundaryY);
1026     simtrackstersCP_tree_->Branch("boundaryZ", &stsCP_boundaryZ);
1027     simtrackstersCP_tree_->Branch("boundaryEta", &stsCP_boundaryEta);
1028     simtrackstersCP_tree_->Branch("boundaryPhi", &stsCP_boundaryPhi);
1029     simtrackstersCP_tree_->Branch("boundaryPx", &stsCP_boundaryPx);
1030     simtrackstersCP_tree_->Branch("boundaryPy", &stsCP_boundaryPy);
1031     simtrackstersCP_tree_->Branch("boundaryPz", &stsCP_boundaryPz);
1032     simtrackstersCP_tree_->Branch("track_boundaryX", &stsCP_track_boundaryX);
1033     simtrackstersCP_tree_->Branch("track_boundaryY", &stsCP_track_boundaryY);
1034     simtrackstersCP_tree_->Branch("track_boundaryZ", &stsCP_track_boundaryZ);
1035     simtrackstersCP_tree_->Branch("track_boundaryEta", &stsCP_track_boundaryEta);
1036     simtrackstersCP_tree_->Branch("track_boundaryPhi", &stsCP_track_boundaryPhi);
1037     simtrackstersCP_tree_->Branch("track_boundaryPx", &stsCP_track_boundaryPx);
1038     simtrackstersCP_tree_->Branch("track_boundaryPy", &stsCP_track_boundaryPy);
1039     simtrackstersCP_tree_->Branch("track_boundaryPz", &stsCP_track_boundaryPz);
1040     simtrackstersCP_tree_->Branch("EV1", &stsCP_trackster_EV1);
1041     simtrackstersCP_tree_->Branch("EV2", &stsCP_trackster_EV2);
1042     simtrackstersCP_tree_->Branch("EV3", &stsCP_trackster_EV3);
1043     simtrackstersCP_tree_->Branch("eVector0_x", &stsCP_trackster_eVector0_x);
1044     simtrackstersCP_tree_->Branch("eVector0_y", &stsCP_trackster_eVector0_y);
1045     simtrackstersCP_tree_->Branch("eVector0_z", &stsCP_trackster_eVector0_z);
1046     simtrackstersCP_tree_->Branch("sigmaPCA1", &stsCP_trackster_sigmaPCA1);
1047     simtrackstersCP_tree_->Branch("sigmaPCA2", &stsCP_trackster_sigmaPCA2);
1048     simtrackstersCP_tree_->Branch("sigmaPCA3", &stsCP_trackster_sigmaPCA3);
1049     simtrackstersCP_tree_->Branch("id_probabilities", &stsCP_trackster_id_probabilities);
1050     simtrackstersCP_tree_->Branch("vertices_indexes", &stsCP_trackster_vertices_indexes);
1051     simtrackstersCP_tree_->Branch("vertices_x", &stsCP_trackster_vertices_x);
1052     simtrackstersCP_tree_->Branch("vertices_y", &stsCP_trackster_vertices_y);
1053     simtrackstersCP_tree_->Branch("vertices_z", &stsCP_trackster_vertices_z);
1054     simtrackstersCP_tree_->Branch("vertices_time", &stsCP_trackster_vertices_time);
1055     simtrackstersCP_tree_->Branch("vertices_timeErr", &stsCP_trackster_vertices_timeErr);
1056     simtrackstersCP_tree_->Branch("vertices_energy", &stsCP_trackster_vertices_energy);
1057     simtrackstersCP_tree_->Branch("vertices_correctedEnergy", &stsCP_trackster_vertices_correctedEnergy);
1058     simtrackstersCP_tree_->Branch("vertices_correctedEnergyUncertainty",
1059                                   &stsCP_trackster_vertices_correctedEnergyUncertainty);
1060     simtrackstersCP_tree_->Branch("vertices_multiplicity", &stsCP_trackster_vertices_multiplicity);
1061   }
1062 
1063   if (saveTracks_) {
1064     tracks_tree_ = fs->make<TTree>("tracks", "Tracks");
1065     tracks_tree_->Branch("event", &ev_event_);
1066     tracks_tree_->Branch("track_id", &track_id);
1067     tracks_tree_->Branch("track_hgcal_pt", &track_hgcal_pt);
1068     tracks_tree_->Branch("track_pt", &track_pt);
1069     tracks_tree_->Branch("track_missing_outer_hits", &track_missing_outer_hits);
1070     tracks_tree_->Branch("track_quality", &track_quality);
1071     tracks_tree_->Branch("track_charge", &track_charge);
1072     tracks_tree_->Branch("track_time", &track_time);
1073     tracks_tree_->Branch("track_time_quality", &track_time_quality);
1074     tracks_tree_->Branch("track_time_err", &track_time_err);
1075     tracks_tree_->Branch("track_nhits", &track_nhits);
1076   }
1077 
1078   if (saveSimTICLCandidate_) {
1079     simTICLCandidate_tree = fs->make<TTree>("simTICLCandidate", "Sim TICL Candidate");
1080     simTICLCandidate_tree->Branch("simTICLCandidate_raw_energy", &simTICLCandidate_raw_energy);
1081     simTICLCandidate_tree->Branch("simTICLCandidate_regressed_energy", &simTICLCandidate_regressed_energy);
1082     simTICLCandidate_tree->Branch("simTICLCandidate_simTracksterCPIndex", &simTICLCandidate_simTracksterCPIndex);
1083     simTICLCandidate_tree->Branch("simTICLCandidate_boundaryX", &simTICLCandidate_boundaryX);
1084     simTICLCandidate_tree->Branch("simTICLCandidate_boundaryY", &simTICLCandidate_boundaryY);
1085     simTICLCandidate_tree->Branch("simTICLCandidate_boundaryZ", &simTICLCandidate_boundaryZ);
1086     simTICLCandidate_tree->Branch("simTICLCandidate_boundaryPx", &simTICLCandidate_boundaryPx);
1087     simTICLCandidate_tree->Branch("simTICLCandidate_boundaryPy", &simTICLCandidate_boundaryPy);
1088     simTICLCandidate_tree->Branch("simTICLCandidate_boundaryPz", &simTICLCandidate_boundaryPz);
1089     simTICLCandidate_tree->Branch("simTICLCandidate_trackTime", &simTICLCandidate_trackTime);
1090     simTICLCandidate_tree->Branch("simTICLCandidate_trackBeta", &simTICLCandidate_trackBeta);
1091     simTICLCandidate_tree->Branch("simTICLCandidate_caloParticleMass", &simTICLCandidate_caloParticleMass);
1092     simTICLCandidate_tree->Branch("simTICLCandidate_pdgId", &simTICLCandidate_pdgId);
1093     simTICLCandidate_tree->Branch("simTICLCandidate_charge", &simTICLCandidate_charge);
1094     simTICLCandidate_tree->Branch("simTICLCandidate_track_in_candidate", &simTICLCandidate_track_in_candidate);
1095   }
1096 }
1097 
1098 void TICLDumper::buildLayers() {
1099   // build disks at HGCal front & EM-Had interface for track propagation
1100 
1101   float zVal = hgcons_->waferZ(1, true);
1102   std::pair<float, float> rMinMax = hgcons_->rangeR(zVal, true);
1103 
1104   float zVal_interface = rhtools_.getPositionLayer(rhtools_.lastLayerEE()).z();
1105   std::pair<float, float> rMinMax_interface = hgcons_->rangeR(zVal_interface, true);
1106 
1107   for (int iSide = 0; iSide < 2; ++iSide) {
1108     float zSide = (iSide == 0) ? (-1. * zVal) : zVal;
1109     firstDisk_[iSide] =
1110         std::make_unique<GeomDet>(Disk::build(Disk::PositionType(0, 0, zSide),
1111                                               Disk::RotationType(),
1112                                               SimpleDiskBounds(rMinMax.first, rMinMax.second, zSide - 0.5, zSide + 0.5))
1113                                       .get());
1114 
1115     zSide = (iSide == 0) ? (-1. * zVal_interface) : zVal_interface;
1116     interfaceDisk_[iSide] = std::make_unique<GeomDet>(
1117         Disk::build(Disk::PositionType(0, 0, zSide),
1118                     Disk::RotationType(),
1119                     SimpleDiskBounds(rMinMax_interface.first, rMinMax_interface.second, zSide - 0.5, zSide + 0.5))
1120             .get());
1121   }
1122 }
1123 
1124 void TICLDumper::initialize(const HGCalDDDConstants* hgcons,
1125                             const hgcal::RecHitTools rhtools,
1126                             const edm::ESHandle<MagneticField> bfieldH,
1127                             const edm::ESHandle<Propagator> propH) {
1128   hgcons_ = hgcons;
1129   rhtools_ = rhtools;
1130   buildLayers();
1131 
1132   bfield_ = bfieldH;
1133   propagator_ = propH;
1134 }
1135 void TICLDumper::analyze(const edm::Event& event, const edm::EventSetup& setup) {
1136   ev_event_ += 1;
1137   clearVariables();
1138   auto bFieldProd = bfield_.product();
1139   const Propagator& prop = (*propagator_);
1140   //get all the tracksters
1141   edm::Handle<std::vector<ticl::Trackster>> tracksters_handle;
1142   event.getByToken(tracksters_token_, tracksters_handle);
1143   const auto& tracksters = *tracksters_handle;
1144 
1145   //get all the layer clusters
1146   edm::Handle<std::vector<reco::CaloCluster>> layer_clusters_h;
1147   event.getByToken(layer_clusters_token_, layer_clusters_h);
1148   const auto& clusters = *layer_clusters_h;
1149 
1150   edm::Handle<edm::ValueMap<std::pair<float, float>>> clustersTime_h;
1151   event.getByToken(clustersTime_token_, clustersTime_h);
1152   const auto& layerClustersTimes = *clustersTime_h;
1153 
1154   //TICL Candidate
1155   edm::Handle<std::vector<TICLCandidate>> candidates_h;
1156   event.getByToken(ticl_candidates_token_, candidates_h);
1157   const auto& ticlcandidates = *candidates_h;
1158 
1159   //Track
1160   edm::Handle<std::vector<reco::Track>> tracks_h;
1161   event.getByToken(tracks_token_, tracks_h);
1162   const auto& tracks = *tracks_h;
1163 
1164   edm::Handle<edm::ValueMap<float>> trackTime_h;
1165   event.getByToken(tracks_time_token_, trackTime_h);
1166   const auto& trackTime = *trackTime_h;
1167 
1168   edm::Handle<edm::ValueMap<float>> trackTimeErr_h;
1169   event.getByToken(tracks_time_err_token_, trackTimeErr_h);
1170   const auto& trackTimeErr = *trackTimeErr_h;
1171 
1172   edm::Handle<edm::ValueMap<float>> trackTimeQual_h;
1173   event.getByToken(tracks_time_quality_token_, trackTimeQual_h);
1174   const auto& trackTimeQual = *trackTimeQual_h;
1175 
1176   //Tracksters merged
1177   edm::Handle<std::vector<ticl::Trackster>> tracksters_merged_h;
1178   event.getByToken(tracksters_merged_token_, tracksters_merged_h);
1179   const auto& trackstersmerged = *tracksters_merged_h;
1180 
1181   // simTracksters from SC
1182   edm::Handle<std::vector<ticl::Trackster>> simTrackstersSC_h;
1183   event.getByToken(simTracksters_SC_token_, simTrackstersSC_h);
1184   const auto& simTrackstersSC = *simTrackstersSC_h;
1185 
1186   // simTracksters from CP
1187   edm::Handle<std::vector<ticl::Trackster>> simTrackstersCP_h;
1188   event.getByToken(simTracksters_CP_token_, simTrackstersCP_h);
1189   const auto& simTrackstersCP = *simTrackstersCP_h;
1190 
1191   // simTracksters from PU
1192   edm::Handle<std::vector<ticl::Trackster>> simTrackstersPU_h;
1193   event.getByToken(simTracksters_PU_token_, simTrackstersPU_h);
1194   const auto& simTrackstersPU = *simTrackstersPU_h;
1195 
1196   edm::Handle<std::vector<TICLCandidate>> simTICLCandidates_h;
1197   event.getByToken(simTICLCandidate_token_, simTICLCandidates_h);
1198   const auto& simTICLCandidates = *simTICLCandidates_h;
1199 
1200   // trackster reco to sim SC
1201   edm::Handle<hgcal::RecoToSimCollectionSimTracksters> tsRecoToSimSC_h;
1202   event.getByToken(tsRecoToSimSC_token_, tsRecoToSimSC_h);
1203   auto const& tsRecoSimSCMap = *tsRecoToSimSC_h;
1204 
1205   // sim simTrackster SC to reco trackster
1206   edm::Handle<hgcal::SimToRecoCollectionSimTracksters> tsSimToRecoSC_h;
1207   event.getByToken(tsSimToRecoSC_token_, tsSimToRecoSC_h);
1208   auto const& tsSimToRecoSCMap = *tsSimToRecoSC_h;
1209 
1210   // trackster reco to sim CP
1211   edm::Handle<hgcal::RecoToSimCollectionSimTracksters> tsRecoToSimCP_h;
1212   event.getByToken(tsRecoToSimCP_token_, tsRecoToSimCP_h);
1213   auto const& tsRecoSimCPMap = *tsRecoToSimCP_h;
1214 
1215   // sim simTrackster CP to reco trackster
1216   edm::Handle<hgcal::SimToRecoCollectionSimTracksters> tsSimToRecoCP_h;
1217   event.getByToken(tsSimToRecoCP_token_, tsSimToRecoCP_h);
1218   auto const& tsSimToRecoCPMap = *tsSimToRecoCP_h;
1219 
1220   edm::Handle<hgcal::RecoToSimCollectionSimTracksters> mergetsRecoToSimSC_h;
1221   event.getByToken(MergeRecoToSimSC_token_, mergetsRecoToSimSC_h);
1222   auto const& MergetsRecoSimSCMap = *mergetsRecoToSimSC_h;
1223 
1224   // sim simTrackster SC to reco trackster
1225   edm::Handle<hgcal::SimToRecoCollectionSimTracksters> mergetsSimToRecoSC_h;
1226   event.getByToken(MergeSimToRecoSC_token_, mergetsSimToRecoSC_h);
1227   auto const& MergetsSimToRecoSCMap = *mergetsSimToRecoSC_h;
1228 
1229   // trackster reco to sim CP
1230   edm::Handle<hgcal::RecoToSimCollectionSimTracksters> mergetsRecoToSimCP_h;
1231   event.getByToken(MergeRecoToSimCP_token_, mergetsRecoToSimCP_h);
1232   auto const& MergetsRecoSimCPMap = *mergetsRecoToSimCP_h;
1233 
1234   // sim simTrackster CP to reco trackster
1235   edm::Handle<hgcal::SimToRecoCollectionSimTracksters> mergetsSimToRecoCP_h;
1236   event.getByToken(MergeSimToRecoCP_token_, mergetsSimToRecoCP_h);
1237   auto const& MergetsSimToRecoCPMap = *mergetsSimToRecoCP_h;
1238 
1239   // trackster reco to sim PU
1240   edm::Handle<hgcal::RecoToSimCollectionSimTracksters> mergetsRecoToSimPU_h;
1241   event.getByToken(MergeRecoToSimPU_token_, mergetsRecoToSimPU_h);
1242   auto const& MergetsRecoSimPUMap = *mergetsRecoToSimPU_h;
1243 
1244   // sim simTrackster PU to reco trackster
1245   edm::Handle<hgcal::SimToRecoCollectionSimTracksters> mergetsSimToRecoPU_h;
1246   event.getByToken(MergeSimToRecoPU_token_, mergetsSimToRecoPU_h);
1247   auto const& MergetsSimToRecoPUMap = *mergetsSimToRecoPU_h;
1248 
1249   edm::Handle<std::vector<CaloParticle>> caloparticles_h;
1250   event.getByToken(caloparticles_token_, caloparticles_h);
1251   const auto& caloparticles = *caloparticles_h;
1252 
1253   const auto& simclusters = event.get(simclusters_token_);
1254 
1255   ntracksters_ = tracksters.size();
1256   nclusters_ = clusters.size();
1257 
1258   for (auto trackster_iterator = tracksters.begin(); trackster_iterator != tracksters.end(); ++trackster_iterator) {
1259     //per-trackster analysis
1260     trackster_time.push_back(trackster_iterator->time());
1261     trackster_timeError.push_back(trackster_iterator->timeError());
1262     trackster_regressed_energy.push_back(trackster_iterator->regressed_energy());
1263     trackster_raw_energy.push_back(trackster_iterator->raw_energy());
1264     trackster_raw_em_energy.push_back(trackster_iterator->raw_em_energy());
1265     trackster_raw_pt.push_back(trackster_iterator->raw_pt());
1266     trackster_raw_em_pt.push_back(trackster_iterator->raw_em_pt());
1267     trackster_barycenter_x.push_back(trackster_iterator->barycenter().x());
1268     trackster_barycenter_y.push_back(trackster_iterator->barycenter().y());
1269     trackster_barycenter_z.push_back(trackster_iterator->barycenter().z());
1270     trackster_barycenter_eta.push_back(trackster_iterator->barycenter().eta());
1271     trackster_barycenter_phi.push_back(trackster_iterator->barycenter().phi());
1272     trackster_EV1.push_back(trackster_iterator->eigenvalues()[0]);
1273     trackster_EV2.push_back(trackster_iterator->eigenvalues()[1]);
1274     trackster_EV3.push_back(trackster_iterator->eigenvalues()[2]);
1275     trackster_eVector0_x.push_back((trackster_iterator->eigenvectors()[0]).x());
1276     trackster_eVector0_y.push_back((trackster_iterator->eigenvectors()[0]).y());
1277     trackster_eVector0_z.push_back((trackster_iterator->eigenvectors()[0]).z());
1278     trackster_sigmaPCA1.push_back(trackster_iterator->sigmasPCA()[0]);
1279     trackster_sigmaPCA2.push_back(trackster_iterator->sigmasPCA()[1]);
1280     trackster_sigmaPCA3.push_back(trackster_iterator->sigmasPCA()[2]);
1281     std::vector<float> id_probs;
1282     for (size_t i = 0; i < 8; i++)
1283       id_probs.push_back(trackster_iterator->id_probabilities(i));
1284     trackster_id_probabilities.push_back(id_probs);
1285 
1286     // Clusters
1287     std::vector<uint32_t> vertices_indexes;
1288     std::vector<float> vertices_x;
1289     std::vector<float> vertices_y;
1290     std::vector<float> vertices_z;
1291     std::vector<float> vertices_time;
1292     std::vector<float> vertices_timeErr;
1293     std::vector<float> vertices_energy;
1294     std::vector<float> vertices_correctedEnergy;
1295     std::vector<float> vertices_correctedEnergyUncertainty;
1296     for (auto idx : trackster_iterator->vertices()) {
1297       vertices_indexes.push_back(idx);
1298       auto associated_cluster = (*layer_clusters_h)[idx];
1299       vertices_x.push_back(associated_cluster.x());
1300       vertices_y.push_back(associated_cluster.y());
1301       vertices_z.push_back(associated_cluster.z());
1302       vertices_energy.push_back(associated_cluster.energy());
1303       vertices_correctedEnergy.push_back(associated_cluster.correctedEnergy());
1304       vertices_correctedEnergyUncertainty.push_back(associated_cluster.correctedEnergyUncertainty());
1305       vertices_time.push_back(layerClustersTimes.get(idx).first);
1306       vertices_timeErr.push_back(layerClustersTimes.get(idx).second);
1307     }
1308     trackster_vertices_indexes.push_back(vertices_indexes);
1309     trackster_vertices_x.push_back(vertices_x);
1310     trackster_vertices_y.push_back(vertices_y);
1311     trackster_vertices_z.push_back(vertices_z);
1312     trackster_vertices_time.push_back(vertices_time);
1313     trackster_vertices_timeErr.push_back(vertices_timeErr);
1314     trackster_vertices_energy.push_back(vertices_energy);
1315     trackster_vertices_correctedEnergy.push_back(vertices_correctedEnergy);
1316     trackster_vertices_correctedEnergyUncertainty.push_back(vertices_correctedEnergyUncertainty);
1317 
1318     // Multiplicity
1319     std::vector<float> vertices_multiplicity;
1320     for (auto multiplicity : trackster_iterator->vertex_multiplicity()) {
1321       vertices_multiplicity.push_back(multiplicity);
1322     }
1323     trackster_vertices_multiplicity.push_back(vertices_multiplicity);
1324   }
1325 
1326   stsSC_ntracksters_ = simTrackstersSC.size();
1327   using CaloObjectVariant = std::variant<CaloParticle, SimCluster>;
1328   for (auto trackster_iterator = simTrackstersSC.begin(); trackster_iterator != simTrackstersSC.end();
1329        ++trackster_iterator) {
1330     //per-trackster analysis
1331     stsSC_trackster_time.push_back(trackster_iterator->time());
1332     stsSC_trackster_timeError.push_back(trackster_iterator->timeError());
1333     stsSC_trackster_regressed_energy.push_back(trackster_iterator->regressed_energy());
1334     stsSC_trackster_raw_energy.push_back(trackster_iterator->raw_energy());
1335     stsSC_trackster_raw_em_energy.push_back(trackster_iterator->raw_em_energy());
1336     stsSC_trackster_raw_pt.push_back(trackster_iterator->raw_pt());
1337     stsSC_trackster_raw_em_pt.push_back(trackster_iterator->raw_em_pt());
1338     stsSC_trackster_barycenter_x.push_back(trackster_iterator->barycenter().x());
1339     stsSC_trackster_barycenter_y.push_back(trackster_iterator->barycenter().y());
1340     stsSC_trackster_barycenter_z.push_back(trackster_iterator->barycenter().z());
1341     stsSC_trackster_barycenter_eta.push_back(trackster_iterator->barycenter().eta());
1342     stsSC_trackster_barycenter_phi.push_back(trackster_iterator->barycenter().phi());
1343     stsSC_trackster_EV1.push_back(trackster_iterator->eigenvalues()[0]);
1344     stsSC_trackster_EV2.push_back(trackster_iterator->eigenvalues()[1]);
1345     stsSC_trackster_EV3.push_back(trackster_iterator->eigenvalues()[2]);
1346     stsSC_trackster_eVector0_x.push_back((trackster_iterator->eigenvectors()[0]).x());
1347     stsSC_trackster_eVector0_y.push_back((trackster_iterator->eigenvectors()[0]).y());
1348     stsSC_trackster_eVector0_z.push_back((trackster_iterator->eigenvectors()[0]).z());
1349     stsSC_trackster_sigmaPCA1.push_back(trackster_iterator->sigmasPCA()[0]);
1350     stsSC_trackster_sigmaPCA2.push_back(trackster_iterator->sigmasPCA()[1]);
1351     stsSC_trackster_sigmaPCA3.push_back(trackster_iterator->sigmasPCA()[2]);
1352     stsSC_pdgID.push_back(simclusters[trackster_iterator->seedIndex()].pdgId());
1353 
1354     CaloObjectVariant caloObj;
1355     if (trackster_iterator->seedID() == caloparticles_h.id()) {
1356       caloObj = caloparticles[trackster_iterator->seedIndex()];
1357     } else {
1358       caloObj = simclusters[trackster_iterator->seedIndex()];
1359     }
1360 
1361     auto const& simTrack = std::visit([](auto&& obj) { return obj.g4Tracks()[0]; }, caloObj);
1362     auto const& caloPt = std::visit([](auto&& obj) { return obj.pt(); }, caloObj);
1363     stsSC_trackster_regressed_pt.push_back(caloPt);
1364     if (simTrack.crossedBoundary()) {
1365       stsSC_boundaryX.push_back(simTrack.getPositionAtBoundary().x());
1366       stsSC_boundaryY.push_back(simTrack.getPositionAtBoundary().y());
1367       stsSC_boundaryZ.push_back(simTrack.getPositionAtBoundary().z());
1368       stsSC_boundaryEta.push_back(simTrack.getPositionAtBoundary().eta());
1369       stsSC_boundaryPhi.push_back(simTrack.getPositionAtBoundary().phi());
1370       stsSC_boundaryPx.push_back(simTrack.getMomentumAtBoundary().x());
1371       stsSC_boundaryPy.push_back(simTrack.getMomentumAtBoundary().y());
1372       stsSC_boundaryPz.push_back(simTrack.getMomentumAtBoundary().z());
1373     } else {
1374       stsSC_boundaryX.push_back(-999);
1375       stsSC_boundaryY.push_back(-999);
1376       stsSC_boundaryZ.push_back(-999);
1377       stsSC_boundaryEta.push_back(-999);
1378       stsSC_boundaryPhi.push_back(-999);
1379       stsSC_boundaryPx.push_back(-999);
1380       stsSC_boundaryPy.push_back(-999);
1381       stsSC_boundaryPz.push_back(-999);
1382     }
1383     auto const trackIdx = trackster_iterator->trackIdx();
1384     stsSC_trackIdx.push_back(trackIdx);
1385     if (trackIdx != -1) {
1386       const auto& track = tracks[trackIdx];
1387 
1388       int iSide = int(track.eta() > 0);
1389 
1390       const auto& fts = trajectoryStateTransform::outerFreeState((track), bFieldProd);
1391       // to the HGCal front
1392       const auto& tsos = prop.propagate(fts, firstDisk_[iSide]->surface());
1393       if (tsos.isValid()) {
1394         const auto& globalPos = tsos.globalPosition();
1395         const auto& globalMom = tsos.globalMomentum();
1396         stsSC_track_boundaryX.push_back(globalPos.x());
1397         stsSC_track_boundaryY.push_back(globalPos.y());
1398         stsSC_track_boundaryZ.push_back(globalPos.z());
1399         stsSC_track_boundaryEta.push_back(globalPos.eta());
1400         stsSC_track_boundaryPhi.push_back(globalPos.phi());
1401         stsSC_track_boundaryPx.push_back(globalMom.x());
1402         stsSC_track_boundaryPy.push_back(globalMom.y());
1403         stsSC_track_boundaryPz.push_back(globalMom.z());
1404         stsSC_trackTime.push_back(track.t0());
1405       } else {
1406         stsSC_track_boundaryX.push_back(-999);
1407         stsSC_track_boundaryY.push_back(-999);
1408         stsSC_track_boundaryZ.push_back(-999);
1409         stsSC_track_boundaryEta.push_back(-999);
1410         stsSC_track_boundaryPhi.push_back(-999);
1411         stsSC_track_boundaryPx.push_back(-999);
1412         stsSC_track_boundaryPy.push_back(-999);
1413         stsSC_track_boundaryPz.push_back(-999);
1414       }
1415     } else {
1416       stsSC_track_boundaryX.push_back(-999);
1417       stsSC_track_boundaryY.push_back(-999);
1418       stsSC_track_boundaryZ.push_back(-999);
1419       stsSC_track_boundaryEta.push_back(-999);
1420       stsSC_track_boundaryPhi.push_back(-999);
1421       stsSC_track_boundaryPx.push_back(-999);
1422       stsSC_track_boundaryPy.push_back(-999);
1423       stsSC_track_boundaryPz.push_back(-999);
1424     }
1425 
1426     std::vector<float> id_probs;
1427     for (size_t i = 0; i < 8; i++)
1428       id_probs.push_back(trackster_iterator->id_probabilities(i));
1429     stsSC_trackster_id_probabilities.push_back(id_probs);
1430 
1431     // Clusters
1432     std::vector<uint32_t> vertices_indexes;
1433     std::vector<float> vertices_x;
1434     std::vector<float> vertices_y;
1435     std::vector<float> vertices_z;
1436     std::vector<float> vertices_time;
1437     std::vector<float> vertices_timeErr;
1438     std::vector<float> vertices_energy;
1439     std::vector<float> vertices_correctedEnergy;
1440     std::vector<float> vertices_correctedEnergyUncertainty;
1441     for (auto idx : trackster_iterator->vertices()) {
1442       vertices_indexes.push_back(idx);
1443       auto associated_cluster = (*layer_clusters_h)[idx];
1444       vertices_x.push_back(associated_cluster.x());
1445       vertices_y.push_back(associated_cluster.y());
1446       vertices_z.push_back(associated_cluster.z());
1447       vertices_energy.push_back(associated_cluster.energy());
1448       vertices_correctedEnergy.push_back(associated_cluster.correctedEnergy());
1449       vertices_correctedEnergyUncertainty.push_back(associated_cluster.correctedEnergyUncertainty());
1450       vertices_time.push_back(layerClustersTimes.get(idx).first);
1451       vertices_timeErr.push_back(layerClustersTimes.get(idx).second);
1452     }
1453     stsSC_trackster_vertices_indexes.push_back(vertices_indexes);
1454     stsSC_trackster_vertices_x.push_back(vertices_x);
1455     stsSC_trackster_vertices_y.push_back(vertices_y);
1456     stsSC_trackster_vertices_z.push_back(vertices_z);
1457     stsSC_trackster_vertices_time.push_back(vertices_time);
1458     stsSC_trackster_vertices_timeErr.push_back(vertices_timeErr);
1459     stsSC_trackster_vertices_energy.push_back(vertices_energy);
1460     stsSC_trackster_vertices_correctedEnergy.push_back(vertices_correctedEnergy);
1461     stsSC_trackster_vertices_correctedEnergyUncertainty.push_back(vertices_correctedEnergyUncertainty);
1462 
1463     // Multiplicity
1464     std::vector<float> vertices_multiplicity;
1465     for (auto multiplicity : trackster_iterator->vertex_multiplicity()) {
1466       vertices_multiplicity.push_back(multiplicity);
1467     }
1468     stsSC_trackster_vertices_multiplicity.push_back(vertices_multiplicity);
1469   }
1470 
1471   stsCP_ntracksters_ = simTrackstersCP.size();
1472 
1473   for (auto trackster_iterator = simTrackstersCP.begin(); trackster_iterator != simTrackstersCP.end();
1474        ++trackster_iterator) {
1475     //per-trackster analysis
1476     stsCP_trackster_time.push_back(trackster_iterator->time());
1477     stsCP_trackster_timeError.push_back(trackster_iterator->timeError());
1478     stsCP_trackster_regressed_energy.push_back(trackster_iterator->regressed_energy());
1479     stsCP_trackster_raw_energy.push_back(trackster_iterator->raw_energy());
1480     stsCP_trackster_raw_em_energy.push_back(trackster_iterator->raw_em_energy());
1481     stsCP_trackster_raw_pt.push_back(trackster_iterator->raw_pt());
1482     stsCP_trackster_raw_em_pt.push_back(trackster_iterator->raw_em_pt());
1483     stsCP_trackster_barycenter_x.push_back(trackster_iterator->barycenter().x());
1484     stsCP_trackster_barycenter_y.push_back(trackster_iterator->barycenter().y());
1485     stsCP_trackster_barycenter_z.push_back(trackster_iterator->barycenter().z());
1486     stsCP_trackster_barycenter_eta.push_back(trackster_iterator->barycenter().eta());
1487     stsCP_trackster_barycenter_phi.push_back(trackster_iterator->barycenter().phi());
1488     stsCP_trackster_EV1.push_back(trackster_iterator->eigenvalues()[0]);
1489     stsCP_trackster_EV2.push_back(trackster_iterator->eigenvalues()[1]);
1490     stsCP_trackster_EV3.push_back(trackster_iterator->eigenvalues()[2]);
1491     stsCP_trackster_eVector0_x.push_back((trackster_iterator->eigenvectors()[0]).x());
1492     stsCP_trackster_eVector0_y.push_back((trackster_iterator->eigenvectors()[0]).y());
1493     stsCP_trackster_eVector0_z.push_back((trackster_iterator->eigenvectors()[0]).z());
1494     stsCP_trackster_sigmaPCA1.push_back(trackster_iterator->sigmasPCA()[0]);
1495     stsCP_trackster_sigmaPCA2.push_back(trackster_iterator->sigmasPCA()[1]);
1496     stsCP_trackster_sigmaPCA3.push_back(trackster_iterator->sigmasPCA()[2]);
1497     stsCP_pdgID.push_back(caloparticles[trackster_iterator->seedIndex()].pdgId());
1498     CaloObjectVariant caloObj;
1499     if (trackster_iterator->seedID() == caloparticles_h.id()) {
1500       caloObj = caloparticles[trackster_iterator->seedIndex()];
1501     } else {
1502       caloObj = simclusters[trackster_iterator->seedIndex()];
1503     }
1504 
1505     auto const& simTrack = std::visit([](auto&& obj) { return obj.g4Tracks()[0]; }, caloObj);
1506     auto const& caloPt = std::visit([](auto&& obj) { return obj.pt(); }, caloObj);
1507     stsCP_trackster_regressed_pt.push_back(caloPt);
1508 
1509     if (simTrack.crossedBoundary()) {
1510       stsCP_boundaryX.push_back(simTrack.getPositionAtBoundary().x());
1511       stsCP_boundaryY.push_back(simTrack.getPositionAtBoundary().y());
1512       stsCP_boundaryZ.push_back(simTrack.getPositionAtBoundary().z());
1513       stsCP_boundaryEta.push_back(simTrack.getPositionAtBoundary().eta());
1514       stsCP_boundaryPhi.push_back(simTrack.getPositionAtBoundary().phi());
1515       stsCP_boundaryPx.push_back(simTrack.getMomentumAtBoundary().x());
1516       stsCP_boundaryPy.push_back(simTrack.getMomentumAtBoundary().y());
1517       stsCP_boundaryPz.push_back(simTrack.getMomentumAtBoundary().z());
1518     } else {
1519       stsCP_boundaryX.push_back(-999);
1520       stsCP_boundaryY.push_back(-999);
1521       stsCP_boundaryZ.push_back(-999);
1522       stsCP_boundaryEta.push_back(-999);
1523       stsCP_boundaryPhi.push_back(-999);
1524       stsCP_boundaryPx.push_back(-999);
1525       stsCP_boundaryPy.push_back(-999);
1526       stsCP_boundaryPz.push_back(-999);
1527     }
1528     auto const trackIdx = trackster_iterator->trackIdx();
1529     stsCP_trackIdx.push_back(trackIdx);
1530     if (trackIdx != -1) {
1531       const auto& track = tracks[trackIdx];
1532 
1533       int iSide = int(track.eta() > 0);
1534 
1535       const auto& fts = trajectoryStateTransform::outerFreeState((track), bFieldProd);
1536       // to the HGCal front
1537       const auto& tsos = prop.propagate(fts, firstDisk_[iSide]->surface());
1538       if (tsos.isValid()) {
1539         const auto& globalPos = tsos.globalPosition();
1540         const auto& globalMom = tsos.globalMomentum();
1541         stsCP_track_boundaryX.push_back(globalPos.x());
1542         stsCP_track_boundaryY.push_back(globalPos.y());
1543         stsCP_track_boundaryZ.push_back(globalPos.z());
1544         stsCP_track_boundaryEta.push_back(globalPos.eta());
1545         stsCP_track_boundaryPhi.push_back(globalPos.phi());
1546         stsCP_track_boundaryPx.push_back(globalMom.x());
1547         stsCP_track_boundaryPy.push_back(globalMom.y());
1548         stsCP_track_boundaryPz.push_back(globalMom.z());
1549         stsCP_trackTime.push_back(track.t0());
1550       } else {
1551         stsCP_track_boundaryX.push_back(-999);
1552         stsCP_track_boundaryY.push_back(-999);
1553         stsCP_track_boundaryZ.push_back(-999);
1554         stsCP_track_boundaryEta.push_back(-999);
1555         stsCP_track_boundaryPhi.push_back(-999);
1556         stsCP_track_boundaryPx.push_back(-999);
1557         stsCP_track_boundaryPy.push_back(-999);
1558         stsCP_track_boundaryPz.push_back(-999);
1559       }
1560     } else {
1561       stsCP_track_boundaryX.push_back(-999);
1562       stsCP_track_boundaryY.push_back(-999);
1563       stsCP_track_boundaryZ.push_back(-999);
1564       stsCP_track_boundaryEta.push_back(-999);
1565       stsCP_track_boundaryPhi.push_back(-999);
1566       stsCP_track_boundaryPx.push_back(-999);
1567       stsCP_track_boundaryPy.push_back(-999);
1568       stsCP_track_boundaryPz.push_back(-999);
1569     }
1570     std::vector<float> id_probs;
1571     for (size_t i = 0; i < 8; i++)
1572       id_probs.push_back(trackster_iterator->id_probabilities(i));
1573     stsCP_trackster_id_probabilities.push_back(id_probs);
1574 
1575     // Clusters
1576     std::vector<uint32_t> vertices_indexes;
1577     std::vector<float> vertices_x;
1578     std::vector<float> vertices_y;
1579     std::vector<float> vertices_z;
1580     std::vector<float> vertices_time;
1581     std::vector<float> vertices_timeErr;
1582     std::vector<float> vertices_energy;
1583     std::vector<float> vertices_correctedEnergy;
1584     std::vector<float> vertices_correctedEnergyUncertainty;
1585     for (auto idx : trackster_iterator->vertices()) {
1586       vertices_indexes.push_back(idx);
1587       auto associated_cluster = (*layer_clusters_h)[idx];
1588       vertices_x.push_back(associated_cluster.x());
1589       vertices_y.push_back(associated_cluster.y());
1590       vertices_z.push_back(associated_cluster.z());
1591       vertices_energy.push_back(associated_cluster.energy());
1592       vertices_correctedEnergy.push_back(associated_cluster.correctedEnergy());
1593       vertices_correctedEnergyUncertainty.push_back(associated_cluster.correctedEnergyUncertainty());
1594       vertices_time.push_back(layerClustersTimes.get(idx).first);
1595       vertices_timeErr.push_back(layerClustersTimes.get(idx).second);
1596     }
1597     stsCP_trackster_vertices_indexes.push_back(vertices_indexes);
1598     stsCP_trackster_vertices_x.push_back(vertices_x);
1599     stsCP_trackster_vertices_y.push_back(vertices_y);
1600     stsCP_trackster_vertices_z.push_back(vertices_z);
1601     stsCP_trackster_vertices_time.push_back(vertices_time);
1602     stsCP_trackster_vertices_timeErr.push_back(vertices_timeErr);
1603     stsCP_trackster_vertices_energy.push_back(vertices_energy);
1604     stsCP_trackster_vertices_correctedEnergy.push_back(vertices_correctedEnergy);
1605     stsCP_trackster_vertices_correctedEnergyUncertainty.push_back(vertices_correctedEnergyUncertainty);
1606 
1607     // Multiplicity
1608     std::vector<float> vertices_multiplicity;
1609     for (auto multiplicity : trackster_iterator->vertex_multiplicity()) {
1610       vertices_multiplicity.push_back(multiplicity);
1611     }
1612     stsCP_trackster_vertices_multiplicity.push_back(vertices_multiplicity);
1613   }
1614 
1615   simTICLCandidate_track_in_candidate.resize(simTICLCandidates.size(), -1);
1616   for (size_t i = 0; i < simTICLCandidates.size(); ++i) {
1617     auto const& cand = simTICLCandidates[i];
1618 
1619     simTICLCandidate_raw_energy.push_back(cand.rawEnergy());
1620     simTICLCandidate_regressed_energy.push_back(cand.p4().energy());
1621     simTICLCandidate_pdgId.push_back(cand.pdgId());
1622     simTICLCandidate_charge.push_back(cand.charge());
1623     std::vector<int> tmpIdxVec;
1624     for (auto const& simTS : cand.tracksters()) {
1625       auto trackster_idx = simTS.get() - (edm::Ptr<ticl::Trackster>(simTrackstersSC_h, 0)).get();
1626       tmpIdxVec.push_back(trackster_idx);
1627     }
1628     simTICLCandidate_simTracksterCPIndex.push_back(tmpIdxVec);
1629     tmpIdxVec.clear();
1630     auto const& trackPtr = cand.trackPtr();
1631     if (!trackPtr.isNull()) {
1632       auto const& track = *trackPtr;
1633       int iSide = int(track.eta() > 0);
1634       int tk_idx = trackPtr.get() - (edm::Ptr<reco::Track>(tracks_h, 0)).get();
1635       simTICLCandidate_track_in_candidate[i] = tk_idx;
1636 
1637       const auto& fts = trajectoryStateTransform::outerFreeState((track), bFieldProd);
1638       // to the HGCal front
1639       const auto& tsos = prop.propagate(fts, firstDisk_[iSide]->surface());
1640       if (tsos.isValid()) {
1641         const auto& globalPos = tsos.globalPosition();
1642         const auto& globalMom = tsos.globalMomentum();
1643         simTICLCandidate_boundaryX.push_back(globalPos.x());
1644         simTICLCandidate_boundaryY.push_back(globalPos.y());
1645         simTICLCandidate_boundaryZ.push_back(globalPos.z());
1646         simTICLCandidate_boundaryPx.push_back(globalMom.x());
1647         simTICLCandidate_boundaryPy.push_back(globalMom.y());
1648         simTICLCandidate_boundaryPz.push_back(globalMom.z());
1649         simTICLCandidate_trackTime.push_back(track.t0());
1650         simTICLCandidate_trackBeta.push_back(track.beta());
1651       } else {
1652         simTICLCandidate_boundaryX.push_back(-999);
1653         simTICLCandidate_boundaryY.push_back(-999);
1654         simTICLCandidate_boundaryZ.push_back(-999);
1655         simTICLCandidate_boundaryPx.push_back(-999);
1656         simTICLCandidate_boundaryPy.push_back(-999);
1657         simTICLCandidate_boundaryPz.push_back(-999);
1658         simTICLCandidate_trackTime.push_back(-999);
1659         simTICLCandidate_trackBeta.push_back(-999);
1660       }
1661     } else {
1662       simTICLCandidate_boundaryX.push_back(-999);
1663       simTICLCandidate_boundaryY.push_back(-999);
1664       simTICLCandidate_boundaryZ.push_back(-999);
1665       simTICLCandidate_boundaryPx.push_back(-999);
1666       simTICLCandidate_boundaryPy.push_back(-999);
1667       simTICLCandidate_boundaryPz.push_back(-999);
1668       simTICLCandidate_trackTime.push_back(-999);
1669       simTICLCandidate_trackBeta.push_back(-999);
1670     }
1671   }
1672 
1673   int c_id = 0;
1674 
1675   for (auto cluster_iterator = clusters.begin(); cluster_iterator != clusters.end(); ++cluster_iterator) {
1676     auto lc_seed = cluster_iterator->seed();
1677     cluster_seedID.push_back(lc_seed);
1678     cluster_energy.push_back(cluster_iterator->energy());
1679     cluster_correctedEnergy.push_back(cluster_iterator->correctedEnergy());
1680     cluster_correctedEnergyUncertainty.push_back(cluster_iterator->correctedEnergyUncertainty());
1681     cluster_position_x.push_back(cluster_iterator->x());
1682     cluster_position_y.push_back(cluster_iterator->y());
1683     cluster_position_z.push_back(cluster_iterator->z());
1684     cluster_position_eta.push_back(cluster_iterator->eta());
1685     cluster_position_phi.push_back(cluster_iterator->phi());
1686     auto haf = cluster_iterator->hitsAndFractions();
1687     auto layerId = rhtools_.getLayerWithOffset(haf[0].first);
1688     cluster_layer_id.push_back(layerId);
1689     uint32_t number_of_hits = cluster_iterator->hitsAndFractions().size();
1690     cluster_number_of_hits.push_back(number_of_hits);
1691     cluster_type.push_back(rhtools_.getCellType(lc_seed));
1692     cluster_timeErr.push_back(layerClustersTimes.get(c_id).second);
1693     cluster_time.push_back(layerClustersTimes.get(c_id).first);
1694     c_id += 1;
1695   }
1696 
1697   tracksters_in_candidate.resize(ticlcandidates.size());
1698   track_in_candidate.resize(ticlcandidates.size(), -1);
1699   nCandidates = ticlcandidates.size();
1700   for (int i = 0; i < static_cast<int>(ticlcandidates.size()); ++i) {
1701     const auto& candidate = ticlcandidates[i];
1702     candidate_charge.push_back(candidate.charge());
1703     candidate_pdgId.push_back(candidate.pdgId());
1704     candidate_energy.push_back(candidate.energy());
1705     candidate_px.push_back(candidate.px());
1706     candidate_py.push_back(candidate.py());
1707     candidate_pz.push_back(candidate.pz());
1708     candidate_time.push_back(candidate.time());
1709     candidate_time_err.push_back(candidate.timeError());
1710     std::vector<float> id_probs;
1711     for (int j = 0; j < 8; j++) {
1712       ticl::Trackster::ParticleType type = static_cast<ticl::Trackster::ParticleType>(j);
1713       id_probs.push_back(candidate.id_probability(type));
1714     }
1715     candidate_id_probabilities.push_back(id_probs);
1716 
1717     auto trackster_ptrs = candidate.tracksters();
1718     auto track_ptr = candidate.trackPtr();
1719     for (const auto& ts_ptr : trackster_ptrs) {
1720       auto ts_idx = ts_ptr.get() - (edm::Ptr<ticl::Trackster>(tracksters_handle, 0)).get();
1721       tracksters_in_candidate[i].push_back(ts_idx);
1722     }
1723     if (track_ptr.isNull())
1724       continue;
1725     int tk_idx = track_ptr.get() - (edm::Ptr<reco::Track>(tracks_h, 0)).get();
1726     track_in_candidate[i] = tk_idx;
1727   }
1728 
1729   nTrackstersMerged = trackstersmerged.size();
1730   for (auto trackster_iterator = trackstersmerged.begin(); trackster_iterator != trackstersmerged.end();
1731        ++trackster_iterator) {
1732     tracksters_merged_time.push_back(trackster_iterator->time());
1733     tracksters_merged_timeError.push_back(trackster_iterator->timeError());
1734     tracksters_merged_regressed_energy.push_back(trackster_iterator->regressed_energy());
1735     tracksters_merged_raw_energy.push_back(trackster_iterator->raw_energy());
1736     tracksters_merged_raw_em_energy.push_back(trackster_iterator->raw_em_energy());
1737     tracksters_merged_raw_pt.push_back(trackster_iterator->raw_pt());
1738     tracksters_merged_raw_em_pt.push_back(trackster_iterator->raw_em_pt());
1739     tracksters_merged_barycenter_x.push_back(trackster_iterator->barycenter().x());
1740     tracksters_merged_barycenter_y.push_back(trackster_iterator->barycenter().y());
1741     tracksters_merged_barycenter_z.push_back(trackster_iterator->barycenter().z());
1742     tracksters_merged_barycenter_eta.push_back(trackster_iterator->barycenter().eta());
1743     tracksters_merged_barycenter_phi.push_back(trackster_iterator->barycenter().phi());
1744     tracksters_merged_EV1.push_back(trackster_iterator->eigenvalues()[0]);
1745     tracksters_merged_EV2.push_back(trackster_iterator->eigenvalues()[1]);
1746     tracksters_merged_EV3.push_back(trackster_iterator->eigenvalues()[2]);
1747     tracksters_merged_eVector0_x.push_back((trackster_iterator->eigenvectors()[0]).x());
1748     tracksters_merged_eVector0_y.push_back((trackster_iterator->eigenvectors()[0]).y());
1749     tracksters_merged_eVector0_z.push_back((trackster_iterator->eigenvectors()[0]).z());
1750     tracksters_merged_sigmaPCA1.push_back(trackster_iterator->sigmasPCA()[0]);
1751     tracksters_merged_sigmaPCA2.push_back(trackster_iterator->sigmasPCA()[1]);
1752     tracksters_merged_sigmaPCA3.push_back(trackster_iterator->sigmasPCA()[2]);
1753 
1754     std::vector<float> id_probs;
1755     for (size_t i = 0; i < 8; i++)
1756       id_probs.push_back(trackster_iterator->id_probabilities(i));
1757     tracksters_merged_id_probabilities.push_back(id_probs);
1758 
1759     std::vector<uint32_t> vertices_indexes;
1760     std::vector<float> vertices_x;
1761     std::vector<float> vertices_y;
1762     std::vector<float> vertices_z;
1763     std::vector<float> vertices_time;
1764     std::vector<float> vertices_timeErr;
1765     std::vector<float> vertices_energy;
1766     std::vector<float> vertices_correctedEnergy;
1767     std::vector<float> vertices_correctedEnergyUncertainty;
1768     for (auto idx : trackster_iterator->vertices()) {
1769       vertices_indexes.push_back(idx);
1770       auto associated_cluster = (*layer_clusters_h)[idx];
1771       vertices_x.push_back(associated_cluster.x());
1772       vertices_y.push_back(associated_cluster.y());
1773       vertices_z.push_back(associated_cluster.z());
1774       vertices_energy.push_back(associated_cluster.energy());
1775       vertices_correctedEnergy.push_back(associated_cluster.correctedEnergy());
1776       vertices_correctedEnergyUncertainty.push_back(associated_cluster.correctedEnergyUncertainty());
1777       vertices_time.push_back(layerClustersTimes.get(idx).first);
1778       vertices_timeErr.push_back(layerClustersTimes.get(idx).second);
1779     }
1780     tracksters_merged_vertices_indexes.push_back(vertices_indexes);
1781     tracksters_merged_vertices_x.push_back(vertices_x);
1782     tracksters_merged_vertices_y.push_back(vertices_y);
1783     tracksters_merged_vertices_z.push_back(vertices_z);
1784     tracksters_merged_vertices_time.push_back(vertices_time);
1785     tracksters_merged_vertices_timeErr.push_back(vertices_timeErr);
1786     tracksters_merged_vertices_energy.push_back(vertices_energy);
1787     tracksters_merged_vertices_correctedEnergy.push_back(vertices_correctedEnergy);
1788     tracksters_merged_vertices_correctedEnergyUncertainty.push_back(vertices_correctedEnergyUncertainty);
1789   }
1790 
1791   // Tackster reco->sim associations
1792   trackstersCLUE3D_recoToSim_SC.resize(tracksters.size());
1793   trackstersCLUE3D_recoToSim_SC_score.resize(tracksters.size());
1794   trackstersCLUE3D_recoToSim_SC_sharedE.resize(tracksters.size());
1795   for (size_t i = 0; i < tracksters.size(); ++i) {
1796     const edm::Ref<ticl::TracksterCollection> tsRef(tracksters_handle, i);
1797 
1798     // CLUE3D -> STS-SC
1799     const auto stsSC_iter = tsRecoSimSCMap.find(tsRef);
1800     if (stsSC_iter != tsRecoSimSCMap.end()) {
1801       const auto& stsSCassociated = stsSC_iter->val;
1802       for (auto& sts : stsSCassociated) {
1803         auto sts_id = (sts.first).get() - (edm::Ref<ticl::TracksterCollection>(simTrackstersSC_h, 0)).get();
1804         trackstersCLUE3D_recoToSim_SC[i].push_back(sts_id);
1805         trackstersCLUE3D_recoToSim_SC_score[i].push_back(sts.second.second);
1806         trackstersCLUE3D_recoToSim_SC_sharedE[i].push_back(sts.second.first);
1807       }
1808     }
1809   }
1810 
1811   // SimTracksters
1812   nsimTrackstersSC = simTrackstersSC.size();
1813   trackstersCLUE3D_simToReco_SC.resize(nsimTrackstersSC);
1814   trackstersCLUE3D_simToReco_SC_score.resize(nsimTrackstersSC);
1815   trackstersCLUE3D_simToReco_SC_sharedE.resize(nsimTrackstersSC);
1816   for (size_t i = 0; i < nsimTrackstersSC; ++i) {
1817     const edm::Ref<ticl::TracksterCollection> stsSCRef(simTrackstersSC_h, i);
1818 
1819     // STS-SC -> CLUE3D
1820     const auto ts_iter = tsSimToRecoSCMap.find(stsSCRef);
1821     if (ts_iter != tsSimToRecoSCMap.end()) {
1822       const auto& tsAssociated = ts_iter->val;
1823       for (auto& ts : tsAssociated) {
1824         auto ts_idx = (ts.first).get() - (edm::Ref<ticl::TracksterCollection>(tracksters_handle, 0)).get();
1825         trackstersCLUE3D_simToReco_SC[i].push_back(ts_idx);
1826         trackstersCLUE3D_simToReco_SC_score[i].push_back(ts.second.second);
1827         trackstersCLUE3D_simToReco_SC_sharedE[i].push_back(ts.second.first);
1828       }
1829     }
1830   }
1831 
1832   // Tackster reco->sim associations
1833   trackstersCLUE3D_recoToSim_CP.resize(tracksters.size());
1834   trackstersCLUE3D_recoToSim_CP_score.resize(tracksters.size());
1835   trackstersCLUE3D_recoToSim_CP_sharedE.resize(tracksters.size());
1836   for (size_t i = 0; i < tracksters.size(); ++i) {
1837     const edm::Ref<ticl::TracksterCollection> tsRef(tracksters_handle, i);
1838 
1839     // CLUE3D -> STS-CP
1840     const auto stsCP_iter = tsRecoSimCPMap.find(tsRef);
1841     if (stsCP_iter != tsRecoSimCPMap.end()) {
1842       const auto& stsCPassociated = stsCP_iter->val;
1843       for (auto& sts : stsCPassociated) {
1844         auto sts_id = (sts.first).get() - (edm::Ref<ticl::TracksterCollection>(simTrackstersCP_h, 0)).get();
1845         trackstersCLUE3D_recoToSim_CP[i].push_back(sts_id);
1846         trackstersCLUE3D_recoToSim_CP_score[i].push_back(sts.second.second);
1847         trackstersCLUE3D_recoToSim_CP_sharedE[i].push_back(sts.second.first);
1848       }
1849     }
1850   }
1851 
1852   // SimTracksters
1853   nsimTrackstersCP = simTrackstersCP.size();
1854   trackstersCLUE3D_simToReco_CP.resize(nsimTrackstersCP);
1855   trackstersCLUE3D_simToReco_CP_score.resize(nsimTrackstersCP);
1856   trackstersCLUE3D_simToReco_CP_sharedE.resize(nsimTrackstersCP);
1857   for (size_t i = 0; i < nsimTrackstersCP; ++i) {
1858     const edm::Ref<ticl::TracksterCollection> stsCPRef(simTrackstersCP_h, i);
1859 
1860     // STS-CP -> CLUE3D
1861     const auto ts_iter = tsSimToRecoCPMap.find(stsCPRef);
1862     if (ts_iter != tsSimToRecoCPMap.end()) {
1863       const auto& tsAssociated = ts_iter->val;
1864       for (auto& ts : tsAssociated) {
1865         auto ts_idx = (ts.first).get() - (edm::Ref<ticl::TracksterCollection>(tracksters_handle, 0)).get();
1866         trackstersCLUE3D_simToReco_CP[i].push_back(ts_idx);
1867         trackstersCLUE3D_simToReco_CP_score[i].push_back(ts.second.second);
1868         trackstersCLUE3D_simToReco_CP_sharedE[i].push_back(ts.second.first);
1869       }
1870     }
1871   }
1872 
1873   // Tackster reco->sim associations
1874   MergeTracksters_recoToSim_SC.resize(trackstersmerged.size());
1875   MergeTracksters_recoToSim_SC_score.resize(trackstersmerged.size());
1876   MergeTracksters_recoToSim_SC_sharedE.resize(trackstersmerged.size());
1877   for (size_t i = 0; i < trackstersmerged.size(); ++i) {
1878     const edm::Ref<ticl::TracksterCollection> tsRef(tracksters_merged_h, i);
1879 
1880     // CLUE3D -> STS-SC
1881     const auto stsSC_iter = MergetsRecoSimSCMap.find(tsRef);
1882     if (stsSC_iter != MergetsRecoSimSCMap.end()) {
1883       const auto& stsSCassociated = stsSC_iter->val;
1884       for (auto& sts : stsSCassociated) {
1885         auto sts_id = (sts.first).get() - (edm::Ref<ticl::TracksterCollection>(simTrackstersSC_h, 0)).get();
1886         MergeTracksters_recoToSim_SC[i].push_back(sts_id);
1887         MergeTracksters_recoToSim_SC_score[i].push_back(sts.second.second);
1888         MergeTracksters_recoToSim_SC_sharedE[i].push_back(sts.second.first);
1889       }
1890     }
1891   }
1892 
1893   // SimTracksters
1894   nsimTrackstersSC = simTrackstersSC.size();
1895   MergeTracksters_simToReco_SC.resize(nsimTrackstersSC);
1896   MergeTracksters_simToReco_SC_score.resize(nsimTrackstersSC);
1897   MergeTracksters_simToReco_SC_sharedE.resize(nsimTrackstersSC);
1898   for (size_t i = 0; i < nsimTrackstersSC; ++i) {
1899     const edm::Ref<ticl::TracksterCollection> stsSCRef(simTrackstersSC_h, i);
1900 
1901     // STS-SC -> CLUE3D
1902     const auto ts_iter = MergetsSimToRecoSCMap.find(stsSCRef);
1903     if (ts_iter != MergetsSimToRecoSCMap.end()) {
1904       const auto& tsAssociated = ts_iter->val;
1905       for (auto& ts : tsAssociated) {
1906         auto ts_idx = (ts.first).get() - (edm::Ref<ticl::TracksterCollection>(tracksters_merged_h, 0)).get();
1907         MergeTracksters_simToReco_SC[i].push_back(ts_idx);
1908         MergeTracksters_simToReco_SC_score[i].push_back(ts.second.second);
1909         MergeTracksters_simToReco_SC_sharedE[i].push_back(ts.second.first);
1910       }
1911     }
1912   }
1913 
1914   // Tackster reco->sim associations
1915   MergeTracksters_recoToSim_CP.resize(trackstersmerged.size());
1916   MergeTracksters_recoToSim_CP_score.resize(trackstersmerged.size());
1917   MergeTracksters_recoToSim_CP_sharedE.resize(trackstersmerged.size());
1918   for (size_t i = 0; i < trackstersmerged.size(); ++i) {
1919     const edm::Ref<ticl::TracksterCollection> tsRef(tracksters_merged_h, i);
1920 
1921     // CLUE3D -> STS-CP
1922     const auto stsCP_iter = MergetsRecoSimCPMap.find(tsRef);
1923     if (stsCP_iter != MergetsRecoSimCPMap.end()) {
1924       const auto& stsCPassociated = stsCP_iter->val;
1925       for (auto& sts : stsCPassociated) {
1926         auto sts_id = (sts.first).get() - (edm::Ref<ticl::TracksterCollection>(simTrackstersCP_h, 0)).get();
1927         MergeTracksters_recoToSim_CP[i].push_back(sts_id);
1928         MergeTracksters_recoToSim_CP_score[i].push_back(sts.second.second);
1929         MergeTracksters_recoToSim_CP_sharedE[i].push_back(sts.second.first);
1930       }
1931     }
1932   }
1933 
1934   // Tackster reco->sim associations
1935   MergeTracksters_recoToSim_PU.resize(trackstersmerged.size());
1936   MergeTracksters_recoToSim_PU_score.resize(trackstersmerged.size());
1937   MergeTracksters_recoToSim_PU_sharedE.resize(trackstersmerged.size());
1938   for (size_t i = 0; i < trackstersmerged.size(); ++i) {
1939     const edm::Ref<ticl::TracksterCollection> tsRef(tracksters_merged_h, i);
1940 
1941     // CLUE3D -> STS-PU
1942     const auto stsPU_iter = MergetsRecoSimPUMap.find(tsRef);
1943     if (stsPU_iter != MergetsRecoSimPUMap.end()) {
1944       const auto& stsPUassociated = stsPU_iter->val;
1945       for (auto& sts : stsPUassociated) {
1946         auto sts_id = (sts.first).get() - (edm::Ref<ticl::TracksterCollection>(simTrackstersPU_h, 0)).get();
1947         MergeTracksters_recoToSim_PU[i].push_back(sts_id);
1948         MergeTracksters_recoToSim_PU_score[i].push_back(sts.second.second);
1949         MergeTracksters_recoToSim_PU_sharedE[i].push_back(sts.second.first);
1950       }
1951     }
1952   }
1953 
1954   // SimTracksters
1955   nsimTrackstersCP = simTrackstersCP.size();
1956   MergeTracksters_simToReco_CP.resize(nsimTrackstersCP);
1957   MergeTracksters_simToReco_CP_score.resize(nsimTrackstersCP);
1958   MergeTracksters_simToReco_CP_sharedE.resize(nsimTrackstersCP);
1959   for (size_t i = 0; i < nsimTrackstersCP; ++i) {
1960     const edm::Ref<ticl::TracksterCollection> stsCPRef(simTrackstersCP_h, i);
1961 
1962     // STS-CP -> TrackstersMerge
1963     const auto ts_iter = MergetsSimToRecoCPMap.find(stsCPRef);
1964     if (ts_iter != MergetsSimToRecoCPMap.end()) {
1965       const auto& tsAssociated = ts_iter->val;
1966       for (auto& ts : tsAssociated) {
1967         auto ts_idx = (ts.first).get() - (edm::Ref<ticl::TracksterCollection>(tracksters_merged_h, 0)).get();
1968         MergeTracksters_simToReco_CP[i].push_back(ts_idx);
1969         MergeTracksters_simToReco_CP_score[i].push_back(ts.second.second);
1970         MergeTracksters_simToReco_CP_sharedE[i].push_back(ts.second.first);
1971       }
1972     }
1973   }
1974 
1975   // SimTracksters
1976   auto nsimTrackstersPU = simTrackstersPU.size();
1977   MergeTracksters_simToReco_PU.resize(nsimTrackstersPU);
1978   MergeTracksters_simToReco_PU_score.resize(nsimTrackstersPU);
1979   MergeTracksters_simToReco_PU_sharedE.resize(nsimTrackstersPU);
1980   for (size_t i = 0; i < nsimTrackstersPU; ++i) {
1981     const edm::Ref<ticl::TracksterCollection> stsPURef(simTrackstersPU_h, i);
1982 
1983     // STS-PU -> Tracksters Merge
1984     const auto ts_iter = MergetsSimToRecoPUMap.find(stsPURef);
1985     if (ts_iter != MergetsSimToRecoPUMap.end()) {
1986       const auto& tsAssociated = ts_iter->val;
1987       for (auto& ts : tsAssociated) {
1988         auto ts_idx = (ts.first).get() - (edm::Ref<ticl::TracksterCollection>(tracksters_merged_h, 0)).get();
1989         MergeTracksters_simToReco_PU[i].push_back(ts_idx);
1990         MergeTracksters_simToReco_PU_score[i].push_back(ts.second.second);
1991         MergeTracksters_simToReco_PU_sharedE[i].push_back(ts.second.first);
1992       }
1993     }
1994   }
1995 
1996   //Tracks
1997   for (size_t i = 0; i < tracks.size(); i++) {
1998     const auto& track = tracks[i];
1999     reco::TrackRef trackref = reco::TrackRef(tracks_h, i);
2000     int iSide = int(track.eta() > 0);
2001     const auto& fts = trajectoryStateTransform::outerFreeState((track), bFieldProd);
2002     // to the HGCal front
2003     const auto& tsos = prop.propagate(fts, firstDisk_[iSide]->surface());
2004     if (tsos.isValid()) {
2005       const auto& globalPos = tsos.globalPosition();
2006       const auto& globalMom = tsos.globalMomentum();
2007       track_id.push_back(i);
2008       track_hgcal_x.push_back(globalPos.x());
2009       track_hgcal_y.push_back(globalPos.y());
2010       track_hgcal_z.push_back(globalPos.z());
2011       track_hgcal_eta.push_back(globalPos.eta());
2012       track_hgcal_phi.push_back(globalPos.phi());
2013       track_hgcal_px.push_back(globalMom.x());
2014       track_hgcal_py.push_back(globalMom.y());
2015       track_hgcal_pz.push_back(globalMom.z());
2016       track_hgcal_pt.push_back(globalMom.perp());
2017       track_pt.push_back(track.pt());
2018       track_quality.push_back(track.quality(reco::TrackBase::highPurity));
2019       track_missing_outer_hits.push_back(track.missingOuterHits());
2020       track_charge.push_back(track.charge());
2021       track_time.push_back(trackTime[trackref]);
2022       track_time_quality.push_back(trackTimeQual[trackref]);
2023       track_time_err.push_back(trackTimeErr[trackref]);
2024       track_nhits.push_back(tracks[i].recHitsSize());
2025     }
2026   }
2027 
2028   if (saveCLUE3DTracksters_)
2029     trackster_tree_->Fill();
2030   if (saveLCs_)
2031     cluster_tree_->Fill();
2032   if (saveTICLCandidate_)
2033     candidate_tree_->Fill();
2034   if (saveTrackstersMerged_)
2035     tracksters_merged_tree_->Fill();
2036   if (saveAssociations_)
2037     associations_tree_->Fill();
2038   if (saveSimTrackstersSC_)
2039     simtrackstersSC_tree_->Fill();
2040   if (saveSimTrackstersCP_)
2041     simtrackstersCP_tree_->Fill();
2042   if (saveTracks_)
2043     tracks_tree_->Fill();
2044   if (saveSimTICLCandidate_)
2045     simTICLCandidate_tree->Fill();
2046 }
2047 
2048 void TICLDumper::endJob() {}
2049 
2050 void TICLDumper::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
2051   edm::ParameterSetDescription desc;
2052   desc.add<edm::InputTag>("trackstersclue3d", edm::InputTag("ticlTrackstersCLUE3DHigh"));
2053   desc.add<edm::InputTag>("layerClusters", edm::InputTag("hgcalMergeLayerClusters"));
2054   desc.add<edm::InputTag>("layer_clustersTime", edm::InputTag("hgcalMergeLayerClusters", "timeLayerCluster"));
2055   desc.add<edm::InputTag>("ticlcandidates", edm::InputTag("ticlTrackstersMerge"));
2056   desc.add<edm::InputTag>("tracks", edm::InputTag("generalTracks"));
2057   desc.add<edm::InputTag>("tracksTime", edm::InputTag("tofPID:t0"));
2058   desc.add<edm::InputTag>("tracksTimeQual", edm::InputTag("mtdTrackQualityMVA:mtdQualMVA"));
2059   desc.add<edm::InputTag>("tracksTimeErr", edm::InputTag("tofPID:sigmat0"));
2060   desc.add<edm::InputTag>("trackstersmerged", edm::InputTag("ticlTrackstersMerge"));
2061   desc.add<edm::InputTag>("simtrackstersSC", edm::InputTag("ticlSimTracksters"));
2062   desc.add<edm::InputTag>("simtrackstersCP", edm::InputTag("ticlSimTracksters", "fromCPs"));
2063   desc.add<edm::InputTag>("simtrackstersPU", edm::InputTag("ticlSimTracksters", "PU"));
2064   desc.add<edm::InputTag>("simTICLCandidates", edm::InputTag("ticlSimTracksters"));
2065   desc.add<edm::InputTag>("recoToSimAssociatorSC",
2066                           edm::InputTag("tracksterSimTracksterAssociationPRbyCLUE3D", "recoToSim"));
2067   desc.add<edm::InputTag>("simToRecoAssociatorSC",
2068                           edm::InputTag("tracksterSimTracksterAssociationPRbyCLUE3D", "simToReco"));
2069   desc.add<edm::InputTag>("recoToSimAssociatorCP",
2070                           edm::InputTag("tracksterSimTracksterAssociationLinkingbyCLUE3D", "recoToSim"));
2071   desc.add<edm::InputTag>("simToRecoAssociatorCP",
2072                           edm::InputTag("tracksterSimTracksterAssociationLinkingbyCLUE3D", "simToReco"));
2073   desc.add<edm::InputTag>("MergerecoToSimAssociatorSC",
2074                           edm::InputTag("tracksterSimTracksterAssociationPR", "recoToSim"));
2075   desc.add<edm::InputTag>("MergesimToRecoAssociatorSC",
2076                           edm::InputTag("tracksterSimTracksterAssociationPR", "simToReco"));
2077   desc.add<edm::InputTag>("MergerecoToSimAssociatorCP",
2078                           edm::InputTag("tracksterSimTracksterAssociationLinking", "recoToSim"));
2079   desc.add<edm::InputTag>("MergesimToRecoAssociatorCP",
2080                           edm::InputTag("tracksterSimTracksterAssociationLinking", "simToReco"));
2081   desc.add<edm::InputTag>("MergerecoToSimAssociatorPU",
2082                           edm::InputTag("tracksterSimTracksterAssociationLinkingPU", "recoToSim"));
2083   desc.add<edm::InputTag>("MergesimToRecoAssociatorPU",
2084                           edm::InputTag("tracksterSimTracksterAssociationLinkingPU", "simToReco"));
2085   desc.add<edm::InputTag>("simclusters", edm::InputTag("mix", "MergedCaloTruth"));
2086   desc.add<edm::InputTag>("caloparticles", edm::InputTag("mix", "MergedCaloTruth"));
2087   desc.add<std::string>("detector", "HGCAL");
2088   desc.add<std::string>("propagator", "PropagatorWithMaterial");
2089 
2090   desc.add<bool>("saveLCs", true);
2091   desc.add<bool>("saveCLUE3DTracksters", true);
2092   desc.add<bool>("saveTrackstersMerged", true);
2093   desc.add<bool>("saveSimTrackstersSC", true);
2094   desc.add<bool>("saveSimTrackstersCP", true);
2095   desc.add<bool>("saveTICLCandidate", true);
2096   desc.add<bool>("saveSimTICLCandidate", true);
2097   desc.add<bool>("saveTracks", true);
2098   desc.add<bool>("saveAssociations", true);
2099   descriptions.add("ticlDumper", desc);
2100 }
2101 
2102 DEFINE_FWK_MODULE(TICLDumper);