FCCAnalyses
Loading...
Searching...
No Matches
JetConstituentsUtils.h
Go to the documentation of this file.
1#ifndef FCCAnalyses_JetConstituentsUtils_h
2#define FCCAnalyses_JetConstituentsUtils_h
3
4#include "ROOT/RVec.hxx"
5#include "edm4hep/ReconstructedParticle.h"
6#include "edm4hep/MCParticle.h"
7#include "edm4hep/Quantity.h"
8#if __has_include("edm4hep/TrackerHit3DData.h")
9#include "edm4hep/TrackerHit3DData.h"
10#else
11#include "edm4hep/TrackerHitData.h"
12namespace edm4hep {
13 using TrackerHit3DData = edm4hep::TrackerHitData;
14}
15#endif
16
17#include "fastjet/JetDefinition.hh"
18
19#include "TMath.h"
20#include "TVector3.h"
21#include "TRotation.h"
22#include "TLorentzVector.h"
23
24namespace FCCAnalyses {
25 namespace JetConstituentsUtils {
26 namespace rv = ROOT::VecOps;
27 using FCCAnalysesJetConstituents = rv::RVec<edm4hep::ReconstructedParticleData>;
28 using FCCAnalysesJetConstituentsData = rv::RVec<float>;
29
31 rv::RVec<FCCAnalysesJetConstituents> build_constituents(const rv::RVec<edm4hep::ReconstructedParticleData>&,
32 const rv::RVec<edm4hep::ReconstructedParticleData>&);
33
34 rv::RVec<FCCAnalysesJetConstituents> build_constituents_cluster(const rv::RVec<edm4hep::ReconstructedParticleData>& rps,
35 const std::vector<std::vector<int>>& indices);
36
38 FCCAnalysesJetConstituents get_jet_constituents(const rv::RVec<FCCAnalysesJetConstituents>&, int);
40 rv::RVec<FCCAnalysesJetConstituents> get_constituents(const rv::RVec<FCCAnalysesJetConstituents>&,
41 const rv::RVec<int>&);
42
43
44 //sorting jets
45 ROOT::VecOps::RVec<edm4hep::ReconstructedParticleData> jets_sorting_on_nconst(const rv::RVec<edm4hep::ReconstructedParticleData>&);
46 ROOT::VecOps::RVec<edm4hep::ReconstructedParticleData> jets_sorting_on_energy(const rv::RVec<edm4hep::ReconstructedParticleData>&);
47
48
49 rv::RVec<FCCAnalysesJetConstituentsData> get_Bz(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
50 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
51
52 rv::RVec<FCCAnalysesJetConstituentsData> get_pt(const rv::RVec<FCCAnalysesJetConstituents>&);
53 rv::RVec<FCCAnalysesJetConstituentsData> get_p(const rv::RVec<FCCAnalysesJetConstituents>&);
54 rv::RVec<FCCAnalysesJetConstituentsData> get_e(const rv::RVec<FCCAnalysesJetConstituents>&);
55 rv::RVec<FCCAnalysesJetConstituentsData> get_theta(const rv::RVec<FCCAnalysesJetConstituents>&);
56 rv::RVec<FCCAnalysesJetConstituentsData> get_phi(const rv::RVec<FCCAnalysesJetConstituents>&);
57 rv::RVec<FCCAnalysesJetConstituentsData> get_type(const rv::RVec<FCCAnalysesJetConstituents>&);
58 rv::RVec<FCCAnalysesJetConstituentsData> get_charge(const rv::RVec<FCCAnalysesJetConstituents>&);
59
60 //displacement
61 rv::RVec<FCCAnalysesJetConstituentsData> get_d0(const rv::RVec<FCCAnalysesJetConstituents>&,
62 const ROOT::VecOps::RVec<edm4hep::TrackState>&);
63
64 rv::RVec<FCCAnalysesJetConstituentsData> get_z0(const rv::RVec<FCCAnalysesJetConstituents>& ,
65 const ROOT::VecOps::RVec<edm4hep::TrackState>&);
66
67 rv::RVec<FCCAnalysesJetConstituentsData> get_phi0(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
68 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
69
70 rv::RVec<FCCAnalysesJetConstituentsData> get_omega(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
71 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
72
73 rv::RVec<FCCAnalysesJetConstituentsData> get_tanLambda(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
74 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
75
76
77 rv::RVec<FCCAnalysesJetConstituentsData> XPtoPar_dxy(const rv::RVec<FCCAnalysesJetConstituents>&,
78 const ROOT::VecOps::RVec<edm4hep::TrackState>&,
79 const TLorentzVector& V, // primary vertex
80 const float&);
81 rv::RVec<FCCAnalysesJetConstituentsData> XPtoPar_dz(const rv::RVec<FCCAnalysesJetConstituents>&,
82 const ROOT::VecOps::RVec<edm4hep::TrackState>&,
83 const TLorentzVector& V, // primary vertex
84 const float&);
85 rv::RVec<FCCAnalysesJetConstituentsData> XPtoPar_phi(const rv::RVec<FCCAnalysesJetConstituents>&,
86 const ROOT::VecOps::RVec<edm4hep::TrackState>&,
87 const TLorentzVector& V, // primary vertex
88 const float&);
89 rv::RVec<FCCAnalysesJetConstituentsData> XPtoPar_C(const rv::RVec<FCCAnalysesJetConstituents>&,
90 const ROOT::VecOps::RVec<edm4hep::TrackState>&,
91 const float&);
92 rv::RVec<FCCAnalysesJetConstituentsData> XPtoPar_ct(const rv::RVec<FCCAnalysesJetConstituents>&,
93 const ROOT::VecOps::RVec<edm4hep::TrackState>&,
94 const float&);
95
96 //covariance matrix
97 //diagonal
98 rv::RVec<FCCAnalysesJetConstituentsData> get_omega_cov(const rv::RVec<FCCAnalysesJetConstituents>&,
99 const ROOT::VecOps::RVec<edm4hep::TrackState>&);
100
101 rv::RVec<FCCAnalysesJetConstituentsData> get_d0_cov(const rv::RVec<FCCAnalysesJetConstituents>&,
102 const ROOT::VecOps::RVec<edm4hep::TrackState>& );
103
104 rv::RVec<FCCAnalysesJetConstituentsData> get_z0_cov(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
105 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
106
107 rv::RVec<FCCAnalysesJetConstituentsData> get_phi0_cov(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
108 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
109
110 rv::RVec<FCCAnalysesJetConstituentsData> get_tanlambda_cov(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
111 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
112 //off-diag
113 rv::RVec<FCCAnalysesJetConstituentsData> get_d0_z0_cov(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
114 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
115
116 rv::RVec<FCCAnalysesJetConstituentsData> get_phi0_d0_cov(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
117 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
118
119 rv::RVec<FCCAnalysesJetConstituentsData> get_phi0_z0_cov(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
120 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
121
122 rv::RVec<FCCAnalysesJetConstituentsData> get_tanlambda_phi0_cov(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
123 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
124
125 rv::RVec<FCCAnalysesJetConstituentsData> get_tanlambda_d0_cov(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
126 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
127
128 rv::RVec<FCCAnalysesJetConstituentsData> get_tanlambda_z0_cov(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
129 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
130
131 rv::RVec<FCCAnalysesJetConstituentsData> get_omega_tanlambda_cov(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
132 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
133
134 rv::RVec<FCCAnalysesJetConstituentsData> get_omega_phi0_cov(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
135 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
136
137 rv::RVec<FCCAnalysesJetConstituentsData> get_omega_d0_cov(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
138 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
139
140 rv::RVec<FCCAnalysesJetConstituentsData> get_omega_z0_cov(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
141 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
142
143
144 rv::RVec<FCCAnalysesJetConstituentsData> get_dndx(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
145 const rv::RVec<edm4hep::Quantity>& dNdx,
146 const rv::RVec<edm4hep::TrackData>& trackdata,
147 const rv::RVec<FCCAnalysesJetConstituentsData> JetsConstituents_isChargedHad);
148
149 rv::RVec<FCCAnalysesJetConstituentsData> get_Sip2dVal(const rv::RVec<edm4hep::ReconstructedParticleData>& jets,
150 const rv::RVec<FCCAnalysesJetConstituents>& jcs,
151 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
152
153 rv::RVec<FCCAnalysesJetConstituentsData> get_Sip2dVal_cluster(const rv::RVec<fastjet::PseudoJet>& jets,
154 const rv::RVec<FCCAnalysesJetConstituents>& jcs,
155 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
156
157
158 rv::RVec<FCCAnalysesJetConstituentsData> get_Sip2dVal_clusterV(const rv::RVec<fastjet::PseudoJet>& jets,
159 const rv::RVec<FCCAnalysesJetConstituentsData>& D0,
160 const rv::RVec<FCCAnalysesJetConstituentsData>& phi0,
161 const float Bz);
162
163
164 rv::RVec<FCCAnalysesJetConstituentsData> get_Sip2dSig(const rv::RVec<FCCAnalysesJetConstituentsData>& Sip2dVals,
165 const rv::RVec<FCCAnalysesJetConstituentsData>& err2_D0);
166
167 rv::RVec<FCCAnalysesJetConstituentsData> get_Sip3dVal(const rv::RVec<edm4hep::ReconstructedParticleData>& jets,
168 const rv::RVec<FCCAnalysesJetConstituents>& jcs,
169 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
170
171
172 rv::RVec<FCCAnalysesJetConstituentsData> get_Sip3dVal_cluster(const rv::RVec<fastjet::PseudoJet>& jets,
173 const rv::RVec<FCCAnalysesJetConstituents>& jcs,
174 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
175
176 rv::RVec<FCCAnalysesJetConstituentsData> get_Sip3dVal_clusterV(const rv::RVec<fastjet::PseudoJet>& jets,
177 const rv::RVec<FCCAnalysesJetConstituentsData>& D0,
178 const rv::RVec<FCCAnalysesJetConstituentsData>& Z0,
179 const rv::RVec<FCCAnalysesJetConstituentsData>& phi0,
180 const float Bz);
181
182 rv::RVec<FCCAnalysesJetConstituentsData> get_Sip3dSig(const rv::RVec<FCCAnalysesJetConstituentsData>& Sip3dVals,
183 const rv::RVec<FCCAnalysesJetConstituentsData>& err2_D0,
184 const rv::RVec<FCCAnalysesJetConstituentsData>& err2_Z0);
185
186 rv::RVec<FCCAnalysesJetConstituentsData> get_JetDistVal(const rv::RVec<edm4hep::ReconstructedParticleData>& jets,
187 const rv::RVec<FCCAnalysesJetConstituents>& jcs,
188 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
189
190 rv::RVec<FCCAnalysesJetConstituentsData> get_JetDistVal_cluster(const rv::RVec<fastjet::PseudoJet>& jets,
191 const rv::RVec<FCCAnalysesJetConstituents>& jcs,
192 const ROOT::VecOps::RVec<edm4hep::TrackState>& tracks);
193
194 rv::RVec<FCCAnalysesJetConstituentsData> get_JetDistVal_clusterV(const rv::RVec<fastjet::PseudoJet>& jets,
195 const rv::RVec<FCCAnalysesJetConstituents>& jcs,
196 const rv::RVec<FCCAnalysesJetConstituentsData>& D0,
197 const rv::RVec<FCCAnalysesJetConstituentsData>& Z0,
198 const rv::RVec<FCCAnalysesJetConstituentsData>& phi0,
199 const float Bz);
200
201 rv::RVec<FCCAnalysesJetConstituentsData> get_JetDistSig(const rv::RVec<FCCAnalysesJetConstituentsData>& JetDistVal,
202 const rv::RVec<FCCAnalysesJetConstituentsData>& err2_D0,
203 const rv::RVec<FCCAnalysesJetConstituentsData>& err2_Z0);
204
205 rv::RVec<FCCAnalysesJetConstituentsData> get_mtof(const rv::RVec<FCCAnalysesJetConstituents>& jcs,
206 const rv::RVec<float>& track_L,
207 const rv::RVec<edm4hep::TrackData>& trackdata,
208 const rv::RVec<edm4hep::TrackerHit3DData>& trackerhits,
209 const rv::RVec<edm4hep::ClusterData>& gammadata,
210 const rv::RVec<edm4hep::ClusterData>& nhdata,
211 const rv::RVec<edm4hep::CalorimeterHitData>& calohits,
212 const TLorentzVector& V // primary vertex
213 );
214
215
216 rv::RVec<FCCAnalysesJetConstituentsData> get_PIDs(const ROOT::VecOps::RVec< int > recin,
217 const ROOT::VecOps::RVec< int > mcin,
218 const rv::RVec<edm4hep::ReconstructedParticleData>& RecPart,
219 const rv::RVec<edm4hep::MCParticleData>& Particle,
220 const rv::RVec<edm4hep::ReconstructedParticleData>& Jets);
221
222 rv::RVec<FCCAnalysesJetConstituentsData> get_PIDs_cluster(const ROOT::VecOps::RVec< int > recin,
223 const ROOT::VecOps::RVec< int > mcin,
224 const rv::RVec<edm4hep::ReconstructedParticleData>& RecPart,
225 const rv::RVec<edm4hep::MCParticleData>& Particle,
226 const std::vector<std::vector<int>>& indices);
227
228 rv::RVec<FCCAnalysesJetConstituentsData> get_isMu(const rv::RVec<FCCAnalysesJetConstituents>& jcs);
229 rv::RVec<FCCAnalysesJetConstituentsData> get_isEl(const rv::RVec<FCCAnalysesJetConstituents>& jcs);
230 rv::RVec<FCCAnalysesJetConstituentsData> get_isChargedHad(const rv::RVec<FCCAnalysesJetConstituents>& jcs);
231 rv::RVec<FCCAnalysesJetConstituentsData> get_isGamma(const rv::RVec<FCCAnalysesJetConstituents>& jcs);
232 rv::RVec<FCCAnalysesJetConstituentsData> get_isNeutralHad(const rv::RVec<FCCAnalysesJetConstituents>& jcs);
233
234 //countings
235 int count_jets(rv::RVec<FCCAnalysesJetConstituents> jets);
236 rv::RVec<int> count_consts(rv::RVec<FCCAnalysesJetConstituents> jets);
237 rv::RVec<int> count_type(const rv::RVec<FCCAnalysesJetConstituentsData>& isType);
238
239
240
241 rv::RVec<FCCAnalysesJetConstituentsData> get_erel(const rv::RVec<edm4hep::ReconstructedParticleData>& jets,
242 const rv::RVec<FCCAnalysesJetConstituents>& jcs);
243 rv::RVec<FCCAnalysesJetConstituentsData> get_erel_cluster(const rv::RVec<fastjet::PseudoJet>& jets,
244 const rv::RVec<FCCAnalysesJetConstituents>& jcs);
245
246 rv::RVec<FCCAnalysesJetConstituentsData> get_erel_log(const rv::RVec<edm4hep::ReconstructedParticleData>& jets,
247 const rv::RVec<FCCAnalysesJetConstituents>& jcs);
248 rv::RVec<FCCAnalysesJetConstituentsData> get_erel_log_cluster(const rv::RVec<fastjet::PseudoJet>& jets,
249 const rv::RVec<FCCAnalysesJetConstituents>& jcs);
250
251 rv::RVec<FCCAnalysesJetConstituentsData> get_thetarel(const rv::RVec<edm4hep::ReconstructedParticleData>& jets,
252 const rv::RVec<FCCAnalysesJetConstituents>& jcs);
253 rv::RVec<FCCAnalysesJetConstituentsData> get_thetarel_cluster(const rv::RVec<fastjet::PseudoJet>& jets,
254 const rv::RVec<FCCAnalysesJetConstituents>& jcs);
255
256 rv::RVec<FCCAnalysesJetConstituentsData> get_phirel(const rv::RVec<edm4hep::ReconstructedParticleData>& jets,
257 const rv::RVec<FCCAnalysesJetConstituents>& jcs);
258 rv::RVec<FCCAnalysesJetConstituentsData> get_phirel_cluster(const rv::RVec<fastjet::PseudoJet>& jets,
259 const rv::RVec<FCCAnalysesJetConstituents>& jcs);
260
261 //residues
262 rv::RVec<TLorentzVector> compute_tlv_jets(const rv::RVec<fastjet::PseudoJet>& jets);
263 rv::RVec<TLorentzVector> sum_tlv_constituents(const rv::RVec<FCCAnalysesJetConstituents>& jets);
264 float InvariantMass(const TLorentzVector& tlv1, const TLorentzVector& tlv2);
265
269 rv::RVec<double> all_invariant_masses(rv::RVec<TLorentzVector> AllJets);
270 rv::RVec<double> compute_residue_energy(const rv::RVec<TLorentzVector>& tlv_jet,
271 const rv::RVec<TLorentzVector>& sum_tlv_jcs);
272 rv::RVec<double> compute_residue_pt(const rv::RVec<TLorentzVector>& tlv_jet,
273 const rv::RVec<TLorentzVector>& sum_tlv_jcs);
274 rv::RVec<double> compute_residue_phi(const rv::RVec<TLorentzVector>& tlv_jet,
275 const rv::RVec<TLorentzVector>& sum_tlv_jcs);
276 rv::RVec<double> compute_residue_theta(const rv::RVec<TLorentzVector>& tlv_jet,
277 const rv::RVec<TLorentzVector>& sum_tlv_jcs);
278 rv::RVec<double> compute_residue_px(const rv::RVec<TLorentzVector>& tlv_jet, const rv::RVec<TLorentzVector>& sum_tlv_jcs);
279 rv::RVec<double> compute_residue_py(const rv::RVec<TLorentzVector>& tlv_jet, const rv::RVec<TLorentzVector>& sum_tlv_jcs);
280 rv::RVec<double> compute_residue_pz(const rv::RVec<TLorentzVector>& tlv_jet, const rv::RVec<TLorentzVector>& sum_tlv_jcs);
281
282 } // namespace JetConstituentsUtils
283} // namespace FCCAnalyses
284
285#endif
rv::RVec< FCCAnalysesJetConstituentsData > get_isMu(const rv::RVec< FCCAnalysesJetConstituents > &jcs)
Definition JetConstituentsUtils.cc:1111
rv::RVec< FCCAnalysesJetConstituentsData > get_isGamma(const rv::RVec< FCCAnalysesJetConstituents > &jcs)
Definition JetConstituentsUtils.cc:1184
rv::RVec< FCCAnalysesJetConstituentsData > get_erel_cluster(const rv::RVec< fastjet::PseudoJet > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs)
Definition JetConstituentsUtils.cc:919
rv::RVec< FCCAnalysesJetConstituentsData > get_phirel(const rv::RVec< edm4hep::ReconstructedParticleData > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs)
Definition JetConstituentsUtils.cc:990
rv::RVec< FCCAnalysesJetConstituentsData > get_thetarel(const rv::RVec< edm4hep::ReconstructedParticleData > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs)
Definition JetConstituentsUtils.cc:938
rv::RVec< FCCAnalysesJetConstituentsData > get_isNeutralHad(const rv::RVec< FCCAnalysesJetConstituents > &jcs)
Definition JetConstituentsUtils.cc:1159
rv::RVec< FCCAnalysesJetConstituentsData > get_z0(const rv::RVec< FCCAnalysesJetConstituents > &, const ROOT::VecOps::RVec< edm4hep::TrackState > &)
Definition JetConstituentsUtils.cc:204
rv::RVec< FCCAnalysesJetConstituentsData > get_JetDistVal(const rv::RVec< edm4hep::ReconstructedParticleData > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:620
rv::RVec< FCCAnalysesJetConstituentsData > get_theta(const rv::RVec< FCCAnalysesJetConstituents > &)
Definition JetConstituentsUtils.cc:144
rv::RVec< FCCAnalysesJetConstituentsData > get_pt(const rv::RVec< FCCAnalysesJetConstituents > &)
Definition JetConstituentsUtils.cc:129
float InvariantMass(const TLorentzVector &tlv1, const TLorentzVector &tlv2)
Definition JetConstituentsUtils.cc:1273
rv::RVec< FCCAnalysesJetConstituentsData > get_phi0_d0_cov(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:309
rv::RVec< TLorentzVector > sum_tlv_constituents(const rv::RVec< FCCAnalysesJetConstituents > &jets)
Definition JetConstituentsUtils.cc:1255
rv::RVec< FCCAnalysesJetConstituentsData > get_JetDistVal_clusterV(const rv::RVec< fastjet::PseudoJet > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs, const rv::RVec< FCCAnalysesJetConstituentsData > &D0, const rv::RVec< FCCAnalysesJetConstituentsData > &Z0, const rv::RVec< FCCAnalysesJetConstituentsData > &phi0, const float Bz)
Definition JetConstituentsUtils.cc:686
rv::RVec< double > all_invariant_masses(rv::RVec< TLorentzVector > AllJets)
all_invariant_masses takes an RVec of TLorentzVectors of jets and computes the invariant masses of al...
Definition JetConstituentsUtils.cc:1283
rv::RVec< TLorentzVector > compute_tlv_jets(const rv::RVec< fastjet::PseudoJet > &jets)
Definition JetConstituentsUtils.cc:1243
rv::RVec< FCCAnalysesJetConstituentsData > get_omega(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:216
rv::RVec< FCCAnalysesJetConstituentsData > get_phi0_cov(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:291
rv::RVec< edm4hep::ReconstructedParticleData > FCCAnalysesJetConstituents
Definition JetConstituentsUtils.h:27
rv::RVec< FCCAnalysesJetConstituentsData > get_tanlambda_z0_cov(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:333
rv::RVec< FCCAnalysesJetConstituents > build_constituents(const rv::RVec< edm4hep::ReconstructedParticleData > &, const rv::RVec< edm4hep::ReconstructedParticleData > &)
Build the collection of constituents (mapping jet -> reconstructed particles) for all jets in event.
Definition JetConstituentsUtils.cc:30
rv::RVec< FCCAnalysesJetConstituents > build_constituents_cluster(const rv::RVec< edm4hep::ReconstructedParticleData > &rps, const std::vector< std::vector< int > > &indices)
Definition JetConstituentsUtils.cc:48
ROOT::VecOps::RVec< edm4hep::ReconstructedParticleData > jets_sorting_on_nconst(const rv::RVec< edm4hep::ReconstructedParticleData > &)
Definition JetConstituentsUtils.cc:165
rv::RVec< FCCAnalysesJetConstituentsData > get_Sip3dSig(const rv::RVec< FCCAnalysesJetConstituentsData > &Sip3dVals, const rv::RVec< FCCAnalysesJetConstituentsData > &err2_D0, const rv::RVec< FCCAnalysesJetConstituentsData > &err2_Z0)
Definition JetConstituentsUtils.cc:596
rv::RVec< FCCAnalysesJetConstituentsData > XPtoPar_dz(const rv::RVec< FCCAnalysesJetConstituents > &, const ROOT::VecOps::RVec< edm4hep::TrackState > &, const TLorentzVector &V, const float &)
Definition JetConstituentsUtils.cc:237
rv::RVec< FCCAnalysesJetConstituentsData > get_PIDs_cluster(const ROOT::VecOps::RVec< int > recin, const ROOT::VecOps::RVec< int > mcin, const rv::RVec< edm4hep::ReconstructedParticleData > &RecPart, const rv::RVec< edm4hep::MCParticleData > &Particle, const std::vector< std::vector< int > > &indices)
Definition JetConstituentsUtils.cc:1065
rv::RVec< FCCAnalysesJetConstituentsData > get_phirel_cluster(const rv::RVec< fastjet::PseudoJet > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs)
Definition JetConstituentsUtils.cc:1016
rv::RVec< FCCAnalysesJetConstituentsData > get_Sip2dVal(const rv::RVec< edm4hep::ReconstructedParticleData > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:396
rv::RVec< FCCAnalysesJetConstituentsData > get_erel_log_cluster(const rv::RVec< fastjet::PseudoJet > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs)
Definition JetConstituentsUtils.cc:881
rv::RVec< FCCAnalysesJetConstituentsData > get_omega_phi0_cov(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:345
rv::RVec< FCCAnalysesJetConstituentsData > get_d0_z0_cov(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:303
rv::RVec< FCCAnalysesJetConstituentsData > get_isEl(const rv::RVec< FCCAnalysesJetConstituents > &jcs)
Definition JetConstituentsUtils.cc:1087
rv::RVec< double > compute_residue_py(const rv::RVec< TLorentzVector > &tlv_jet, const rv::RVec< TLorentzVector > &sum_tlv_jcs)
Definition JetConstituentsUtils.cc:1337
rv::RVec< FCCAnalysesJetConstituentsData > get_tanlambda_phi0_cov(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:321
rv::RVec< FCCAnalysesJetConstituentsData > get_Sip2dSig(const rv::RVec< FCCAnalysesJetConstituentsData > &Sip2dVals, const rv::RVec< FCCAnalysesJetConstituentsData > &err2_D0)
The functions get_Sip2dSig and get_Sip2dVal can be made independent; I passed to the former the resul...
Definition JetConstituentsUtils.cc:484
rv::RVec< FCCAnalysesJetConstituentsData > get_erel_log(const rv::RVec< edm4hep::ReconstructedParticleData > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs)
Definition JetConstituentsUtils.cc:862
rv::RVec< FCCAnalysesJetConstituentsData > get_Sip3dVal(const rv::RVec< edm4hep::ReconstructedParticleData > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:507
rv::RVec< FCCAnalysesJetConstituentsData > get_phi0(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:210
rv::RVec< FCCAnalysesJetConstituentsData > get_tanlambda_d0_cov(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:327
rv::RVec< int > count_type(const rv::RVec< FCCAnalysesJetConstituentsData > &isType)
Definition JetConstituentsUtils.cc:1225
rv::RVec< double > compute_residue_phi(const rv::RVec< TLorentzVector > &tlv_jet, const rv::RVec< TLorentzVector > &sum_tlv_jcs)
Definition JetConstituentsUtils.cc:1372
rv::RVec< FCCAnalysesJetConstituentsData > get_JetDistVal_cluster(const rv::RVec< fastjet::PseudoJet > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:653
rv::RVec< FCCAnalysesJetConstituentsData > get_omega_d0_cov(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:351
rv::RVec< FCCAnalysesJetConstituentsData > get_phi0_z0_cov(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:315
rv::RVec< FCCAnalysesJetConstituentsData > get_JetDistSig(const rv::RVec< FCCAnalysesJetConstituentsData > &JetDistVal, const rv::RVec< FCCAnalysesJetConstituentsData > &err2_D0, const rv::RVec< FCCAnalysesJetConstituentsData > &err2_Z0)
Definition JetConstituentsUtils.cc:720
rv::RVec< FCCAnalysesJetConstituentsData > XPtoPar_dxy(const rv::RVec< FCCAnalysesJetConstituents > &, const ROOT::VecOps::RVec< edm4hep::TrackState > &, const TLorentzVector &V, const float &)
Definition JetConstituentsUtils.cc:228
rv::RVec< int > count_consts(rv::RVec< FCCAnalysesJetConstituents > jets)
Definition JetConstituentsUtils.cc:1215
rv::RVec< float > FCCAnalysesJetConstituentsData
Definition JetConstituentsUtils.h:28
rv::RVec< FCCAnalysesJetConstituentsData > get_tanlambda_cov(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:297
rv::RVec< FCCAnalysesJetConstituentsData > get_d0_cov(const rv::RVec< FCCAnalysesJetConstituents > &, const ROOT::VecOps::RVec< edm4hep::TrackState > &)
Definition JetConstituentsUtils.cc:279
rv::RVec< FCCAnalysesJetConstituentsData > get_Sip2dVal_cluster(const rv::RVec< fastjet::PseudoJet > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:425
rv::RVec< FCCAnalysesJetConstituentsData > get_thetarel_cluster(const rv::RVec< fastjet::PseudoJet > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs)
Definition JetConstituentsUtils.cc:964
rv::RVec< FCCAnalysesJetConstituentsData > get_d0(const rv::RVec< FCCAnalysesJetConstituents > &, const ROOT::VecOps::RVec< edm4hep::TrackState > &)
Definition JetConstituentsUtils.cc:198
rv::RVec< FCCAnalysesJetConstituentsData > get_omega_z0_cov(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:357
rv::RVec< FCCAnalysesJetConstituentsData > get_tanLambda(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:222
rv::RVec< FCCAnalysesJetConstituentsData > get_omega_cov(const rv::RVec< FCCAnalysesJetConstituents > &, const ROOT::VecOps::RVec< edm4hep::TrackState > &)
Definition JetConstituentsUtils.cc:273
rv::RVec< FCCAnalysesJetConstituentsData > get_Sip3dVal_clusterV(const rv::RVec< fastjet::PseudoJet > &jets, const rv::RVec< FCCAnalysesJetConstituentsData > &D0, const rv::RVec< FCCAnalysesJetConstituentsData > &Z0, const rv::RVec< FCCAnalysesJetConstituentsData > &phi0, const float Bz)
Definition JetConstituentsUtils.cc:567
rv::RVec< double > compute_residue_pz(const rv::RVec< TLorentzVector > &tlv_jet, const rv::RVec< TLorentzVector > &sum_tlv_jcs)
Definition JetConstituentsUtils.cc:1348
rv::RVec< FCCAnalysesJetConstituentsData > XPtoPar_C(const rv::RVec< FCCAnalysesJetConstituents > &, const ROOT::VecOps::RVec< edm4hep::TrackState > &, const float &)
Definition JetConstituentsUtils.cc:255
rv::RVec< FCCAnalysesJetConstituentsData > get_Bz(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:123
rv::RVec< FCCAnalysesJetConstituentsData > get_e(const rv::RVec< FCCAnalysesJetConstituents > &)
Definition JetConstituentsUtils.cc:139
rv::RVec< FCCAnalysesJetConstituentsData > get_isChargedHad(const rv::RVec< FCCAnalysesJetConstituents > &jcs)
Definition JetConstituentsUtils.cc:1135
rv::RVec< double > compute_residue_energy(const rv::RVec< TLorentzVector > &tlv_jet, const rv::RVec< TLorentzVector > &sum_tlv_jcs)
Definition JetConstituentsUtils.cc:1315
rv::RVec< double > compute_residue_theta(const rv::RVec< TLorentzVector > &tlv_jet, const rv::RVec< TLorentzVector > &sum_tlv_jcs)
Definition JetConstituentsUtils.cc:1385
rv::RVec< FCCAnalysesJetConstituents > get_constituents(const rv::RVec< FCCAnalysesJetConstituents > &, const rv::RVec< int > &)
Retrieve the constituents of a collection of indexed jets in event.
Definition JetConstituentsUtils.cc:71
rv::RVec< FCCAnalysesJetConstituentsData > get_erel(const rv::RVec< edm4hep::ReconstructedParticleData > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs)
Definition JetConstituentsUtils.cc:900
ROOT::VecOps::RVec< edm4hep::ReconstructedParticleData > jets_sorting_on_energy(const rv::RVec< edm4hep::ReconstructedParticleData > &)
Definition JetConstituentsUtils.cc:181
int count_jets(rv::RVec< FCCAnalysesJetConstituents > jets)
Definition JetConstituentsUtils.cc:1210
rv::RVec< double > compute_residue_px(const rv::RVec< TLorentzVector > &tlv_jet, const rv::RVec< TLorentzVector > &sum_tlv_jcs)
Definition JetConstituentsUtils.cc:1326
rv::RVec< FCCAnalysesJetConstituentsData > get_Sip3dVal_cluster(const rv::RVec< fastjet::PseudoJet > &jets, const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:537
rv::RVec< FCCAnalysesJetConstituentsData > get_dndx(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const rv::RVec< edm4hep::Quantity > &dNdx, const rv::RVec< edm4hep::TrackData > &trackdata, const rv::RVec< FCCAnalysesJetConstituentsData > JetsConstituents_isChargedHad)
Definition JetConstituentsUtils.cc:365
rv::RVec< FCCAnalysesJetConstituentsData > get_z0_cov(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:285
rv::RVec< FCCAnalysesJetConstituentsData > get_omega_tanlambda_cov(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const ROOT::VecOps::RVec< edm4hep::TrackState > &tracks)
Definition JetConstituentsUtils.cc:339
rv::RVec< FCCAnalysesJetConstituentsData > XPtoPar_phi(const rv::RVec< FCCAnalysesJetConstituents > &, const ROOT::VecOps::RVec< edm4hep::TrackState > &, const TLorentzVector &V, const float &)
Definition JetConstituentsUtils.cc:246
rv::RVec< FCCAnalysesJetConstituentsData > get_charge(const rv::RVec< FCCAnalysesJetConstituents > &)
Definition JetConstituentsUtils.cc:159
FCCAnalysesJetConstituents get_jet_constituents(const rv::RVec< FCCAnalysesJetConstituents > &, int)
Retrieve the constituents of an indexed jet in event.
Definition JetConstituentsUtils.cc:64
rv::RVec< FCCAnalysesJetConstituentsData > get_p(const rv::RVec< FCCAnalysesJetConstituents > &)
Definition JetConstituentsUtils.cc:134
rv::RVec< FCCAnalysesJetConstituentsData > get_type(const rv::RVec< FCCAnalysesJetConstituents > &)
Definition JetConstituentsUtils.cc:154
rv::RVec< FCCAnalysesJetConstituentsData > get_phi(const rv::RVec< FCCAnalysesJetConstituents > &)
Definition JetConstituentsUtils.cc:149
rv::RVec< FCCAnalysesJetConstituentsData > get_mtof(const rv::RVec< FCCAnalysesJetConstituents > &jcs, const rv::RVec< float > &track_L, const rv::RVec< edm4hep::TrackData > &trackdata, const rv::RVec< edm4hep::TrackerHit3DData > &trackerhits, const rv::RVec< edm4hep::ClusterData > &gammadata, const rv::RVec< edm4hep::ClusterData > &nhdata, const rv::RVec< edm4hep::CalorimeterHitData > &calohits, const TLorentzVector &V)
Definition JetConstituentsUtils.cc:758
rv::RVec< FCCAnalysesJetConstituentsData > XPtoPar_ct(const rv::RVec< FCCAnalysesJetConstituents > &, const ROOT::VecOps::RVec< edm4hep::TrackState > &, const float &)
Definition JetConstituentsUtils.cc:263
rv::RVec< FCCAnalysesJetConstituentsData > get_PIDs(const ROOT::VecOps::RVec< int > recin, const ROOT::VecOps::RVec< int > mcin, const rv::RVec< edm4hep::ReconstructedParticleData > &RecPart, const rv::RVec< edm4hep::MCParticleData > &Particle, const rv::RVec< edm4hep::ReconstructedParticleData > &Jets)
Definition JetConstituentsUtils.cc:1044
rv::RVec< double > compute_residue_pt(const rv::RVec< TLorentzVector > &tlv_jet, const rv::RVec< TLorentzVector > &sum_tlv_jcs)
Definition JetConstituentsUtils.cc:1359
rv::RVec< FCCAnalysesJetConstituentsData > get_Sip2dVal_clusterV(const rv::RVec< fastjet::PseudoJet > &jets, const rv::RVec< FCCAnalysesJetConstituentsData > &D0, const rv::RVec< FCCAnalysesJetConstituentsData > &phi0, const float Bz)
Definition JetConstituentsUtils.cc:454
FCC analyzers collection.
Definition Algorithms.h:15
Definition JetConstituentsUtils.h:12
edm4hep::TrackerHitData TrackerHit3DData
Definition JetConstituentsUtils.h:13