diff --git a/DeepNtuplizer/BuildFile.xml b/DeepNtuplizer/BuildFile.xml index 3759283a42b..67d94077ab3 100644 --- a/DeepNtuplizer/BuildFile.xml +++ b/DeepNtuplizer/BuildFile.xml @@ -4,7 +4,7 @@ - + diff --git a/DeepNtuplizer/bin/mergeSamples.cc b/DeepNtuplizer/bin/mergeSamples.cc index 574d07a05ad..10ddaf96447 100644 --- a/DeepNtuplizer/bin/mergeSamples.cc +++ b/DeepNtuplizer/bin/mergeSamples.cc @@ -40,6 +40,8 @@ #include "DeepNTuples/DeepNtuplizer/interface/ntuple_JetInfo.h" #include "DeepNTuples/DeepNtuplizer/interface/ntuple_pfCands.h" #include "DeepNTuples/DeepNtuplizer/interface/ntuple_SV.h" +#include "DeepNTuples/DeepNtuplizer/interface/ntuple_LT.h" +#include "DeepNTuples/DeepNtuplizer/interface/ntuple_pairwise.h" #include #include @@ -97,9 +99,10 @@ int main(int argc, char *argv[]){ std::vector branchinfos; branchinfos.push_back(new ntuple_JetInfo()); branchinfos.push_back(new ntuple_SV()); + branchinfos.push_back(new ntuple_LT()); branchinfos.push_back(new ntuple_bTagVars()); branchinfos.push_back(new ntuple_pfCands()); - + branchinfos.push_back(new ntuple_pairwise()); //simple opt parsing diff --git a/DeepNtuplizer/bin/mergeSamples_parallel.cc b/DeepNtuplizer/bin/mergeSamples_parallel.cc index d8d400b005a..8c01cd0b524 100644 --- a/DeepNtuplizer/bin/mergeSamples_parallel.cc +++ b/DeepNtuplizer/bin/mergeSamples_parallel.cc @@ -40,6 +40,8 @@ #include "DeepNTuples/DeepNtuplizer/interface/ntuple_JetInfo.h" #include "DeepNTuples/DeepNtuplizer/interface/ntuple_pfCands.h" #include "DeepNTuples/DeepNtuplizer/interface/ntuple_SV.h" +#include "DeepNTuples/DeepNtuplizer/interface/ntuple_LT.h" +#include "DeepNTuples/DeepNtuplizer/interface/ntuple_pairwise.h" #include "DeepNTuples/DeepNtuplizer/interface/ntuple_FatJetInfo.h" #include @@ -121,8 +123,10 @@ std::vector createChains(const std::vector >& inf branchinfos.push_back(new ntuple_JetInfo()); branchinfos.push_back(new ntuple_SV()); + branchinfos.push_back(new ntuple_LT()); branchinfos.push_back(new ntuple_bTagVars()); branchinfos.push_back(new ntuple_pfCands()); + branchinfos.push_back(new ntuple_pairwise()); branchinfos.push_back(new ntuple_FatJetInfo()); std::vector chains; diff --git a/DeepNtuplizer/bin/printBranches.cc b/DeepNtuplizer/bin/printBranches.cc index 1b5a72072a2..c3467c794f4 100644 --- a/DeepNtuplizer/bin/printBranches.cc +++ b/DeepNtuplizer/bin/printBranches.cc @@ -2,10 +2,12 @@ #include "../interface/ntuple_JetInfo.h" #include "../interface/ntuple_SV.h" +#include "../interface/ntuple_LT.h" #include "../interface/ntuple_bTagVars.h" #include "../interface/ntuple_pfCands.h" -#include "../interface/ntuple_FatJetInfo.h" -#include "../interface/ntuple_DeepVertex.h" +//#include "../interface/ntuple_FatJetInfo.h" +//#include "../interface/ntuple_DeepVertex.h" +#include "../interface/ntuple_pairwise.h" #include "TFile.h" #include #include "TH1F.h" @@ -21,10 +23,12 @@ int main(int argc, char *argv[]) { std::vector branchinfos; branchinfos.push_back(new ntuple_JetInfo()); branchinfos.push_back(new ntuple_SV()); + branchinfos.push_back(new ntuple_LT()); branchinfos.push_back(new ntuple_bTagVars()); branchinfos.push_back(new ntuple_pfCands()); - // branchinfos.push_back(new ntuple_FatJetInfo()); + branchinfos.push_back(new ntuple_pairwise()); // branchinfos.push_back(new ntuple_DeepVertex()); + // branchinfos.push_back(new ntuple_GraphB()); if (argc < 3) return -1; diff --git a/DeepNtuplizer/interface/TrackPairInfoBuilder.h b/DeepNtuplizer/interface/TrackPairInfoBuilder.h new file mode 100644 index 00000000000..0094793b665 --- /dev/null +++ b/DeepNtuplizer/interface/TrackPairInfoBuilder.h @@ -0,0 +1,102 @@ +#ifndef BTAGHELPERS_INTERFACE_TRACKPARINFOBUILDER_H_ +#define BTAGHELPERS_INTERFACE_TRACKPAIRINFOBUILDER_H_ + +#include "DataFormats/GeometrySurface/interface/Line.h" +#include "TrackingTools/Records/interface/TransientTrackRecord.h" +#include "TrackingTools/IPTools/interface/IPTools.h" +#include "TrackingTools/PatternTools/interface/TwoTrackMinimumDistance.h" +#include "DataFormats/PatCandidates/interface/Jet.h" + +namespace deepntuples { + + class TrackPairInfoBuilder { + public: + TrackPairInfoBuilder(); + + void buildTrackPairInfo(const reco::TransientTrack it, + const reco::TransientTrack tt, + const reco::Vertex& pv, + const pat::Jet &jet + ); + + const float track_i_pt() const { return track_i_pt_; } + const float track_t_pt() const { return track_t_pt_; } + const float track_eta() const { return track_eta_; } + const float track_phi() const { return track_phi_; } + const float track_dz() const { return track_dz_; } + const float track_dxy() const { return track_dxy_; } + const float pca_distance() const { return pca_distance_; } + const float pca_significance() const { return pca_significance_; } + const float pcaSeed_x() const { return pcaSeed_x_; } + const float pcaSeed_y() const { return pcaSeed_y_; } + const float pcaSeed_z() const { return pcaSeed_z_; } + const float pcaSeed_xerr() const { return pcaSeed_xerr_; } + const float pcaSeed_yerr() const { return pcaSeed_yerr_; } + const float pcaSeed_zerr() const { return pcaSeed_zerr_; } + const float pcaTrack_x() const { return pcaTrack_x_; } + const float pcaTrack_y() const { return pcaTrack_y_; } + const float pcaTrack_z() const { return pcaTrack_z_; } + const float pcaTrack_xerr() const { return pcaTrack_xerr_; } + const float pcaTrack_yerr() const { return pcaTrack_yerr_; } + const float pcaTrack_zerr() const { return pcaTrack_zerr_; } + const float dotprodTrack() const { return dotprodTrack_; } + const float dotprodSeed() const { return dotprodSeed_; } + const float pcaSeed_dist() const { return pcaSeed_dist_; } + const float pcaTrack_dist() const { return pcaTrack_dist_; } + const float track_candMass() const { return track_candMass_; } + const float track_ip2d() const { return track_ip2d_; } + const float track_ip2dSig() const { return track_ip2dSig_; } + const float track_ip3d() const { return track_ip3d_; } + const float track_ip3dSig() const { return track_ip3dSig_; } + const float dotprodTrackSeed2D() const { return dotprodTrackSeed2D_; } + const float dotprodTrackSeed2DV() const { return dotprodTrackSeed2DV_; } + const float dotprodTrackSeed3D() const { return dotprodTrackSeed3D_; } + const float dotprodTrackSeed3DV() const { return dotprodTrackSeed3DV_; } + const float pca_jetAxis_dist() const { return pca_jetAxis_dist_; } + const float pca_jetAxis_dotprod() const { return pca_jetAxis_dotprod_; } + const float pca_jetAxis_dEta() const { return pca_jetAxis_dEta_; } + const float pca_jetAxis_dPhi() const { return pca_jetAxis_dPhi_; } + + private: + float track_i_pt_; + float track_t_pt_; + float track_eta_; + float track_phi_; + float track_dz_; + float track_dxy_; + float pca_distance_; + float pca_significance_; + float pcaSeed_x_; + float pcaSeed_y_; + float pcaSeed_z_; + float pcaSeed_xerr_; + float pcaSeed_yerr_; + float pcaSeed_zerr_; + float pcaTrack_x_; + float pcaTrack_y_; + float pcaTrack_z_; + float pcaTrack_xerr_; + float pcaTrack_yerr_; + float pcaTrack_zerr_; + float dotprodTrack_; + float dotprodSeed_; + float pcaSeed_dist_; + float pcaTrack_dist_; + float track_candMass_; + float track_ip2d_; + float track_ip2dSig_; + float track_ip3d_; + float track_ip3dSig_; + float dotprodTrackSeed2D_; + float dotprodTrackSeed2DV_; + float dotprodTrackSeed3D_; + float dotprodTrackSeed3DV_; + float pca_jetAxis_dist_; + float pca_jetAxis_dotprod_; + float pca_jetAxis_dEta_; + float pca_jetAxis_dPhi_; + }; + +} // namespace btagbtvdeep + +#endif //RecoBTag_FeatureTools_TrackPairInfoBuilder_h diff --git a/DeepNtuplizer/interface/ntuple_DeepVertex.h b/DeepNtuplizer/interface/ntuple_DeepVertex.h deleted file mode 100644 index 2726b671f38..00000000000 --- a/DeepNtuplizer/interface/ntuple_DeepVertex.h +++ /dev/null @@ -1,148 +0,0 @@ -/* - * ntuple_DeepVertex.h - * - * Created on: 23 June 2017 - * Author: Seth Moortgat - */ - -#ifndef DEEPNTUPLES_DEEPNTUPLIZER_INTERFACE_NTUPLE_DEEPVERTEX_H_ -#define DEEPNTUPLES_DEEPNTUPLIZER_INTERFACE_NTUPLE_DEEPVERTEX_H_ - -#include "ntuple_content.h" -#include "neighbourTrackVars.h" -#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h" -#include "RecoBTag/TrackProbability/interface/HistogramProbabilityEstimator.h" - -class HistogramProbabilityEstimator; -#include - -class ntuple_DeepVertex: public ntuple_content{ -public: - - ntuple_DeepVertex(double jetR = 0.4); - ~ntuple_DeepVertex(); - - void getInput(const edm::ParameterSet& iConfig); - void initBranches(TTree* ); - void readEvent(const edm::Event& iEvent); - void readSetup(const edm::EventSetup& iSetup); - void checkEventSetup(const edm::EventSetup & iSetup); - - //use either of these functions - - bool fillBranches(const pat::Jet &, const size_t& jetidx, const edm::View * coll=0); - - - void setCandidatesToken(const edm::EDGetTokenT > & t){ - CandidateToken=t; - } - -private: - - // seed candidates - static constexpr size_t max_seeds=10; - - unsigned int n_seeds=0; - float nSeeds=0; - - float seed_pt[max_seeds]; - float seed_eta[max_seeds]; - float seed_phi[max_seeds]; - float seed_mass[max_seeds]; - - float seed_dz[max_seeds]; - float seed_dxy[max_seeds]; - float seed_3D_ip[max_seeds]; - float seed_3D_sip[max_seeds]; - float seed_2D_ip[max_seeds]; - float seed_2D_sip[max_seeds]; - float seed_3D_signedIp[max_seeds]; - float seed_3D_signedSip[max_seeds]; - float seed_2D_signedIp[max_seeds]; - float seed_2D_signedSip[max_seeds]; - float seed_3D_TrackProbability[max_seeds]; - float seed_2D_TrackProbability[max_seeds]; - - float seed_chi2reduced[max_seeds]; - float seed_nPixelHits[max_seeds]; - float seed_nHits[max_seeds]; - float seed_jetAxisDistance[max_seeds]; - float seed_jetAxisDlength[max_seeds]; - - unsigned int seed_n_NearTracks[max_seeds]; - - //nearest track candidates - static constexpr size_t max_nearestTrk=200; // 20 per seed - - unsigned int n_NearTracksTotal=0; - - float nearTracks_pt[max_nearestTrk]; - float nearTracks_eta[max_nearestTrk]; - float nearTracks_phi[max_nearestTrk]; - float nearTracks_mass[max_nearestTrk]; - float nearTracks_dz[max_nearestTrk]; - float nearTracks_dxy[max_nearestTrk]; - float nearTracks_3D_ip[max_nearestTrk]; - float nearTracks_3D_sip[max_nearestTrk]; - float nearTracks_2D_ip[max_nearestTrk]; - float nearTracks_2D_sip[max_nearestTrk]; - float nearTracks_PCAdist[max_nearestTrk]; - float nearTracks_PCAdsig[max_nearestTrk]; - float nearTracks_PCAonSeed_x[max_nearestTrk]; - float nearTracks_PCAonSeed_y[max_nearestTrk]; - float nearTracks_PCAonSeed_z[max_nearestTrk]; - float nearTracks_PCAonSeed_xerr[max_nearestTrk]; - float nearTracks_PCAonSeed_yerr[max_nearestTrk]; - float nearTracks_PCAonSeed_zerr[max_nearestTrk]; - float nearTracks_PCAonTrack_x[max_nearestTrk]; - float nearTracks_PCAonTrack_y[max_nearestTrk]; - float nearTracks_PCAonTrack_z[max_nearestTrk]; - float nearTracks_PCAonTrack_xerr[max_nearestTrk]; - float nearTracks_PCAonTrack_yerr[max_nearestTrk]; - float nearTracks_PCAonTrack_zerr[max_nearestTrk]; - float nearTracks_dotprodTrack[max_nearestTrk]; - float nearTracks_dotprodSeed[max_nearestTrk]; - float nearTracks_dotprodTrackSeed2D[max_nearestTrk]; - float nearTracks_dotprodTrackSeed3D[max_nearestTrk]; - float nearTracks_dotprodTrackSeedVectors2D[max_nearestTrk]; - float nearTracks_dotprodTrackSeedVectors3D[max_nearestTrk]; - float nearTracks_PCAonSeed_pvd[max_nearestTrk]; - float nearTracks_PCAonTrack_pvd[max_nearestTrk]; - float nearTracks_PCAjetAxis_dist[max_nearestTrk]; - float nearTracks_PCAjetMomenta_dotprod[max_nearestTrk]; - float nearTracks_PCAjetDirs_DEta[max_nearestTrk]; - float nearTracks_PCAjetDirs_DPhi[max_nearestTrk]; - - - // IVF cut parameters (HARDCODED?? OR CONFIGURABLE IN PYTHON CONFIG) - float min3DIPValue=0.005; - float min3DIPSignificance=1.2; - int max3DIPValue=9999.; - int max3DIPSignificance=9999.; - - - //tokens to be defined from main analyzer - edm::EDGetTokenT > CandidateToken; - - //helper: - edm::Handle > tracks; - - // builder - edm::ESHandle builder; - - std::auto_ptr m_probabilityEstimator; - bool m_computeProbabilities=1; - unsigned long long m_calibrationCacheId2D; - unsigned long long m_calibrationCacheId3D; - - // temporary containers - neighbourTrackVars myTrack; - std::vector nearTracks; - std::multimap > > SortedSeedsMap; - - -}; - - - -#endif /* DEEPNTUPLES_DEEPNTUPLIZER_INTERFACE_NTUPLE_DEEPVERTEX_H_ */ diff --git a/DeepNtuplizer/interface/ntuple_GraphB.h b/DeepNtuplizer/interface/ntuple_GraphB.h deleted file mode 100644 index 36a916b662d..00000000000 --- a/DeepNtuplizer/interface/ntuple_GraphB.h +++ /dev/null @@ -1,96 +0,0 @@ -/* - * ntuple_GraphB.h - * - * Created on: 23 June 2017 - * Author: Seth Moortgat - */ - -#ifndef DEEPNTUPLES_DEEPNTUPLIZER_INTERFACE_NTUPLE_GRAPHB_H_ -#define DEEPNTUPLES_DEEPNTUPLIZER_INTERFACE_NTUPLE_GRAPHB_H_ - -#include "ntuple_content.h" -#include "neighbourTrackVars.h" -#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h" -#include "RecoBTag/TrackProbability/interface/HistogramProbabilityEstimator.h" - -class HistogramProbabilityEstimator; -#include - -class ntuple_GraphB: public ntuple_content{ -public: - - ntuple_GraphB(double jetR = 0.4); - ~ntuple_GraphB(); - - void getInput(const edm::ParameterSet& iConfig); - void initBranches(TTree* ); - void readEvent(const edm::Event& iEvent); - void readSetup(const edm::EventSetup& iSetup); - void checkEventSetup(const edm::EventSetup & iSetup); - - //use either of these functions - - bool fillBranches(const pat::Jet &, const size_t& jetidx, const edm::View * coll=0); - - - void setCandidatesToken(const edm::EDGetTokenT > & t){ - CandidateToken=t; - } - -private: - - // seed candidates - static constexpr size_t max_gtracks=100; - - unsigned int n_gtracks=0; - float nGtracks=0; - - float gtrack_pt[max_gtracks]; - float gtrack_eta[max_gtracks]; - float gtrack_phi[max_gtracks]; - float gtrack_mass[max_gtracks]; - float gtrack_dR[max_gtracks]; - float gtrack_dz[max_gtracks]; - float gtrack_dxy[max_gtracks]; - float gtrack_3D_ip[max_gtracks]; - float gtrack_3D_sip[max_gtracks]; - float gtrack_2D_ip[max_gtracks]; - float gtrack_2D_sip[max_gtracks]; - float gtrack_3D_TrackProbability[max_gtracks]; - float gtrack_2D_TrackProbability[max_gtracks]; - float gtrack_chi2reduced[max_gtracks]; - float gtrack_nPixelHits[max_gtracks]; - float gtrack_nHits[max_gtracks]; - float gtrack_jetAxisDistance[max_gtracks]; - float gtrack_jetAxisDlength[max_gtracks]; - float gtrack_dotProdTrack[max_gtracks]; - float gtrack_dotProdTrack2D[max_gtracks]; - float gtrack_PCAtrackFromPV[max_gtracks]; - - - // IVF cut parameters (HARDCODED?? OR CONFIGURABLE IN PYTHON CONFIG) - float min3DIPValue=0.005; - float min3DIPSignificance=1.2; - int max3DIPValue=9999.; - int max3DIPSignificance=9999.; - - - //tokens to be defined from main analyzer - edm::EDGetTokenT > CandidateToken; - - //helper: - edm::Handle > tracks; - - // builder - edm::ESHandle builder; - - std::auto_ptr m_probabilityEstimator; - bool m_computeProbabilities=1; - unsigned long long m_calibrationCacheId2D; - unsigned long long m_calibrationCacheId3D; - -}; - - - -#endif /* DEEPNTUPLES_DEEPNTUPLIZER_INTERFACE_NTUPLE_GRAPHB_H_ */ diff --git a/DeepNtuplizer/interface/ntuple_LT.h b/DeepNtuplizer/interface/ntuple_LT.h new file mode 100644 index 00000000000..958564f74b5 --- /dev/null +++ b/DeepNtuplizer/interface/ntuple_LT.h @@ -0,0 +1,105 @@ +/* + * ntuple_LT.h + * + * Created on: 28 Mai 2023 + * Author: Alexandre De Moor + */ + +#ifndef DEEPNTUPLES_DEEPNTUPLIZER_INTERFACE_NTUPLE_LT_H_ +#define DEEPNTUPLES_DEEPNTUPLIZER_INTERFACE_NTUPLE_LT_H_ + +#include "ntuple_content.h" +#include "DataFormats/PatCandidates/interface/PackedCandidate.h" +#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h" + +class ntuple_LT: public ntuple_content{ + public: + + ntuple_LT():ntuple_content(),jetradius_(0.4), + n_Cpfcand_(0),n_Npfcand_(0){} + + void setJetRadius(const float& radius){jetradius_=radius;} + void getInput(const edm::ParameterSet& iConfig); + void initBranches(TTree* ); + virtual void readConfig(const edm::ParameterSet& iConfig, edm::ConsumesCollector && cc) override; + void readEvent(const edm::Event& iEvent); + void readSetup(const edm::EventSetup& iSetup); + + + //use either of these functions + + bool fillBranches(const pat::Jet &, const size_t& jetidx, const edm::View * coll=0); + + private: + + float jetradius_; + float min_candidate_pt_ = -1; + + edm::ESHandle builder; + + edm::EDGetTokenT> ltToken_; + edm::Handle> LTs; + + int n_LT_; + int n_Cpfcand_; + int n_Npfcand_; + + static constexpr size_t max_pfcand_=50; + + float LT_pt_[max_pfcand_]; + float LT_eta_[max_pfcand_]; + float LT_phi_[max_pfcand_]; + float LT_e_[max_pfcand_]; + + float LT_puppiw_[max_pfcand_]; + float LT_VTX_ass_[max_pfcand_]; + float LT_dz_[max_pfcand_]; + + float LT_BtagPf_trackEtaRel_[max_pfcand_]; + float LT_BtagPf_trackPtRel_[max_pfcand_]; + float LT_BtagPf_trackPPar_[max_pfcand_]; + float LT_BtagPf_trackDeltaR_[max_pfcand_]; + float LT_BtagPf_trackPParRatio_[max_pfcand_]; + float LT_BtagPf_trackSip3dVal_[max_pfcand_]; + float LT_BtagPf_trackSip3dSig_[max_pfcand_]; + float LT_BtagPf_trackSip2dVal_[max_pfcand_]; + float LT_BtagPf_trackSip2dSig_[max_pfcand_]; + float LT_BtagPf_trackDecayLen_[max_pfcand_]; + float LT_BtagPf_trackJetDistVal_[max_pfcand_]; + + float LT_charge_[max_pfcand_]; + float LT_chi2_[max_pfcand_]; + float LT_quality_[max_pfcand_]; + float LT_drminsv_[max_pfcand_]; + float LT_distminsvold_[max_pfcand_]; + float LT_distminsv_[max_pfcand_]; + float LT_distminsv2_[max_pfcand_]; + + float LT_nhitpixelBarrelLayer1_[max_pfcand_]; + float LT_nhitpixelBarrelLayer2_[max_pfcand_]; + float LT_nhitpixelBarrelLayer3_[max_pfcand_]; + float LT_nhitpixelBarrelLayer4_[max_pfcand_]; + float LT_nhitpixelEndcapLayer1_[max_pfcand_]; + float LT_nhitpixelEndcapLayer2_[max_pfcand_]; + float LT_numberOfValidHits_[max_pfcand_]; + float LT_numberOfValidPixelHits_[max_pfcand_]; + float LT_numberOfValidStripHits_[max_pfcand_]; + float LT_numberOfValidStripTIBHits_[max_pfcand_]; + float LT_numberOfValidStripTIDHits_[max_pfcand_]; + float LT_numberOfValidStripTOBHits_[max_pfcand_]; + float LT_numberOfValidStripTECHits_[max_pfcand_]; + + float LT_pdgID_[max_pfcand_]; + float LT_HadFrac_[max_pfcand_]; + float LT_CaloFrac_[max_pfcand_]; + + float mindrsvpfcand(const pat::PackedCandidate* pfcand); + float mindistsvpfcandold(const reco::TransientTrack track); + float mindistsvpfcand(const reco::TransientTrack track); + GlobalPoint mingpsvpfcand(const reco::TransientTrack track); + GlobalPoint gppvpfcand(const reco::TransientTrack track, GlobalVector direction, const reco::Vertex vertex); + +}; + + +#endif /* DEEPNTUPLES_DEEPNTUPLIZER_INTERFACE_NTUPLE_LT_H_ */ diff --git a/DeepNtuplizer/interface/ntuple_SV.h b/DeepNtuplizer/interface/ntuple_SV.h index 7c3ca0a00d2..28df27a3c31 100644 --- a/DeepNtuplizer/interface/ntuple_SV.h +++ b/DeepNtuplizer/interface/ntuple_SV.h @@ -9,6 +9,8 @@ #define DEEPNTUPLES_DEEPNTUPLIZER_INTERFACE_NTUPLE_SV_H_ #include "ntuple_content.h" +#include "DataFormats/PatCandidates/interface/PackedCandidate.h" +#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h" class ntuple_SV: public ntuple_content{ public: @@ -19,6 +21,7 @@ class ntuple_SV: public ntuple_content{ void getInput(const edm::ParameterSet& iConfig); void initBranches(TTree* ); void readEvent(const edm::Event& iEvent); + void readSetup(const edm::EventSetup& iSetup); //use either of these functions @@ -33,11 +36,14 @@ class ntuple_SV: public ntuple_content{ float nsv_; std::string prefix_; + edm::ESHandle builder; + static constexpr size_t max_sv=10; float sv_pt_[max_sv]; float sv_eta_[max_sv]; float sv_phi_[max_sv]; + float sv_e_[max_sv]; float sv_etarel_[max_sv]; float sv_phirel_[max_sv]; float sv_deltaR_[max_sv]; @@ -57,7 +63,15 @@ class ntuple_SV: public ntuple_content{ float sv_costhetasvpv_[max_sv]; float sv_enratio_[max_sv]; - + float sv_hcal_frac_[max_sv]; + float sv_calo_frac_[max_sv]; + float sv_dz_[max_sv]; + float sv_pfd2dval_[max_sv]; + float sv_pfd2dsig_[max_sv]; + float sv_pfd3dval_[max_sv]; + float sv_pfd3dsig_[max_sv]; + float sv_puppiw_[max_sv]; + float sv_charge_sum_[max_sv]; static const reco::Vertex * spvp_; diff --git a/DeepNtuplizer/interface/ntuple_content.h b/DeepNtuplizer/interface/ntuple_content.h index b8c864ad7f8..3979c607c3c 100644 --- a/DeepNtuplizer/interface/ntuple_content.h +++ b/DeepNtuplizer/interface/ntuple_content.h @@ -19,6 +19,8 @@ #include "DataFormats/VertexReco/interface/Vertex.h" #include "SimDataFormats/PileupSummaryInfo/interface/PileupSummaryInfo.h" #include "DataFormats/Candidate/interface/VertexCompositePtrCandidate.h" +#include "FWCore/Framework/interface/ConsumesCollector.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" #include #include #include @@ -34,6 +36,7 @@ class ntuple_content{ virtual void getInput(const edm::ParameterSet& iConfig){} virtual void initBranches(TTree* )=0; + virtual void readConfig(const edm::ParameterSet& iConfig, edm::ConsumesCollector && cc) {} virtual void readEvent(const edm::Event& iEvent)=0; virtual void readSetup(const edm::EventSetup& iSetup){} //use either of these functions diff --git a/DeepNtuplizer/interface/ntuple_pairwise.h b/DeepNtuplizer/interface/ntuple_pairwise.h new file mode 100644 index 00000000000..522c5d78ce7 --- /dev/null +++ b/DeepNtuplizer/interface/ntuple_pairwise.h @@ -0,0 +1,88 @@ +/* + * ntuple_pairwise.h + * + * Created on: 01 sep 2022 + * Author: Alexandre De Moor + */ + +#ifndef DEEPNTUPLES_DEEPNTUPLIZER_INTERFACE_NTUPLE_PAIRWISE_H_ +#define DEEPNTUPLES_DEEPNTUPLIZER_INTERFACE_NTUPLE_PAIRWISE_H_ + +#include "ntuple_content.h" +#include "DataFormats/PatCandidates/interface/PackedCandidate.h" +#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h" +#include "TrackPairInfoBuilder.h" + +class ntuple_pairwise: public ntuple_content{ +public: + + ntuple_pairwise():ntuple_content(),jetradius_(0.4){}//, + //n_Cpfcand2_(0),n_Npfcand2_(0){} + + void setJetRadius(const float& radius){jetradius_=radius;} + void getInput(const edm::ParameterSet& iConfig); + void initBranches(TTree* ); + void readEvent(const edm::Event& iEvent); + void readSetup(const edm::EventSetup& iSetup); + + + //use either of these functions + + bool fillBranches(const pat::Jet &, const size_t& jetidx, const edm::View * coll=0); + +private: + + float jetradius_; + float min_candidate_pt_ = -1; + + edm::ESHandle builder; + + int n_Npfcand2_; + int n_Cpfcand2_; + int n_Cpfpairs_; + float nCpfpairs_; + + static constexpr size_t max_pfcand_=800; + + float pair_pca_distance_[max_pfcand_]; + float pair_pca_significance_[max_pfcand_]; + + float pair_pcaSeed_x1_[max_pfcand_]; + float pair_pcaSeed_y1_[max_pfcand_]; + float pair_pcaSeed_z1_[max_pfcand_]; + + float pair_pcaSeed_x2_[max_pfcand_]; + float pair_pcaSeed_y2_[max_pfcand_]; + float pair_pcaSeed_z2_[max_pfcand_]; + + float pair_pcaSeed_xerr1_[max_pfcand_]; + float pair_pcaSeed_yerr1_[max_pfcand_]; + float pair_pcaSeed_zerr1_[max_pfcand_]; + + float pair_pcaSeed_xerr2_[max_pfcand_]; + float pair_pcaSeed_yerr2_[max_pfcand_]; + float pair_pcaSeed_zerr2_[max_pfcand_]; + + float pair_dotprod1_[max_pfcand_]; + float pair_dotprod2_[max_pfcand_]; + + float pair_pca_dist1_[max_pfcand_]; + float pair_pca_dist2_[max_pfcand_]; + + float pair_dotprod12_2D_[max_pfcand_]; + float pair_dotprod12_2DV_[max_pfcand_]; + float pair_dotprod12_3D_[max_pfcand_]; + float pair_dotprod12_3DV_[max_pfcand_]; + + float pair_pca_jetAxis_dist_[max_pfcand_]; + float pair_pca_jetAxis_dotprod_[max_pfcand_]; + float pair_pca_jetAxis_dEta_[max_pfcand_]; + float pair_pca_jetAxis_dPhi_[max_pfcand_]; + + float pfcand_dist_vtx_12_[max_pfcand_]; + + float mindrsvpfcand(const pat::PackedCandidate* pfcand); + +}; + +#endif /* DEEPNTUPLES_DEEPNTUPLIZER_INTERFACE_NTUPLE_PAIRWISE_H_ */ diff --git a/DeepNtuplizer/interface/ntuple_pfCands.h b/DeepNtuplizer/interface/ntuple_pfCands.h index 60e1665c699..6c21bde72f8 100644 --- a/DeepNtuplizer/interface/ntuple_pfCands.h +++ b/DeepNtuplizer/interface/ntuple_pfCands.h @@ -12,125 +12,240 @@ #include "DataFormats/PatCandidates/interface/PackedCandidate.h" #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h" +#include "DataFormats/PatCandidates/interface/PackedGenParticle.h" +#include "DataFormats/Candidate/interface/Candidate.h" +#include "DataFormats/HepMCCandidate/interface/GenParticle.h" + class ntuple_pfCands: public ntuple_content{ -public: + public: - ntuple_pfCands():ntuple_content(),jetradius_(0.4), + ntuple_pfCands():ntuple_content(),jetradius_(0.4), n_Cpfcand_(0),n_Npfcand_(0){} - void setJetRadius(const float& radius){jetradius_=radius;} - void getInput(const edm::ParameterSet& iConfig); - void initBranches(TTree* ); - void readEvent(const edm::Event& iEvent); - void readSetup(const edm::EventSetup& iSetup); - - - //use either of these functions - - bool fillBranches(const pat::Jet &, const size_t& jetidx, const edm::View * coll=0); - -private: - - float jetradius_; - float min_candidate_pt_ = -1; - - edm::ESHandle builder; - - unsigned int n_Cpfcand_; - float nCpfcand_; - - static constexpr size_t max_pfcand_=50; - - float Cpfcan_pt_[max_pfcand_]; - float Cpfcan_eta_[max_pfcand_]; - float Cpfcan_phi_[max_pfcand_]; - float Cpfcan_ptrel_[max_pfcand_]; - float Cpfcan_erel_[max_pfcand_]; - float Cpfcan_phirel_[max_pfcand_]; - float Cpfcan_etarel_[max_pfcand_]; - float Cpfcan_deltaR_[max_pfcand_]; - float Cpfcan_puppiw_[max_pfcand_]; - float Cpfcan_VTX_ass_[max_pfcand_]; - - float Cpfcan_fromPV_[max_pfcand_]; - - float Cpfcan_vertexChi2_[max_pfcand_]; - float Cpfcan_vertexNdof_[max_pfcand_]; - float Cpfcan_vertexNormalizedChi2_[max_pfcand_]; - float Cpfcan_vertex_rho_[max_pfcand_]; - float Cpfcan_vertex_phirel_[max_pfcand_]; - float Cpfcan_vertex_etarel_[max_pfcand_]; - float Cpfcan_vertexRef_mass_[max_pfcand_]; - - // covariance - float Cpfcan_dz_[max_pfcand_]; - float Cpfcan_dxy_[max_pfcand_]; - - float Cpfcan_dxyerrinv_[max_pfcand_]; - float Cpfcan_dxysig_[max_pfcand_]; - - float Cpfcan_dptdpt_[max_pfcand_]; - float Cpfcan_detadeta_[max_pfcand_]; - float Cpfcan_dphidphi_[max_pfcand_]; - float Cpfcan_dxydxy_[max_pfcand_]; - float Cpfcan_dzdz_[max_pfcand_]; - float Cpfcan_dxydz_[max_pfcand_]; - float Cpfcan_dphidxy_[max_pfcand_]; - float Cpfcan_dlambdadz_[max_pfcand_]; - - - - float Cpfcan_BtagPf_trackMomentum_[max_pfcand_]; - float Cpfcan_BtagPf_trackEta_[max_pfcand_]; - float Cpfcan_BtagPf_trackEtaRel_[max_pfcand_]; - float Cpfcan_BtagPf_trackPtRel_[max_pfcand_]; - float Cpfcan_BtagPf_trackPPar_[max_pfcand_]; - float Cpfcan_BtagPf_trackDeltaR_[max_pfcand_]; - float Cpfcan_BtagPf_trackPtRatio_[max_pfcand_]; - float Cpfcan_BtagPf_trackPParRatio_[max_pfcand_]; - float Cpfcan_BtagPf_trackSip3dVal_[max_pfcand_]; - float Cpfcan_BtagPf_trackSip3dSig_[max_pfcand_]; - float Cpfcan_BtagPf_trackSip2dVal_[max_pfcand_]; - float Cpfcan_BtagPf_trackSip2dSig_[max_pfcand_]; - - float Cpfcan_BtagPf_trackDecayLen_[max_pfcand_]; - - float Cpfcan_BtagPf_trackJetDistVal_[max_pfcand_]; - float Cpfcan_BtagPf_trackJetDistSig_[max_pfcand_]; - - // ID, skipped "charged hadron" as that is true if now the other - // TODO (comment of Markus Stoye) add reco information - float Cpfcan_isMu_[max_pfcand_]; // pitty that the quality is missing - float Cpfcan_isEl_[max_pfcand_]; // pitty that the quality is missing - float Cpfcan_charge_[max_pfcand_]; - - // track quality - float Cpfcan_lostInnerHits_[max_pfcand_]; - float Cpfcan_numberOfPixelHits_[max_pfcand_]; - - float Cpfcan_chi2_[max_pfcand_]; - float Cpfcan_quality_[max_pfcand_]; - - float Cpfcan_drminsv_[max_pfcand_]; - - //Neutral Pf candidates - unsigned int n_Npfcand_; - float nNpfcand_; - float Npfcan_pt_[max_pfcand_]; - float Npfcan_eta_[max_pfcand_]; - float Npfcan_phi_[max_pfcand_]; - float Npfcan_ptrel_[max_pfcand_]; - float Npfcan_erel_[max_pfcand_]; - float Npfcan_puppiw_[max_pfcand_]; - float Npfcan_phirel_[max_pfcand_]; - float Npfcan_etarel_[max_pfcand_]; - float Npfcan_deltaR_[max_pfcand_]; - float Npfcan_isGamma_[max_pfcand_]; - float Npfcan_HadFrac_[max_pfcand_]; - float Npfcan_drminsv_[max_pfcand_]; - - - float mindrsvpfcand(const pat::PackedCandidate* pfcand); + void setJetRadius(const float& radius){jetradius_=radius;} + void getInput(const edm::ParameterSet& iConfig); + void initBranches(TTree* ); + virtual void readConfig(const edm::ParameterSet& iConfig, edm::ConsumesCollector && cc) override; + void readEvent(const edm::Event& iEvent); + void readSetup(const edm::EventSetup& iSetup); + + + //use either of these functions + + bool fillBranches(const pat::Jet &, const size_t& jetidx, const edm::View * coll=0); + + private: + + float jetradius_; + float min_candidate_pt_ = -1; + + edm::ESHandle builder; + + edm::EDGetTokenT> ltToken_; + edm::Handle> LTs; + + edm::EDGetTokenT> packedToken_; + edm::Handle> packed; + + edm::EDGetTokenT> prunedToken_; + edm::Handle> pruned; + + int n_Cpfcand_; + float nCpfcand_; + + static constexpr size_t max_pfcand_=50; + + float Cpfcan_pt_[max_pfcand_]; + float Cpfcan_eta_[max_pfcand_]; + float Cpfcan_phi_[max_pfcand_]; + float Cpfcan_ptrel_[max_pfcand_]; + float Cpfcan_e_[max_pfcand_]; + float Cpfcan_erel_[max_pfcand_]; + float Cpfcan_phirel_[max_pfcand_]; + float Cpfcan_etarel_[max_pfcand_]; + float Cpfcan_deltaR_[max_pfcand_]; + float Cpfcan_puppiw_[max_pfcand_]; + float Cpfcan_VTX_ass_[max_pfcand_]; + + float Cpfcan_fromPV_[max_pfcand_]; + + float Cpfcan_vertexChi2_[max_pfcand_]; + float Cpfcan_vertexNdof_[max_pfcand_]; + float Cpfcan_vertexNormalizedChi2_[max_pfcand_]; + float Cpfcan_vertex_rho_[max_pfcand_]; + float Cpfcan_vertex_phirel_[max_pfcand_]; + float Cpfcan_vertex_etarel_[max_pfcand_]; + float Cpfcan_vertexRef_mass_[max_pfcand_]; + + // covariance + float Cpfcan_dz_[max_pfcand_]; + float Cpfcan_dxy_[max_pfcand_]; + + float Cpfcan_dxyerrinv_[max_pfcand_]; + float Cpfcan_dxysig_[max_pfcand_]; + + float Cpfcan_dptdpt_[max_pfcand_]; + float Cpfcan_detadeta_[max_pfcand_]; + float Cpfcan_dphidphi_[max_pfcand_]; + float Cpfcan_dxydxy_[max_pfcand_]; + float Cpfcan_dzdz_[max_pfcand_]; + float Cpfcan_dxydz_[max_pfcand_]; + float Cpfcan_dphidxy_[max_pfcand_]; + float Cpfcan_dlambdadz_[max_pfcand_]; + + + + float Cpfcan_BtagPf_trackMomentum_[max_pfcand_]; + float Cpfcan_BtagPf_trackEta_[max_pfcand_]; + float Cpfcan_BtagPf_trackEtaRel_[max_pfcand_]; + float Cpfcan_BtagPf_trackPtRel_[max_pfcand_]; + float Cpfcan_BtagPf_trackPPar_[max_pfcand_]; + float Cpfcan_BtagPf_trackDeltaR_[max_pfcand_]; + float Cpfcan_BtagPf_trackPtRatio_[max_pfcand_]; + float Cpfcan_BtagPf_trackPParRatio_[max_pfcand_]; + float Cpfcan_BtagPf_trackSip3dVal_[max_pfcand_]; + float Cpfcan_BtagPf_trackSip3dSig_[max_pfcand_]; + float Cpfcan_BtagPf_trackSip2dVal_[max_pfcand_]; + float Cpfcan_BtagPf_trackSip2dSig_[max_pfcand_]; + + float Cpfcan_BtagPf_trackDecayLen_[max_pfcand_]; + + float Cpfcan_BtagPf_trackJetDistVal_[max_pfcand_]; + float Cpfcan_BtagPf_trackJetDistSig_[max_pfcand_]; + + // ID, skipped "charged hadron" as that is true if now the other + // TODO (comment of Markus Stoye) add reco information + float Cpfcan_isMu_[max_pfcand_]; // pitty that the quality is missing + float Cpfcan_isEl_[max_pfcand_]; // pitty that the quality is missing + float Cpfcan_charge_[max_pfcand_]; + + // track quality + float Cpfcan_lostInnerHits_[max_pfcand_]; + float Cpfcan_numberOfPixelHits_[max_pfcand_]; + float Cpfcan_chi2_[max_pfcand_]; + float Cpfcan_quality_[max_pfcand_]; + float Cpfcan_drminsv_[max_pfcand_]; + float Cpfcan_distminsvold_[max_pfcand_]; + float Cpfcan_distminsv_[max_pfcand_]; + float Cpfcan_distminsv2_[max_pfcand_]; + float Cpfcan_dxminsv_[max_pfcand_]; + float Cpfcan_dyminsv_[max_pfcand_]; + float Cpfcan_dzminsv_[max_pfcand_]; + float Cpfcan_dxpv_[max_pfcand_]; + float Cpfcan_dypv_[max_pfcand_]; + float Cpfcan_dzpv_[max_pfcand_]; + + //hit pattern variables, as defined here https://github.com/cms-sw/cmssw/blob/master/DataFormats/TrackReco/interface/HitPattern.h + //Tracker per layer + //Pixel barrel + float Cpfcan_nhitpixelBarrelLayer1_[max_pfcand_]; + float Cpfcan_nhitpixelBarrelLayer2_[max_pfcand_]; + float Cpfcan_nhitpixelBarrelLayer3_[max_pfcand_]; + float Cpfcan_nhitpixelBarrelLayer4_[max_pfcand_]; + //Pixel Endcap + float Cpfcan_nhitpixelEndcapLayer1_[max_pfcand_]; + float Cpfcan_nhitpixelEndcapLayer2_[max_pfcand_]; + //Strip TIB + float Cpfcan_nhitstripTIBLayer1_[max_pfcand_]; + float Cpfcan_nhitstripTIBLayer2_[max_pfcand_]; + float Cpfcan_nhitstripTIBLayer3_[max_pfcand_]; + float Cpfcan_nhitstripTIBLayer4_[max_pfcand_]; + //Strip TID + float Cpfcan_nhitstripTIDLayer1_[max_pfcand_]; + float Cpfcan_nhitstripTIDLayer2_[max_pfcand_]; + float Cpfcan_nhitstripTIDLayer3_[max_pfcand_]; + //Strip TOB + float Cpfcan_nhitstripTOBLayer1_[max_pfcand_]; + float Cpfcan_nhitstripTOBLayer2_[max_pfcand_]; + float Cpfcan_nhitstripTOBLayer3_[max_pfcand_]; + float Cpfcan_nhitstripTOBLayer4_[max_pfcand_]; + float Cpfcan_nhitstripTOBLayer5_[max_pfcand_]; + float Cpfcan_nhitstripTOBLayer6_[max_pfcand_]; + //Strip TEC + float Cpfcan_nhitstripTECLayer1_[max_pfcand_]; + float Cpfcan_nhitstripTECLayer2_[max_pfcand_]; + float Cpfcan_nhitstripTECLayer3_[max_pfcand_]; + float Cpfcan_nhitstripTECLayer4_[max_pfcand_]; + float Cpfcan_nhitstripTECLayer5_[max_pfcand_]; + float Cpfcan_nhitstripTECLayer6_[max_pfcand_]; + float Cpfcan_nhitstripTECLayer7_[max_pfcand_]; + float Cpfcan_nhitstripTECLayer8_[max_pfcand_]; + float Cpfcan_nhitstripTECLayer9_[max_pfcand_]; + //Tracker all layers together + //Valid hits + float Cpfcan_numberOfValidHits_[max_pfcand_]; + float Cpfcan_numberOfValidTrackerHits_[max_pfcand_]; + float Cpfcan_numberOfValidPixelHits_[max_pfcand_]; + float Cpfcan_numberOfValidPixelBarrelHits_[max_pfcand_]; + float Cpfcan_numberOfValidPixelEndcapHits_[max_pfcand_]; + float Cpfcan_numberOfValidStripHits_[max_pfcand_]; + float Cpfcan_numberOfValidStripTIBHits_[max_pfcand_]; + float Cpfcan_numberOfValidStripTIDHits_[max_pfcand_]; + float Cpfcan_numberOfValidStripTOBHits_[max_pfcand_]; + float Cpfcan_numberOfValidStripTECHits_[max_pfcand_]; + //LayersWithMeasurement + float Cpfcan_trackerLayersWithMeasurementOld_[max_pfcand_]; + float Cpfcan_trackerLayersWithMeasurement_[max_pfcand_]; + float Cpfcan_pixelLayersWithMeasurementOld_[max_pfcand_]; + float Cpfcan_pixelLayersWithMeasurement_[max_pfcand_]; + float Cpfcan_stripLayersWithMeasurement_[max_pfcand_]; + float Cpfcan_pixelBarrelLayersWithMeasurement_[max_pfcand_]; + float Cpfcan_pixelEndcapLayersWithMeasurement_[max_pfcand_]; + float Cpfcan_stripTIBLayersWithMeasurement_[max_pfcand_]; + float Cpfcan_stripTIDLayersWithMeasurement_[max_pfcand_]; + float Cpfcan_stripTOBLayersWithMeasurement_[max_pfcand_]; + float Cpfcan_stripTECLayersWithMeasurement_[max_pfcand_]; + //Null + float Cpfcan_trackerLayersNull_[max_pfcand_]; + float Cpfcan_pixelLayersNull_[max_pfcand_]; + float Cpfcan_stripLayersNull_[max_pfcand_]; + float Cpfcan_pixelBarrelLayersNull_[max_pfcand_]; + float Cpfcan_pixelEndcapLayersNull_[max_pfcand_]; + float Cpfcan_stripTIBLayersNull_[max_pfcand_]; + float Cpfcan_stripTIDLayersNull_[max_pfcand_]; + float Cpfcan_stripTOBLayersNull_[max_pfcand_]; + float Cpfcan_stripTECLayersNull_[max_pfcand_]; + + //Neutral Pf candidates + int n_Npfcand_; + float nNpfcand_; + float Npfcan_pt_[max_pfcand_]; + float Npfcan_eta_[max_pfcand_]; + float Npfcan_phi_[max_pfcand_]; + float Npfcan_ptrel_[max_pfcand_]; + float Npfcan_e_[max_pfcand_]; + float Npfcan_erel_[max_pfcand_]; + float Npfcan_puppiw_[max_pfcand_]; + float Npfcan_phirel_[max_pfcand_]; + float Npfcan_etarel_[max_pfcand_]; + float Npfcan_deltaR_[max_pfcand_]; + float Npfcan_isGamma_[max_pfcand_]; + float Npfcan_HadFrac_[max_pfcand_]; + //float Npfcan_CaloFrac_[max_pfcand_]; + float Npfcan_drminsv_[max_pfcand_]; + + float Npfcan_pdgID_[max_pfcand_]; + float Cpfcan_pdgID_[max_pfcand_]; + + float Cpfcan_HadFrac_[max_pfcand_]; + float Cpfcan_CaloFrac_[max_pfcand_]; + + int Cpfcan_c_tag_[max_pfcand_]; + int Cpfcan_b_tag_[max_pfcand_]; + int Cpfcan_g_tag_[max_pfcand_]; + + float Cpfcan_vtx_x_[max_pfcand_]; + float Cpfcan_vtx_y_[max_pfcand_]; + float Cpfcan_vtx_z_[max_pfcand_]; + + float Cpfcan_dist_from_pv_[max_pfcand_]; + + float mindrsvpfcand(const pat::PackedCandidate* pfcand); + float mindistsvpfcandold(const reco::TransientTrack track); + float mindistsvpfcand(const reco::TransientTrack track); + GlobalPoint mingpsvpfcand(const reco::TransientTrack track); + GlobalPoint gppvpfcand(const reco::TransientTrack track, GlobalVector direction, const reco::Vertex vertex); + bool containParton(const reco::Candidate *pruned_part, int pdgid); }; diff --git a/DeepNtuplizer/plugins/DeepNtuplizer.cc b/DeepNtuplizer/plugins/DeepNtuplizer.cc index 761bd99e7ca..a8201c701ff 100644 --- a/DeepNtuplizer/plugins/DeepNtuplizer.cc +++ b/DeepNtuplizer/plugins/DeepNtuplizer.cc @@ -8,12 +8,13 @@ #include "../interface/ntuple_content.h" #include "../interface/ntuple_SV.h" +#include "../interface/ntuple_LT.h" #include "../interface/ntuple_JetInfo.h" #include "../interface/ntuple_pfCands.h" #include "../interface/ntuple_bTagVars.h" -#include "../interface/ntuple_FatJetInfo.h" -#include "../interface/ntuple_DeepVertex.h" -#include "../interface/ntuple_GraphB.h" +//#include "../interface/ntuple_FatJetInfo.h" +//#include "../interface/ntuple_DeepVertex.h" +#include "../interface/ntuple_pairwise.h" //ROOT includes #include "TTree.h" #include @@ -156,9 +157,10 @@ DeepNtuplizer::DeepNtuplizer(const edm::ParameterSet& iConfig): } */ - //ntuple_GraphB* deepvertexmodule=new ntuple_GraphB(jetR); - //deepvertexmodule->setCandidatesToken(consumes >(iConfig.getParameter("candidates"))); - //addModule(deepvertexmodule); + + /*ntuple_GraphB* deepvertexmodule=new ntuple_GraphB(jetR); + deepvertexmodule->setCandidatesToken(consumes >(iConfig.getParameter("candidates"))); + addModule(deepvertexmodule, "DeepVertextuple");*/ ntuple_JetInfo* jetinfo=new ntuple_JetInfo(); @@ -198,26 +200,34 @@ DeepNtuplizer::DeepNtuplizer(const edm::ParameterSet& iConfig): ntuple_pfCands * pfcands = new ntuple_pfCands(); pfcands->setJetRadius(jetR); - addModule(pfcands, "pfcands"); + ntuple_LT * LT = new ntuple_LT(); + LT->setJetRadius(jetR); + addModule(LT, "LT"); + + ntuple_pairwise * pairwise = new ntuple_pairwise(); + pairwise->setJetRadius(jetR); + addModule(pairwise, "pairwise"); + addModule(new ntuple_bTagVars(), "bTagVars"); - if(runFatJets_){ + /*if(runFatJets_){ auto *fatjetinfo = new ntuple_FatJetInfo(jetR); fatjetinfo->setGenParticleToken(consumes(iConfig.getParameter("pruned"))); fatjetinfo->setFatJetToken(consumes(iConfig.getParameter("fatjets"))); addModule(fatjetinfo, "fatJets"); - } + }*/ /* * * Modules initialized * * parse the input parameters (if any) */ - for(auto& m: modules_) - m->getInput(iConfig); - + for(auto& m: modules_){ + m->readConfig(iConfig, consumesCollector()); + m->getInput(iConfig); + } } diff --git a/DeepNtuplizer/production/DV_small.cfg b/DeepNtuplizer/production/DV_small.cfg new file mode 100644 index 00000000000..4ca95109d34 --- /dev/null +++ b/DeepNtuplizer/production/DV_small.cfg @@ -0,0 +1,12 @@ +500 /BulkGravitonToHHTo4Q_MX-600to6000_MH-15to250_part1_TuneCP5_13TeV-madgraph_pythia8/RunIISummer19UL17MiniAODv2-multigridpack_106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_bulk_grav_1 +500 /BulkGravitonToHHTo4Q_MX-600to6000_MH-15to250_part2_TuneCP5_13TeV-madgraph_pythia8/RunIISummer19UL17MiniAODv2-multigridpack_106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_bulk_grav_2 +500 /BulkGravitonToHHTo4Q_MX-600to6000_MH-15to250_part3_TuneCP5_13TeV-madgraph_pythia8/RunIISummer19UL17MiniAODv2-multigridpack_106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_bulk_grav_3 + +500 /QCD_HT1000to1500_TuneCP5_PSWeights_13TeV-madgraphMLM-pythia8/RunIISummer19UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1000_1500 +500 /QCD_HT100to200_TuneCP5_PSWeights_13TeV-madgraphMLM-pythia8/RunIISummer19UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_100_200 +500 /QCD_HT1500to2000_TuneCP5_PSWeights_13TeV-madgraphMLM-pythia8/RunIISummer19UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1500_2000 +500 /QCD_HT2000toInf_TuneCP5_PSWeights_13TeV-madgraphMLM-pythia8/RunIISummer19UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_2000_Inf +500 /QCD_HT200to300_TuneCP5_PSWeights_13TeV-madgraphMLM-pythia8/RunIISummer19UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_200_300 +500 /QCD_HT300to500_TuneCP5_PSWeights_13TeV-madgraphMLM-pythia8/RunIISummer19UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_300_500 +500 /QCD_HT500to700_TuneCP5_PSWeights_13TeV-madgraphMLM-pythia8/RunIISummer19UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_500_700 +500 /QCD_HT700to1000_TuneCP5_PSWeights_13TeV-madgraphMLM-pythia8/RunIISummer19UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_700_1000 diff --git a/DeepNtuplizer/production/DeepNtuplizer.py b/DeepNtuplizer/production/DeepNtuplizer.py index a4729d6b610..e3ebf6c5afb 100644 --- a/DeepNtuplizer/production/DeepNtuplizer.py +++ b/DeepNtuplizer/production/DeepNtuplizer.py @@ -9,11 +9,11 @@ options.register('inputScript','',VarParsing.VarParsing.multiplicity.singleton,VarParsing.VarParsing.varType.string,"input Script") options.register('outputFile','output',VarParsing.VarParsing.multiplicity.singleton,VarParsing.VarParsing.varType.string,"output File (w/o .root)") -options.register('maxEvents',-1,VarParsing.VarParsing.multiplicity.singleton,VarParsing.VarParsing.varType.int,"maximum events") +options.register('maxEvents',50001,VarParsing.VarParsing.multiplicity.singleton,VarParsing.VarParsing.varType.int,"maximum events") options.register('skipEvents', 0, VarParsing.VarParsing.multiplicity.singleton, VarParsing.VarParsing.varType.int, "skip N events") options.register('job', 0, VarParsing.VarParsing.multiplicity.singleton, VarParsing.VarParsing.varType.int, "job number") options.register('nJobs', 1, VarParsing.VarParsing.multiplicity.singleton, VarParsing.VarParsing.varType.int, "total jobs") -options.register('reportEvery', 100, VarParsing.VarParsing.multiplicity.singleton, VarParsing.VarParsing.varType.int, "report every") +options.register('reportEvery', 1000, VarParsing.VarParsing.multiplicity.singleton, VarParsing.VarParsing.varType.int, "report every") options.register('gluonReduction', 0.0, VarParsing.VarParsing.multiplicity.singleton, VarParsing.VarParsing.varType.float, "gluon reduction") options.register('selectJets', True, VarParsing.VarParsing.multiplicity.singleton, VarParsing.VarParsing.varType.bool, "select jets with good gen match") options.register('phase2', False, VarParsing.VarParsing.multiplicity.singleton, VarParsing.VarParsing.varType.bool, "apply jet selection for phase 2. Currently sets JetEtaMax to 3.0 and picks slimmedJetsPuppi as jet collection.") @@ -94,12 +94,15 @@ bTagInfos = [ 'pfImpactParameterTagInfos', 'pfInclusiveSecondaryVertexFinderTagInfos', - 'pfDeepCSVTagInfos' ] + 'pfDeepCSVTagInfos', + 'pfParticleNetAK4TagInfos', +] else : bTagInfos = [ 'pfImpactParameterTagInfos', 'pfInclusiveSecondaryVertexFinderTagInfos', 'deepNNTagInfos', + 'pfParticleNetAK4TagInfos', ] @@ -114,6 +117,18 @@ 'pfDeepCSVJetTags:probb', 'pfDeepCSVJetTags:probc', 'pfDeepCSVJetTags:probbb', + 'pfDeepFlavourJetTags:probb', + 'pfDeepFlavourJetTags:probbb', + 'pfDeepFlavourJetTags:problepb', + 'pfDeepFlavourJetTags:probc', + 'pfDeepFlavourJetTags:probuds', + 'pfDeepFlavourJetTags:probg', + 'pfParticleNetAK4JetTags:probb', + 'pfParticleNetAK4JetTags:probbb', + 'pfParticleNetAK4JetTags:probc', + 'pfParticleNetAK4JetTags:probcc', + 'pfParticleNetAK4JetTags:probuds', + 'pfParticleNetAK4JetTags:probg', ] else : bTagDiscriminators = [ @@ -127,6 +142,18 @@ 'pfDeepCSVJetTags:probc', 'pfDeepCSVJetTags:probbb', 'pfDeepCSVJetTags:probcc', + 'pfDeepFlavourJetTags:probb', + 'pfDeepFlavourJetTags:probbb', + 'pfDeepFlavourJetTags:problepb', + 'pfDeepFlavourJetTags:probc', + 'pfDeepFlavourJetTags:probuds', + 'pfDeepFlavourJetTags:probg', + 'pfParticleNetAK4JetTags:probb', + 'pfParticleNetAK4JetTags:probbb', + 'pfParticleNetAK4JetTags:probc', + 'pfParticleNetAK4JetTags:probcc', + 'pfParticleNetAK4JetTags:probuds', + 'pfParticleNetAK4JetTags:probg', ] jetCorrectionsAK4 = ('AK4PFchs', ['L1FastJet', 'L2Relative', 'L3Absolute'], 'None') @@ -242,7 +269,7 @@ # DeepNtuplizer process.load("DeepNTuples.DeepNtuplizer.DeepNtuplizer_cfi") -process.deepntuplizer.jets = cms.InputTag('selectedUpdatedPatJetsDeepFlavour'); +process.deepntuplizer.jets = cms.InputTag('selectedUpdatedPatJetsDeepFlavour') process.deepntuplizer.bDiscriminators = bTagDiscriminators process.deepntuplizer.bDiscriminators.append('pfCombinedMVAV2BJetTags') process.deepntuplizer.LooseSVs = cms.InputTag("looseIVFinclusiveCandidateSecondaryVertices") @@ -274,9 +301,9 @@ #Trick to make it work in 9_1_X process.tsk = cms.Task() -for mod in process.producers_().itervalues(): +for mod in process.producers_().values(): #.itervalues(): process.tsk.add(mod) -for mod in process.filters_().itervalues(): +for mod in process.filters_().values(): #.itervalues(): process.tsk.add(mod) process.p = cms.Path( diff --git a/DeepNtuplizer/production/DeepVertexing_campaign.cfg b/DeepNtuplizer/production/DeepVertexing_campaign.cfg new file mode 100644 index 00000000000..e8cb9253563 --- /dev/null +++ b/DeepNtuplizer/production/DeepVertexing_campaign.cfg @@ -0,0 +1,18 @@ + + +14 /TTToHadronic_TuneCP5_13TeV-powheg-pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM ntuple_ttbar_had_puppi_102X + +2 /QCD_Pt_15to30_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15_ext1-v2/MINIAODSIM ntuple_qcd_15_30_puppi_102X +2 /QCD_Pt_30to50_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM ntuple_qcd_30_50_puppi_102X +2 /QCD_Pt_50to80_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM ntuple_qcd_50_80_puppi_102X +2 /QCD_Pt_80to120_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM ntuple_qcd_80_120_puppi_102X +2 /QCD_Pt_120to170_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM ntuple_qcd_120_170_puppi_102X +2 /QCD_Pt_170to300_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM ntuple_qcd_170_300_puppi_102X +2 /QCD_Pt_300to470_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM ntuple_qcd_300_470_puppi_102X +2 /QCD_Pt_470to600_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM ntuple_qcd_470_600_puppi_102X +2 /QCD_Pt_600to800_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM ntuple_qcd_600_800_puppi_102X +2 /QCD_Pt_800to1000_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15_ext1-v2/MINIAODSIM ntuple_qcd_800_1000_puppi_102X +2 /QCD_Pt_1000to1400_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM ntuple_qcd_1000_1400_puppi_102X +2 /QCD_Pt_1400to1800_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM ntuple_qcd_1400_1800_puppi_102X +2 /QCD_Pt_1800to2400_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM ntuple_qcd_1800_2400_puppi_102X +2 /QCD_Pt_2400to3200_TuneCP5_13TeV_pythia8/RunIIAutumn18MiniAOD-102X_upgrade2018_realistic_v15-v1/MINIAODSIM ntuple_qcd_2400_3200_puppi_102X \ No newline at end of file diff --git a/DeepNtuplizer/production/New_UL17.cfg b/DeepNtuplizer/production/New_UL17.cfg new file mode 100644 index 00000000000..59b65ada24b --- /dev/null +++ b/DeepNtuplizer/production/New_UL17.cfg @@ -0,0 +1,14 @@ +500 /TTToHadronic_TuneCP5_13TeV-powheg-pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v2/MINIAODSIM ntuple_ttbar_had +500 /TTToSemiLeptonic_TuneCP5_13TeV-powheg-pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_ttbar_semilep + +500 /QCD_Pt_120to170_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1000_1400 gluonReduction=0.5 +500 /QCD_Pt_170to300_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1000_1400 gluonReduction=0.5 +500 /QCD_Pt_300to470_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1000_1400 gluonReduction=0.5 +500 /QCD_Pt_470to600_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1000_1400 gluonReduction=0.5 +500 /QCD_Pt_600to800_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1000_1400 gluonReduction=0.5 +500 /QCD_Pt_800to1000_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1000_1400 gluonReduction=0.5 +500 /QCD_Pt_1000to1400_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1000_1400 gluonReduction=0.5 +500 /QCD_Pt_1400to1800_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1400_1800 gluonReduction=0.5 +500 /QCD_Pt_1800to2400_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1800_2400 gluonReduction=0.5 +500 /QCD_Pt_2400to3200_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_2400_3200 gluonReduction=0.5 +500 /QCD_Pt_3200toInf_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_3200_inf gluonReduction=0.5 \ No newline at end of file diff --git a/DeepNtuplizer/production/New_UL17_HT.cfg b/DeepNtuplizer/production/New_UL17_HT.cfg new file mode 100644 index 00000000000..82ba7538ce3 --- /dev/null +++ b/DeepNtuplizer/production/New_UL17_HT.cfg @@ -0,0 +1,4 @@ +400 /QCD_HT300to500_TuneCP5_13TeV-madgraphMLM-pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v2/MINIAODSIM ntuple_qcd_HT_300_500 gluonReduction=0.5 +400 /QCD_HT500to700_TuneCP5_13TeV-madgraphMLM-pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v2/MINIAODSIM ntuple_qcd_HT_500_700 gluonReduction=0.5 +400 /QCD_HT700to1000_TuneCP5_13TeV-madgraphMLM-pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v2/MINIAODSIM ntuple_qcd_HT_700_1000 gluonReduction=0.5 +400 /QCD_HT1000to1500_TuneCP5_13TeV-madgraphMLM-pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v2/MINIAODSIM ntuple_qcd_HT_1000_1500 gluonReduction=0.5 \ No newline at end of file diff --git a/DeepNtuplizer/production/New_UL17_QCD.cfg b/DeepNtuplizer/production/New_UL17_QCD.cfg new file mode 100644 index 00000000000..1f4c35eb682 --- /dev/null +++ b/DeepNtuplizer/production/New_UL17_QCD.cfg @@ -0,0 +1,5 @@ +850 /QCD_Pt_1000to1400_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1000_1400 gluonReduction=0.5 +850 /QCD_Pt_1400to1800_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1400_1800 gluonReduction=0.5 +850 /QCD_Pt_1800to2400_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_1800_2400 gluonReduction=0.5 +850 /QCD_Pt_2400to3200_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_2400_3200 gluonReduction=0.5 +850 /QCD_Pt_3200toInf_TuneCP5_13TeV_pythia8/RunIISummer20UL17MiniAODv2-106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_qcd_3200_inf gluonReduction=0.5 diff --git a/DeepNtuplizer/production/New_UL17_grav.cfg b/DeepNtuplizer/production/New_UL17_grav.cfg new file mode 100644 index 00000000000..3b577e9413f --- /dev/null +++ b/DeepNtuplizer/production/New_UL17_grav.cfg @@ -0,0 +1,3 @@ +750 /BulkGravitonToHHTo4Q_MX-600to6000_MH-15to250_part1_TuneCP5_13TeV-madgraph_pythia8/RunIISummer19UL17MiniAODv2-multigridpack_106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_blukgrav_part1 gluonReduction=0.5 +750 /BulkGravitonToHHTo4Q_MX-600to6000_MH-15to250_part2_TuneCP5_13TeV-madgraph_pythia8/RunIISummer19UL17MiniAODv2-multigridpack_106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_blukgrav_part2 gluonReduction=0.5 +750 /BulkGravitonToHHTo4Q_MX-600to6000_MH-15to250_part3_TuneCP5_13TeV-madgraph_pythia8/RunIISummer19UL17MiniAODv2-multigridpack_106X_mc2017_realistic_v9-v1/MINIAODSIM ntuple_blukgrav_part3 gluonReduction=0.5 \ No newline at end of file diff --git a/DeepNtuplizer/production/ParT_test.cfg b/DeepNtuplizer/production/ParT_test.cfg new file mode 100644 index 00000000000..0ed89f9f5bf --- /dev/null +++ b/DeepNtuplizer/production/ParT_test.cfg @@ -0,0 +1 @@ +15 /TTToHadronic_TuneCP5_13TeV-powheg-pythia8/RunIISummer20UL17MiniAOD-106X_mc2017_realistic_v6-v2/MINIAODSIM ntuple_ttbar_had diff --git a/DeepNtuplizer/production/samples_phase2.cfg b/DeepNtuplizer/production/samples_phase2.cfg index 10c7708e18d..bc9ba01e103 100644 --- a/DeepNtuplizer/production/samples_phase2.cfg +++ b/DeepNtuplizer/production/samples_phase2.cfg @@ -1,7 +1,7 @@ -10000 /TT_TuneCUETP8M2T4_14TeV-powheg-pythia8/PhaseIITDRFall17MiniAOD-PU200_93X_upgrade2023_realistic_v2_ext1-v1/MINIAODSIM ntuple_tt_phase2_v2_ext1-v1 phase2=True -10000 /TT_TuneCUETP8M2T4_14TeV-powheg-pythia8/PhaseIITDRFall17MiniAOD-PU200_93X_upgrade2023_realistic_v2-v3/MINIAODSIM ntuple_tt_phase2_v2-v3 phase2=True -10000 /TT_TuneCUETP8M2T4_14TeV-powheg-pythia8/PhaseIITDRFall17MiniAOD-PU200_93X_upgrade2023_realistic_v2-v1/MINIAODSIM ntuple_tt_phase2_v2-v1 phase2=True -10000 /QCD_Flat_Pt-15to7000_TuneCUETP8M1_14TeV_pythia8/PhaseIITDRFall17MiniAOD-PU200_93X_upgrade2023_realistic_v2-v1/MINIAODSIM ntuple_qcd_15_7000_phase2_v2-v1 phase2=True -10000 /QCD_Flat_Pt-15to7000_TuneCUETP8M1_14TeV_pythia8/PhaseIITDRFall17MiniAOD-PU200_93X_upgrade2023_realistic_v2-v2/MINIAODSIM ntuple_qcd_15_7000_phase2_v2-v2 phase2=True -10000 /QCD_Flat_Pt-15to7000_TuneCUETP8M1_14TeV_pythia8/PhaseIITDRFall17MiniAOD-PU200FEVT_93X_upgrade2023_realistic_v2_ext1-v1/MINIAODSIM ntuple_qcd_15_7000_phase2_v2_ext-v1 phase2=True +3 /TT_TuneCUETP8M2T4_14TeV-powheg-pythia8/PhaseIITDRFall17MiniAOD-PU200_93X_upgrade2023_realistic_v2_ext1-v1/MINIAODSIM ntuple_tt_phase2_v2_ext1-v1 phase2=True +3 /TT_TuneCUETP8M2T4_14TeV-powheg-pythia8/PhaseIITDRFall17MiniAOD-PU200_93X_upgrade2023_realistic_v2-v3/MINIAODSIM ntuple_tt_phase2_v2-v3 phase2=True +3 /TT_TuneCUETP8M2T4_14TeV-powheg-pythia8/PhaseIITDRFall17MiniAOD-PU200_93X_upgrade2023_realistic_v2-v1/MINIAODSIM ntuple_tt_phase2_v2-v1 phase2=True +3 /QCD_Flat_Pt-15to7000_TuneCUETP8M1_14TeV_pythia8/PhaseIITDRFall17MiniAOD-PU200_93X_upgrade2023_realistic_v2-v1/MINIAODSIM ntuple_qcd_15_7000_phase2_v2-v1 phase2=True +3 /QCD_Flat_Pt-15to7000_TuneCUETP8M1_14TeV_pythia8/PhaseIITDRFall17MiniAOD-PU200_93X_upgrade2023_realistic_v2-v2/MINIAODSIM ntuple_qcd_15_7000_phase2_v2-v2 phase2=True +3 /QCD_Flat_Pt-15to7000_TuneCUETP8M1_14TeV_pythia8/PhaseIITDRFall17MiniAOD-PU200FEVT_93X_upgrade2023_realistic_v2_ext1-v1/MINIAODSIM ntuple_qcd_15_7000_phase2_v2_ext-v1 phase2=True diff --git a/DeepNtuplizer/python/DeepNtuplizer_cfi.py b/DeepNtuplizer/python/DeepNtuplizer_cfi.py index d9280f1ec50..fdd9c8daa15 100644 --- a/DeepNtuplizer/python/DeepNtuplizer_cfi.py +++ b/DeepNtuplizer/python/DeepNtuplizer_cfi.py @@ -1,39 +1,41 @@ import FWCore.ParameterSet.Config as cms deepntuplizer = cms.EDAnalyzer('DeepNtuplizer', - vertices = cms.InputTag("offlineSlimmedPrimaryVertices"), - secVertices = cms.InputTag("slimmedSecondaryVertices"), - jets = cms.InputTag("slimmedJetsPuppi"), - jetR = cms.double(0.4), + vertices = cms.InputTag("offlineSlimmedPrimaryVertices"), + secVertices = cms.InputTag("slimmedSecondaryVertices"), + jets = cms.InputTag("slimmedJetsPuppi"), + losttracks = cms.InputTag("lostTracks"), + packed = cms.InputTag("packedGenParticles"), + jetR = cms.double(0.4), runFatJet = cms.bool(False), eta = cms.bool(True), puppi = cms.bool(True), runDeepVertex = cms.bool(False), - pupInfo = cms.InputTag("slimmedAddPileupInfo"), - rhoInfo = cms.InputTag("fixedGridRhoFastjetAll"), - SVs = cms.InputTag("slimmedSecondaryVertices"), - LooseSVs = cms.InputTag("inclusiveCandidateSecondaryVertices"), - genJetMatchWithNu = cms.InputTag("patGenJetMatchWithNu"), - genJetMatchRecluster = cms.InputTag("patGenJetMatchRecluster"), - genJetMatchAllowDuplicates = cms.InputTag("patGenJetMatchAllowDuplicates"), - pruned = cms.InputTag("prunedGenParticles"), - fatjets = cms.InputTag('slimmedJetsAK8'), - muons = cms.InputTag("slimmedMuons"), - electrons = cms.InputTag("slimmedElectrons"), - jetPtMin = cms.double(10.0), - jetPtMax = cms.double(2000), - jetAbsEtaMin = cms.double(0.0), - jetAbsEtaMax = cms.double(5.0), - gluonReduction = cms.double(0.0), - tagInfoName = cms.string('deepNN'), - tagInfoFName = cms.string('pfBoostedDoubleSVAK8'), - bDiscriminators = cms.vstring(), - qgtagger = cms.string("QGTagger"), - candidates = cms.InputTag("packedPFCandidates"), - minCandidatePt = cms.double(0.95), - - useHerwigCompatible=cms.bool(False), - isHerwig=cms.bool(False), - useOffsets=cms.bool(True), - applySelection=cms.bool(True) - ) + pupInfo = cms.InputTag("slimmedAddPileupInfo"), + rhoInfo = cms.InputTag("fixedGridRhoFastjetAll"), + SVs = cms.InputTag("slimmedSecondaryVertices"), + LooseSVs = cms.InputTag("inclusiveCandidateSecondaryVertices"), + genJetMatchWithNu = cms.InputTag("patGenJetMatchWithNu"), + genJetMatchRecluster = cms.InputTag("patGenJetMatchRecluster"), + genJetMatchAllowDuplicates = cms.InputTag("patGenJetMatchAllowDuplicates"), + pruned = cms.InputTag("prunedGenParticles"), + fatjets = cms.InputTag('slimmedJetsAK8'), + muons = cms.InputTag("slimmedMuons"), + electrons = cms.InputTag("slimmedElectrons"), + jetPtMin = cms.double(10.0), + jetPtMax = cms.double(2000), + jetAbsEtaMin = cms.double(0.0), + jetAbsEtaMax = cms.double(5.0), + gluonReduction = cms.double(0.0), + tagInfoName = cms.string('deepNN'), + tagInfoFName = cms.string('pfBoostedDoubleSVAK8'), + bDiscriminators = cms.vstring(), + qgtagger = cms.string("QGTagger"), + candidates = cms.InputTag("packedPFCandidates"), + minCandidatePt = cms.double(0.95), + + useHerwigCompatible=cms.bool(False), + isHerwig=cms.bool(False), + useOffsets=cms.bool(True), + applySelection=cms.bool(True) +) diff --git a/DeepNtuplizer/scripts/jobSub.py b/DeepNtuplizer/scripts/jobSub.py index abedd1dd942..7c3054890ba 100755 --- a/DeepNtuplizer/scripts/jobSub.py +++ b/DeepNtuplizer/scripts/jobSub.py @@ -29,8 +29,8 @@ def doSub(): parser.add_argument('--file',default='samples.cfg',help='file containing a sample list') parser.add_argument('--nosubmit',default=False,help='no submission') parser.add_argument('--outpath',default='',help='set path to store the .root output') - parser.add_argument('--walltime',default='21600',help='set job wall time in seconds') - parser.add_argument('--maxsize',default='2000',help='set maximum allowed size of output ntuple') + parser.add_argument('--walltime',default='28800',help='set job wall time in seconds') + parser.add_argument('--maxsize',default='15000',help='set maximum allowed size of output ntuple') args = parser.parse_args() @@ -38,7 +38,9 @@ def doSub(): maxSize=args.maxsize #eosGlobalOutDir='/eos/user/'+os.environ['USER'][0]+'/'+os.environ['USER']+'/DeepNtuples' - eosGlobalOutDir='/eos/cms/store/cmst3/group/dehep/DeepJet/NTuples/' + #eosGlobalOutDir='/eos/cms/store/cmst3/group/dehep/DeepJet/NTuples/' + eosGlobalOutDir='/afs/cern.ch/work/a/ademoor/NTuples_DV/' + if len(args.outpath): eosGlobalOutDir=args.outpath if not os.path.isdir(eosGlobalOutDir): diff --git a/DeepNtuplizer/scripts/mergeSamples.py b/DeepNtuplizer/scripts/mergeSamples.py index b8a0cf34f59..3d08b8b9e98 100755 --- a/DeepNtuplizer/scripts/mergeSamples.py +++ b/DeepNtuplizer/scripts/mergeSamples.py @@ -1,11 +1,11 @@ -#!/usr/bin/env python2 +#!/usr/bin/env python3 from argparse import ArgumentParser import os import subprocess def syscall(cmd): - print 'Executing: %s' % cmd + print('Executing: %s' % cmd) retval = os.system(cmd) if retval != 0: raise RuntimeError('Command failed!') @@ -25,7 +25,10 @@ def syscall(cmd): allins='' for l in args.infiles: allins+=' '+l - + + print(str(args.nsamples)) + print(args.outdir) + print(allins) syscall('createMergeList '+str(args.nsamples)+' '+args.outdir+' '+allins) @@ -87,7 +90,7 @@ def syscall(cmd): condor.write(''.join(jobs)) os.chdir(dname) syscall('condor_submit condor.sub') - print 'Once all the jobs are run please run again this command to ensure everything worked' + print('Once all the jobs are run please run again this command to ensure everything worked') else: import multiprocessing as mp diff --git a/DeepNtuplizer/src/TrackPairInfoBuilder.cc b/DeepNtuplizer/src/TrackPairInfoBuilder.cc new file mode 100644 index 00000000000..5bf800a8108 --- /dev/null +++ b/DeepNtuplizer/src/TrackPairInfoBuilder.cc @@ -0,0 +1,145 @@ +#include "DataFormats/Candidate/interface/Candidate.h" +#include "DataFormats/Candidate/interface/VertexCompositePtrCandidate.h" +#include "DataFormats/VertexReco/interface/Vertex.h" +#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h" +#include "TrackingTools/Records/interface/TransientTrackRecord.h" +#include "TVector3.h" +#include "DataFormats/PatCandidates/interface/Jet.h" + +#include "../interface/TrackPairInfoBuilder.h" + +namespace deepntuples { + + TrackPairInfoBuilder::TrackPairInfoBuilder() + : + + track_i_pt_(0), + track_t_pt_(0), + track_eta_(0), + track_phi_(0), + track_dz_(0), + track_dxy_(0), + + pca_distance_(0), + pca_significance_(0), + + pcaSeed_x_(0), + pcaSeed_y_(0), + pcaSeed_z_(0), + pcaSeed_xerr_(0), + pcaSeed_yerr_(0), + pcaSeed_zerr_(0), + pcaTrack_x_(0), + pcaTrack_y_(0), + pcaTrack_z_(0), + pcaTrack_xerr_(0), + pcaTrack_yerr_(0), + pcaTrack_zerr_(0), + + dotprodTrack_(0), + dotprodSeed_(0), + pcaSeed_dist_(0), + pcaTrack_dist_(0), + + track_candMass_(0), + track_ip2d_(0), + track_ip2dSig_(0), + track_ip3d_(0), + track_ip3dSig_(0), + + dotprodTrackSeed2D_(0), + dotprodTrackSeed2DV_(0), + dotprodTrackSeed3D_(0), + dotprodTrackSeed3DV_(0), + + pca_jetAxis_dist_(0), + pca_jetAxis_dotprod_(0), + pca_jetAxis_dEta_(0), + pca_jetAxis_dPhi_(0) + + {} + + void TrackPairInfoBuilder::buildTrackPairInfo(const reco::TransientTrack it, + const reco::TransientTrack tt, + const reco::Vertex& pv, + const pat::Jet &jet + ) { + + GlobalVector jetdirection(jet.px(), jet.py(), jet.pz()); + GlobalPoint pvp(pv.x(), pv.y(), pv.z()); + + VertexDistance3D distanceComputer; + TwoTrackMinimumDistance dist; + + auto const& iImpactState = it.impactPointState(); + auto const& tImpactState = tt.impactPointState(); + + if (dist.calculate(tImpactState, iImpactState)) { + GlobalPoint ttPoint = dist.points().first; + GlobalError ttPointErr = tImpactState.cartesianError().position(); + GlobalPoint seedPosition = dist.points().second; + GlobalError seedPositionErr = iImpactState.cartesianError().position(); + + Measurement1D m = + distanceComputer.distance(VertexState(seedPosition, seedPositionErr), VertexState(ttPoint, ttPointErr)); + + GlobalPoint cp(dist.crossingPoint()); + + GlobalVector pairMomentum((Basic3DVector)(it.track().momentum() + tt.track().momentum())); + GlobalVector pvToPCA(cp - pvp); + + float pvToPCAseed = (seedPosition - pvp).mag(); + float pvToPCAtrack = (ttPoint - pvp).mag(); + float distance = dist.distance(); + + GlobalVector trackDir2D(tImpactState.globalDirection().x(), tImpactState.globalDirection().y(), 0.); + GlobalVector seedDir2D(iImpactState.globalDirection().x(), iImpactState.globalDirection().y(), 0.); + GlobalVector trackPCADir2D(ttPoint.x() - pvp.x(), ttPoint.y() - pvp.y(), 0.); + GlobalVector seedPCADir2D(seedPosition.x() - pvp.x(), seedPosition.y() - pvp.y(), 0.); + + float dotprodTrack = (ttPoint - pvp).unit().dot(tImpactState.globalDirection().unit()); + float dotprodSeed = (seedPosition - pvp).unit().dot(iImpactState.globalDirection().unit()); + + Line::PositionType pos(pvp); + Line::DirectionType dir(jetdirection); + Line::DirectionType pairMomentumDir(pairMomentum); + Line jetLine(pos, dir); + Line PCAMomentumLine(cp, pairMomentumDir); + + track_t_pt_ = tt.track().pt(); + track_i_pt_ = it.track().pt(); + + pca_distance_ = distance; + pca_significance_ = m.significance(); + + pcaSeed_x_ = seedPosition.x(); + pcaSeed_y_ = seedPosition.y(); + pcaSeed_z_ = seedPosition.z(); + pcaSeed_xerr_ = seedPositionErr.cxx(); + pcaSeed_yerr_ = seedPositionErr.cyy(); + pcaSeed_zerr_ = seedPositionErr.czz(); + pcaTrack_x_ = ttPoint.x(); + pcaTrack_y_ = ttPoint.y(); + pcaTrack_z_ = ttPoint.z(); + pcaTrack_xerr_ = ttPointErr.cxx(); + pcaTrack_yerr_ = ttPointErr.cyy(); + pcaTrack_zerr_ = ttPointErr.czz(); + + dotprodTrack_ = dotprodTrack; + dotprodSeed_ = dotprodSeed; + pcaSeed_dist_ = pvToPCAseed; + pcaTrack_dist_ = pvToPCAtrack; + + dotprodTrackSeed2D_ = trackDir2D.unit().dot(seedDir2D.unit()); + dotprodTrackSeed3D_ = iImpactState.globalDirection().unit().dot(tImpactState.globalDirection().unit()); + dotprodTrackSeed2DV_ = trackPCADir2D.unit().dot(seedPCADir2D.unit()); + dotprodTrackSeed3DV_ = (seedPosition - pvp).unit().dot((ttPoint - pvp).unit()); + + pca_jetAxis_dist_ = jetLine.distance(cp).mag(); + pca_jetAxis_dotprod_ = pairMomentum.unit().dot(jetdirection.unit()); + pca_jetAxis_dEta_ = std::fabs(pvToPCA.eta() - jetdirection.eta()); + pca_jetAxis_dPhi_ = std::fabs(pvToPCA.phi() - jetdirection.phi()); + } + } + +} // namespace btagbtvdeep diff --git a/DeepNtuplizer/src/mergeDescriptor.cc b/DeepNtuplizer/src/mergeDescriptor.cc index 61f47edbcb5..374c53e94fb 100644 --- a/DeepNtuplizer/src/mergeDescriptor.cc +++ b/DeepNtuplizer/src/mergeDescriptor.cc @@ -24,8 +24,9 @@ #include "DeepNTuples/DeepNtuplizer/interface/ntuple_JetInfo.h" #include "DeepNTuples/DeepNtuplizer/interface/ntuple_pfCands.h" #include "DeepNTuples/DeepNtuplizer/interface/ntuple_SV.h" -#include "DeepNTuples/DeepNtuplizer/interface/ntuple_DeepVertex.h" -#include "DeepNTuples/DeepNtuplizer/interface/ntuple_GraphB.h" +#include "DeepNTuples/DeepNtuplizer/interface/ntuple_LT.h" +#include "DeepNTuples/DeepNtuplizer/interface/ntuple_pairwise.h" +//#include "DeepNTuples/DeepNtuplizer/interface/ntuple_GraphB.h" static bool debug=true; @@ -131,9 +132,10 @@ std::vector mergeDescriptor::createChains( branchinfos.push_back(new ntuple_JetInfo()); // branchinfos.push_back(new ntuple_FatJetInfo()); branchinfos.push_back(new ntuple_SV()); + branchinfos.push_back(new ntuple_LT()); branchinfos.push_back(new ntuple_bTagVars()); branchinfos.push_back(new ntuple_pfCands()); - //branchinfos.push_back(new ntuple_DeepVertex()); + branchinfos.push_back(new ntuple_pairwise()); //branchinfos.push_back(new ntuple_GraphB()); std::vector chains; diff --git a/DeepNtuplizer/src/ntuple_DeepVertex.cc b/DeepNtuplizer/src/ntuple_DeepVertex.cc deleted file mode 100644 index e6bd3ec74d1..00000000000 --- a/DeepNtuplizer/src/ntuple_DeepVertex.cc +++ /dev/null @@ -1,442 +0,0 @@ -/* - * ntuple_DeepVertex.cc - * - * Created on: 23 June 2017 - * Author: Seth Moortgat - */ - - -#include "../interface/ntuple_DeepVertex.h" - -#include "DataFormats/GeometrySurface/interface/Line.h" - -#include "TrackingTools/Records/interface/TransientTrackRecord.h" -#include "TrackingTools/IPTools/interface/IPTools.h" -#include "TrackingTools/PatternTools/interface/TwoTrackMinimumDistance.h" - -#include "CondFormats/BTauObjects/interface/TrackProbabilityCalibration.h" -#include "CondFormats/DataRecord/interface/BTagTrackProbability2DRcd.h" -#include "CondFormats/DataRecord/interface/BTagTrackProbability3DRcd.h" -#include "FWCore/Framework/interface/EventSetupRecord.h" -#include "FWCore/Framework/interface/EventSetupRecordImplementation.h" -#include "FWCore/Framework/interface/EventSetupRecordKey.h" - -ntuple_DeepVertex::ntuple_DeepVertex(double jetR):ntuple_content(jetR){} - -ntuple_DeepVertex::~ntuple_DeepVertex(){} - -void ntuple_DeepVertex::getInput(const edm::ParameterSet& iConfig){} - -void ntuple_DeepVertex::initBranches(TTree* tree){ - - addBranch(tree,"n_seeds",&n_seeds, "n_seeds/i"); - addBranch(tree,"nSeeds",&nSeeds, "nSeeds/f"); - - addBranch(tree,"seed_pt",&seed_pt, "seed_pt[n_seeds]/f"); - addBranch(tree,"seed_eta",&seed_eta, "seed_eta[n_seeds]/f"); - addBranch(tree,"seed_phi",&seed_phi, "seed_phi[n_seeds]/f"); - addBranch(tree,"seed_mass",&seed_mass, "seed_mass[n_seeds]/f"); - - addBranch(tree,"seed_dz", &seed_dz, "seed_dz[n_seeds]/f"); - addBranch(tree,"seed_dxy", &seed_dxy, "seed_dxy[n_seeds]/f"); - addBranch(tree,"seed_3D_ip", &seed_3D_ip, "seed_3D_ip[n_seeds]/f"); - addBranch(tree,"seed_3D_sip", &seed_3D_sip, "seed_3D_sip[n_seeds]/f"); - addBranch(tree,"seed_2D_ip", &seed_2D_ip, "seed_2D_ip[n_seeds]/f"); - addBranch(tree,"seed_2D_sip", &seed_2D_sip, "seed_2D_sip[n_seeds]/f"); - - addBranch(tree,"seed_3D_signedIp", &seed_3D_signedIp, "seed_3D_signedIp[n_seeds]/f"); - addBranch(tree,"seed_3D_signedSip", &seed_3D_signedSip, "seed_3D_signedSip[n_seeds]/f"); - addBranch(tree,"seed_2D_signedIp", &seed_2D_signedIp, "seed_2D_signedIp[n_seeds]/f"); - addBranch(tree,"seed_2D_signedSip", &seed_2D_signedSip, "seed_2D_signedSip[n_seeds]/f"); - addBranch(tree,"seed_3D_TrackProbability", &seed_3D_TrackProbability, "seed_3D_TrackProbability[n_seeds]/f"); - addBranch(tree,"seed_2D_TrackProbability", &seed_2D_TrackProbability, "seed_2D_TrackProbability[n_seeds]/f"); - - addBranch(tree,"seed_chi2reduced",&seed_chi2reduced, "seed_chi2reduced[n_seeds]/f"); - addBranch(tree,"seed_nPixelHits",&seed_nPixelHits, "seed_nPixelHits[n_seeds]/f"); - addBranch(tree,"seed_nHits",&seed_nHits, "seed_nHits[n_seeds]/f"); - addBranch(tree,"seed_jetAxisDistance",&seed_jetAxisDistance, "seed_jetAxisDistance[n_seeds]/f"); - addBranch(tree,"seed_jetAxisDlength",&seed_jetAxisDlength, "seed_jetAxisDlength[n_seeds]/f"); - - addBranch(tree,"seed_n_NearTracks",&seed_n_NearTracks, "seed_n_NearTracks[n_seeds]/i"); - - // near Tracks - addBranch(tree,"n_NearTracksTotal",&n_NearTracksTotal, "n_NearTracksTotal/i"); - - addBranch(tree,"nearTracks_pt", &nearTracks_pt, "nearTracks_pt[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_eta", &nearTracks_eta, "nearTracks_eta[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_phi", &nearTracks_phi, "nearTracks_phi[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_mass", &nearTracks_mass, "nearTracks_mass[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_dz", &nearTracks_dz, "nearTracks_dz[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_dxy", &nearTracks_dxy, "nearTracks_dxy[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_3D_ip", &nearTracks_3D_ip, "nearTracks_3D_ip[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_3D_sip", &nearTracks_3D_sip, "nearTracks_3D_sip[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_2D_ip", &nearTracks_2D_ip, "nearTracks_2D_ip[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_2D_sip", &nearTracks_2D_sip, "nearTracks_2D_sip[n_NearTracksTotal]/f"); - - addBranch(tree,"nearTracks_PCAdist", &nearTracks_PCAdist, "nearTracks_PCAdist[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAdsig", &nearTracks_PCAdsig, "nearTracks_PCAdsig[n_NearTracksTotal]/f"); - - addBranch(tree,"nearTracks_PCAonSeed_x", &nearTracks_PCAonSeed_x, "nearTracks_PCAonSeed_x[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAonSeed_y", &nearTracks_PCAonSeed_y, "nearTracks_PCAonSeed_y[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAonSeed_z", &nearTracks_PCAonSeed_z, "nearTracks_PCAonSeed_z[n_NearTracksTotal]/f"); - - addBranch(tree,"nearTracks_PCAonSeed_xerr", &nearTracks_PCAonSeed_xerr, "nearTracks_PCAonSeed_xerr[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAonSeed_yerr", &nearTracks_PCAonSeed_yerr, "nearTracks_PCAonSeed_yerr[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAonSeed_zerr", &nearTracks_PCAonSeed_zerr, "nearTracks_PCAonSeed_zerr[n_NearTracksTotal]/f"); - - addBranch(tree,"nearTracks_PCAonTrack_x", &nearTracks_PCAonTrack_x, "nearTracks_PCAonTrack_x[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAonTrack_y", &nearTracks_PCAonTrack_y, "nearTracks_PCAonTrack_y[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAonTrack_z", &nearTracks_PCAonTrack_z, "nearTracks_PCAonTrack_z[n_NearTracksTotal]/f"); - - addBranch(tree,"nearTracks_PCAonTrack_xerr", &nearTracks_PCAonTrack_xerr, "nearTracks_PCAonTrack_xerr[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAonTrack_yerr", &nearTracks_PCAonTrack_yerr, "nearTracks_PCAonTrack_yerr[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAonTrack_zerr", &nearTracks_PCAonTrack_zerr, "nearTracks_PCAonTrack_zerr[n_NearTracksTotal]/f"); - - addBranch(tree,"nearTracks_dotprodTrack", &nearTracks_dotprodTrack, "nearTracks_dotprodTrack[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_dotprodSeed", &nearTracks_dotprodSeed, "nearTracks_dotprodSeed[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_dotprodTrackSeed2D", &nearTracks_dotprodTrackSeed2D, "nearTracks_dotprodTrackSeed2D[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_dotprodTrackSeed3D", &nearTracks_dotprodTrackSeed3D, "nearTracks_dotprodTrackSeed3D[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_dotprodTrackSeedVectors2D", &nearTracks_dotprodTrackSeedVectors2D, "nearTracks_dotprodTrackSeedVectors2D[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_dotprodTrackSeedVectors3D", &nearTracks_dotprodTrackSeedVectors3D, "nearTracks_dotprodTrackSeedVectors3D[n_NearTracksTotal]/f"); - - addBranch(tree,"nearTracks_PCAonSeed_pvd", &nearTracks_PCAonSeed_pvd, "nearTracks_PCAonSeed_pvd[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAonTrack_pvd", &nearTracks_PCAonTrack_pvd, "nearTracks_PCAonTrack_pvd[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAjetAxis_dist",&nearTracks_PCAjetAxis_dist,"nearTracks_PCAjetAxis_dist[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAjetMomenta_dotprod",&nearTracks_PCAjetMomenta_dotprod,"nearTracks_PCAjetMomenta_dotprod[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAjetDirs_DEta",&nearTracks_PCAjetDirs_DEta,"nearTracks_PCAjetDirs_DEta[n_NearTracksTotal]/f"); - addBranch(tree,"nearTracks_PCAjetDirs_DPhi",&nearTracks_PCAjetDirs_DPhi,"nearTracks_PCAjetDirs_DPhi[n_NearTracksTotal]/f"); - - -} - - -void ntuple_DeepVertex::readEvent(const edm::Event& iEvent) -{ - iEvent.getByToken(CandidateToken, tracks); -} - - -void ntuple_DeepVertex::readSetup(const edm::EventSetup& iSetup){ - - //this part was to be in checkEventSetup, but idk how to call it - using namespace edm; - using namespace edm::eventsetup; - - const EventSetupRecord & re2D= iSetup.get(); - const EventSetupRecord & re3D= iSetup.get(); - unsigned long long cacheId2D= re2D.cacheIdentifier(); - unsigned long long cacheId3D= re3D.cacheIdentifier(); - - if(cacheId2D!=m_calibrationCacheId2D || cacheId3D!=m_calibrationCacheId3D ) //Calibration changed - { - //iSetup.get().get(calib); - ESHandle calib2DHandle; - iSetup.get().get(calib2DHandle); - ESHandle calib3DHandle; - iSetup.get().get(calib3DHandle); - - const TrackProbabilityCalibration * ca2D= calib2DHandle.product(); - const TrackProbabilityCalibration * ca3D= calib3DHandle.product(); - - m_probabilityEstimator.reset(new HistogramProbabilityEstimator(ca3D,ca2D)); - - } - - m_calibrationCacheId3D=cacheId3D; - m_calibrationCacheId2D=cacheId2D; - - //readEvent only line - iSetup.get().get("TransientTrackBuilder", builder); - -} - -void ntuple_DeepVertex::checkEventSetup(const edm::EventSetup & iSetup) { - - using namespace edm; - using namespace edm::eventsetup; - - const EventSetupRecord & re2D= iSetup.get(); - const EventSetupRecord & re3D= iSetup.get(); - unsigned long long cacheId2D= re2D.cacheIdentifier(); - unsigned long long cacheId3D= re3D.cacheIdentifier(); - - if(cacheId2D!=m_calibrationCacheId2D || cacheId3D!=m_calibrationCacheId3D ) //Calibration changed - { - //iSetup.get().get(calib); - ESHandle calib2DHandle; - iSetup.get().get(calib2DHandle); - ESHandle calib3DHandle; - iSetup.get().get(calib3DHandle); - - const TrackProbabilityCalibration * ca2D= calib2DHandle.product(); - const TrackProbabilityCalibration * ca3D= calib3DHandle.product(); - - m_probabilityEstimator.reset(new HistogramProbabilityEstimator(ca3D,ca2D)); - - } - - m_calibrationCacheId3D=cacheId3D; - m_calibrationCacheId2D=cacheId2D; - -} - - -bool ntuple_DeepVertex::fillBranches(const pat::Jet & jet, const size_t& jetidx, const edm::View * coll){ - - // pv info - const reco::Vertex &pv = vertices()->at(0); - GlobalPoint pvp(pv.x(),pv.y(),pv.z()); - - - std::vector selectedTracks; - std::vector masses; - - - for(size_t k = 0; ksize(); ++k) { - if((*tracks)[k].bestTrack() != 0 && (*tracks)[k].pt()>0.5 && std::fabs(pvp.z()-builder->build(tracks->ptrAt(k)).track().vz())<0.5) { - selectedTracks.push_back(builder->build(tracks->ptrAt(k))); - masses.push_back(tracks->ptrAt(k)->mass()); - } - } - - double jet_radius = jetR(); - GlobalVector direction(jet.px(), jet.py(), jet.pz()); - - for(std::vector::const_iterator it = selectedTracks.begin(); it != selectedTracks.end(); it++){ - - //is the track in the jet cone? - float angular_distance=reco::deltaR(jet,it->track());//std::sqrt(std::pow(jet.eta()-it->track().eta(),2) + std::pow(jet.phi()-it->track().phi(),2) ); - if (angular_distance>jet_radius) { continue; } - - // is it a seed track? - std::pair ip = IPTools::absoluteImpactParameter3D(*it, pv); - std::pair ip2d = IPTools::absoluteTransverseImpactParameter(*it, pv); - std::pair jet_dist =IPTools::jetTrackDistance(*it, direction, pv); - TrajectoryStateOnSurface closest = IPTools::closestApproachToJet(it->impactPointState(),pv, direction,it->field()); - float length=999; - if (closest.isValid()) length=(closest.globalPosition() - pvp).mag(); - - // shouldn't it be like this, including the minimal 3DIP cuts? more conform with IVF! - // https://github.com/cms-sw/cmssw/blob/09c3fce6626f70fd04223e7dacebf0b485f73f54/RecoVertex/AdaptiveVertexFinder/src/TracksClusteringFromDisplacedSeed.cc#L96 - //maybe... but qulity cuts for now - - bool is_seed_candidate = (ip.first && ip.second.value() >= 0.0 && ip.second.significance() >= 1.0 && - ip.second.value() <= max3DIPValue && ip.second.significance() <= max3DIPSignificance && - it->track().normalizedChi2()<5. && std::fabs(it->track().dxy(pv.position())) < 2 && - std::fabs(it->track().dz(pv.position())) < 17 && jet_dist.second.value()<0.07 && length<5. ); - - //if (!is_seed_candidate){continue;} - - std::pair ipSigned = IPTools::signedImpactParameter3D(*it,direction, pv); - //n_seeds++; - - nearTracks.clear(); - //now that we found a seed, loop over all other tracks and look for neighbours - for(std::vector::const_iterator tt = selectedTracks.begin();tt!=selectedTracks.end(); ++tt ) { - float near_angular_distance=reco::deltaR(jet,tt->track()); - if(near_angular_distancetrack().vz())>0.1) continue; - - VertexDistance3D distanceComputer; - TwoTrackMinimumDistance dist; - - if(dist.calculate(tt->impactPointState(),it->impactPointState())) { - GlobalPoint ttPoint = dist.points().first; - GlobalError ttPointErr = tt->impactPointState().cartesianError().position(); - GlobalPoint seedPosition = dist.points().second; - GlobalError seedPositionErr = it->impactPointState().cartesianError().position(); - - Measurement1D m = distanceComputer.distance(VertexState(seedPosition,seedPositionErr), VertexState(ttPoint, ttPointErr)); - GlobalPoint cp(dist.crossingPoint()); - - GlobalVector PairMomentum(it->track().px()+tt->track().px(), it->track().py()+tt->track().py(), it->track().pz()+tt->track().pz()); - GlobalVector PCA_pv(cp-pvp); - - float PCAseedFromPV = (dist.points().second-pvp).mag(); - float PCAtrackFromPV = (dist.points().first-pvp).mag(); - float distance = dist.distance(); - - GlobalVector trackDir2D(tt->impactPointState().globalDirection().x(),tt->impactPointState().globalDirection().y(),0.); - GlobalVector seedDir2D(it->impactPointState().globalDirection().x(),it->impactPointState().globalDirection().y(),0.); - GlobalVector trackPCADir2D(dist.points().first.x()-pvp.x(),dist.points().first.y()-pvp.y(),0.); - GlobalVector seedPCADir2D(dist.points().second.x()-pvp.x(),dist.points().second.y()-pvp.y(),0.); - - float dotprodTrack = (dist.points().first-pvp).unit().dot(tt->impactPointState().globalDirection().unit()); - float dotprodSeed = (dist.points().second-pvp).unit().dot(it->impactPointState().globalDirection().unit()); - float dotprodTrackSeed2D = trackDir2D.unit().dot(seedDir2D.unit()); - float dotprodTrackSeed3D = it->impactPointState().globalDirection().unit().dot(tt->impactPointState().globalDirection().unit()); - float dotprodTrackSeed2DV = trackPCADir2D.unit().dot(seedPCADir2D.unit()); - float dotprodTrackSeed3DV = (dist.points().second-pvp).unit().dot((dist.points().first-pvp).unit()); - - std::pair t_ip = IPTools::absoluteImpactParameter3D(*tt,pv); - std::pair t_ip2d = IPTools::absoluteTransverseImpactParameter(*tt,pv); - - myTrack.set_PtEtaPhiMassDzDxy(tt->track().pt(), tt->track().eta(), tt->track().phi(), masses[tt-selectedTracks.begin()], - tt->track().dz(pv.position()), tt->track().dxy(pv.position())); - - myTrack.set_IPs(t_ip2d.second.value() , t_ip2d.second.significance(), t_ip.second.value() , t_ip.second.significance()); - - myTrack.set_PCAdistance (distance, m.significance()); - myTrack.set_PCAonSeedXYZ(seedPosition.x(), seedPosition.y(), seedPosition.z(), seedPositionErr.cxx(), seedPositionErr.cyy(), seedPositionErr.czz()); - myTrack.set_PCAonTrackXYZ(ttPoint.x(), ttPoint.y(), ttPoint.z(), ttPointErr.cxx(), ttPointErr.cyy(), ttPointErr.czz()); - - myTrack.set_dotProds(dotprodTrack, dotprodSeed, dotprodTrackSeed2D, dotprodTrackSeed3D, dotprodTrackSeed2DV, dotprodTrackSeed3DV); - myTrack.set_PVdistance(PCAseedFromPV, PCAtrackFromPV); - - Line::PositionType pos(pvp); - Line::DirectionType dir(direction); - Line::DirectionType pairMomentumDir(PairMomentum); - Line jetLine(pos,dir); - Line PCAMomentumLine(cp,pairMomentumDir); - float PCA_JetAxis_dist=jetLine.distance(cp).mag(); - float dotprodMomenta=PairMomentum.unit().dot(direction.unit()); - float dEta=std::fabs(PCA_pv.eta()-jet.eta()); - float dPhi=std::fabs(PCA_pv.phi()-jet.phi()); - - myTrack.setSeedMass(masses[it-selectedTracks.begin()]); - myTrack.set_JetAxisVars(PCA_JetAxis_dist,dotprodMomenta,dEta,dPhi); - nearTracks.push_back(myTrack); - - } - } - - std::sort (nearTracks.begin(), nearTracks.end(), sortfunctionNTracks()); - if (nearTracks.size() > 20){nearTracks.resize(20);} - //for(int n = 0; n<20; n++){ - //std::cout << "n = " << n << std::endl; - //std::cout << "dR = " << reco::deltaR( jet.eta(), jet.phi(), nearTracks.at(n).eta, nearTracks.at(n).phi) << std::endl; - //std::cout << "dist = " << nearTracks.at(n).dist << std::endl; - - //} - SortedSeedsMap.insert(std::make_pair(-ipSigned.second.significance(), std::make_pair(&(*it), nearTracks))); - - } - // std::cout << "jet tracks = " << counter << std::endl; - //std::cout << "ang dist 4 = " << big_counter << std::endl; - unsigned int seeds_max_counter=0; - unsigned int neartracks_max_counter=0; - - for(std::multimap > >::const_iterator im = SortedSeedsMap.begin(); im != SortedSeedsMap.end(); im++){ - - if(seeds_max_counter>=max_seeds) break; - - std::pair ipSigned = IPTools::signedImpactParameter3D(*im->second.first,direction, pv); - std::pair ip2dSigned = IPTools::signedTransverseImpactParameter(*im->second.first,direction, pv); - std::pair ip = IPTools::absoluteImpactParameter3D(*im->second.first, pv); - std::pair ip2d = IPTools::absoluteTransverseImpactParameter(*im->second.first, pv); - - seed_pt[seeds_max_counter]=im->second.first->track().pt(); - seed_eta[seeds_max_counter]=im->second.first->track().eta(); - seed_phi[seeds_max_counter]=im->second.first->track().phi(); - if (im->second.second.size() != 0){ seed_mass[seeds_max_counter]=im->second.second.at(0).seedMass;} - else{seed_mass[seeds_max_counter]=-1;} - seed_dz[seeds_max_counter]=im->second.first->track().dz(pv.position()); - seed_dxy[seeds_max_counter]=im->second.first->track().dxy(pv.position()); - seed_3D_ip[seeds_max_counter]=ip.second.value(); - seed_3D_sip[seeds_max_counter]=ip.second.significance(); - seed_2D_ip[seeds_max_counter]=ip2d.second.value(); - seed_2D_sip[seeds_max_counter]=ip2d.second.significance(); - seed_3D_signedIp[seeds_max_counter]=ipSigned.second.value(); - seed_3D_signedSip[seeds_max_counter]=ipSigned.second.significance(); - seed_2D_signedIp[seeds_max_counter]=ip2dSigned.second.value(); - seed_2D_signedSip[seeds_max_counter]=ip2dSigned.second.significance(); - - seed_3D_TrackProbability[seeds_max_counter]=1; - seed_2D_TrackProbability[seeds_max_counter]=1; - - if (m_computeProbabilities) { - - std::pair probability; - - //probability with 3D ip - probability = m_probabilityEstimator->probability(0, 0,ip.second.significance(),im->second.first->track(),jet,pv); - double prob3D=(probability.first ? probability.second : -1.); - - //probability with 2D ip - probability = m_probabilityEstimator->probability(0, 1,ip2d.second.significance(),im->second.first->track(),jet,pv); - double prob2D=(probability.first ? probability.second : -1.); - - seed_3D_TrackProbability[seeds_max_counter]=prob3D; - seed_2D_TrackProbability[seeds_max_counter]=prob2D; - - } - - seed_chi2reduced[seeds_max_counter]=im->second.first->track().normalizedChi2(); - seed_nPixelHits[seeds_max_counter]=im->second.first->track().hitPattern().numberOfValidPixelHits(); - seed_nHits[seeds_max_counter]=im->second.first->track().hitPattern().numberOfValidHits(); - - std::pair jet_distance =IPTools::jetTrackDistance(*im->second.first, direction, pv); - seed_jetAxisDistance[seeds_max_counter]=std::fabs(jet_distance.second.value()); - - TrajectoryStateOnSurface closest = IPTools::closestApproachToJet(im->second.first->impactPointState(),pv, direction,im->second.first->field()); - if (closest.isValid()) seed_jetAxisDlength[seeds_max_counter]=(closest.globalPosition() - pvp).mag(); - else seed_jetAxisDlength[seeds_max_counter]= -99; - - seed_n_NearTracks[seeds_max_counter]=im->second.second.size(); - - // FILL NEAREAST VARIABLES - for(unsigned int i=0; i< im->second.second.size(); i++) { - - if((neartracks_max_counter+i)>=max_nearestTrk) break; - - nearTracks_pt[neartracks_max_counter+i]=im->second.second.at(i).pt; - nearTracks_eta[neartracks_max_counter+i]=im->second.second.at(i).eta; - nearTracks_phi[neartracks_max_counter+i]=im->second.second.at(i).phi; - nearTracks_dz[neartracks_max_counter+i]=im->second.second.at(i).dz; - nearTracks_dxy[neartracks_max_counter+i]=im->second.second.at(i).dxy; - nearTracks_mass[neartracks_max_counter+i]=im->second.second.at(i).mass; - nearTracks_3D_ip[neartracks_max_counter+i]=im->second.second.at(i).t3Dip; - nearTracks_3D_sip[neartracks_max_counter+i]=im->second.second.at(i).t3Dsip; - nearTracks_2D_ip[neartracks_max_counter+i]=im->second.second.at(i).t2Dip; - nearTracks_2D_sip[neartracks_max_counter+i]=im->second.second.at(i).t2Dsip; - nearTracks_PCAdist[neartracks_max_counter+i]=im->second.second.at(i).dist; - nearTracks_PCAdsig[neartracks_max_counter+i]=im->second.second.at(i).dsig; - nearTracks_PCAonSeed_x[neartracks_max_counter+i]=im->second.second.at(i).PCA_sx; - nearTracks_PCAonSeed_y[neartracks_max_counter+i]=im->second.second.at(i).PCA_sy; - nearTracks_PCAonSeed_z[neartracks_max_counter+i]=im->second.second.at(i).PCA_sz; - nearTracks_PCAonSeed_xerr[neartracks_max_counter+i]=im->second.second.at(i).PCA_sxerr; - nearTracks_PCAonSeed_yerr[neartracks_max_counter+i]=im->second.second.at(i).PCA_syerr; - nearTracks_PCAonSeed_zerr[neartracks_max_counter+i]=im->second.second.at(i).PCA_szerr; - nearTracks_PCAonTrack_x[neartracks_max_counter+i]=im->second.second.at(i).PCA_tx; - nearTracks_PCAonTrack_y[neartracks_max_counter+i]=im->second.second.at(i).PCA_ty; - nearTracks_PCAonTrack_z[neartracks_max_counter+i]=im->second.second.at(i).PCA_tz; - nearTracks_PCAonTrack_xerr[neartracks_max_counter+i]=im->second.second.at(i).PCA_txerr; - nearTracks_PCAonTrack_yerr[neartracks_max_counter+i]=im->second.second.at(i).PCA_tyerr; - nearTracks_PCAonTrack_zerr[neartracks_max_counter+i]=im->second.second.at(i).PCA_tzerr; - nearTracks_dotprodTrack[neartracks_max_counter+i]=im->second.second.at(i).dotprodTrack; - nearTracks_dotprodSeed[neartracks_max_counter+i]=im->second.second.at(i).dotprodSeed; - nearTracks_dotprodTrackSeed2D[neartracks_max_counter+i]=im->second.second.at(i).dotprodTrackSeed2D; - nearTracks_dotprodTrackSeed3D[neartracks_max_counter+i]=im->second.second.at(i).dotprodTrackSeed3D; - nearTracks_dotprodTrackSeedVectors2D[neartracks_max_counter+i]=im->second.second.at(i).dotprodTrackSeedVectors2D; - nearTracks_dotprodTrackSeedVectors3D[neartracks_max_counter+i]=im->second.second.at(i).dotprodTrackSeedVectors3D; - - nearTracks_PCAonSeed_pvd[neartracks_max_counter+i]=im->second.second.at(i).seedPCA_pv; - nearTracks_PCAonTrack_pvd[neartracks_max_counter+i]=im->second.second.at(i).trackPCA_pv; - - nearTracks_PCAjetAxis_dist[neartracks_max_counter+i]=im->second.second.at(i).PCA_JetAxis_distance; - nearTracks_PCAjetMomenta_dotprod[neartracks_max_counter+i]=im->second.second.at(i).PCAPair_Jet_dotprod; - - nearTracks_PCAjetDirs_DEta[neartracks_max_counter+i]=im->second.second.at(i).PCAAxis_JetAxis_DEta; - nearTracks_PCAjetDirs_DPhi[neartracks_max_counter+i]=im->second.second.at(i).PCAAxis_JetAxis_DPhi; - - } - - // ********* - neartracks_max_counter += im->second.second.size(); - seeds_max_counter++; - } - - n_NearTracksTotal = neartracks_max_counter; - n_seeds = seeds_max_counter; - - SortedSeedsMap.clear(); - nearTracks.clear(); - masses.clear(); - - return true; -} - diff --git a/DeepNtuplizer/src/ntuple_GraphB.cc b/DeepNtuplizer/src/ntuple_GraphB.cc deleted file mode 100644 index 0464b50ea4d..00000000000 --- a/DeepNtuplizer/src/ntuple_GraphB.cc +++ /dev/null @@ -1,260 +0,0 @@ -/* - * ntuple_GraphB.cc - * - * Created on: 23 June 2017 - * Author: Seth Moortgat - - */ - - -#include "../interface/ntuple_GraphB.h" - -#include "DataFormats/GeometrySurface/interface/Line.h" - -#include "TrackingTools/Records/interface/TransientTrackRecord.h" -#include "TrackingTools/IPTools/interface/IPTools.h" -#include "TrackingTools/PatternTools/interface/TwoTrackMinimumDistance.h" -#include "../interface/sorting_modules.h" -#include "CondFormats/BTauObjects/interface/TrackProbabilityCalibration.h" -#include "CondFormats/DataRecord/interface/BTagTrackProbability2DRcd.h" -#include "CondFormats/DataRecord/interface/BTagTrackProbability3DRcd.h" -#include "FWCore/Framework/interface/EventSetupRecord.h" -#include "FWCore/Framework/interface/EventSetupRecordImplementation.h" -#include "FWCore/Framework/interface/EventSetupRecordKey.h" - -ntuple_GraphB::ntuple_GraphB(double jetR):ntuple_content(jetR){} - -ntuple_GraphB::~ntuple_GraphB(){} - -void ntuple_GraphB::getInput(const edm::ParameterSet& iConfig){} - -void ntuple_GraphB::initBranches(TTree* tree){ - - addBranch(tree,"n_gtracks",&n_gtracks, "n_gtracks/i"); - addBranch(tree,"nGtracks",&nGtracks, "nGtracks/f"); - - addBranch(tree,"gtrack_pt",>rack_pt, "gtrack_pt[n_gtracks]/f"); - addBranch(tree,"gtrack_eta",>rack_eta, "gtrack_eta[n_gtracks]/f"); - addBranch(tree,"gtrack_phi",>rack_phi, "gtrack_phi[n_gtracks]/f"); - addBranch(tree,"gtrack_mass",>rack_mass, "gtrack_mass[n_gtracks]/f"); - - addBranch(tree,"gtrack_dz", >rack_dz, "gtrack_dz[n_gtracks]/f"); - addBranch(tree,"gtrack_dxy", >rack_dxy, "gtrack_dxy[n_gtracks]/f"); - addBranch(tree,"gtrack_3D_ip", >rack_3D_ip, "gtrack_3D_ip[n_gtracks]/f"); - addBranch(tree,"gtrack_3D_sip", >rack_3D_sip, "gtrack_3D_sip[n_gtracks]/f"); - addBranch(tree,"gtrack_2D_ip", >rack_2D_ip, "gtrack_2D_ip[n_gtracks]/f"); - addBranch(tree,"gtrack_2D_sip", >rack_2D_sip, "gtrack_2D_sip[n_gtracks]/f"); - addBranch(tree,"gtrack_dR", >rack_dR, "gtrack_dR[n_gtracks]/f"); - - addBranch(tree,"gtrack_3D_TrackProbability", >rack_3D_TrackProbability, "gtrack_3D_TrackProbability[n_gtracks]/f"); - addBranch(tree,"gtrack_2D_TrackProbability", >rack_2D_TrackProbability, "gtrack_2D_TrackProbability[n_gtracks]/f"); - - addBranch(tree,"gtrack_chi2reduced",>rack_chi2reduced, "gtrack_chi2reduced[n_gtracks]/f"); - addBranch(tree,"gtrack_nPixelHits",>rack_nPixelHits, "gtrack_nPixelHits[n_gtracks]/f"); - addBranch(tree,"gtrack_nHits",>rack_nHits, "gtrack_nHits[n_gtracks]/f"); - addBranch(tree,"gtrack_jetAxisDistance",>rack_jetAxisDistance, "gtrack_jetAxisDistance[n_gtracks]/f"); - addBranch(tree,"gtrack_jetAxisDlength",>rack_jetAxisDlength, "gtrack_jetAxisDlength[n_gtracks]/f"); - addBranch(tree,"gtrack_PCAtrackFromPV",>rack_PCAtrackFromPV, "gtrack_PCAtrackFromPV[n_gtracks]/f"); - addBranch(tree,"gtrack_dotProdTrack",>rack_dotProdTrack, "gtrack_dotProdTrack[n_gtracks]/f"); - addBranch(tree,"gtrack_dotProdTrack2D",>rack_dotProdTrack2D, "gtrack_dotProdTrack2D[n_gtracks]/f"); - - -} - - -void ntuple_GraphB::readEvent(const edm::Event& iEvent) -{ - iEvent.getByToken(CandidateToken, tracks); -} - - -void ntuple_GraphB::readSetup(const edm::EventSetup& iSetup){ - - //this part was to be in checkEventSetup, but idk how to call it - using namespace edm; - using namespace edm::eventsetup; - - const EventSetupRecord & re2D= iSetup.get(); - const EventSetupRecord & re3D= iSetup.get(); - unsigned long long cacheId2D= re2D.cacheIdentifier(); - unsigned long long cacheId3D= re3D.cacheIdentifier(); - - if(cacheId2D!=m_calibrationCacheId2D || cacheId3D!=m_calibrationCacheId3D ) //Calibration changed - { - //iSetup.get().get(calib); - ESHandle calib2DHandle; - iSetup.get().get(calib2DHandle); - ESHandle calib3DHandle; - iSetup.get().get(calib3DHandle); - - const TrackProbabilityCalibration * ca2D= calib2DHandle.product(); - const TrackProbabilityCalibration * ca3D= calib3DHandle.product(); - - m_probabilityEstimator.reset(new HistogramProbabilityEstimator(ca3D,ca2D)); - - } - - m_calibrationCacheId3D=cacheId3D; - m_calibrationCacheId2D=cacheId2D; - - //readEvent only line - iSetup.get().get("TransientTrackBuilder", builder); - -} - -void ntuple_GraphB::checkEventSetup(const edm::EventSetup & iSetup) { - - using namespace edm; - using namespace edm::eventsetup; - - const EventSetupRecord & re2D= iSetup.get(); - const EventSetupRecord & re3D= iSetup.get(); - unsigned long long cacheId2D= re2D.cacheIdentifier(); - unsigned long long cacheId3D= re3D.cacheIdentifier(); - - if(cacheId2D!=m_calibrationCacheId2D || cacheId3D!=m_calibrationCacheId3D ) //Calibration changed - { - //iSetup.get().get(calib); - ESHandle calib2DHandle; - iSetup.get().get(calib2DHandle); - ESHandle calib3DHandle; - iSetup.get().get(calib3DHandle); - - const TrackProbabilityCalibration * ca2D= calib2DHandle.product(); - const TrackProbabilityCalibration * ca3D= calib3DHandle.product(); - - m_probabilityEstimator.reset(new HistogramProbabilityEstimator(ca3D,ca2D)); - - } - - m_calibrationCacheId3D=cacheId3D; - m_calibrationCacheId2D=cacheId2D; - -} - - -bool ntuple_GraphB::fillBranches(const pat::Jet & jet, const size_t& jetidx, const edm::View * coll){ - - // pv info - const reco::Vertex &pv = vertices()->at(0); - GlobalPoint pvp(pv.x(),pv.y(),pv.z()); - - - std::vector selectedTracks; - std::vector masses; - - - for(size_t k = 0; ksize(); ++k) { - if((*tracks)[k].bestTrack() != 0 && (*tracks)[k].pt()>0.5 && std::fabs(pvp.z()-builder->build(tracks->ptrAt(k)).track().vz())<0.5) { - selectedTracks.push_back(builder->build(tracks->ptrAt(k))); - masses.push_back(tracks->ptrAt(k)->mass()); - } - } - - double jet_radius = jetR(); - GlobalVector direction(jet.px(), jet.py(), jet.pz()); - std::vector > sorted_tracks; - for(std::vector::const_iterator it = selectedTracks.begin(); it != selectedTracks.end(); it++){ - float angular_distance=reco::deltaR(jet,it->track()); - sorted_tracks.push_back(sorting::sortingClass - (it-selectedTracks.begin(), -angular_distance, - -1, -1)); - } - std::sort(sorted_tracks.begin(),sorted_tracks.end(),sorting::sortingClass::compareByABCInv); - std::vector sorted_track_indices; - sorted_track_indices=sorting::invertSortingVector(sorted_tracks); - //n_gtracks = std::min(sorted_tracks.size(),max_gtracks); - //nGtracks = n_gtracks; - - size_t counter= 0; - for(std::vector::const_iterator it = selectedTracks.begin(); it != selectedTracks.end(); it++){ - //is the track in the jet cone? - float angular_distance=reco::deltaR(jet,it->track()); - //std::sqrt(std::pow(jet.eta()-it->track().eta(),2) + std::pow(jet.phi()-it->track().phi(),2) ); - bool hasNeighbour = false; - bool include = false; - if (angular_distance>jet_radius) { - if(std::fabs(pvp.z()-it->track().vz())>0.1) continue; - for(std::vector::const_iterator tt = selectedTracks.begin();tt!=selectedTracks.end(); ++tt ) { - TwoTrackMinimumDistance dist; - float near_angular_distance=reco::deltaR(jet,tt->track()); - if(*tt==*it) continue; - if(near_angular_distance ip = IPTools::absoluteImpactParameter3D(*tt, pv); - if(ip.second.significance() < 1.0) continue; - if(dist.calculate(tt->impactPointState(),it->impactPointState())) { - float distance = dist.distance(); - if(distance < 0.02){ - hasNeighbour = true; - } - } - } - } - if(hasNeighbour){ - int matched_jets = 0; - for (std::size_t jet_n = 0; jet_n < coll->size(); jet_n++) { - const auto& test_jet = coll->at(jet_n); - if(test_jet.pt() < 5.0){continue;} - float new_angular_distance=reco::deltaR(test_jet,it->track()); - if(new_angular_distance < jet_radius){ - matched_jets = matched_jets + 1; - } - } - if(matched_jets != 0){continue;} - include = true; - } - } - if(angular_distance ip = IPTools::signedImpactParameter3D(*it, direction, pv); - std::pair ip2d = IPTools::signedTransverseImpactParameter(*it, direction, pv); - std::pair jet_dist =IPTools::jetTrackDistance(*it, direction, pv); - TrajectoryStateOnSurface closest = IPTools::closestApproachToJet(it->impactPointState(),pv, direction,it->field()); - if(counter >= max_gtracks){continue;} - gtrack_dR[counter] = catchInfsAndBound( reco::deltaR(jet,it->track()), -1.0,0.0,10.0); - float length=999; - if (closest.isValid()) length=(closest.globalPosition() - pvp).mag(); - gtrack_jetAxisDlength[counter] = catchInfsAndBound(length,-1.0,-100.0,100.0); - gtrack_jetAxisDistance[counter] = catchInfsAndBound(jet_dist.second.value(),-1.0,-100.0,100.0); - gtrack_chi2reduced[counter] = catchInfsAndBound(it->track().normalizedChi2(),-1.0,-100.0,300.0); - gtrack_3D_ip[counter] = catchInfsAndBound(ip.second.value(),0.0,-500.0,500.0); - gtrack_3D_sip[counter] = catchInfsAndBound(ip.second.significance(),0.0,-500.0,500.0); - gtrack_2D_ip[counter] = catchInfsAndBound(ip2d.second.value(),0.0,-500.0,500.0); - gtrack_2D_sip[counter] = catchInfsAndBound(ip2d.second.significance(),0.0,-500.0,500.0); - gtrack_pt[counter] = catchInfsAndBound(it->track().pt(),0.0,-100.0,1000.0); - gtrack_eta[counter] = catchInfsAndBound(it->track().eta(),0.0,-2.5,2.5); - gtrack_phi[counter] = catchInfsAndBound(it->track().phi(),0.0,-5,5); - gtrack_mass[counter] = catchInfsAndBound(masses[it-selectedTracks.begin()],0.0,-1.0,500.0); - gtrack_dxy[counter] = catchInfsAndBound(it->track().dxy(pv.position()),0.0,-100.0,100.0); - gtrack_dz[counter] = catchInfsAndBound(it->track().dz(pv.position()),0.0,-100.0,100.0); - gtrack_PCAtrackFromPV[counter] = catchInfsAndBound( (it->impactPointState().globalPosition()-pvp).mag() , -1.0,-100.0,100.0); - gtrack_dotProdTrack[counter] = catchInfsAndBound( (it->impactPointState().globalPosition()-pvp).unit().dot(it->impactPointState().globalDirection().unit()), -3.0,-3.0,3.0); - GlobalVector trackDir2D(it->impactPointState().globalDirection().x(),it->impactPointState().globalDirection().y(),0.); - GlobalPoint pvp2d(pv.x(),pv.y(),0.0); - GlobalPoint trackPos2D(it->impactPointState().globalPosition().x(),it->impactPointState().globalPosition().y(),0.0); - gtrack_dotProdTrack2D[counter] = catchInfsAndBound( (trackPos2D-pvp2d).unit().dot(trackDir2D.unit()), -3.0,-3.0,3.0); - std::pair probability; - //probability with 3D ip - probability = m_probabilityEstimator->probability(0, 0,ip.second.significance(),it->track(),jet,pv); - double prob3D=(probability.first ? probability.second : -1.); - gtrack_3D_TrackProbability[counter] = catchInfsAndBound(prob3D,-2.0,-100.0,100.0); - //probability with 2D ip - probability = m_probabilityEstimator->probability(0, 1,ip2d.second.significance(),it->track(),jet,pv); - double prob2D=(probability.first ? probability.second : -1.); - gtrack_2D_TrackProbability[counter] = catchInfsAndBound(prob2D,-2.0,-100.0,100.0); - - gtrack_nPixelHits[counter] = catchInfsAndBound(it->hitPattern().numberOfValidPixelHits(),-1.0,0.0,100.0); - gtrack_nHits[counter] = catchInfsAndBound(it->hitPattern().numberOfValidHits(),-1.0,0.0,100.0); - counter++; - } - } - n_gtracks = counter; - nGtracks = n_gtracks; - - masses.clear(); - - return true; -} - diff --git a/DeepNtuplizer/src/ntuple_LT.cc b/DeepNtuplizer/src/ntuple_LT.cc new file mode 100644 index 00000000000..44fe63c7ad0 --- /dev/null +++ b/DeepNtuplizer/src/ntuple_LT.cc @@ -0,0 +1,617 @@ +/* + * ntuple_LT.cc + * + * Created on: 28 Mai 2023 + * Author: Alexandre De Moor + */ + + +#include "../interface/ntuple_LT.h" +#include "DataFormats/PatCandidates/interface/PackedCandidate.h" +#include "DataFormats/Candidate/interface/VertexCompositePtrCandidate.h" +#include "../interface/sorting_modules.h" + + +#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h" +#include "TrackingTools/IPTools/interface/IPTools.h" +#include "TrackingTools/Records/interface/TransientTrackRecord.h" +#include "RecoVertex/VertexTools/interface/VertexDistance3D.h" +#include "TVector3.h" + +#include "DataFormats/GeometrySurface/interface/Line.h" +#include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h" +#include "TrackingTools/GeomPropagators/interface/AnalyticalTrajectoryExtrapolatorToLine.h" +#include "TrackingTools/GeomPropagators/interface/AnalyticalImpactPointExtrapolator.h" +#include "DataFormats/GeometryCommonDetAlgo/interface/Measurement1D.h" +#include "RecoVertex/VertexTools/interface/VertexDistanceXY.h" +#include "RecoVertex/VertexPrimitives/interface/ConvertToFromReco.h" + +class TrackInfoBuilder{ +public: + TrackInfoBuilder(edm::ESHandle & build): + builder(build), + trackMomentum_(0), + trackEta_(0), + trackEtaRel_(0), + trackPtRel_(0), + trackPPar_(0), + trackDeltaR_(0), + trackPtRatio_(0), + trackPParRatio_(0), + trackSip2dVal_(0), + trackSip2dSig_(0), + trackSip3dVal_(0), + trackSip3dSig_(0), + + trackJetDecayLen_(0), + trackJetDistVal_(0), + trackJetDistSig_(0), + ttrack_(0) + { + + + } + + void buildTrackInfo(const pat::PackedCandidate* PackedCandidate_ ,const math::XYZVector& jetDir, GlobalVector refjetdirection, const reco::Vertex & pv){ + TVector3 jetDir3(jetDir.x(),jetDir.y(),jetDir.z()); + if(!PackedCandidate_->hasTrackDetails()) { + TVector3 trackMom3( + PackedCandidate_->momentum().x(), + PackedCandidate_->momentum().y(), + PackedCandidate_->momentum().z() + ); + trackMomentum_=PackedCandidate_->p(); + trackEta_= PackedCandidate_->eta(); + trackEtaRel_=reco::btau::etaRel(jetDir, PackedCandidate_->momentum()); + trackPtRel_=trackMom3.Perp(jetDir3); + trackPPar_=jetDir.Dot(PackedCandidate_->momentum()); + trackDeltaR_=reco::deltaR(PackedCandidate_->momentum(), jetDir); + trackPtRatio_=trackMom3.Perp(jetDir3) / PackedCandidate_->p(); + trackPParRatio_=jetDir.Dot(PackedCandidate_->momentum()) / PackedCandidate_->p(); + trackSip2dVal_=0.; + trackSip2dSig_=0.; + trackSip3dVal_=0.; + trackSip3dSig_=0.; + trackJetDecayLen_=0.; + trackJetDistVal_=0.; + trackJetDistSig_=0.; + return; + } + + const reco::Track & PseudoTrack = PackedCandidate_->pseudoTrack(); + + reco::TransientTrack transientTrack; + transientTrack=builder->build(PseudoTrack); + Measurement1D meas_ip2d=IPTools::signedTransverseImpactParameter(transientTrack, refjetdirection, pv).second; + Measurement1D meas_ip3d=IPTools::signedImpactParameter3D(transientTrack, refjetdirection, pv).second; + Measurement1D jetdist=IPTools::jetTrackDistance(transientTrack, refjetdirection, pv).second; + Measurement1D decayl = IPTools::signedDecayLength3D(transientTrack, refjetdirection, pv).second; + math::XYZVector trackMom = PseudoTrack.momentum(); + double trackMag = std::sqrt(trackMom.Mag2()); + TVector3 trackMom3(trackMom.x(),trackMom.y(),trackMom.z()); + + + trackMomentum_=std::sqrt(trackMom.Mag2()); + trackEta_= trackMom.Eta(); + trackEtaRel_=reco::btau::etaRel(jetDir, trackMom); + trackPtRel_=trackMom3.Perp(jetDir3); + trackPPar_=jetDir.Dot(trackMom); + trackDeltaR_=reco::deltaR(trackMom, jetDir); + trackPtRatio_=trackMom3.Perp(jetDir3) / trackMag; + trackPParRatio_=jetDir.Dot(trackMom) / trackMag; + + trackSip2dVal_=(meas_ip2d.value()); + trackSip2dSig_=(meas_ip2d.significance()); + trackSip3dVal_=(meas_ip3d.value()); + trackSip3dSig_=meas_ip3d.significance(); + + trackJetDecayLen_= decayl.value(); + trackJetDistVal_= jetdist.value(); + trackJetDistSig_= jetdist.significance(); + + ttrack_ = transientTrack; + + } + + const float& getTrackDeltaR() const {return trackDeltaR_;} + const float& getTrackEta() const {return trackEta_;} + const float& getTrackEtaRel() const {return trackEtaRel_;} + const float& getTrackJetDecayLen() const {return trackJetDecayLen_;} + const float& getTrackJetDistSig() const {return trackJetDistSig_;} + const float& getTrackJetDistVal() const {return trackJetDistVal_;} + const float& getTrackMomentum() const {return trackMomentum_;} + const float& getTrackPPar() const {return trackPPar_;} + const float& getTrackPParRatio() const {return trackPParRatio_;} + const float& getTrackPtRatio() const {return trackPtRatio_;} + const float& getTrackPtRel() const {return trackPtRel_;} + const float& getTrackSip2dSig() const {return trackSip2dSig_;} + const float& getTrackSip2dVal() const {return trackSip2dVal_;} + const float& getTrackSip3dSig() const {return trackSip3dSig_;} + const float& getTrackSip3dVal() const {return trackSip3dVal_;} + const reco::TransientTrack getTTrack() const {return ttrack_;} + +private: + + edm::ESHandle& builder; + + float trackMomentum_; + float trackEta_; + float trackEtaRel_; + float trackPtRel_; + float trackPPar_; + float trackDeltaR_; + float trackPtRatio_; + float trackPParRatio_; + float trackSip2dVal_; + float trackSip2dSig_; + float trackSip3dVal_; + float trackSip3dSig_; + + float trackJetDecayLen_; + float trackJetDistVal_; + float trackJetDistSig_; + reco::TransientTrack ttrack_; + +}; + +void ntuple_LT::readConfig(const edm::ParameterSet& iConfig, edm::ConsumesCollector&& cc) { + ltToken_ = cc.consumes>(iConfig.getParameter("losttracks")); +} + +void ntuple_LT::readSetup(const edm::EventSetup& iSetup){ + iSetup.get().get("TransientTrackBuilder", builder); +} + +void ntuple_LT::getInput(const edm::ParameterSet& iConfig){ + min_candidate_pt_ = (iConfig.getParameter("minCandidatePt")); +} + +void ntuple_LT::initBranches(TTree* tree){ + + addBranch(tree,"n_LT", &n_LT_,"n_LT_/I"); + + addBranch(tree,"LT_pt", <_pt_,"LT_pt_[n_Cpfcand_]/F"); + addBranch(tree,"LT_eta", <_eta_,"LT_eta_[n_Cpfcand_]/F"); + addBranch(tree,"LT_phi", <_phi_,"LT_phi_[n_Cpfcand_]/F"); + addBranch(tree,"LT_e", <_e_,"LT_e_[n_Cpfcand_]/F"); + + addBranch(tree,"LT_puppiw",<_puppiw_,"LT_puppiw_[n_Cpfcand_]/F"); + addBranch(tree,"LT_dz",<_dz_,"LT_dz_[n_Cpfcand_]/F"); + + addBranch(tree,"LT_drminsv",<_drminsv_,"LT_drminsv_[n_Cpfcand_]/F"); + addBranch(tree,"LT_distminsvold",<_distminsvold_,"LT_distminsvold_[n_Cpfcand_]/F"); + addBranch(tree,"LT_distminsv",<_distminsv_,"LT_distminsv_[n_Cpfcand_]/F"); + addBranch(tree,"LT_distminsv2",<_distminsv2_,"LT_distminsv2_[n_Cpfcand_]/F"); + + addBranch(tree,"LT_BtagPf_trackEtaRel",<_BtagPf_trackEtaRel_,"LT_BtagPf_trackEtaRel_[n_Cpfcand_]/F"); + addBranch(tree,"LT_BtagPf_trackPtRel",<_BtagPf_trackPtRel_,"LT_BtagPf_trackPtRel_[n_Cpfcand_]/F"); + addBranch(tree,"LT_BtagPf_trackPPar",<_BtagPf_trackPPar_,"LT_BtagPf_trackPPar_[n_Cpfcand_]/F"); + addBranch(tree,"LT_BtagPf_trackDeltaR",<_BtagPf_trackDeltaR_,"LT_BtagPf_trackDeltaR_[n_Cpfcand_]/F"); + addBranch(tree,"LT_BtagPf_trackPParRatio",<_BtagPf_trackPParRatio_,"LT_BtagPf_trackPParRatio[n_Cpfcand_]/F"); + addBranch(tree,"LT_BtagPf_trackSip3dVal",<_BtagPf_trackSip3dVal_,"LT_BtagPf_trackSip3dVal_[n_Cpfcand_]/F"); + addBranch(tree,"LT_BtagPf_trackSip3dSig",<_BtagPf_trackSip3dSig_,"LT_BtagPf_trackSip3dSig_[n_Cpfcand_]/F"); + addBranch(tree,"LT_BtagPf_trackSip2dVal",<_BtagPf_trackSip2dVal_,"LT_BtagPf_trackSip2dVal_[n_Cpfcand_]/F"); + addBranch(tree,"LT_BtagPf_trackSip2dSig",<_BtagPf_trackSip2dSig_,"LT_BtagPf_trackSip2dSig_[n_Cpfcand_]/F"); + addBranch(tree,"LT_BtagPf_trackDecayLen",<_BtagPf_trackDecayLen_,"LT_BtagPf_trackDecayLen_[n_Cpfcand_]/F"); + addBranch(tree,"LT_BtagPf_trackJetDistVal",<_BtagPf_trackJetDistVal_,"LT_BtagPf_trackJetDistVal_[n_Cpfcand_]/F"); + + addBranch(tree,"LT_charge",<_charge_,"LT_charge_[n_Cpfcand_]/F"); + addBranch(tree,"LT_chi2",<_chi2_,"LT_chi2_[n_Cpfcand_]/F"); + addBranch(tree,"LT_quality",<_quality_,"LT_quality_[n_Cpfcand_]/F"); + + addBranch(tree,"LT_nhitpixelBarrelLayer1",<_nhitpixelBarrelLayer1_,"LT_nhitpixelBarrelLayer1_[n_Cpfcand_]/F"); + addBranch(tree,"LT_nhitpixelBarrelLayer2",<_nhitpixelBarrelLayer2_,"LT_nhitpixelBarrelLayer2_[n_Cpfcand_]/F"); + addBranch(tree,"LT_nhitpixelBarrelLayer3",<_nhitpixelBarrelLayer3_,"LT_nhitpixelBarrelLayer3_[n_Cpfcand_]/F"); + addBranch(tree,"LT_nhitpixelBarrelLayer4",<_nhitpixelBarrelLayer4_,"LT_nhitpixelBarrelLayer4_[n_Cpfcand_]/F"); + + addBranch(tree,"LT_nhitpixelEndcapLayer1",<_nhitpixelEndcapLayer1_,"LT_nhitpixelEndcapLayer1_[n_Cpfcand_]/F"); + addBranch(tree,"LT_nhitpixelEndcapLayer2",<_nhitpixelEndcapLayer2_,"LT_nhitpixelEndcapLayer2_[n_Cpfcand_]/F"); + + addBranch(tree,"LT_numberOfValidHits",<_numberOfValidHits_,"LT_numberOfValidHits_[n_Cpfcand_]/F"); + addBranch(tree,"LT_numberOfValidPixelHits",<_numberOfValidPixelHits_,"LT_numberOfValidPixelHits_[n_Cpfcand_]/F"); + addBranch(tree,"LT_numberOfValidStripHits",<_numberOfValidStripHits_,"LT_numberOfValidStripHits_[n_Cpfcand_]/F"); + addBranch(tree,"LT_numberOfValidStripTIBHits",<_numberOfValidStripTIBHits_,"LT_numberOfValidStripTIBHits_[n_Cpfcand_]/F"); + addBranch(tree,"LT_numberOfValidStripTIDHits",<_numberOfValidStripTIDHits_,"LT_numberOfValidStripTIDHits_[n_Cpfcand_]/F"); + addBranch(tree,"LT_numberOfValidStripTOBHits",<_numberOfValidStripTOBHits_,"LT_numberOfValidStripTOBHits_[n_Cpfcand_]/F"); + addBranch(tree,"LT_numberOfValidStripTECHits",<_numberOfValidStripTECHits_,"LT_numberOfValidStripTECHits_[n_Cpfcand_]/F"); + + addBranch(tree,"LT_pdgID",<_pdgID_,"LT_pdgID_[n_Cpfcand_]/F"); + addBranch(tree,"LT_HadFrac",<_HadFrac_,"LT_HadFrac_[n_Cpfcand_]/F"); + addBranch(tree,"LT_CaloFrac",<_CaloFrac_,"LT_CaloFrac_[n_Cpfcand_]/F"); + +} + +void ntuple_LT::readEvent(const edm::Event& iEvent){ + + iEvent.getByToken(ltToken_, LTs); + n_Npfcand_=0; + n_Cpfcand_=0; + +} + + + +//use either of these functions + +bool ntuple_LT::fillBranches(const pat::Jet & jet, const size_t& jetidx, const edm::View * coll){ + + math::XYZVector jetDir = jet.momentum().Unit(); + GlobalVector jetRefTrackDir(jet.px(),jet.py(),jet.pz()); + const reco::Vertex & pv = vertices()->at(0); + + std::vector > sortedcharged; + + const float jet_uncorr_pt=jet.correctedJet("Uncorrected").pt(); + const float jet_uncorr_e=jet.correctedJet("Uncorrected").energy(); + + TrackInfoBuilder trackinfo(builder); + int n_lts = 0; + + std::vector cpfPtrs; + + for (size_t i = 0; i < LTs->size(); ++i) { + auto cand = LTs->ptrAt(i); + if (reco::deltaR(*cand, jet) < jetradius_) { + const auto *PackedCandidate = dynamic_cast(&(*cand)); + if(PackedCandidate){ + if(PackedCandidate->pt() < min_candidate_pt_) continue; + trackinfo.buildTrackInfo(PackedCandidate,jetDir,jetRefTrackDir,pv); + sortedcharged.push_back(sorting::sortingClass + (i, trackinfo.getTrackSip2dSig(), + -mindrsvpfcand(PackedCandidate), PackedCandidate->pt()/jet_uncorr_pt)); + cpfPtrs.push_back(cand); + n_lts++; + } + } + } + + std::sort(sortedcharged.begin(),sortedcharged.end(),sorting::sortingClass::compareByABCInv); + n_Cpfcand_ = std::min(sortedcharged.size(),max_pfcand_); + + std::vector sortedchargedindices; + sortedchargedindices=sorting::invertSortingVector(sortedcharged); + + for (unsigned int i = 0; i < (unsigned int)n_Cpfcand_; i++){ + const auto *PackedCandidate_ = dynamic_cast(&(*cpfPtrs.at(i))); + //const auto& PackedCandidate_=s.get(); + if(!PackedCandidate_) continue; + if(PackedCandidate_->pt() < min_candidate_pt_) continue; + + // get the dr with the closest sv + float drminpfcandsv_ = mindrsvpfcand(PackedCandidate_); + + float pdgid_; + if (abs(PackedCandidate_->pdgId()) == 11 and PackedCandidate_->charge() != 0){ + pdgid_ = 0.0; + } + else if (abs(PackedCandidate_->pdgId()) == 13 and PackedCandidate_->charge() != 0){ + pdgid_ = 1.0; + } + else if (abs(PackedCandidate_->pdgId()) == 22 and PackedCandidate_->charge() == 0){ + pdgid_ = 2.0; + } + else if (abs(PackedCandidate_->pdgId()) != 22 and PackedCandidate_->charge() == 0 and abs(PackedCandidate_->pdgId()) != 1 and abs(PackedCandidate_->pdgId()) != 2){ + pdgid_ = 3.0; + } + else if (abs(PackedCandidate_->pdgId()) != 11 and abs(PackedCandidate_->pdgId()) != 13 and PackedCandidate_->charge() != 0){ + pdgid_ = 4.0; + } + else if (PackedCandidate_->charge() == 0 and abs(PackedCandidate_->pdgId()) == 1){ + pdgid_ = 5.0; + } + else if (PackedCandidate_->charge() == 0 and abs(PackedCandidate_->pdgId()) == 2){ + pdgid_ = 6.0; + } + else{ + pdgid_ = 7.0; + } + + /// This might include more than PF candidates, e.g. Reco muons and could + /// be double counting. Needs to be checked.!!!! + /// + /// Split to charged and neutral candidates + if(PackedCandidate_->charge()!=0 ){ + + size_t fillntupleentry= sortedchargedindices.at(i); + if(fillntupleentry>=max_pfcand_) continue; + + LT_pt_[fillntupleentry] = PackedCandidate_->pt(); + LT_eta_[fillntupleentry] = PackedCandidate_->eta(); + LT_phi_[fillntupleentry] = PackedCandidate_->phi(); + LT_e_[fillntupleentry] = PackedCandidate_->energy(); + + LT_puppiw_[fillntupleentry] = PackedCandidate_->puppiWeight(); + LT_dz_[fillntupleentry] = PackedCandidate_->dz(); + LT_VTX_ass_[fillntupleentry] = PackedCandidate_->pvAssociationQuality(); + + float tempdontopt=PackedCandidate_->vx(); + tempdontopt++; + + LT_pdgID_[fillntupleentry] = pdgid_; + LT_CaloFrac_[fillntupleentry] = PackedCandidate_->caloFraction(); + LT_HadFrac_[fillntupleentry] = PackedCandidate_->hcalFraction(); + + + trackinfo.buildTrackInfo(PackedCandidate_,jetDir,jetRefTrackDir,pv); + + const reco::TransientTrack ttrack = trackinfo.getTTrack(); + float mindistsvold = mindistsvpfcandold(ttrack); + + LT_distminsvold_[fillntupleentry] = mindistsvold; + + float mindistsv = mindistsvpfcand(ttrack); + float eng_mindistsv = std::log(std::fabs(mindistsv)+1.0); + + LT_distminsv_[fillntupleentry] = mindistsv; + LT_distminsv2_[fillntupleentry] = eng_mindistsv; + + LT_BtagPf_trackEtaRel_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackEtaRel(), 0,-5,15); + LT_BtagPf_trackPtRel_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackPtRel(), 0,-1,4); + LT_BtagPf_trackPPar_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackPPar(), 0,-1e5,1e5 ); + LT_BtagPf_trackDeltaR_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackDeltaR(), 0,-5,5 ); + LT_BtagPf_trackPParRatio_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackPParRatio(),0,-10,100); + LT_BtagPf_trackSip3dVal_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackSip3dVal(), 0, -1,1e5 ); + LT_BtagPf_trackSip3dSig_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackSip3dSig(), 0, -1,4e4 ); + LT_BtagPf_trackSip2dVal_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackSip2dVal(), 0, -1,70 ); + LT_BtagPf_trackSip2dSig_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackSip2dSig(), 0, -1,4e4 ); + LT_BtagPf_trackDecayLen_[fillntupleentry] =trackinfo.getTrackJetDecayLen(); + LT_BtagPf_trackJetDistVal_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackJetDistVal(),0,-20,1 ); + + float cand_charge_ = PackedCandidate_->charge(); + LT_charge_[fillntupleentry] = cand_charge_; + + LT_drminsv_[fillntupleentry] = catchInfsAndBound(drminpfcandsv_,0,-0.4,0,-0.4); + + //hit pattern variables, as defined here https://github.com/cms-sw/cmssw/blob/master/DataFormats/TrackReco/interface/HitPattern.h + //get track associated to a jet constituent + const reco::Track *track_ptr = nullptr; + auto pf_candidate = dynamic_cast(PackedCandidate_); + auto packed_candidate = dynamic_cast(PackedCandidate_); + if(pf_candidate){ + track_ptr = pf_candidate->bestTrack(); //trackRef was sometimes null + }else if(packed_candidate && packed_candidate->hasTrackDetails()){//if PackedCandidate does not have TrackDetails this gives an Exception because unpackCovariance might be called for pseudoTrack/bestTrack + track_ptr = &(packed_candidate->pseudoTrack()); + } + //get hit pattern information + if(track_ptr){ + const reco::HitPattern &p = track_ptr->hitPattern(); + //Tracker per layer + //Pixel barrel + int Cpfcan_nhitpixelBarrelLayer1 = 0; + int Cpfcan_nhitpixelBarrelLayer2 = 0; + int Cpfcan_nhitpixelBarrelLayer3 = 0; + int Cpfcan_nhitpixelBarrelLayer4 = 0; + //Pixel Endcap + int Cpfcan_nhitpixelEndcapLayer1 = 0; + int Cpfcan_nhitpixelEndcapLayer2 = 0; + // loop over the hits of the track. + //const static unsigned short LayerOffset = 3; + //const static unsigned short LayerMask = 0xF; + for(int nh = 0; nh < p.numberOfAllHits(reco::HitPattern::TRACK_HITS); nh++){ + uint32_t hit = p.getHitPattern(reco::HitPattern::TRACK_HITS, nh); + if(p.validHitFilter(hit)){// if the hit is valid + //Pixel Barrel // it is in pixel barrel + if(p.pixelBarrelHitFilter(hit)){ + //std::cout << "valid hit found in pixel Barrel layer " << p.getLayer(hit) << std::endl; + //if(p.getLayer(hit)==1){ + // std::cout<< (hit >> LayerOffset) << " " << ((hit >> LayerOffset) & LayerMask) << std::endl; + //} + if(p.getLayer(hit)==1) Cpfcan_nhitpixelBarrelLayer1 = Cpfcan_nhitpixelBarrelLayer1+1; + if(p.getLayer(hit)==2) Cpfcan_nhitpixelBarrelLayer2 = Cpfcan_nhitpixelBarrelLayer2+1; + if(p.getLayer(hit)==3) Cpfcan_nhitpixelBarrelLayer3 = Cpfcan_nhitpixelBarrelLayer3+1; + if(p.getLayer(hit)==4) Cpfcan_nhitpixelBarrelLayer4 = Cpfcan_nhitpixelBarrelLayer4+1; + } + //Pixel Endcap + if(p.pixelEndcapHitFilter(hit)){ + //std::cout << "valid hit found in pixel Endcap layer " << p.getLayer(hit) << std::endl; + if(p.getLayer(hit)==1) Cpfcan_nhitpixelEndcapLayer1 = Cpfcan_nhitpixelEndcapLayer1+1; + if(p.getLayer(hit)==2) Cpfcan_nhitpixelEndcapLayer2 = Cpfcan_nhitpixelEndcapLayer2+1; + } + } + } + //Pixel Barrel + LT_nhitpixelBarrelLayer1_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelBarrelHits()) ? catchInfsAndBound(Cpfcan_nhitpixelBarrelLayer1,-1,0,100,0) : -1; + LT_nhitpixelBarrelLayer2_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelBarrelHits()) ? catchInfsAndBound(Cpfcan_nhitpixelBarrelLayer2,-1,0,100,0) : -1; + LT_nhitpixelBarrelLayer3_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelBarrelHits()) ? catchInfsAndBound(Cpfcan_nhitpixelBarrelLayer3,-1,0,100,0) : -1; + LT_nhitpixelBarrelLayer4_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelBarrelHits()) ? catchInfsAndBound(Cpfcan_nhitpixelBarrelLayer4,-1,0,100,0) : -1; + + LT_nhitpixelEndcapLayer1_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelEndcapHits()) ? catchInfsAndBound(Cpfcan_nhitpixelEndcapLayer1,-1,0,100,0) : -1; + LT_nhitpixelEndcapLayer2_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelEndcapHits()) ? catchInfsAndBound(Cpfcan_nhitpixelEndcapLayer2,-1,0,100,0) : -1; + //Tracker all layers together + //Valid hits + LT_numberOfValidHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidHits(),-1,0,100,0) : -1; + LT_numberOfValidPixelHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidPixelHits(),-1,0,100,0) : -1; + LT_numberOfValidStripHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidStripHits(),-1,0,100,0) : -1; + LT_numberOfValidStripTIBHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTIBHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidStripTIBHits(),-1,0,100,0) : -1; + LT_numberOfValidStripTIDHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTIDHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidStripTIDHits(),-1,0,100,0) : -1; + LT_numberOfValidStripTOBHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTOBHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidStripTOBHits(),-1,0,100,0) : -1; + LT_numberOfValidStripTECHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTECHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidStripTECHits(),-1,0,100,0) : -1; + } + else{ + //Tracker per layer + //Pixel barrel + LT_nhitpixelBarrelLayer1_[fillntupleentry] = -1; + LT_nhitpixelBarrelLayer2_[fillntupleentry] = -1; + LT_nhitpixelBarrelLayer3_[fillntupleentry] = -1; + LT_nhitpixelBarrelLayer4_[fillntupleentry] = -1; + //Pixel Endcap + LT_nhitpixelEndcapLayer1_[fillntupleentry] = -1; + LT_nhitpixelEndcapLayer2_[fillntupleentry] = -1; + + LT_numberOfValidHits_[fillntupleentry] = -1; + LT_numberOfValidPixelHits_[fillntupleentry] = -1; + LT_numberOfValidStripHits_[fillntupleentry] = -1; + LT_numberOfValidStripTIBHits_[fillntupleentry] = -1; + LT_numberOfValidStripTIDHits_[fillntupleentry] = -1; + LT_numberOfValidStripTOBHits_[fillntupleentry] = -1; + LT_numberOfValidStripTECHits_[fillntupleentry] = -1; + } + } + } + + n_LT_ = n_lts; + + return true; +} + + +float ntuple_LT::mindrsvpfcand(const pat::PackedCandidate* pfcand) { + + float mindr_ = jetradius_; + for (unsigned int i=0; isize(); ++i) { + if(!pfcand) continue; + //if(!svs.at(i)) continue; + float tempdr_ = reco::deltaR(secVertices()->at(i),*pfcand); + if (tempdr_size(); ++i) { + if (!track.isValid()) {continue;} + reco::Vertex::CovarianceMatrix csv; secVertices()->at(i).fillVertexCovariance(csv); + reco::Vertex vertex(secVertices()->at(i).vertex(), csv); + if (!vertex.isValid()) {continue;} + + GlobalVector direction(secVertices()->at(i).px(),secVertices()->at(i).py(),secVertices()->at(i).pz()); + + + AnalyticalImpactPointExtrapolator extrapolator(track.field()); + TrajectoryStateOnSurface tsos = extrapolator.extrapolate(track.impactPointState(), RecoVertex::convertPos(vertex.position())); + + VertexDistance3D dist; + + if (!tsos.isValid()) {continue;} + GlobalPoint refPoint = tsos.globalPosition(); + GlobalError refPointErr = tsos.cartesianError().position(); + GlobalPoint vertexPosition = RecoVertex::convertPos(vertex.position()); + GlobalError vertexPositionErr = RecoVertex::convertError(vertex.error()); + + std::pair result(true, dist.distance(VertexState(vertexPosition, vertexPositionErr), VertexState(refPoint, refPointErr))); + if (!result.first) {continue;} + + GlobalPoint impactPoint = tsos.globalPosition(); + GlobalVector IPVec(impactPoint.x() - vertex.x(), impactPoint.y() - vertex.y(), impactPoint.z() - vertex.z()); + double prod = IPVec.dot(direction); + double sign = (prod >= 0) ? 1. : -1.; + + if(result.second.value() < mindist_){ + out_dist = sign * result.second.value(); + mindist_ = result.second.value(); + } + } + return out_dist; +} + +GlobalPoint ntuple_LT::mingpsvpfcand(const reco::TransientTrack track) { + + float mindist_ = 999.999; + GlobalPoint out_dist(0.0,0.0,0.0); + GlobalPoint pca; + for (unsigned int i=0; isize(); ++i) { + if (!track.isValid()) {continue;} + + reco::Vertex::CovarianceMatrix csv; secVertices()->at(i).fillVertexCovariance(csv); + reco::Vertex vertex(secVertices()->at(i).vertex(), csv); + if (!vertex.isValid()) {continue;} + + GlobalVector direction(secVertices()->at(i).px(),secVertices()->at(i).py(),secVertices()->at(i).pz()); + + AnalyticalImpactPointExtrapolator extrapolator(track.field()); + TrajectoryStateOnSurface tsos = extrapolator.extrapolate(track.impactPointState(), RecoVertex::convertPos(vertex.position())); + + VertexDistance3D dist; + + if (!tsos.isValid()) {continue;} + GlobalPoint refPoint = tsos.globalPosition(); + GlobalError refPointErr = tsos.cartesianError().position(); + GlobalPoint vertexPosition = RecoVertex::convertPos(vertex.position()); + GlobalError vertexPositionErr = RecoVertex::convertError(vertex.error()); + + std::pair result(true, dist.distance(VertexState(vertexPosition, vertexPositionErr), VertexState(refPoint, refPointErr))); + if (!result.first) {continue;} + + pca = tsos.globalPosition(); + + if(result.second.value() < mindist_){ + out_dist = pca; + mindist_ = result.second.value(); + } + } + return out_dist; +} + +GlobalPoint ntuple_LT::gppvpfcand(const reco::TransientTrack track, const GlobalVector direction, const reco::Vertex vertex) { + + float mindist_ = 999.999; + float dist = 0.; + GlobalPoint out_dist(0.0,0.0,0.0); + GlobalPoint pca; + if ((track.isValid()) && (vertex.isValid())){ + + AnalyticalImpactPointExtrapolator extrapolator(track.field()); + TrajectoryStateOnSurface tsos = extrapolator.extrapolate(track.impactPointState(), RecoVertex::convertPos(vertex.position())); + + VertexDistance3D dist; + + if (tsos.isValid()) { + GlobalPoint refPoint = tsos.globalPosition(); + GlobalError refPointErr = tsos.cartesianError().position(); + GlobalPoint vertexPosition = RecoVertex::convertPos(vertex.position()); + GlobalError vertexPositionErr = RecoVertex::convertError(vertex.error()); + + std::pair result(true, dist.distance(VertexState(vertexPosition, vertexPositionErr), VertexState(refPoint, refPointErr))); + if (result.first) { + pca = tsos.globalPosition(); + } + } + } + out_dist = pca; + return out_dist; +} + +float ntuple_LT::mindistsvpfcand(const reco::TransientTrack track) { + + float mindist_ = 999.999; + float out_dist = 0.0; + for (unsigned int i=0; isize(); ++i) { + if (!track.isValid()) {continue;} + reco::Vertex::CovarianceMatrix csv; secVertices()->at(i).fillVertexCovariance(csv); + reco::Vertex vertex(secVertices()->at(i).vertex(), csv); + if (!vertex.isValid()) {continue;} + + GlobalVector direction(secVertices()->at(i).px(),secVertices()->at(i).py(),secVertices()->at(i).pz()); + + + AnalyticalImpactPointExtrapolator extrapolator(track.field()); + TrajectoryStateOnSurface tsos = extrapolator.extrapolate(track.impactPointState(), RecoVertex::convertPos(vertex.position())); + + VertexDistance3D dist; + + if (!tsos.isValid()) {continue;} + GlobalPoint refPoint = tsos.globalPosition(); + GlobalError refPointErr = tsos.cartesianError().position(); + GlobalPoint vertexPosition = RecoVertex::convertPos(vertex.position()); + GlobalError vertexPositionErr = RecoVertex::convertError(vertex.error()); + + std::pair result(true, dist.distance(VertexState(vertexPosition, vertexPositionErr), VertexState(refPoint, refPointErr))); + if (!result.first) {continue;} + + GlobalPoint impactPoint = tsos.globalPosition(); + GlobalVector IPVec(impactPoint.x() - vertex.x(), impactPoint.y() - vertex.y(), impactPoint.z() - vertex.z()); + double prod = IPVec.dot(direction); + double sign = (prod >= 0) ? 1. : -1.; + + if(result.second.value() < mindist_){ + float inv_dist = 1.0 / (sign * result.second.value() + 0.001); + out_dist = inv_dist; + mindist_ = result.second.value(); + } + } + return out_dist; +} diff --git a/DeepNtuplizer/src/ntuple_SV.cc b/DeepNtuplizer/src/ntuple_SV.cc index a32fc6757f8..ac43425cd3a 100644 --- a/DeepNtuplizer/src/ntuple_SV.cc +++ b/DeepNtuplizer/src/ntuple_SV.cc @@ -16,7 +16,139 @@ #include "DataFormats/Math/interface/deltaR.h" #include "TrackingTools/IPTools/interface/IPTools.h" #include "DataFormats/GeometryCommonDetAlgo/interface/Measurement1D.h" +#include "DataFormats/PatCandidates/interface/PackedCandidate.h" +#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h" +#include "TrackingTools/IPTools/interface/IPTools.h" +#include "TrackingTools/Records/interface/TransientTrackRecord.h" +#include "RecoVertex/VertexTools/interface/VertexDistance3D.h" +#include "TVector3.h" + +class TrackInfoBuilder{ +public: + TrackInfoBuilder(edm::ESHandle & build): + builder(build), + trackMomentum_(0), + trackEta_(0), + trackEtaRel_(0), + trackPtRel_(0), + trackPPar_(0), + trackDeltaR_(0), + trackPtRatio_(0), + trackPParRatio_(0), + trackSip2dVal_(0), + trackSip2dSig_(0), + trackSip3dVal_(0), + trackSip3dSig_(0), + + trackJetDecayLen_(0), + trackJetDistVal_(0), + trackJetDistSig_(0), + ttrack_(0) + { + + } + + void buildTrackInfo(const pat::PackedCandidate* PackedCandidate_ ,const math::XYZVector& jetDir, GlobalVector refjetdirection, const reco::Vertex & pv){ + TVector3 jetDir3(jetDir.x(),jetDir.y(),jetDir.z()); + if(!PackedCandidate_->hasTrackDetails()) { + TVector3 trackMom3( + PackedCandidate_->momentum().x(), + PackedCandidate_->momentum().y(), + PackedCandidate_->momentum().z() + ); + trackMomentum_=PackedCandidate_->p(); + trackEta_= PackedCandidate_->eta(); + trackEtaRel_=reco::btau::etaRel(jetDir, PackedCandidate_->momentum()); + trackPtRel_=trackMom3.Perp(jetDir3); + trackPPar_=jetDir.Dot(PackedCandidate_->momentum()); + trackDeltaR_=reco::deltaR(PackedCandidate_->momentum(), jetDir); + trackPtRatio_=trackMom3.Perp(jetDir3) / PackedCandidate_->p(); + trackPParRatio_=jetDir.Dot(PackedCandidate_->momentum()) / PackedCandidate_->p(); + trackSip2dVal_=0.; + trackSip2dSig_=0.; + trackSip3dVal_=0.; + trackSip3dSig_=0.; + trackJetDecayLen_=0.; + trackJetDistVal_=0.; + trackJetDistSig_=0.; + return; + } + + const reco::Track & PseudoTrack = PackedCandidate_->pseudoTrack(); + + reco::TransientTrack transientTrack; + transientTrack=builder->build(PseudoTrack); + Measurement1D meas_ip2d=IPTools::signedTransverseImpactParameter(transientTrack, refjetdirection, pv).second; + Measurement1D meas_ip3d=IPTools::signedImpactParameter3D(transientTrack, refjetdirection, pv).second; + Measurement1D jetdist=IPTools::jetTrackDistance(transientTrack, refjetdirection, pv).second; + Measurement1D decayl = IPTools::signedDecayLength3D(transientTrack, refjetdirection, pv).second; + math::XYZVector trackMom = PseudoTrack.momentum(); + double trackMag = std::sqrt(trackMom.Mag2()); + TVector3 trackMom3(trackMom.x(),trackMom.y(),trackMom.z()); + + trackMomentum_=std::sqrt(trackMom.Mag2()); + trackEta_= trackMom.Eta(); + trackEtaRel_=reco::btau::etaRel(jetDir, trackMom); + trackPtRel_=trackMom3.Perp(jetDir3); + trackPPar_=jetDir.Dot(trackMom); + trackDeltaR_=reco::deltaR(trackMom, jetDir); + trackPtRatio_=trackMom3.Perp(jetDir3) / trackMag; + trackPParRatio_=jetDir.Dot(trackMom) / trackMag; + + trackSip2dVal_=(meas_ip2d.value()); + trackSip2dSig_=(meas_ip2d.significance()); + trackSip3dVal_=(meas_ip3d.value()); + trackSip3dSig_=meas_ip3d.significance(); + + trackJetDecayLen_= decayl.value(); + trackJetDistVal_= jetdist.value(); + trackJetDistSig_= jetdist.significance(); + + ttrack_ = transientTrack; + + } + + const float& getTrackDeltaR() const {return trackDeltaR_;} + const float& getTrackEta() const {return trackEta_;} + const float& getTrackEtaRel() const {return trackEtaRel_;} + const float& getTrackJetDecayLen() const {return trackJetDecayLen_;} + const float& getTrackJetDistSig() const {return trackJetDistSig_;} + const float& getTrackJetDistVal() const {return trackJetDistVal_;} + const float& getTrackMomentum() const {return trackMomentum_;} + const float& getTrackPPar() const {return trackPPar_;} + const float& getTrackPParRatio() const {return trackPParRatio_;} + const float& getTrackPtRatio() const {return trackPtRatio_;} + const float& getTrackPtRel() const {return trackPtRel_;} + const float& getTrackSip2dSig() const {return trackSip2dSig_;} + const float& getTrackSip2dVal() const {return trackSip2dVal_;} + const float& getTrackSip3dSig() const {return trackSip3dSig_;} + const float& getTrackSip3dVal() const {return trackSip3dVal_;} + const reco::TransientTrack getTTrack() const {return ttrack_;} + +private: + + edm::ESHandle& builder; + + float trackMomentum_; + float trackEta_; + float trackEtaRel_; + float trackPtRel_; + float trackPPar_; + float trackDeltaR_; + float trackPtRatio_; + float trackPParRatio_; + float trackSip2dVal_; + float trackSip2dSig_; + float trackSip3dVal_; + float trackSip3dSig_; + + float trackJetDecayLen_; + float trackJetDistVal_; + float trackJetDistSig_; + reco::TransientTrack ttrack_; + +}; const reco::Vertex * ntuple_SV::spvp_; @@ -31,32 +163,48 @@ void ntuple_SV::getInput(const edm::ParameterSet& iConfig){ } void ntuple_SV::initBranches(TTree* tree){ - // SV candidates - addBranch(tree,(prefix_+"n_sv").c_str() ,&sv_num_ ,(prefix_+"sv_num_/i").c_str() ); - addBranch(tree,(prefix_+"nsv").c_str() ,&nsv_ ,(prefix_+"nsv_/f").c_str() ); - addBranch(tree,(prefix_+"sv_pt").c_str() ,&sv_pt_ ,(prefix_+"sv_pt_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_eta").c_str() ,&sv_eta_ ,(prefix_+"sv_eta_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_phi").c_str() ,&sv_phi_ ,(prefix_+"sv_phi_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_etarel").c_str() ,&sv_etarel_ ,(prefix_+"sv_etarel_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_phirel").c_str() ,&sv_phirel_ ,(prefix_+"sv_phirel_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_deltaR").c_str() ,&sv_deltaR_ ,(prefix_+"sv_deltaR_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_mass").c_str() ,&sv_mass_ ,(prefix_+"sv_mass_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_ntracks").c_str() ,&sv_ntracks_ ,(prefix_+"sv_ntracks_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_chi2").c_str() ,&sv_chi2_ ,(prefix_+"sv_chi2_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_ndf").c_str() ,&sv_ndf_ ,(prefix_+"sv_ndf_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_normchi2").c_str() ,&sv_normchi2_ ,(prefix_+"sv_normchi2_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_dxy").c_str() ,&sv_dxy_ ,(prefix_+"sv_dxy_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_dxyerr").c_str() ,&sv_dxyerr_ ,(prefix_+"sv_dxyerr_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_dxysig").c_str() ,&sv_dxysig_ ,(prefix_+"sv_dxysig_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_d3d").c_str() ,&sv_d3d_ ,(prefix_+"sv_d3d_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_d3derr").c_str() ,&sv_d3derr_ ,(prefix_+"sv_d3err_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_d3dsig").c_str() ,&sv_d3dsig_ ,(prefix_+"sv_d3dsig_["+prefix_+"sv_num_]/f").c_str() ); - addBranch(tree,(prefix_+"sv_costhetasvpv").c_str(),&sv_costhetasvpv_,(prefix_+"sv_costhetasvpv_["+prefix_+"sv_num_]/f").c_str()); - addBranch(tree,(prefix_+"sv_enratio").c_str() ,&sv_enratio_ ,(prefix_+"sv_enratio_["+prefix_+"sv_num_]/f").c_str()); - + // SV candidates + addBranch(tree,(prefix_+"n_sv").c_str() ,&sv_num_ ,(prefix_+"sv_num_/I").c_str() ); + addBranch(tree,(prefix_+"nsv").c_str() ,&nsv_ ,(prefix_+"nsv_/F").c_str() ); + addBranch(tree,(prefix_+"sv_pt").c_str() ,&sv_pt_ ,(prefix_+"sv_pt_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_eta").c_str() ,&sv_eta_ ,(prefix_+"sv_eta_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_phi").c_str() ,&sv_phi_ ,(prefix_+"sv_phi_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_e").c_str() ,&sv_e_ ,(prefix_+"sv_e_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_etarel").c_str() ,&sv_etarel_ ,(prefix_+"sv_etarel_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_phirel").c_str() ,&sv_phirel_ ,(prefix_+"sv_phirel_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_deltaR").c_str() ,&sv_deltaR_ ,(prefix_+"sv_deltaR_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_mass").c_str() ,&sv_mass_ ,(prefix_+"sv_mass_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_ntracks").c_str() ,&sv_ntracks_ ,(prefix_+"sv_ntracks_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_chi2").c_str() ,&sv_chi2_ ,(prefix_+"sv_chi2_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_ndf").c_str() ,&sv_ndf_ ,(prefix_+"sv_ndf_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_normchi2").c_str() ,&sv_normchi2_ ,(prefix_+"sv_normchi2_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_dxy").c_str() ,&sv_dxy_ ,(prefix_+"sv_dxy_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_dxyerr").c_str() ,&sv_dxyerr_ ,(prefix_+"sv_dxyerr_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_dxysig").c_str() ,&sv_dxysig_ ,(prefix_+"sv_dxysig_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_d3d").c_str() ,&sv_d3d_ ,(prefix_+"sv_d3d_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_d3derr").c_str() ,&sv_d3derr_ ,(prefix_+"sv_d3err_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_d3dsig").c_str() ,&sv_d3dsig_ ,(prefix_+"sv_d3dsig_["+prefix_+"sv_num_]/F").c_str() ); + addBranch(tree,(prefix_+"sv_costhetasvpv").c_str(),&sv_costhetasvpv_,(prefix_+"sv_costhetasvpv_["+prefix_+"sv_num_]/F").c_str()); + addBranch(tree,(prefix_+"sv_enratio").c_str() ,&sv_enratio_ ,(prefix_+"sv_enratio_["+prefix_+"sv_num_]/F").c_str()); + + addBranch(tree,(prefix_+"sv_hcal_frac").c_str() ,&sv_hcal_frac_ ,(prefix_+"sv_hcal_frac_["+prefix_+"sv_num_]/F").c_str()); + addBranch(tree,(prefix_+"sv_calo_frac").c_str() ,&sv_calo_frac_ ,(prefix_+"sv_calo_frac_["+prefix_+"sv_num_]/F").c_str()); + + addBranch(tree,(prefix_+"sv_dz").c_str() ,&sv_dz_ ,(prefix_+"sv_dz_["+prefix_+"sv_num_]/F").c_str()); + addBranch(tree,(prefix_+"sv_pfd2dval").c_str() ,&sv_pfd2dval_ ,(prefix_+"sv_pfd2dval_["+prefix_+"sv_num_]/F").c_str()); + addBranch(tree,(prefix_+"sv_pfd2dsig").c_str() ,&sv_pfd2dsig_ ,(prefix_+"sv_pfd2dsig_["+prefix_+"sv_num_]/F").c_str()); + addBranch(tree,(prefix_+"sv_pfd3dval").c_str() ,&sv_pfd3dval_ ,(prefix_+"sv_pfd3dval_["+prefix_+"sv_num_]/F").c_str()); + addBranch(tree,(prefix_+"sv_pfd3dsig").c_str() ,&sv_pfd3dsig_ ,(prefix_+"sv_pfd3dsig_["+prefix_+"sv_num_]/F").c_str()); + addBranch(tree,(prefix_+"sv_puppiw").c_str() ,&sv_puppiw_ ,(prefix_+"sv_puppiw_["+prefix_+"sv_num_]/F").c_str()); + addBranch(tree,(prefix_+"sv_charge_sum").c_str() ,&sv_charge_sum_ ,(prefix_+"sv_charge_sum_["+prefix_+"sv_num_]/F").c_str()); } +void ntuple_SV::readSetup(const edm::EventSetup& iSetup){ + + iSetup.get().get("TransientTrackBuilder", builder); + +} void ntuple_SV::readEvent(const edm::Event& iEvent){ @@ -78,18 +226,19 @@ bool ntuple_SV::compareDxyDxyErr(const reco::VertexCompositePtrCandidate &sva,co bool ntuple_SV::fillBranches(const pat::Jet & jet, const size_t& jetidx, const edm::View * coll){ - const float jet_uncorr_e=jet.correctedJet("Uncorrected").energy(); - const reco::Vertex & pv = vertices()->at(0); + math::XYZVector jetDir = jet.momentum().Unit(); + GlobalVector jetRefTrackDir(jet.px(),jet.py(),jet.pz()); sv_num_ = 0; reco::VertexCompositePtrCandidateCollection cpvtx=*secVertices(); - spvp_ = & vertices()->at(0); std::sort(cpvtx.begin(),cpvtx.end(),ntuple_SV::compareDxyDxyErr); + TrackInfoBuilder trackinfo(builder); + float etasign=1; etasign++; //avoid unused warning if(jet.eta()<0)etasign=-1; @@ -131,8 +280,48 @@ bool ntuple_SV::fillBranches(const pat::Jet & jet, const size_t& jetidx, const // of the tracks in the SV and the flight direction betwen PV and SV) sv_enratio_[sv_num_]=sv.energy()/jet_uncorr_e; - - + sv_e_[sv_num_]=sv.energy(); + + float calo_frac = 0.0; + float hcal_frac = 0.0; + float puppiw = 0.0; + float charge = 0.0; + float dz = 0.0; + + float pfd3dval = 0.0; + float pfd3dsig = 0.0; + float pfd2dval = 0.0; + float pfd2dsig = 0.0; + float pfcount = 0.0; + + for (unsigned idx=0; idx(sv.daughter(idx)); + + calo_frac = calo_frac + PackedCandidate_->caloFraction(); + hcal_frac = hcal_frac + PackedCandidate_->hcalFraction(); + puppiw = puppiw + PackedCandidate_->puppiWeight(); + charge = charge + PackedCandidate_->charge(); + dz = dz + PackedCandidate_->dz(); + if(PackedCandidate_->charge() != 0 and PackedCandidate_->pt() > 0.95){ + trackinfo.buildTrackInfo(PackedCandidate_,jetDir,jetRefTrackDir,pv); + pfd3dval = pfd3dval + catchInfsAndBound(trackinfo.getTrackSip3dVal(), 0, -1,1e5 ); + pfd3dsig = pfd3dsig + catchInfsAndBound(trackinfo.getTrackSip3dSig(), 0, -1,4e4 ); + pfd2dval = pfd2dval + catchInfsAndBound(trackinfo.getTrackSip2dVal(), 0, -1,70 ); + pfd2dsig = pfd2dsig + catchInfsAndBound(trackinfo.getTrackSip2dSig(), 0, -1,4e4 ); + pfcount = pfcount + 1.0; + } + } + + sv_calo_frac_[sv_num_] = calo_frac / sv.numberOfDaughters(); + sv_hcal_frac_[sv_num_] = hcal_frac / sv.numberOfDaughters(); + sv_puppiw_[sv_num_] = puppiw / sv.numberOfDaughters(); + sv_dz_[sv_num_] = dz / sv.numberOfDaughters(); + sv_charge_sum_[sv_num_] = charge; + + sv_pfd3dval_[sv_num_] = pfd3dval / pfcount; + sv_pfd3dsig_[sv_num_] = pfd3dsig / pfcount; + sv_pfd2dval_[sv_num_] = pfd2dval / pfcount; + sv_pfd2dsig_[sv_num_] = pfd2dsig / pfcount; sv_num_++; } @@ -142,18 +331,7 @@ bool ntuple_SV::fillBranches(const pat::Jet & jet, const size_t& jetidx, const return true; } - - - - - - - - ///helpers seldomly touched - - - Measurement1D ntuple_SV::vertexDxy(const reco::VertexCompositePtrCandidate &svcand, const reco::Vertex &pv) { VertexDistanceXY dist; reco::Vertex::CovarianceMatrix csv; svcand.fillVertexCovariance(csv); diff --git a/DeepNtuplizer/src/ntuple_pairwise.cc b/DeepNtuplizer/src/ntuple_pairwise.cc new file mode 100644 index 00000000000..beb2f68ba87 --- /dev/null +++ b/DeepNtuplizer/src/ntuple_pairwise.cc @@ -0,0 +1,246 @@ +/* + * ntuple_pairwise.cc + * + * Created on: 01 Sep 2022 + * Authors: Matteo Malucchi & Alexandre De Moor + */ + +#include "../interface/ntuple_pairwise.h" +#include "DataFormats/PatCandidates/interface/PackedCandidate.h" +#include "DataFormats/Candidate/interface/VertexCompositePtrCandidate.h" +#include "../interface/sorting_modules.h" + +#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h" +#include "TrackingTools/IPTools/interface/IPTools.h" +#include "TrackingTools/Records/interface/TransientTrackRecord.h" +#include "RecoVertex/VertexTools/interface/VertexDistance3D.h" +#include "TrackingTools/PatternTools/interface/TwoTrackMinimumDistance.h" +#include "TVector3.h" + + +class TrackInfoBuilder{ +public: + TrackInfoBuilder(edm::ESHandle & build): + builder(build), + trackSip2dSig_(0), + ttrack_(0) + +{ + +} + + void buildTrackInfo(const pat::PackedCandidate* PackedCandidate_ ,const math::XYZVector& jetDir, GlobalVector refjetdirection, const reco::Vertex & pv){ + TVector3 jetDir3(jetDir.x(),jetDir.y(),jetDir.z()); + if(!PackedCandidate_->hasTrackDetails()) { + trackSip2dSig_=0.; + return; + } + const reco::Track & PseudoTrack = PackedCandidate_->pseudoTrack(); + reco::TransientTrack transientTrack; + transientTrack=builder->build(PseudoTrack); + ttrack_ = transientTrack; + + Measurement1D meas_ip2d=IPTools::signedTransverseImpactParameter(transientTrack, refjetdirection, pv).second; + trackSip2dSig_=(meas_ip2d.significance()); + } + + const reco::TransientTrack getTTrack() const {return ttrack_;} + const float& getTrackSip2dSig() const {return trackSip2dSig_;} + +private: + + edm::ESHandle& builder; + float trackSip2dSig_; + reco::TransientTrack ttrack_; + +}; + +void ntuple_pairwise::readSetup(const edm::EventSetup& iSetup){ + + iSetup.get().get("TransientTrackBuilder", builder); + +} + +void ntuple_pairwise::getInput(const edm::ParameterSet& iConfig){ + min_candidate_pt_ = (iConfig.getParameter("minCandidatePt")); +} + +void ntuple_pairwise::initBranches(TTree* tree){ + + addBranch(tree,"n_Cpfpairs", &n_Cpfpairs_,"n_Cpfpairs_/I"); + addBranch(tree,"nCpfpairs", &nCpfpairs_,"nCpfpairs_/F"); + + addBranch(tree,"pair_pca_distance", &pair_pca_distance_,"pair_pca_distance_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_pca_significance", &pair_pca_significance_,"pair_pca_significance_[n_Cpfpairs_]/F"); + + addBranch(tree,"pair_pcaSeed_x1", &pair_pcaSeed_x1_,"pair_pcaSeed_x1_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_pcaSeed_y1", &pair_pcaSeed_y1_,"pair_pcaSeed_y1_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_pcaSeed_z1", &pair_pcaSeed_z1_,"pair_pcaSeed_z1_[n_Cpfpairs_]/F"); + + addBranch(tree,"pair_pcaSeed_x2", &pair_pcaSeed_x2_,"pair_pcaSeed_x2_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_pcaSeed_y2", &pair_pcaSeed_y2_,"pair_pcaSeed_y2_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_pcaSeed_z2", &pair_pcaSeed_z2_,"pair_pcaSeed_z2_[n_Cpfpairs_]/F"); + + addBranch(tree,"pair_pcaSeed_xerr1", &pair_pcaSeed_xerr1_,"pair_pcaSeed_xerr1_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_pcaSeed_yerr1", &pair_pcaSeed_yerr1_,"pair_pcaSeed_yerr1_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_pcaSeed_zerr1", &pair_pcaSeed_zerr1_,"pair_pcaSeed_zerr1_[n_Cpfpairs_]/F"); + + addBranch(tree,"pair_pcaSeed_xerr2", &pair_pcaSeed_xerr2_,"pair_pcaSeed_xerr2_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_pcaSeed_yerr2", &pair_pcaSeed_yerr2_,"pair_pcaSeed_yerr2_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_pcaSeed_zerr2", &pair_pcaSeed_zerr2_,"pair_pcaSeed_zerr2_[n_Cpfpairs_]/F"); + + addBranch(tree,"pair_dotprod1", &pair_dotprod1_,"pair_dotprod1_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_dotprod2", &pair_dotprod2_,"pair_dotprod2_[n_Cpfpairs_]/F"); + + addBranch(tree,"pair_pca_dist1", &pair_pca_dist1_,"pair_pca_dist1_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_pca_dist2", &pair_pca_dist2_,"pair_pca_dist2_[n_Cpfpairs_]/F"); + + addBranch(tree,"pair_dotprod12_2D", &pair_dotprod12_2D_,"pair_dotprod12_2D_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_dotprod12_2DV", &pair_dotprod12_2DV_,"pair_dotprod12_2DV_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_dotprod12_3D", &pair_dotprod12_3D_,"pair_dotprod12_3D_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_dotprod12_3DV", &pair_dotprod12_3DV_,"pair_dotprod12_3DV_[n_Cpfpairs_]/F"); + + addBranch(tree,"pair_pca_jetAxis_dist", &pair_pca_jetAxis_dist_,"pair_pca_jetAxis_dist_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_pca_jetAxis_dotprod", &pair_pca_jetAxis_dotprod_,"pair_pca_jetAxis_dotprod_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_pca_jetAxis_dEta", &pair_pca_jetAxis_dEta_,"pair_pca_jetAxis_dEta_[n_Cpfpairs_]/F"); + addBranch(tree,"pair_pca_jetAxis_dPhi", &pair_pca_jetAxis_dPhi_,"pair_pca_jetAxis_dPhi_[n_Cpfpairs_]/F"); + + addBranch(tree,"pfcand_dist_vtx_12", &pfcand_dist_vtx_12_,"pfcand_dist_vtx_12_[n_Cpfpairs_]/F"); + +} + +void ntuple_pairwise::readEvent(const edm::Event& iEvent){ + + n_Npfcand2_=0; + n_Cpfcand2_=0; + +} + +//use either of these functions + +bool ntuple_pairwise::fillBranches(const pat::Jet & jet, const size_t& jetidx, const edm::View * coll){ + + math::XYZVector jetDir = jet.momentum().Unit(); + GlobalVector jetRefTrackDir(jet.px(),jet.py(),jet.pz()); + const reco::Vertex & pv = vertices()->at(0); + + std::vector > sortedcharged, sortedneutrals; + + const float jet_uncorr_pt=jet.correctedJet("Uncorrected").pt(); + const float jet_uncorr_e=jet.correctedJet("Uncorrected").energy(); + + TrackInfoBuilder trackinfo(builder); + //create collection first, to be able to do some sorting + for (unsigned int i = 0; i < jet.numberOfDaughters(); i++){ + const pat::PackedCandidate* PackedCandidate = dynamic_cast(jet.daughter(i)); + if(PackedCandidate){ + if(PackedCandidate->pt() < min_candidate_pt_) continue; + if(PackedCandidate->charge()!=0){ + trackinfo.buildTrackInfo(PackedCandidate,jetDir,jetRefTrackDir,pv); + sortedcharged.push_back(sorting::sortingClass + (i, trackinfo.getTrackSip2dSig(), + -mindrsvpfcand(PackedCandidate), PackedCandidate->pt()/jet_uncorr_pt)); + } + /* else{ + sortedneutrals.push_back(sorting::sortingClass + (i, -1, -mindrsvpfcand(PackedCandidate), PackedCandidate->pt()/jet_uncorr_pt)); + }*/ + } + } + std::sort(sortedcharged.begin(),sortedcharged.end(),sorting::sortingClass::compareByABCInv); + n_Cpfcand2_ = std::min(sortedcharged.size(),max_pfcand_); + + //std::sort(sortedneutrals.begin(),sortedneutrals.end(),sorting::sortingClass::compareByABCInv); + std::vector sortedchargedindices; //,sortedneutralsindices; + //n_Npfcand2_ = std::min(sortedneutrals.size(),max_pfcand_); + sortedchargedindices=sorting::invertSortingVector(sortedcharged); + //sortedneutralsindices=sorting::invertSortingVector(sortedneutrals); + + size_t counter = 0; + int n_cpf_ = std::min((int)25, n_Cpfcand2_); + + for (int i = 0; i < n_cpf_; i++){ + for (int j = 0; j < i; j++){ + + deepntuples::TrackPairInfoBuilder trkpairinfo; + int ind_i = sortedcharged.at(i).get(); + int ind_j = sortedcharged.at(j).get(); + const pat::PackedCandidate* Part_i_ = dynamic_cast(jet.daughter(ind_i)); + + if(!Part_i_){ + std::cout << i << " Bug PackedCandidate " << j << std::endl; + } + const pat::PackedCandidate* Part_j_ = dynamic_cast(jet.daughter(ind_j)); + if(!Part_j_){ + std::cout << i << " Bug PackedCandidate " << j << std::endl; + } + + trackinfo.buildTrackInfo(Part_i_,jetDir,jetRefTrackDir,pv); + const reco::TransientTrack it = trackinfo.getTTrack(); + trackinfo.buildTrackInfo(Part_j_,jetDir,jetRefTrackDir,pv); + const reco::TransientTrack tt = trackinfo.getTTrack(); + + trkpairinfo.buildTrackPairInfo(it,tt,vertices()->at(0),jet); + + //const reco::Candidate * pruned_part_match1 = Part_i_.lastPrunedRef().get(); + //const reco::Candidate * pruned_part_match2 = Part_j_.lastPrunedRef().get(); + float dist_vtx_12 = -1.0; //sqrt((pruned_part_match1->vertex()- pruned_part_match2->vertex()).mag2()); + + pair_pca_distance_[counter] = trkpairinfo.pca_distance(); + pair_pca_significance_[counter] = trkpairinfo.pca_significance(); + + pair_pcaSeed_x1_[counter] = trkpairinfo.pcaSeed_x(); + pair_pcaSeed_y1_[counter] = trkpairinfo.pcaSeed_y(); + pair_pcaSeed_z1_[counter] = trkpairinfo.pcaSeed_z(); + + pair_pcaSeed_x2_[counter] = trkpairinfo.pcaTrack_x(); + pair_pcaSeed_y2_[counter] = trkpairinfo.pcaTrack_y(); + pair_pcaSeed_z2_[counter] = trkpairinfo.pcaTrack_z(); + + pair_pcaSeed_xerr1_[counter] = trkpairinfo.pcaSeed_xerr(); + pair_pcaSeed_yerr1_[counter] = trkpairinfo.pcaSeed_yerr(); + pair_pcaSeed_zerr1_[counter] = trkpairinfo.pcaSeed_zerr(); + + pair_pcaSeed_xerr2_[counter] = trkpairinfo.pcaTrack_xerr(); + pair_pcaSeed_yerr2_[counter] = trkpairinfo.pcaTrack_yerr(); + pair_pcaSeed_zerr2_[counter] = trkpairinfo.pcaTrack_zerr(); + + pair_dotprod1_[counter] = trkpairinfo.dotprodTrack(); + pair_dotprod2_[counter] = trkpairinfo.dotprodSeed(); + + pair_pca_dist1_[counter] = trkpairinfo.pcaSeed_dist(); + pair_pca_dist2_[counter] = trkpairinfo.pcaTrack_dist(); + + pair_dotprod12_2D_[counter] = trkpairinfo.dotprodTrackSeed2D(); + pair_dotprod12_2DV_[counter] = trkpairinfo.dotprodTrackSeed2DV(); + pair_dotprod12_3D_[counter] = trkpairinfo.dotprodTrackSeed3D(); + pair_dotprod12_3DV_[counter] = trkpairinfo.dotprodTrackSeed3DV(); + + pair_pca_jetAxis_dist_[counter] = trkpairinfo.pca_jetAxis_dist(); + pair_pca_jetAxis_dotprod_[counter] = trkpairinfo.pca_jetAxis_dotprod(); + pair_pca_jetAxis_dEta_[counter] = trkpairinfo.pca_jetAxis_dEta(); + pair_pca_jetAxis_dPhi_[counter] = trkpairinfo.pca_jetAxis_dPhi(); + + pfcand_dist_vtx_12_[counter] = dist_vtx_12; + + counter++; + } + } + nCpfpairs_ = counter; + n_Cpfpairs_ = counter; + + return true; //for making cuts +} + +float ntuple_pairwise::mindrsvpfcand(const pat::PackedCandidate* pfcand) { + + float mindr_ = jetradius_; + for (unsigned int i=0; isize(); ++i) { + if(!pfcand) continue; + //if(!svs.at(i)) continue; + float tempdr_ = reco::deltaR(secVertices()->at(i),*pfcand); + if (tempdr_ & build): - builder(build), - trackMomentum_(0), - trackEta_(0), - trackEtaRel_(0), - trackPtRel_(0), - trackPPar_(0), - trackDeltaR_(0), - trackPtRatio_(0), - trackPParRatio_(0), - trackSip2dVal_(0), - trackSip2dSig_(0), - trackSip3dVal_(0), - trackSip3dSig_(0), - - trackJetDistVal_(0), - trackJetDistSig_(0) -{ - - -} - - void buildTrackInfo(const pat::PackedCandidate* PackedCandidate_ ,const math::XYZVector& jetDir, GlobalVector refjetdirection, const reco::Vertex & pv){ - TVector3 jetDir3(jetDir.x(),jetDir.y(),jetDir.z()); - if(!PackedCandidate_->hasTrackDetails()) { - TVector3 trackMom3( - PackedCandidate_->momentum().x(), - PackedCandidate_->momentum().y(), - PackedCandidate_->momentum().z() - ); - trackMomentum_=PackedCandidate_->p(); - trackEta_= PackedCandidate_->eta(); - trackEtaRel_=reco::btau::etaRel(jetDir, PackedCandidate_->momentum()); - trackPtRel_=trackMom3.Perp(jetDir3); - trackPPar_=jetDir.Dot(PackedCandidate_->momentum()); - trackDeltaR_=reco::deltaR(PackedCandidate_->momentum(), jetDir); - trackPtRatio_=trackMom3.Perp(jetDir3) / PackedCandidate_->p(); - trackPParRatio_=jetDir.Dot(PackedCandidate_->momentum()) / PackedCandidate_->p(); - trackSip2dVal_=0.; - trackSip2dSig_=0.; - trackSip3dVal_=0.; - trackSip3dSig_=0.; - trackJetDistVal_=0.; - trackJetDistSig_=0.; - return; - } - - const reco::Track & PseudoTrack = PackedCandidate_->pseudoTrack(); - - reco::TransientTrack transientTrack; - transientTrack=builder->build(PseudoTrack); - Measurement1D meas_ip2d=IPTools::signedTransverseImpactParameter(transientTrack, refjetdirection, pv).second; - Measurement1D meas_ip3d=IPTools::signedImpactParameter3D(transientTrack, refjetdirection, pv).second; - Measurement1D jetdist=IPTools::jetTrackDistance(transientTrack, refjetdirection, pv).second; - math::XYZVector trackMom = PseudoTrack.momentum(); - double trackMag = std::sqrt(trackMom.Mag2()); - TVector3 trackMom3(trackMom.x(),trackMom.y(),trackMom.z()); - - - trackMomentum_=std::sqrt(trackMom.Mag2()); - trackEta_= trackMom.Eta(); - trackEtaRel_=reco::btau::etaRel(jetDir, trackMom); - trackPtRel_=trackMom3.Perp(jetDir3); - trackPPar_=jetDir.Dot(trackMom); - trackDeltaR_=reco::deltaR(trackMom, jetDir); - trackPtRatio_=trackMom3.Perp(jetDir3) / trackMag; - trackPParRatio_=jetDir.Dot(trackMom) / trackMag; - trackSip2dVal_=(meas_ip2d.value()); - - trackSip2dSig_=(meas_ip2d.significance()); - trackSip3dVal_=(meas_ip3d.value()); - - - trackSip3dSig_=meas_ip3d.significance(); - trackJetDistVal_= jetdist.value(); - trackJetDistSig_= jetdist.significance(); - + TrackInfoBuilder(edm::ESHandle & build): + builder(build), + trackMomentum_(0), + trackEta_(0), + trackEtaRel_(0), + trackPtRel_(0), + trackPPar_(0), + trackDeltaR_(0), + trackPtRatio_(0), + trackPParRatio_(0), + trackSip2dVal_(0), + trackSip2dSig_(0), + trackSip3dVal_(0), + trackSip3dSig_(0), + + trackJetDecayLen_(0), + trackJetDistVal_(0), + trackJetDistSig_(0), + ttrack_(0) + { + + + } + + void buildTrackInfo(const pat::PackedCandidate* PackedCandidate_ ,const math::XYZVector& jetDir, GlobalVector refjetdirection, const reco::Vertex & pv){ + TVector3 jetDir3(jetDir.x(),jetDir.y(),jetDir.z()); + if(!PackedCandidate_->hasTrackDetails()) { + TVector3 trackMom3( + PackedCandidate_->momentum().x(), + PackedCandidate_->momentum().y(), + PackedCandidate_->momentum().z() + ); + trackMomentum_=PackedCandidate_->p(); + trackEta_= PackedCandidate_->eta(); + trackEtaRel_=reco::btau::etaRel(jetDir, PackedCandidate_->momentum()); + trackPtRel_=trackMom3.Perp(jetDir3); + trackPPar_=jetDir.Dot(PackedCandidate_->momentum()); + trackDeltaR_=reco::deltaR(PackedCandidate_->momentum(), jetDir); + trackPtRatio_=trackMom3.Perp(jetDir3) / PackedCandidate_->p(); + trackPParRatio_=jetDir.Dot(PackedCandidate_->momentum()) / PackedCandidate_->p(); + trackSip2dVal_=0.; + trackSip2dSig_=0.; + trackSip3dVal_=0.; + trackSip3dSig_=0.; + trackJetDecayLen_=0.; + trackJetDistVal_=0.; + trackJetDistSig_=0.; + return; } - const float& getTrackDeltaR() const {return trackDeltaR_;} - const float& getTrackEta() const {return trackEta_;} - const float& getTrackEtaRel() const {return trackEtaRel_;} - const float& getTrackJetDistSig() const {return trackJetDistSig_;} - const float& getTrackJetDistVal() const {return trackJetDistVal_;} - const float& getTrackMomentum() const {return trackMomentum_;} - const float& getTrackPPar() const {return trackPPar_;} - const float& getTrackPParRatio() const {return trackPParRatio_;} - const float& getTrackPtRatio() const {return trackPtRatio_;} - const float& getTrackPtRel() const {return trackPtRel_;} - const float& getTrackSip2dSig() const {return trackSip2dSig_;} - const float& getTrackSip2dVal() const {return trackSip2dVal_;} - const float& getTrackSip3dSig() const {return trackSip3dSig_;} - const float& getTrackSip3dVal() const {return trackSip3dVal_;} + const reco::Track & PseudoTrack = PackedCandidate_->pseudoTrack(); + + reco::TransientTrack transientTrack; + transientTrack=builder->build(PseudoTrack); + Measurement1D meas_ip2d=IPTools::signedTransverseImpactParameter(transientTrack, refjetdirection, pv).second; + Measurement1D meas_ip3d=IPTools::signedImpactParameter3D(transientTrack, refjetdirection, pv).second; + Measurement1D jetdist=IPTools::jetTrackDistance(transientTrack, refjetdirection, pv).second; + Measurement1D decayl = IPTools::signedDecayLength3D(transientTrack, refjetdirection, pv).second; + math::XYZVector trackMom = PseudoTrack.momentum(); + double trackMag = std::sqrt(trackMom.Mag2()); + TVector3 trackMom3(trackMom.x(),trackMom.y(),trackMom.z()); + + + trackMomentum_=std::sqrt(trackMom.Mag2()); + trackEta_= trackMom.Eta(); + trackEtaRel_=reco::btau::etaRel(jetDir, trackMom); + trackPtRel_=trackMom3.Perp(jetDir3); + trackPPar_=jetDir.Dot(trackMom); + trackDeltaR_=reco::deltaR(trackMom, jetDir); + trackPtRatio_=trackMom3.Perp(jetDir3) / trackMag; + trackPParRatio_=jetDir.Dot(trackMom) / trackMag; + + trackSip2dVal_=(meas_ip2d.value()); + trackSip2dSig_=(meas_ip2d.significance()); + trackSip3dVal_=(meas_ip3d.value()); + trackSip3dSig_=meas_ip3d.significance(); + + trackJetDecayLen_= decayl.value(); + trackJetDistVal_= jetdist.value(); + trackJetDistSig_= jetdist.significance(); + + ttrack_ = transientTrack; + + } + + const float& getTrackDeltaR() const {return trackDeltaR_;} + const float& getTrackEta() const {return trackEta_;} + const float& getTrackEtaRel() const {return trackEtaRel_;} + const float& getTrackJetDecayLen() const {return trackJetDecayLen_;} + const float& getTrackJetDistSig() const {return trackJetDistSig_;} + const float& getTrackJetDistVal() const {return trackJetDistVal_;} + const float& getTrackMomentum() const {return trackMomentum_;} + const float& getTrackPPar() const {return trackPPar_;} + const float& getTrackPParRatio() const {return trackPParRatio_;} + const float& getTrackPtRatio() const {return trackPtRatio_;} + const float& getTrackPtRel() const {return trackPtRel_;} + const float& getTrackSip2dSig() const {return trackSip2dSig_;} + const float& getTrackSip2dVal() const {return trackSip2dVal_;} + const float& getTrackSip3dSig() const {return trackSip3dSig_;} + const float& getTrackSip3dVal() const {return trackSip3dVal_;} + const reco::TransientTrack getTTrack() const {return ttrack_;} private: edm::ESHandle& builder; - float trackMomentum_; - float trackEta_; - float trackEtaRel_; - float trackPtRel_; - float trackPPar_; - float trackDeltaR_; - float trackPtRatio_; - float trackPParRatio_; - float trackSip2dVal_; - float trackSip2dSig_; - float trackSip3dVal_; - float trackSip3dSig_; - - float trackJetDistVal_; - float trackJetDistSig_; + float trackMomentum_; + float trackEta_; + float trackEtaRel_; + float trackPtRel_; + float trackPPar_; + float trackDeltaR_; + float trackPtRatio_; + float trackPParRatio_; + float trackSip2dVal_; + float trackSip2dSig_; + float trackSip3dVal_; + float trackSip3dSig_; + + float trackJetDecayLen_; + float trackJetDistVal_; + float trackJetDistSig_; + reco::TransientTrack ttrack_; }; - - +void ntuple_pfCands::readConfig(const edm::ParameterSet& iConfig, edm::ConsumesCollector&& cc) { + ltToken_ = cc.consumes>(iConfig.getParameter("losttracks")); + packedToken_ = cc.consumes>(iConfig.getParameter("packed")); + prunedToken_ = cc.consumes>(iConfig.getParameter("pruned")); +} void ntuple_pfCands::readSetup(const edm::EventSetup& iSetup){ - iSetup.get().get("TransientTrackBuilder", builder); - } void ntuple_pfCands::getInput(const edm::ParameterSet& iConfig){ @@ -151,123 +170,190 @@ void ntuple_pfCands::getInput(const edm::ParameterSet& iConfig){ void ntuple_pfCands::initBranches(TTree* tree){ - addBranch(tree,"n_Cpfcand", &n_Cpfcand_,"n_Cpfcand_/i"); - - addBranch(tree,"nCpfcand", &nCpfcand_,"nCpfcand_/f"); - - addBranch(tree,"Cpfcan_pt", &Cpfcan_pt_,"Cpfcan_pt_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_eta", &Cpfcan_eta_,"Cpfcan_eta_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_phi", &Cpfcan_phi_,"Cpfcan_phi_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_ptrel", &Cpfcan_ptrel_,"Cpfcan_ptrel_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_erel", &Cpfcan_erel_,"Cpfcan_erel_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_phirel",&Cpfcan_phirel_,"Cpfcan_phirel_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_etarel",&Cpfcan_etarel_,"Cpfcan_etarel_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_deltaR",&Cpfcan_deltaR_,"Cpfcan_deltaR_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_puppiw",&Cpfcan_puppiw_,"Cpfcan_puppiw_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_dxy",&Cpfcan_dxy_,"Cpfcan_dxy_[n_Cpfcand_]/f"); - - addBranch(tree,"Cpfcan_dxyerrinv",&Cpfcan_dxyerrinv_,"Cpfcan_dxyerrinv_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_dxysig",&Cpfcan_dxysig_,"Cpfcan_dxysig_[n_Cpfcand_]/f"); - - addBranch(tree,"Cpfcan_dz",&Cpfcan_dz_,"Cpfcan_dz_[n_Cpfcand_]/f"); - - addBranch(tree,"Cpfcan_VTX_ass",&Cpfcan_VTX_ass_,"Cpfcan_VTX_ass_[n_Cpfcand_]/f"); - - addBranch(tree,"Cpfcan_fromPV",&Cpfcan_fromPV_,"Cpfcan_fromPV_[n_Cpfcand_]/f"); - - addBranch(tree,"Cpfcan_drminsv",&Cpfcan_drminsv_,"Cpfcan_drminsv_[n_Cpfcand_]/f"); - - //commented ones don't work - /**///addBranch(tree,"Cpfcan_vertexChi2",&Cpfcan_vertexChi2_,"Cpfcan_vertexChi2_[n_Cpfcand_]/f"); - /**///addBranch(tree,"Cpfcan_vertexNdof",&Cpfcan_vertexNdof_,"Cpfcan_vertexNdof_[n_Cpfcand_]/f"); - /**///addBranch(tree,"Cpfcan_vertexNormalizedChi2",&Cpfcan_vertexNormalizedChi2_,"Cpfcan_vertexNormalizedChi2_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_vertex_rho",&Cpfcan_vertex_rho_,"Cpfcan_vertex_rho_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_vertex_phirel",&Cpfcan_vertex_phirel_,"Cpfcan_vertex_phirel_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_vertex_etarel",&Cpfcan_vertex_etarel_,"Cpfcan_vertex_etarel_[n_Cpfcand_]/f"); - /**///addBranch(tree,"Cpfcan_vertexRef_mass",&Cpfcan_vertexRef_mass_,"Cpfcan_vertexRef_mass_[n_Cpfcand_]/f"); - - /* - addBranch(tree,"Cpfcan_dptdpt",&Cpfcan_dptdpt_,"Cpfcan_dptdpt_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_detadeta",&Cpfcan_detadeta_,"Cpfcan_detadeta_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_dphidphi",&Cpfcan_dphidphi_,"Cpfcan_dphidphi_[n_Cpfcand_]/f"); - - - addBranch(tree,"Cpfcan_dxydxy",&Cpfcan_dxydxy_,"Cpfcan_dxydxy_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_dzdz",&Cpfcan_dzdz_,"Cpfcan_dzdz_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_dxydz",&Cpfcan_dxydz_,"Cpfcan_dxydz_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_dphidxy",&Cpfcan_dphidxy_,"Cpfcan_dphidxy_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_dlambdadz",&Cpfcan_dlambdadz_,"Cpfcan_dlambdadz_[n_Cpfcand_]/f"); - */ - - - - addBranch(tree,"Cpfcan_BtagPf_trackMomentum",&Cpfcan_BtagPf_trackMomentum_,"Cpfcan_BtagPf_trackMomentum_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackEta",&Cpfcan_BtagPf_trackEta_,"Cpfcan_BtagPf_trackEta_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackEtaRel",&Cpfcan_BtagPf_trackEtaRel_,"Cpfcan_BtagPf_trackEtaRel_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackPtRel",&Cpfcan_BtagPf_trackPtRel_,"Cpfcan_BtagPf_trackPtRel_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackPPar",&Cpfcan_BtagPf_trackPPar_,"Cpfcan_BtagPf_trackPPar_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackDeltaR",&Cpfcan_BtagPf_trackDeltaR_,"Cpfcan_BtagPf_trackDeltaR_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackPtRatio",&Cpfcan_BtagPf_trackPtRatio_,"Cpfcan_BtagPf_trackPtRatio_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackPParRatio",&Cpfcan_BtagPf_trackPParRatio_,"Cpfcan_BtagPf_trackPParRatio[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackSip3dVal",&Cpfcan_BtagPf_trackSip3dVal_,"Cpfcan_BtagPf_trackSip3dVal_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackSip3dSig",&Cpfcan_BtagPf_trackSip3dSig_,"Cpfcan_BtagPf_trackSip3dSig_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackSip2dVal",&Cpfcan_BtagPf_trackSip2dVal_,"Cpfcan_BtagPf_trackSip2dVal_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackSip2dSig",&Cpfcan_BtagPf_trackSip2dSig_,"Cpfcan_BtagPf_trackSip2dSig_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackDecayLen",&Cpfcan_BtagPf_trackDecayLen_,"Cpfcan_BtagPf_trackDecayLen_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackJetDistVal",&Cpfcan_BtagPf_trackJetDistVal_,"Cpfcan_BtagPf_trackJetDistVal_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_BtagPf_trackJetDistSig",&Cpfcan_BtagPf_trackJetDistSig_,"Cpfcan_BtagPf_trackJetDistSig_[n_Cpfcand_]/f"); - - - - - addBranch(tree,"Cpfcan_isMu",&Cpfcan_isMu_,"Cpfcan_isMu_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_isEl",&Cpfcan_isEl_,"Cpfcan_isEl_[n_Cpfcand_]/f"); - - //in16 conversion broken - addBranch(tree,"Cpfcan_lostInnerHits",&Cpfcan_lostInnerHits_,"Cpfcan_lostInnerHits_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_numberOfPixelHits",&Cpfcan_numberOfPixelHits_,"Cpfcan_numberOfPixelHits_[n_Cpfcand_]/f"); - - addBranch(tree,"Cpfcan_chi2",&Cpfcan_chi2_,"Cpfcan_chi2_[n_Cpfcand_]/f"); - addBranch(tree,"Cpfcan_quality",&Cpfcan_quality_,"Cpfcan_quality_[n_Cpfcand_]/f"); - - // did not give integers !! - addBranch(tree,"Cpfcan_charge",&Cpfcan_charge_,"Cpfcan_charge_[n_Cpfcand_]/f"); - - //Neutral Pf candidates - addBranch(tree,"n_Npfcand", &n_Npfcand_,"n_Npfcand_/i"); - addBranch(tree,"nNpfcand", &nNpfcand_,"nNpfcand/f"); - - addBranch(tree,"Npfcan_pt", &Npfcan_pt_,"Npfcan_pt_[n_Npfcand_]/f"); - addBranch(tree,"Npfcan_eta", &Npfcan_eta_,"Npfcan_eta_[n_Npfcand_]/f"); - addBranch(tree,"Npfcan_phi", &Npfcan_phi_,"Npfcan_phi_[n_Npfcand_]/f"); - addBranch(tree,"Npfcan_ptrel", &Npfcan_ptrel_,"Npfcan_ptrel_[n_Npfcand_]/f"); - addBranch(tree,"Npfcan_erel", &Npfcan_erel_,"Npfcan_erel_[n_Npfcand_]/f"); - - addBranch(tree,"Npfcan_puppiw", &Npfcan_puppiw_,"Npfcan_puppiw_[n_Npfcand_]/f"); - - - addBranch(tree,"Npfcan_phirel",&Npfcan_phirel_,"Npfcan_phirel_[n_Npfcand_]/f"); - addBranch(tree,"Npfcan_etarel",&Npfcan_etarel_,"Npfcan_etarel_[n_Npfcand_]/f"); - addBranch(tree,"Npfcan_deltaR",&Npfcan_deltaR_,"Npfcan_deltaR_[n_Npfcand_]/f"); - addBranch(tree,"Npfcan_isGamma",&Npfcan_isGamma_,"Npfcan_isGamma_[n_Npfcand_]/f"); - addBranch(tree,"Npfcan_HadFrac",&Npfcan_HadFrac_,"Npfcan_HadFrac_[n_Npfcand_]/f"); - addBranch(tree,"Npfcan_drminsv",&Npfcan_drminsv_,"Npfcan_drminsv_[n_Npfcand_]/f"); - + addBranch(tree,"n_Cpfcand", &n_Cpfcand_,"n_Cpfcand_/I"); + + addBranch(tree,"nCpfcand", &nCpfcand_,"nCpfcand_/F"); + + addBranch(tree,"Cpfcan_pt", &Cpfcan_pt_,"Cpfcan_pt_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_eta", &Cpfcan_eta_,"Cpfcan_eta_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_phi", &Cpfcan_phi_,"Cpfcan_phi_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_ptrel", &Cpfcan_ptrel_,"Cpfcan_ptrel_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_e", &Cpfcan_e_,"Cpfcan_e_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_erel", &Cpfcan_erel_,"Cpfcan_erel_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_phirel",&Cpfcan_phirel_,"Cpfcan_phirel_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_etarel",&Cpfcan_etarel_,"Cpfcan_etarel_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_deltaR",&Cpfcan_deltaR_,"Cpfcan_deltaR_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_puppiw",&Cpfcan_puppiw_,"Cpfcan_puppiw_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_dxy",&Cpfcan_dxy_,"Cpfcan_dxy_[n_Cpfcand_]/F"); + + addBranch(tree,"Cpfcan_dxyerrinv",&Cpfcan_dxyerrinv_,"Cpfcan_dxyerrinv_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_dxysig",&Cpfcan_dxysig_,"Cpfcan_dxysig_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_dz",&Cpfcan_dz_,"Cpfcan_dz_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_VTX_ass",&Cpfcan_VTX_ass_,"Cpfcan_VTX_ass_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_fromPV",&Cpfcan_fromPV_,"Cpfcan_fromPV_[n_Cpfcand_]/F"); + + addBranch(tree,"Cpfcan_drminsv",&Cpfcan_drminsv_,"Cpfcan_drminsv_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_distminsvold",&Cpfcan_distminsvold_,"Cpfcan_distminsvold_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_distminsv",&Cpfcan_distminsv_,"Cpfcan_distminsv_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_distminsv2",&Cpfcan_distminsv2_,"Cpfcan_distminsv2_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_dxminsv",&Cpfcan_dxminsv_,"Cpfcan_dxminsv_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_dyminsv",&Cpfcan_dyminsv_,"Cpfcan_dyminsv_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_dzminsv",&Cpfcan_dzminsv_,"Cpfcan_dzminsv_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_dxpv",&Cpfcan_dxpv_,"Cpfcan_dxpv_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_dypv",&Cpfcan_dypv_,"Cpfcan_dypv_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_dzpv",&Cpfcan_dzpv_,"Cpfcan_dzpv_[n_Cpfcand_]/F"); + + //commented ones don't work + addBranch(tree,"Cpfcan_vertex_rho",&Cpfcan_vertex_rho_,"Cpfcan_vertex_rho_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_vertex_phirel",&Cpfcan_vertex_phirel_,"Cpfcan_vertex_phirel_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_vertex_etarel",&Cpfcan_vertex_etarel_,"Cpfcan_vertex_etarel_[n_Cpfcand_]/F"); + + addBranch(tree,"Cpfcan_BtagPf_trackMomentum",&Cpfcan_BtagPf_trackMomentum_,"Cpfcan_BtagPf_trackMomentum_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackEta",&Cpfcan_BtagPf_trackEta_,"Cpfcan_BtagPf_trackEta_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackEtaRel",&Cpfcan_BtagPf_trackEtaRel_,"Cpfcan_BtagPf_trackEtaRel_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackPtRel",&Cpfcan_BtagPf_trackPtRel_,"Cpfcan_BtagPf_trackPtRel_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackPPar",&Cpfcan_BtagPf_trackPPar_,"Cpfcan_BtagPf_trackPPar_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackDeltaR",&Cpfcan_BtagPf_trackDeltaR_,"Cpfcan_BtagPf_trackDeltaR_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackPtRatio",&Cpfcan_BtagPf_trackPtRatio_,"Cpfcan_BtagPf_trackPtRatio_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackPParRatio",&Cpfcan_BtagPf_trackPParRatio_,"Cpfcan_BtagPf_trackPParRatio[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackSip3dVal",&Cpfcan_BtagPf_trackSip3dVal_,"Cpfcan_BtagPf_trackSip3dVal_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackSip3dSig",&Cpfcan_BtagPf_trackSip3dSig_,"Cpfcan_BtagPf_trackSip3dSig_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackSip2dVal",&Cpfcan_BtagPf_trackSip2dVal_,"Cpfcan_BtagPf_trackSip2dVal_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackSip2dSig",&Cpfcan_BtagPf_trackSip2dSig_,"Cpfcan_BtagPf_trackSip2dSig_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackDecayLen",&Cpfcan_BtagPf_trackDecayLen_,"Cpfcan_BtagPf_trackDecayLen_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackJetDistVal",&Cpfcan_BtagPf_trackJetDistVal_,"Cpfcan_BtagPf_trackJetDistVal_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_BtagPf_trackJetDistSig",&Cpfcan_BtagPf_trackJetDistSig_,"Cpfcan_BtagPf_trackJetDistSig_[n_Cpfcand_]/F"); + + addBranch(tree,"Cpfcan_isMu",&Cpfcan_isMu_,"Cpfcan_isMu_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_isEl",&Cpfcan_isEl_,"Cpfcan_isEl_[n_Cpfcand_]/F"); + // did not give integers !! + addBranch(tree,"Cpfcan_charge",&Cpfcan_charge_,"Cpfcan_charge_[n_Cpfcand_]/F"); + + //in16 conversion broken + addBranch(tree,"Cpfcan_lostInnerHits",&Cpfcan_lostInnerHits_,"Cpfcan_lostInnerHits_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_numberOfPixelHits",&Cpfcan_numberOfPixelHits_,"Cpfcan_numberOfPixelHits_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_chi2",&Cpfcan_chi2_,"Cpfcan_chi2_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_quality",&Cpfcan_quality_,"Cpfcan_quality_[n_Cpfcand_]/F"); + //hit pattern variables, as defined here https://github.com/cms-sw/cmssw/blob/master/DataFormats/TrackReco/interface/HitPattern.h + //Tracker per layer + //Pixel barrel + addBranch(tree,"Cpfcan_nhitpixelBarrelLayer1",&Cpfcan_nhitpixelBarrelLayer1_,"Cpfcan_nhitpixelBarrelLayer1_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitpixelBarrelLayer2",&Cpfcan_nhitpixelBarrelLayer2_,"Cpfcan_nhitpixelBarrelLayer2_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitpixelBarrelLayer3",&Cpfcan_nhitpixelBarrelLayer3_,"Cpfcan_nhitpixelBarrelLayer3_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitpixelBarrelLayer4",&Cpfcan_nhitpixelBarrelLayer4_,"Cpfcan_nhitpixelBarrelLayer4_[n_Cpfcand_]/F"); + //Pixel Endcap + addBranch(tree,"Cpfcan_nhitpixelEndcapLayer1",&Cpfcan_nhitpixelEndcapLayer1_,"Cpfcan_nhitpixelEndcapLayer1_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitpixelEndcapLayer2",&Cpfcan_nhitpixelEndcapLayer2_,"Cpfcan_nhitpixelEndcapLayer2_[n_Cpfcand_]/F"); + //Strip TIB + addBranch(tree,"Cpfcan_nhitstripTIBLayer1",&Cpfcan_nhitstripTIBLayer1_,"Cpfcan_nhitstripTIBLayer1_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTIBLayer2",&Cpfcan_nhitstripTIBLayer2_,"Cpfcan_nhitstripTIBLayer2_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTIBLayer3",&Cpfcan_nhitstripTIBLayer3_,"Cpfcan_nhitstripTIBLayer3_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTIBLayer4",&Cpfcan_nhitstripTIBLayer4_,"Cpfcan_nhitstripTIBLayer4_[n_Cpfcand_]/F"); + //Strip TID + addBranch(tree,"Cpfcan_nhitstripTIDLayer1",&Cpfcan_nhitstripTIDLayer1_,"Cpfcan_nhitstripTIDLayer1_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTIDLayer2",&Cpfcan_nhitstripTIDLayer2_,"Cpfcan_nhitstripTIDLayer2_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTIDLayer3",&Cpfcan_nhitstripTIDLayer3_,"Cpfcan_nhitstripTIDLayer3_[n_Cpfcand_]/F"); + //Strip TOB + addBranch(tree,"Cpfcan_nhitstripTOBLayer1",&Cpfcan_nhitstripTOBLayer1_,"Cpfcan_nhitstripTOBLayer1_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTOBLayer2",&Cpfcan_nhitstripTOBLayer2_,"Cpfcan_nhitstripTOBLayer2_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTOBLayer3",&Cpfcan_nhitstripTOBLayer3_,"Cpfcan_nhitstripTOBLayer3_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTOBLayer4",&Cpfcan_nhitstripTOBLayer4_,"Cpfcan_nhitstripTOBLayer4_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTOBLayer5",&Cpfcan_nhitstripTOBLayer5_,"Cpfcan_nhitstripTOBLayer5_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTOBLayer6",&Cpfcan_nhitstripTOBLayer6_,"Cpfcan_nhitstripTOBLayer6_[n_Cpfcand_]/F"); + //Strip TEC + addBranch(tree,"Cpfcan_nhitstripTECLayer1",&Cpfcan_nhitstripTECLayer1_,"Cpfcan_nhitstripTECLayer1_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTECLayer2",&Cpfcan_nhitstripTECLayer2_,"Cpfcan_nhitstripTECLayer2_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTECLayer3",&Cpfcan_nhitstripTECLayer3_,"Cpfcan_nhitstripTECLayer3_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTECLayer4",&Cpfcan_nhitstripTECLayer4_,"Cpfcan_nhitstripTECLayer4_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTECLayer5",&Cpfcan_nhitstripTECLayer5_,"Cpfcan_nhitstripTECLayer5_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTECLayer6",&Cpfcan_nhitstripTECLayer6_,"Cpfcan_nhitstripTECLayer6_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTECLayer7",&Cpfcan_nhitstripTECLayer7_,"Cpfcan_nhitstripTECLayer7_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTECLayer8",&Cpfcan_nhitstripTECLayer8_,"Cpfcan_nhitstripTECLayer8_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_nhitstripTECLayer9",&Cpfcan_nhitstripTECLayer9_,"Cpfcan_nhitstripTECLayer9_[n_Cpfcand_]/F"); + //Tracker all layers together + //Valid hits + addBranch(tree,"Cpfcan_numberOfValidHits",&Cpfcan_numberOfValidHits_,"Cpfcan_numberOfValidHits_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_numberOfValidTrackerHits",&Cpfcan_numberOfValidTrackerHits_,"Cpfcan_numberOfValidTrackerHits_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_numberOfValidPixelHits",&Cpfcan_numberOfValidPixelHits_,"Cpfcan_numberOfValidPixelHits_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_numberOfValidPixelBarrelHits",&Cpfcan_numberOfValidPixelBarrelHits_,"Cpfcan_numberOfValidPixelBarrelHits_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_numberOfValidPixelEndcapHits",&Cpfcan_numberOfValidPixelEndcapHits_,"Cpfcan_numberOfValidPixelEndcapHits_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_numberOfValidStripHits",&Cpfcan_numberOfValidStripHits_,"Cpfcan_numberOfValidStripHits_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_numberOfValidStripTIBHits",&Cpfcan_numberOfValidStripTIBHits_,"Cpfcan_numberOfValidStripTIBHits_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_numberOfValidStripTIDHits",&Cpfcan_numberOfValidStripTIDHits_,"Cpfcan_numberOfValidStripTIDHits_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_numberOfValidStripTOBHits",&Cpfcan_numberOfValidStripTOBHits_,"Cpfcan_numberOfValidStripTOBHits_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_numberOfValidStripTECHits",&Cpfcan_numberOfValidStripTECHits_,"Cpfcan_numberOfValidStripTECHits_[n_Cpfcand_]/F"); + //LayersWithMeasuremen + addBranch(tree,"Cpfcan_trackerLayersWithMeasurementOld",&Cpfcan_trackerLayersWithMeasurementOld_,"Cpfcan_trackerLayersWithMeasurementOld_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_trackerLayersWithMeasurement",&Cpfcan_trackerLayersWithMeasurement_,"Cpfcan_trackerLayersWithMeasurement_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_pixelLayersWithMeasurementOld",&Cpfcan_pixelLayersWithMeasurementOld_,"Cpfcan_pixelLayersWithMeasurementOld_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_pixelLayersWithMeasurement",&Cpfcan_pixelLayersWithMeasurement_,"Cpfcan_pixelLayersWithMeasurement_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_stripLayersWithMeasurement",&Cpfcan_stripLayersWithMeasurement_,"Cpfcan_stripLayersWithMeasurement_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_pixelBarrelLayersWithMeasurement",&Cpfcan_pixelBarrelLayersWithMeasurement_,"Cpfcan_pixelBarrelLayersWithMeasurement_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_pixelEndcapLayersWithMeasurement",&Cpfcan_pixelEndcapLayersWithMeasurement_,"Cpfcan_pixelEndcapLayersWithMeasurement_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_stripTIBLayersWithMeasurement",&Cpfcan_stripTIBLayersWithMeasurement_,"Cpfcan_stripTIBLayersWithMeasurement_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_stripTIDLayersWithMeasurement",&Cpfcan_stripTIDLayersWithMeasurement_,"Cpfcan_stripTIDLayersWithMeasurement_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_stripTOBLayersWithMeasurement",&Cpfcan_stripTOBLayersWithMeasurement_,"Cpfcan_stripTOBLayersWithMeasurement_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_stripTECLayersWithMeasurement",&Cpfcan_stripTECLayersWithMeasurement_,"Cpfcan_stripTECLayersWithMeasurement_[n_Cpfcand_]/F"); + //Null + addBranch(tree,"Cpfcan_trackerLayersNull",&Cpfcan_trackerLayersNull_,"Cpfcan_trackerLayersNull_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_pixelLayersNull",&Cpfcan_pixelLayersNull_,"Cpfcan_pixelLayersNull_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_stripLayersNull",&Cpfcan_stripLayersNull_,"Cpfcan_stripLayersNull_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_pixelBarrelLayersNull",&Cpfcan_pixelBarrelLayersNull_,"Cpfcan_pixelBarrelLayersNull_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_pixelEndcapLayersNull",&Cpfcan_pixelEndcapLayersNull_,"Cpfcan_pixelEndcapLayersNull_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_stripTIBLayersNull",&Cpfcan_stripTIBLayersNull_,"Cpfcan_stripTIBLayersNull_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_stripTIDLayersNull",&Cpfcan_stripTIDLayersNull_,"Cpfcan_stripTIDLayersNull_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_stripTOBLayersNull",&Cpfcan_stripTOBLayersNull_,"Cpfcan_stripTOBLayersNull_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_stripTECLayersNull",&Cpfcan_stripTECLayersNull_,"Cpfcan_stripTECLayersNull_[n_Cpfcand_]/F"); + + //Neutral Pf candidates + addBranch(tree,"n_Npfcand", &n_Npfcand_,"n_Npfcand_/I"); + addBranch(tree,"nNpfcand", &nNpfcand_,"nNpfcand/F"); + + addBranch(tree,"Npfcan_pt", &Npfcan_pt_,"Npfcan_pt_[n_Npfcand_]/F"); + addBranch(tree,"Npfcan_eta", &Npfcan_eta_,"Npfcan_eta_[n_Npfcand_]/F"); + addBranch(tree,"Npfcan_phi", &Npfcan_phi_,"Npfcan_phi_[n_Npfcand_]/F"); + addBranch(tree,"Npfcan_ptrel", &Npfcan_ptrel_,"Npfcan_ptrel_[n_Npfcand_]/F"); + addBranch(tree,"Npfcan_e", &Npfcan_e_,"Npfcan_e_[n_Npfcand_]/F"); + addBranch(tree,"Npfcan_erel", &Npfcan_erel_,"Npfcan_erel_[n_Npfcand_]/F"); + + addBranch(tree,"Npfcan_puppiw", &Npfcan_puppiw_,"Npfcan_puppiw_[n_Npfcand_]/F"); + + addBranch(tree,"Npfcan_phirel",&Npfcan_phirel_,"Npfcan_phirel_[n_Npfcand_]/F"); + addBranch(tree,"Npfcan_etarel",&Npfcan_etarel_,"Npfcan_etarel_[n_Npfcand_]/F"); + addBranch(tree,"Npfcan_deltaR",&Npfcan_deltaR_,"Npfcan_deltaR_[n_Npfcand_]/F"); + addBranch(tree,"Npfcan_isGamma",&Npfcan_isGamma_,"Npfcan_isGamma_[n_Npfcand_]/F"); + addBranch(tree,"Npfcan_HadFrac",&Npfcan_HadFrac_,"Npfcan_HadFrac_[n_Npfcand_]/F"); + addBranch(tree,"Npfcan_drminsv",&Npfcan_drminsv_,"Npfcan_drminsv_[n_Npfcand_]/F"); + + addBranch(tree,"Npfcan_pdgID",&Npfcan_pdgID_,"Npfcan_pdgID_[n_Npfcand_]/F"); + addBranch(tree,"Cpfcan_pdgID",&Cpfcan_pdgID_,"Cpfcan_pdgID_[n_Cpfcand_]/F"); + + addBranch(tree,"Cpfcan_HadFrac",&Cpfcan_HadFrac_,"Cpfcan_HadFrac_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_CaloFrac",&Cpfcan_CaloFrac_,"Cpfcan_CaloFrac_[n_Cpfcand_]/F"); + + addBranch(tree,"Cpfcan_b_tag",&Cpfcan_b_tag_,"Cpfcan_b_tag_[n_Cpfcand_]/I"); + addBranch(tree,"Cpfcan_c_tag",&Cpfcan_c_tag_,"Cpfcan_c_tag_[n_Cpfcand_]/I"); + addBranch(tree,"Cpfcan_g_tag",&Cpfcan_g_tag_,"Cpfcan_g_tag_[n_Cpfcand_]/I"); + + addBranch(tree,"Cpfcan_vtx_x",&Cpfcan_vtx_x_,"Cpfcan_vtx_x_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_vtx_y",&Cpfcan_vtx_y_,"Cpfcan_vtx_y_[n_Cpfcand_]/F"); + addBranch(tree,"Cpfcan_vtx_z",&Cpfcan_vtx_z_,"Cpfcan_vtx_z_[n_Cpfcand_]/F"); + + addBranch(tree,"Cpfcan_dist_from_pv",&Cpfcan_dist_from_pv_,"Cpfcan_dist_from_pv_[n_Cpfcand_]/F"); } void ntuple_pfCands::readEvent(const edm::Event& iEvent){ - + iEvent.getByToken(ltToken_, LTs); + iEvent.getByToken(packedToken_, packed); + iEvent.getByToken(prunedToken_, pruned); n_Npfcand_=0; n_Cpfcand_=0; } - - //use either of these functions - bool ntuple_pfCands::fillBranches(const pat::Jet & jet, const size_t& jetidx, const edm::View * coll){ float etasign = 1.; if (jet.eta()<0) etasign =-1.; @@ -275,6 +361,9 @@ bool ntuple_pfCands::fillBranches(const pat::Jet & jet, const size_t& jetidx, co GlobalVector jetRefTrackDir(jet.px(),jet.py(),jet.pz()); const reco::Vertex & pv = vertices()->at(0); + for(const auto &pruned_part : *pruned){ + if(pruned_part.pdgId()!=2212) { + const auto pv_0 = pruned_part.vertex(); std::vector > sortedcharged, sortedneutrals; @@ -299,14 +388,15 @@ bool ntuple_pfCands::fillBranches(const pat::Jet & jet, const size_t& jetidx, co } } } - std::sort(sortedcharged.begin(),sortedcharged.end(),sorting::sortingClass::compareByABCInv); + + std::sort(sortedcharged.begin(),sortedcharged.end(),sorting::sortingClass::compareByABCInv); n_Cpfcand_ = std::min(sortedcharged.size(),max_pfcand_); std::sort(sortedneutrals.begin(),sortedneutrals.end(),sorting::sortingClass::compareByABCInv); std::vector sortedchargedindices,sortedneutralsindices; n_Npfcand_ = std::min(sortedneutrals.size(),max_pfcand_); - sortedchargedindices=sorting::invertSortingVector(sortedcharged); - sortedneutralsindices=sorting::invertSortingVector(sortedneutrals); + sortedchargedindices=sorting::invertSortingVector(sortedcharged); + sortedneutralsindices=sorting::invertSortingVector(sortedneutrals); for (unsigned int i = 0; i < jet.numberOfDaughters(); i++){ const pat::PackedCandidate* PackedCandidate_ = dynamic_cast(jet.daughter(i)); @@ -317,6 +407,31 @@ bool ntuple_pfCands::fillBranches(const pat::Jet & jet, const size_t& jetidx, co // get the dr with the closest sv float drminpfcandsv_ = mindrsvpfcand(PackedCandidate_); + float pdgid_; + if (abs(PackedCandidate_->pdgId()) == 11 and PackedCandidate_->charge() != 0){ + pdgid_ = 0.0; + } + else if (abs(PackedCandidate_->pdgId()) == 13 and PackedCandidate_->charge() != 0){ + pdgid_ = 1.0; + } + else if (abs(PackedCandidate_->pdgId()) == 22 and PackedCandidate_->charge() == 0){ + pdgid_ = 2.0; + } + else if (abs(PackedCandidate_->pdgId()) != 22 and PackedCandidate_->charge() == 0 and abs(PackedCandidate_->pdgId()) != 1 and abs(PackedCandidate_->pdgId()) != 2){ + pdgid_ = 3.0; + } + else if (abs(PackedCandidate_->pdgId()) != 11 and abs(PackedCandidate_->pdgId()) != 13 and PackedCandidate_->charge() != 0){ + pdgid_ = 4.0; + } + else if (PackedCandidate_->charge() == 0 and abs(PackedCandidate_->pdgId()) == 1){ + pdgid_ = 5.0; + } + else if (PackedCandidate_->charge() == 0 and abs(PackedCandidate_->pdgId()) == 2){ + pdgid_ = 6.0; + } + else{ + pdgid_ = 7.0; + } /// This might include more than PF candidates, e.g. Reco muons and could /// be double counting. Needs to be checked.!!!! @@ -327,12 +442,53 @@ bool ntuple_pfCands::fillBranches(const pat::Jet & jet, const size_t& jetidx, co size_t fillntupleentry= sortedchargedindices.at(i); if(fillntupleentry>=max_pfcand_) continue; + int b_tag=-1, c_tag=-1, g_tag=-1; + float dist_from_pv=-1; + float vtx_x=0, vtx_y=0, vtx_z=0; + + double dR_min=pow10(6); + + const reco::Candidate * pruned_part_match=nullptr; + for (const auto &packed_part : *packed){ + double dR = reco::deltaR(*PackedCandidate_, packed_part); + double dpt = std::abs((PackedCandidate_->pt()- packed_part.pt())/PackedCandidate_->pt()); + + if(dR<0.01 && dpt<0.1 && PackedCandidate_->charge()==packed_part.charge()){ + if (dRvertex()).mag2()); + + vtx_x=pruned_part_match->vertex().x(); + vtx_y=pruned_part_match->vertex().y(); + vtx_z=pruned_part_match->vertex().z(); + } + Cpfcan_c_tag_[fillntupleentry] = c_tag; + Cpfcan_b_tag_[fillntupleentry] = b_tag; + Cpfcan_g_tag_[fillntupleentry] = g_tag; + + Cpfcan_vtx_x_[fillntupleentry] = vtx_x; + Cpfcan_vtx_y_[fillntupleentry] = vtx_y; + Cpfcan_vtx_z_[fillntupleentry] = vtx_z; + + Cpfcan_dist_from_pv_[fillntupleentry] = dist_from_pv; + + Cpfcan_pdgID_[fillntupleentry] = pdgid_; Cpfcan_pt_[fillntupleentry] = PackedCandidate_->pt(); Cpfcan_eta_[fillntupleentry] = PackedCandidate_->eta(); Cpfcan_phi_[fillntupleentry] = PackedCandidate_->phi(); Cpfcan_ptrel_[fillntupleentry] = catchInfsAndBound(PackedCandidate_->pt()/jet_uncorr_pt,0,-1,0,-1); Cpfcan_erel_[fillntupleentry] = catchInfsAndBound(PackedCandidate_->energy()/jet_uncorr_e,0,-1,0,-1); + Cpfcan_e_[fillntupleentry] = PackedCandidate_->energy(); Cpfcan_phirel_[fillntupleentry] = catchInfsAndBound(fabs(reco::deltaPhi(PackedCandidate_->phi(),jet.phi())),0,-2,0,-0.5); Cpfcan_etarel_[fillntupleentry] = catchInfsAndBound(fabs(PackedCandidate_->eta()-jet.eta()),0,-2,0,-0.5); Cpfcan_deltaR_[fillntupleentry] =catchInfsAndBound(reco::deltaR(*PackedCandidate_,jet),0,-0.6,0,-0.6); @@ -353,33 +509,38 @@ bool ntuple_pfCands::fillBranches(const pat::Jet & jet, const size_t& jetidx, co Cpfcan_vertexChi2_[fillntupleentry]=PackedCandidate_->vertexChi2(); Cpfcan_vertexNdof_[fillntupleentry]=PackedCandidate_->vertexNdof(); + + Cpfcan_CaloFrac_[fillntupleentry] = PackedCandidate_->caloFraction(); + Cpfcan_HadFrac_[fillntupleentry] = PackedCandidate_->hcalFraction(); + //divided Cpfcan_vertexNormalizedChi2_[fillntupleentry]=PackedCandidate_->vertexNormalizedChi2(); Cpfcan_vertex_rho_[fillntupleentry]=catchInfsAndBound(PackedCandidate_->vertex().rho(),0,-1,50); Cpfcan_vertex_phirel_[fillntupleentry]=reco::deltaPhi(PackedCandidate_->vertex().phi(),jet.phi()); Cpfcan_vertex_etarel_[fillntupleentry]=etasign*(PackedCandidate_->vertex().eta()-jet.eta()); Cpfcan_vertexRef_mass_[fillntupleentry]=PackedCandidate_->vertexRef()->p4().M(); - - Cpfcan_puppiw_[fillntupleentry] = PackedCandidate_->puppiWeight(); + trackinfo.buildTrackInfo(PackedCandidate_,jetDir,jetRefTrackDir,pv); - /* - reco::Track::CovarianceMatrix myCov = PseudoTrack.covariance (); - //https://github.com/cms-sw/cmssw/blob/CMSSW_9_0_X/DataFormats/PatCandidates/interface/PackedCandidate.h#L394 + const reco::TransientTrack ttrack = trackinfo.getTTrack(); + float mindistsvold = mindistsvpfcandold(ttrack); + GlobalPoint mindistgpsv = mingpsvpfcand(ttrack); + GlobalPoint gppv = gppvpfcand(ttrack, jetRefTrackDir, pv); - Cpfcan_dptdpt_[fillntupleentry] = catchInfsAndBound(myCov[0][0],0,-1,1); - Cpfcan_detadeta_[fillntupleentry]= catchInfsAndBound(myCov[1][1],0,-1,0.01); - Cpfcan_dphidphi_[fillntupleentry]= catchInfsAndBound(myCov[2][2],0,-1,0.1); + Cpfcan_distminsvold_[fillntupleentry] = mindistsvold; + Cpfcan_dxminsv_[fillntupleentry] = mindistgpsv.x(); + Cpfcan_dyminsv_[fillntupleentry] = mindistgpsv.y(); + Cpfcan_dzminsv_[fillntupleentry] = mindistgpsv.z(); + Cpfcan_dxpv_[fillntupleentry] = gppv.x(); + Cpfcan_dypv_[fillntupleentry] = gppv.y(); + Cpfcan_dzpv_[fillntupleentry] = gppv.z(); - Cpfcan_dxydxy_[fillntupleentry] = catchInfsAndBound(myCov[3][3],7.,-1,7); //zero if pvAssociationQuality ==7 ? - Cpfcan_dzdz_[fillntupleentry] = catchInfsAndBound(myCov[4][4],6.5,-1,6.5); //zero if pvAssociationQuality ==7 ? - Cpfcan_dxydz_[fillntupleentry] = catchInfsAndBound(myCov[3][4],6.,-6,6); //zero if pvAssociationQuality ==7 ? - Cpfcan_dphidxy_[fillntupleentry] = catchInfs(myCov[2][3],-0.03); //zero if pvAssociationQuality ==7 ? - Cpfcan_dlambdadz_[fillntupleentry]= catchInfs(myCov[1][4],-0.03); //zero if pvAssociationQuality ==7 ? - */ + float mindistsv = mindistsvpfcand(ttrack); + float eng_mindistsv = std::log(std::fabs(mindistsv)+1.0); - trackinfo.buildTrackInfo(PackedCandidate_,jetDir,jetRefTrackDir,pv); + Cpfcan_distminsv_[fillntupleentry] = mindistsv; + Cpfcan_distminsv2_[fillntupleentry] = eng_mindistsv; Cpfcan_BtagPf_trackMomentum_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackMomentum(),0,0 ,1000); Cpfcan_BtagPf_trackEta_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackEta() , 0,-5,5); @@ -393,7 +554,7 @@ bool ntuple_pfCands::fillBranches(const pat::Jet & jet, const size_t& jetidx, co Cpfcan_BtagPf_trackSip3dSig_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackSip3dSig(), 0, -1,4e4 ); Cpfcan_BtagPf_trackSip2dVal_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackSip2dVal(), 0, -1,70 ); Cpfcan_BtagPf_trackSip2dSig_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackSip2dSig(), 0, -1,4e4 ); - Cpfcan_BtagPf_trackDecayLen_[fillntupleentry] =0; + Cpfcan_BtagPf_trackDecayLen_[fillntupleentry] =trackinfo.getTrackJetDecayLen(); Cpfcan_BtagPf_trackJetDistVal_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackJetDistVal(),0,-20,1 ); Cpfcan_BtagPf_trackJetDistSig_[fillntupleentry] =catchInfsAndBound(trackinfo.getTrackJetDistSig(),0,-1,1e5 ); @@ -412,19 +573,262 @@ bool ntuple_pfCands::fillBranches(const pat::Jet & jet, const size_t& jetidx, co Cpfcan_charge_[fillntupleentry] = cand_charge_; Cpfcan_lostInnerHits_[fillntupleentry] = catchInfs(PackedCandidate_->lostInnerHits(),2); Cpfcan_numberOfPixelHits_[fillntupleentry] = catchInfs(PackedCandidate_->numberOfPixelHits(),-1); - - //std::cout << PackedCandidate_->lostInnerHits()<< " inner hits " <numberOfPixelHits()<< " Pixel hits + masked " <pixelLayersWithMeasurement()<< " Pixel hits " <hasTrackDetails() ? \ - catchInfsAndBound(PackedCandidate_->pseudoTrack().normalizedChi2(),300,-1,300) : -1; - //for some reason this returns the quality enum not a mask. - Cpfcan_quality_[fillntupleentry] = PackedCandidate_->hasTrackDetails() ? - PackedCandidate_->pseudoTrack().qualityMask() : (1 << reco::TrackBase::loose); - + Cpfcan_chi2_[fillntupleentry] = PackedCandidate_->hasTrackDetails() ? \ + catchInfsAndBound(PackedCandidate_->pseudoTrack().normalizedChi2(),300,-1,300) : -1; + //for some reason this returns the quality enum not a mask. + Cpfcan_quality_[fillntupleentry] = PackedCandidate_->hasTrackDetails() ? + PackedCandidate_->pseudoTrack().qualityMask() : (1 << reco::TrackBase::loose); Cpfcan_drminsv_[fillntupleentry] = catchInfsAndBound(drminpfcandsv_,0,-0.4,0,-0.4); - + //hit pattern variables, as defined here https://github.com/cms-sw/cmssw/blob/master/DataFormats/TrackReco/interface/HitPattern.h + //get track associated to a jet constituent + const reco::Track *track_ptr = nullptr; + auto pf_candidate = dynamic_cast(PackedCandidate_); + auto packed_candidate = dynamic_cast(PackedCandidate_); + if(pf_candidate){ + track_ptr = pf_candidate->bestTrack(); //trackRef was sometimes null + }else if(packed_candidate && packed_candidate->hasTrackDetails()){//if PackedCandidate does not have TrackDetails this gives an Exception because unpackCovariance might be called for pseudoTrack/bestTrack + track_ptr = &(packed_candidate->pseudoTrack()); + } + //get hit pattern information + if(track_ptr){ + const reco::HitPattern &p = track_ptr->hitPattern(); + //Tracker per layer + //Pixel barrel + int Cpfcan_nhitpixelBarrelLayer1 = 0; + int Cpfcan_nhitpixelBarrelLayer2 = 0; + int Cpfcan_nhitpixelBarrelLayer3 = 0; + int Cpfcan_nhitpixelBarrelLayer4 = 0; + //Pixel Endcap + int Cpfcan_nhitpixelEndcapLayer1 = 0; + int Cpfcan_nhitpixelEndcapLayer2 = 0; + //Strip TIB + int Cpfcan_nhitstripTIBLayer1 = 0; + int Cpfcan_nhitstripTIBLayer2 = 0; + int Cpfcan_nhitstripTIBLayer3 = 0; + int Cpfcan_nhitstripTIBLayer4 = 0; + //Strip TID + int Cpfcan_nhitstripTIDLayer1 = 0; + int Cpfcan_nhitstripTIDLayer2 = 0; + int Cpfcan_nhitstripTIDLayer3 = 0; + //Strip TOB + int Cpfcan_nhitstripTOBLayer1 = 0; + int Cpfcan_nhitstripTOBLayer2 = 0; + int Cpfcan_nhitstripTOBLayer3 = 0; + int Cpfcan_nhitstripTOBLayer4 = 0; + int Cpfcan_nhitstripTOBLayer5 = 0; + int Cpfcan_nhitstripTOBLayer6 = 0; + //Strip TEC + int Cpfcan_nhitstripTECLayer1 = 0; + int Cpfcan_nhitstripTECLayer2 = 0; + int Cpfcan_nhitstripTECLayer3 = 0; + int Cpfcan_nhitstripTECLayer4 = 0; + int Cpfcan_nhitstripTECLayer5 = 0; + int Cpfcan_nhitstripTECLayer6 = 0; + int Cpfcan_nhitstripTECLayer7 = 0; + int Cpfcan_nhitstripTECLayer8 = 0; + int Cpfcan_nhitstripTECLayer9 = 0; + // loop over the hits of the track. + //const static unsigned short LayerOffset = 3; + //const static unsigned short LayerMask = 0xF; + for(int nh = 0; nh < p.numberOfAllHits(reco::HitPattern::TRACK_HITS); nh++){ + uint32_t hit = p.getHitPattern(reco::HitPattern::TRACK_HITS, nh); + if(p.validHitFilter(hit)){// if the hit is valid + //Pixel Barrel // it is in pixel barrel + if(p.pixelBarrelHitFilter(hit)){ + //std::cout << "valid hit found in pixel Barrel layer " << p.getLayer(hit) << std::endl; + //if(p.getLayer(hit)==1){ + // std::cout<< (hit >> LayerOffset) << " " << ((hit >> LayerOffset) & LayerMask) << std::endl; + //} + if(p.getLayer(hit)==1) Cpfcan_nhitpixelBarrelLayer1 = Cpfcan_nhitpixelBarrelLayer1+1; + if(p.getLayer(hit)==2) Cpfcan_nhitpixelBarrelLayer2 = Cpfcan_nhitpixelBarrelLayer2+1; + if(p.getLayer(hit)==3) Cpfcan_nhitpixelBarrelLayer3 = Cpfcan_nhitpixelBarrelLayer3+1; + if(p.getLayer(hit)==4) Cpfcan_nhitpixelBarrelLayer4 = Cpfcan_nhitpixelBarrelLayer4+1; + } + //Pixel Endcap + if(p.pixelEndcapHitFilter(hit)){ + //std::cout << "valid hit found in pixel Endcap layer " << p.getLayer(hit) << std::endl; + if(p.getLayer(hit)==1) Cpfcan_nhitpixelEndcapLayer1 = Cpfcan_nhitpixelEndcapLayer1+1; + if(p.getLayer(hit)==2) Cpfcan_nhitpixelEndcapLayer2 = Cpfcan_nhitpixelEndcapLayer2+1; + } + //Strip TIB + if(p.stripTIBHitFilter(hit)){ + //std::cout << "valid hit found in TIB layer " << p.getLayer(hit) << std::endl; + if(p.getLayer(hit)==1) Cpfcan_nhitstripTIBLayer1 = Cpfcan_nhitstripTIBLayer1+1; + if(p.getLayer(hit)==2) Cpfcan_nhitstripTIBLayer2 = Cpfcan_nhitstripTIBLayer2+1; + if(p.getLayer(hit)==3) Cpfcan_nhitstripTIBLayer3 = Cpfcan_nhitstripTIBLayer3+1; + if(p.getLayer(hit)==4) Cpfcan_nhitstripTIBLayer4 = Cpfcan_nhitstripTIBLayer4+1; + } + //Strip TID + if(p.stripTIDHitFilter(hit)){ + //std::cout << "valid hit found in TID layer " << p.getLayer(hit) << std::endl; + if(p.getLayer(hit)==1) Cpfcan_nhitstripTIDLayer1 = Cpfcan_nhitstripTIDLayer1+1; + if(p.getLayer(hit)==2) Cpfcan_nhitstripTIDLayer2 = Cpfcan_nhitstripTIDLayer2+1; + if(p.getLayer(hit)==3) Cpfcan_nhitstripTIDLayer3 = Cpfcan_nhitstripTIDLayer3+1; + } + //Strip TOB + if(p.stripTOBHitFilter(hit)){ + //std::cout << "valid hit found in TOB layer " << p.getLayer(hit) << std::endl; + if(p.getLayer(hit)==1) Cpfcan_nhitstripTOBLayer1 = Cpfcan_nhitstripTOBLayer1+1; + if(p.getLayer(hit)==2) Cpfcan_nhitstripTOBLayer2 = Cpfcan_nhitstripTOBLayer2+1; + if(p.getLayer(hit)==3) Cpfcan_nhitstripTOBLayer3 = Cpfcan_nhitstripTOBLayer3+1; + if(p.getLayer(hit)==4) Cpfcan_nhitstripTOBLayer4 = Cpfcan_nhitstripTOBLayer4+1; + if(p.getLayer(hit)==5) Cpfcan_nhitstripTOBLayer5 = Cpfcan_nhitstripTOBLayer5+1; + if(p.getLayer(hit)==6) Cpfcan_nhitstripTOBLayer6 = Cpfcan_nhitstripTOBLayer6+1; + } + //Strip TEC + if(p.stripTECHitFilter(hit)){ + //std::cout << "valid hit found in TEC layer " << p.getLayer(hit) << std::endl; + if(p.getLayer(hit)==1) Cpfcan_nhitstripTECLayer1 = Cpfcan_nhitstripTECLayer1+1; + if(p.getLayer(hit)==2) Cpfcan_nhitstripTECLayer2 = Cpfcan_nhitstripTECLayer2+1; + if(p.getLayer(hit)==3) Cpfcan_nhitstripTECLayer3 = Cpfcan_nhitstripTECLayer3+1; + if(p.getLayer(hit)==4) Cpfcan_nhitstripTECLayer4 = Cpfcan_nhitstripTECLayer4+1; + if(p.getLayer(hit)==5) Cpfcan_nhitstripTECLayer5 = Cpfcan_nhitstripTECLayer5+1; + if(p.getLayer(hit)==6) Cpfcan_nhitstripTECLayer6 = Cpfcan_nhitstripTECLayer6+1; + if(p.getLayer(hit)==7) Cpfcan_nhitstripTECLayer7 = Cpfcan_nhitstripTECLayer7+1; + if(p.getLayer(hit)==8) Cpfcan_nhitstripTECLayer8 = Cpfcan_nhitstripTECLayer8+1; + if(p.getLayer(hit)==9) Cpfcan_nhitstripTECLayer9 = Cpfcan_nhitstripTECLayer9+1; + } + } + } + //Pixel Barrel + Cpfcan_nhitpixelBarrelLayer1_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelBarrelHits()) ? catchInfsAndBound(Cpfcan_nhitpixelBarrelLayer1,-1,0,100,0) : -1; + Cpfcan_nhitpixelBarrelLayer2_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelBarrelHits()) ? catchInfsAndBound(Cpfcan_nhitpixelBarrelLayer2,-1,0,100,0) : -1; + Cpfcan_nhitpixelBarrelLayer3_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelBarrelHits()) ? catchInfsAndBound(Cpfcan_nhitpixelBarrelLayer3,-1,0,100,0) : -1; + Cpfcan_nhitpixelBarrelLayer4_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelBarrelHits()) ? catchInfsAndBound(Cpfcan_nhitpixelBarrelLayer4,-1,0,100,0) : -1; + //Pixel Endcap + Cpfcan_nhitpixelEndcapLayer1_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelEndcapHits()) ? catchInfsAndBound(Cpfcan_nhitpixelEndcapLayer1,-1,0,100,0) : -1; + Cpfcan_nhitpixelEndcapLayer2_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelEndcapHits()) ? catchInfsAndBound(Cpfcan_nhitpixelEndcapLayer2,-1,0,100,0) : -1; + //Strip TIB + Cpfcan_nhitstripTIBLayer1_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTIBHits()) ? catchInfsAndBound(Cpfcan_nhitstripTIBLayer1,-1,0,100,0) : -1; + Cpfcan_nhitstripTIBLayer2_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTIBHits()) ? catchInfsAndBound(Cpfcan_nhitstripTIBLayer2,-1,0,100,0) : -1; + Cpfcan_nhitstripTIBLayer3_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTIBHits()) ? catchInfsAndBound(Cpfcan_nhitstripTIBLayer3,-1,0,100,0) : -1; + Cpfcan_nhitstripTIBLayer4_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTIBHits()) ? catchInfsAndBound(Cpfcan_nhitstripTIBLayer4,-1,0,100,0) : -1; + //Strip TID + Cpfcan_nhitstripTIDLayer1_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTIDHits()) ? catchInfsAndBound(Cpfcan_nhitstripTIDLayer1,-1,0,100,0) : -1; + Cpfcan_nhitstripTIDLayer2_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTIDHits()) ? catchInfsAndBound(Cpfcan_nhitstripTIDLayer2,-1,0,100,0) : -1; + Cpfcan_nhitstripTIDLayer3_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTIDHits()) ? catchInfsAndBound(Cpfcan_nhitstripTIDLayer3,-1,0,100,0) : -1; + //Strip TOB + Cpfcan_nhitstripTOBLayer1_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTOBHits()) ? catchInfsAndBound(Cpfcan_nhitstripTOBLayer1,-1,0,100,0) : -1; + Cpfcan_nhitstripTOBLayer2_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTOBHits()) ? catchInfsAndBound(Cpfcan_nhitstripTOBLayer2,-1,0,100,0) : -1; + Cpfcan_nhitstripTOBLayer3_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTOBHits()) ? catchInfsAndBound(Cpfcan_nhitstripTOBLayer3,-1,0,100,0) : -1; + Cpfcan_nhitstripTOBLayer4_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTOBHits()) ? catchInfsAndBound(Cpfcan_nhitstripTOBLayer4,-1,0,100,0) : -1; + Cpfcan_nhitstripTOBLayer5_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTOBHits()) ? catchInfsAndBound(Cpfcan_nhitstripTOBLayer5,-1,0,100,0) : -1; + Cpfcan_nhitstripTOBLayer6_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTOBHits()) ? catchInfsAndBound(Cpfcan_nhitstripTOBLayer6,-1,0,100,0) : -1; + //Strip TEC + Cpfcan_nhitstripTECLayer1_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTECHits()) ? catchInfsAndBound(Cpfcan_nhitstripTECLayer1,-1,0,100,0) : -1; + Cpfcan_nhitstripTECLayer2_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTECHits()) ? catchInfsAndBound(Cpfcan_nhitstripTECLayer2,-1,0,100,0) : -1; + Cpfcan_nhitstripTECLayer3_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTECHits()) ? catchInfsAndBound(Cpfcan_nhitstripTECLayer3,-1,0,100,0) : -1; + Cpfcan_nhitstripTECLayer4_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTECHits()) ? catchInfsAndBound(Cpfcan_nhitstripTECLayer4,-1,0,100,0) : -1; + Cpfcan_nhitstripTECLayer5_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTECHits()) ? catchInfsAndBound(Cpfcan_nhitstripTECLayer5,-1,0,100,0) : -1; + Cpfcan_nhitstripTECLayer6_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTECHits()) ? catchInfsAndBound(Cpfcan_nhitstripTECLayer6,-1,0,100,0) : -1; + Cpfcan_nhitstripTECLayer7_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTECHits()) ? catchInfsAndBound(Cpfcan_nhitstripTECLayer7,-1,0,100,0) : -1; + Cpfcan_nhitstripTECLayer8_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTECHits()) ? catchInfsAndBound(Cpfcan_nhitstripTECLayer8,-1,0,100,0) : -1; + Cpfcan_nhitstripTECLayer9_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTECHits()) ? catchInfsAndBound(Cpfcan_nhitstripTECLayer9,-1,0,100,0) : -1; + //Tracker all layers together + //Valid hits + Cpfcan_numberOfValidHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidHits(),-1,0,100,0) : -1; + Cpfcan_numberOfValidTrackerHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidTrackerHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidTrackerHits(),-1,0,100,0) : -1; + Cpfcan_numberOfValidPixelHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidPixelHits(),-1,0,100,0) : -1; + Cpfcan_numberOfValidPixelBarrelHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelBarrelHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidPixelBarrelHits(),-1,0,100,0) : -1; + Cpfcan_numberOfValidPixelEndcapHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidPixelEndcapHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidPixelEndcapHits(),-1,0,100,0) : -1; + Cpfcan_numberOfValidStripHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidStripHits(),-1,0,100,0) : -1; + Cpfcan_numberOfValidStripTIBHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTIBHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidStripTIBHits(),-1,0,100,0) : -1; + Cpfcan_numberOfValidStripTIDHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTIDHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidStripTIDHits(),-1,0,100,0) : -1; + Cpfcan_numberOfValidStripTOBHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTOBHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidStripTOBHits(),-1,0,100,0) : -1; + Cpfcan_numberOfValidStripTECHits_[fillntupleentry] = (track_ptr->hitPattern().numberOfValidStripTECHits()) ? catchInfsAndBound(track_ptr->hitPattern().numberOfValidStripTECHits(),-1,0,100,0) : -1; + //LayersWithMeasurement + Cpfcan_trackerLayersWithMeasurementOld_[fillntupleentry] = (track_ptr->hitPattern().trackerLayersWithMeasurementOld()) ? catchInfsAndBound(track_ptr->hitPattern().trackerLayersWithMeasurementOld(),-1,0,100,0) : -1; + Cpfcan_trackerLayersWithMeasurement_[fillntupleentry] = (track_ptr->hitPattern().trackerLayersWithMeasurement()) ? catchInfsAndBound(track_ptr->hitPattern().trackerLayersWithMeasurement(),-1,0,100,0) : -1; + Cpfcan_pixelLayersWithMeasurementOld_[fillntupleentry] = (track_ptr->hitPattern().pixelLayersWithMeasurementOld()) ? catchInfsAndBound(track_ptr->hitPattern().pixelLayersWithMeasurementOld(),-1,0,100,0) : -1; + Cpfcan_pixelLayersWithMeasurement_[fillntupleentry] = (track_ptr->hitPattern().pixelLayersWithMeasurement()) ? catchInfsAndBound(track_ptr->hitPattern().pixelLayersWithMeasurement(),-1,0,100,0) : -1; + Cpfcan_stripLayersWithMeasurement_[fillntupleentry] = (track_ptr->hitPattern().stripLayersWithMeasurement()) ? catchInfsAndBound(track_ptr->hitPattern().stripLayersWithMeasurement(),-1,0,100,0) : -1; + Cpfcan_pixelBarrelLayersWithMeasurement_[fillntupleentry] = (track_ptr->hitPattern().pixelBarrelLayersWithMeasurement()) ? catchInfsAndBound(track_ptr->hitPattern().pixelBarrelLayersWithMeasurement(),-1,0,100,0) : -1; + Cpfcan_pixelEndcapLayersWithMeasurement_[fillntupleentry] = (track_ptr->hitPattern().pixelEndcapLayersWithMeasurement()) ? catchInfsAndBound(track_ptr->hitPattern().pixelEndcapLayersWithMeasurement(),-1,0,100,0) : -1; + Cpfcan_stripTIBLayersWithMeasurement_[fillntupleentry] = (track_ptr->hitPattern().stripTIBLayersWithMeasurement()) ? catchInfsAndBound(track_ptr->hitPattern().stripTIBLayersWithMeasurement(),-1,0,100,0) : -1; + Cpfcan_stripTIDLayersWithMeasurement_[fillntupleentry] = (track_ptr->hitPattern().stripTIDLayersWithMeasurement()) ? catchInfsAndBound(track_ptr->hitPattern().stripTIDLayersWithMeasurement(),-1,0,100,0) : -1; + Cpfcan_stripTOBLayersWithMeasurement_[fillntupleentry] = (track_ptr->hitPattern().stripTOBLayersWithMeasurement()) ? catchInfsAndBound(track_ptr->hitPattern().stripTOBLayersWithMeasurement(),-1,0,100,0) : -1; + Cpfcan_stripTECLayersWithMeasurement_[fillntupleentry] = (track_ptr->hitPattern().stripTECLayersWithMeasurement()) ? catchInfsAndBound(track_ptr->hitPattern().stripTECLayersWithMeasurement(),-1,0,100,0) : -1; + //Null + Cpfcan_trackerLayersNull_[fillntupleentry] = (track_ptr->hitPattern().trackerLayersNull()) ? catchInfsAndBound(track_ptr->hitPattern().trackerLayersNull(),-1,0,100,0) : -1; + Cpfcan_pixelLayersNull_[fillntupleentry] = (track_ptr->hitPattern().pixelLayersNull()) ? catchInfsAndBound(track_ptr->hitPattern().pixelLayersNull(),-1,0,100,0) : -1; + Cpfcan_stripLayersNull_[fillntupleentry] = (track_ptr->hitPattern().stripLayersNull()) ? catchInfsAndBound(track_ptr->hitPattern().stripLayersNull(),-1,0,100,0) : -1; + Cpfcan_pixelBarrelLayersNull_[fillntupleentry] = (track_ptr->hitPattern().pixelBarrelLayersNull()) ? catchInfsAndBound(track_ptr->hitPattern().pixelBarrelLayersNull(),-1,0,100,0) : -1; + Cpfcan_pixelEndcapLayersNull_[fillntupleentry] = (track_ptr->hitPattern().pixelEndcapLayersNull()) ? catchInfsAndBound(track_ptr->hitPattern().pixelEndcapLayersNull(),-1,0,100,0) : -1; + Cpfcan_stripTIBLayersNull_[fillntupleentry] = (track_ptr->hitPattern().stripTIBLayersNull()) ? catchInfsAndBound(track_ptr->hitPattern().stripTIBLayersNull(),-1,0,100,0) : -1; + Cpfcan_stripTIDLayersNull_[fillntupleentry] = (track_ptr->hitPattern().stripTIDLayersNull()) ? catchInfsAndBound(track_ptr->hitPattern().stripTIDLayersNull(),-1,0,100,0) : -1; + Cpfcan_stripTOBLayersNull_[fillntupleentry] = (track_ptr->hitPattern().stripTOBLayersNull()) ? catchInfsAndBound(track_ptr->hitPattern().stripTOBLayersNull(),-1,0,100,0) : -1; + Cpfcan_stripTECLayersNull_[fillntupleentry] = (track_ptr->hitPattern().stripTECLayersNull()) ? catchInfsAndBound(track_ptr->hitPattern().stripTECLayersNull(),-1,0,100,0) : -1; + }else{ + //Tracker per layer + //Pixel barrel + Cpfcan_nhitpixelBarrelLayer1_[fillntupleentry] = -1; + Cpfcan_nhitpixelBarrelLayer2_[fillntupleentry] = -1; + Cpfcan_nhitpixelBarrelLayer3_[fillntupleentry] = -1; + Cpfcan_nhitpixelBarrelLayer4_[fillntupleentry] = -1; + //Pixel Endcap + Cpfcan_nhitpixelEndcapLayer1_[fillntupleentry] = -1; + Cpfcan_nhitpixelEndcapLayer2_[fillntupleentry] = -1; + //Strip TIB + Cpfcan_nhitstripTIBLayer1_[fillntupleentry] = -1; + Cpfcan_nhitstripTIBLayer2_[fillntupleentry] = -1; + Cpfcan_nhitstripTIBLayer3_[fillntupleentry] = -1; + Cpfcan_nhitstripTIBLayer4_[fillntupleentry] = -1; + //Strip TID + Cpfcan_nhitstripTIDLayer1_[fillntupleentry] = -1; + Cpfcan_nhitstripTIDLayer2_[fillntupleentry] = -1; + Cpfcan_nhitstripTIDLayer3_[fillntupleentry] = -1; + //Strip TOB + Cpfcan_nhitstripTOBLayer1_[fillntupleentry] = -1; + Cpfcan_nhitstripTOBLayer2_[fillntupleentry] = -1; + Cpfcan_nhitstripTOBLayer3_[fillntupleentry] = -1; + Cpfcan_nhitstripTOBLayer4_[fillntupleentry] = -1; + Cpfcan_nhitstripTOBLayer5_[fillntupleentry] = -1; + Cpfcan_nhitstripTOBLayer6_[fillntupleentry] = -1; + //Strip TEC + Cpfcan_nhitstripTECLayer1_[fillntupleentry] = -1; + Cpfcan_nhitstripTECLayer2_[fillntupleentry] = -1; + Cpfcan_nhitstripTECLayer3_[fillntupleentry] = -1; + Cpfcan_nhitstripTECLayer4_[fillntupleentry] = -1; + Cpfcan_nhitstripTECLayer5_[fillntupleentry] = -1; + Cpfcan_nhitstripTECLayer6_[fillntupleentry] = -1; + Cpfcan_nhitstripTECLayer7_[fillntupleentry] = -1; + Cpfcan_nhitstripTECLayer8_[fillntupleentry] = -1; + Cpfcan_nhitstripTECLayer9_[fillntupleentry] = -1; + //Tracker all layers together + //Valid hits + Cpfcan_numberOfValidHits_[fillntupleentry] = -1; + Cpfcan_numberOfValidTrackerHits_[fillntupleentry] = -1; + Cpfcan_numberOfValidPixelHits_[fillntupleentry] = -1; + Cpfcan_numberOfValidPixelBarrelHits_[fillntupleentry] = -1; + Cpfcan_numberOfValidPixelEndcapHits_[fillntupleentry] = -1; + Cpfcan_numberOfValidStripHits_[fillntupleentry] = -1; + Cpfcan_numberOfValidStripTIBHits_[fillntupleentry] = -1; + Cpfcan_numberOfValidStripTIDHits_[fillntupleentry] = -1; + Cpfcan_numberOfValidStripTOBHits_[fillntupleentry] = -1; + Cpfcan_numberOfValidStripTECHits_[fillntupleentry] = -1; + //LayersWithMeasurement + Cpfcan_trackerLayersWithMeasurementOld_[fillntupleentry] = -1; + Cpfcan_trackerLayersWithMeasurement_[fillntupleentry] = -1; + Cpfcan_pixelLayersWithMeasurementOld_[fillntupleentry] = -1; + Cpfcan_pixelLayersWithMeasurement_[fillntupleentry] = -1; + Cpfcan_stripLayersWithMeasurement_[fillntupleentry] = -1; + Cpfcan_pixelBarrelLayersWithMeasurement_[fillntupleentry] = -1; + Cpfcan_pixelEndcapLayersWithMeasurement_[fillntupleentry] = -1; + Cpfcan_stripTIBLayersWithMeasurement_[fillntupleentry] = -1; + Cpfcan_stripTIDLayersWithMeasurement_[fillntupleentry] = -1; + Cpfcan_stripTOBLayersWithMeasurement_[fillntupleentry] = -1; + Cpfcan_stripTECLayersWithMeasurement_[fillntupleentry] = -1; + //Null + Cpfcan_trackerLayersNull_[fillntupleentry] = -1; + Cpfcan_pixelLayersNull_[fillntupleentry] = -1; + Cpfcan_stripLayersNull_[fillntupleentry] = -1; + Cpfcan_pixelBarrelLayersNull_[fillntupleentry] = -1; + Cpfcan_pixelEndcapLayersNull_[fillntupleentry] = -1; + Cpfcan_stripTIBLayersNull_[fillntupleentry] = -1; + Cpfcan_stripTIDLayersNull_[fillntupleentry] = -1; + Cpfcan_stripTOBLayersNull_[fillntupleentry] = -1; + Cpfcan_stripTECLayersNull_[fillntupleentry] = -1; + } } else{// neutral candidates @@ -437,6 +841,7 @@ bool ntuple_pfCands::fillBranches(const pat::Jet & jet, const size_t& jetidx, co Npfcan_phi_[fillntupleentry] = PackedCandidate_->phi(); Npfcan_ptrel_[fillntupleentry] = catchInfsAndBound(PackedCandidate_->pt()/jet_uncorr_pt,0,-1,0,-1); Npfcan_erel_[fillntupleentry] = catchInfsAndBound(PackedCandidate_->energy()/jet_uncorr_e,0,-1,0,-1); + Npfcan_e_[fillntupleentry] = PackedCandidate_->energy(); Npfcan_puppiw_[fillntupleentry] = PackedCandidate_->puppiWeight(); Npfcan_phirel_[fillntupleentry] = catchInfsAndBound(fabs(reco::deltaPhi(PackedCandidate_->phi(),jet.phi())),0,-2,0,-0.5); Npfcan_etarel_[fillntupleentry] = catchInfsAndBound(fabs(PackedCandidate_->eta()-jet.eta()),0,-2,0,-0.5); @@ -444,26 +849,15 @@ bool ntuple_pfCands::fillBranches(const pat::Jet & jet, const size_t& jetidx, co Npfcan_isGamma_[fillntupleentry] = 0; if(fabs(PackedCandidate_->pdgId())==22) Npfcan_isGamma_[fillntupleentry] = 1; Npfcan_HadFrac_[fillntupleentry] = PackedCandidate_->hcalFraction(); + Npfcan_pdgID_[fillntupleentry] = pdgid_; Npfcan_drminsv_[fillntupleentry] = catchInfsAndBound(drminpfcandsv_,0,-0.4,0,-0.4); } - - } // end loop over jet.numberOfDaughters() - - /* - std::cout <<"numbers charged/neutrals"<size(); ++i) { + if (!track.isValid()) {continue;} + reco::Vertex::CovarianceMatrix csv; secVertices()->at(i).fillVertexCovariance(csv); + reco::Vertex vertex(secVertices()->at(i).vertex(), csv); + if (!vertex.isValid()) {continue;} + + GlobalVector direction(secVertices()->at(i).px(),secVertices()->at(i).py(),secVertices()->at(i).pz()); + + + AnalyticalImpactPointExtrapolator extrapolator(track.field()); + TrajectoryStateOnSurface tsos = extrapolator.extrapolate(track.impactPointState(), RecoVertex::convertPos(vertex.position())); + + VertexDistance3D dist; + + if (!tsos.isValid()) {continue;} + GlobalPoint refPoint = tsos.globalPosition(); + GlobalError refPointErr = tsos.cartesianError().position(); + GlobalPoint vertexPosition = RecoVertex::convertPos(vertex.position()); + GlobalError vertexPositionErr = RecoVertex::convertError(vertex.error()); + + std::pair result(true, dist.distance(VertexState(vertexPosition, vertexPositionErr), VertexState(refPoint, refPointErr))); + if (!result.first) {continue;} + + GlobalPoint impactPoint = tsos.globalPosition(); + GlobalVector IPVec(impactPoint.x() - vertex.x(), impactPoint.y() - vertex.y(), impactPoint.z() - vertex.z()); + double prod = IPVec.dot(direction); + double sign = (prod >= 0) ? 1. : -1.; + + if(result.second.value() < mindist_){ + out_dist = sign * result.second.value(); + mindist_ = result.second.value(); + } + } + return out_dist; +} + +GlobalPoint ntuple_pfCands::mingpsvpfcand(const reco::TransientTrack track) { + + float mindist_ = 999.999; + GlobalPoint out_dist(0.0,0.0,0.0); + GlobalPoint pca; + for (unsigned int i=0; isize(); ++i) { + if (!track.isValid()) {continue;} + + reco::Vertex::CovarianceMatrix csv; secVertices()->at(i).fillVertexCovariance(csv); + reco::Vertex vertex(secVertices()->at(i).vertex(), csv); + if (!vertex.isValid()) {continue;} + + GlobalVector direction(secVertices()->at(i).px(),secVertices()->at(i).py(),secVertices()->at(i).pz()); + + AnalyticalImpactPointExtrapolator extrapolator(track.field()); + TrajectoryStateOnSurface tsos = extrapolator.extrapolate(track.impactPointState(), RecoVertex::convertPos(vertex.position())); + + VertexDistance3D dist; + + if (!tsos.isValid()) {continue;} + GlobalPoint refPoint = tsos.globalPosition(); + GlobalError refPointErr = tsos.cartesianError().position(); + GlobalPoint vertexPosition = RecoVertex::convertPos(vertex.position()); + GlobalError vertexPositionErr = RecoVertex::convertError(vertex.error()); + + std::pair result(true, dist.distance(VertexState(vertexPosition, vertexPositionErr), VertexState(refPoint, refPointErr))); + if (!result.first) {continue;} + + pca = tsos.globalPosition(); + + if(result.second.value() < mindist_){ + out_dist = pca; + mindist_ = result.second.value(); + } + } + return out_dist; +} + +GlobalPoint ntuple_pfCands::gppvpfcand(const reco::TransientTrack track, const GlobalVector direction, const reco::Vertex vertex) { + + float mindist_ = 999.999; + float dist = 0.; + GlobalPoint out_dist(0.0,0.0,0.0); + GlobalPoint pca; + if ((track.isValid()) && (vertex.isValid())){ + + AnalyticalImpactPointExtrapolator extrapolator(track.field()); + TrajectoryStateOnSurface tsos = extrapolator.extrapolate(track.impactPointState(), RecoVertex::convertPos(vertex.position())); + + VertexDistance3D dist; + + if (tsos.isValid()) { + GlobalPoint refPoint = tsos.globalPosition(); + GlobalError refPointErr = tsos.cartesianError().position(); + GlobalPoint vertexPosition = RecoVertex::convertPos(vertex.position()); + GlobalError vertexPositionErr = RecoVertex::convertError(vertex.error()); + + std::pair result(true, dist.distance(VertexState(vertexPosition, vertexPositionErr), VertexState(refPoint, refPointErr))); + if (result.first) { + pca = tsos.globalPosition(); + } + } + } + out_dist = pca; + return out_dist; +} + +float ntuple_pfCands::mindistsvpfcand(const reco::TransientTrack track) { + + float mindist_ = 999.999; + float out_dist = 0.0; + for (unsigned int i=0; isize(); ++i) { + if (!track.isValid()) {continue;} + reco::Vertex::CovarianceMatrix csv; secVertices()->at(i).fillVertexCovariance(csv); + reco::Vertex vertex(secVertices()->at(i).vertex(), csv); + if (!vertex.isValid()) {continue;} + + GlobalVector direction(secVertices()->at(i).px(),secVertices()->at(i).py(),secVertices()->at(i).pz()); + + + AnalyticalImpactPointExtrapolator extrapolator(track.field()); + TrajectoryStateOnSurface tsos = extrapolator.extrapolate(track.impactPointState(), RecoVertex::convertPos(vertex.position())); + + VertexDistance3D dist; + + if (!tsos.isValid()) {continue;} + GlobalPoint refPoint = tsos.globalPosition(); + GlobalError refPointErr = tsos.cartesianError().position(); + GlobalPoint vertexPosition = RecoVertex::convertPos(vertex.position()); + GlobalError vertexPositionErr = RecoVertex::convertError(vertex.error()); + + std::pair result(true, dist.distance(VertexState(vertexPosition, vertexPositionErr), VertexState(refPoint, refPointErr))); + if (!result.first) {continue;} + + GlobalPoint impactPoint = tsos.globalPosition(); + GlobalVector IPVec(impactPoint.x() - vertex.x(), impactPoint.y() - vertex.y(), impactPoint.z() - vertex.z()); + double prod = IPVec.dot(direction); + double sign = (prod >= 0) ? 1. : -1.; + + if(result.second.value() < mindist_){ + float inv_dist = 1.0 / (sign * result.second.value() + 0.001); + out_dist = inv_dist; + mindist_ = result.second.value(); + } + } + return out_dist; +} + +bool ntuple_pfCands::containParton(const reco::Candidate * pruned_part, int pdgid) { + if (abs(pruned_part->pdgId())==pdgid) return true; + for(size_t i=0;i< pruned_part->numberOfMothers();i++){ + if (containParton(pruned_part->mother(i), pdgid)) return true; + } + return false; +} diff --git a/README.md b/README.md index 3644d0c8b39..c1ffceeb304 100644 --- a/README.md +++ b/README.md @@ -6,18 +6,18 @@ Installation (CMSSW 10_6_X) ============ ``` -cmsrel CMSSW_10_6_0 -cd CMSSW_10_6_0/src/ +cmsrel CMSSW_10_6_30 +cd CMSSW_10_6_30/src/ cmsenv git cms-init -git clone https://github.com/emilbols/DeepNTuples +git clone https://github.com/AlexDeMoor/DeepNTuples cd DeepNTuples git checkout 106X # Add JetToolBox git submodule init git submodule update -scram b -j 4 +scram b -j 8 ``` Further settings