ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/CmsHi/JetAnalysis/modifiedFiles/PFBlockAlgo.h
Revision: 1.1
Committed: Thu Sep 15 22:04:26 2011 UTC (13 years, 7 months ago) by yilmaz
Content type: text/plain
Branch: MAIN
CVS Tags: HiForest_V02_85, HiForest_V02_84, HiForest_V02_83, HiForest_V02_82, HiForest_V02_81, HiForest_V02_80, HiForest_V02_79, HiForest_V02_78, HiForest_V02_77, HiForest_V02_76, HiForest_V02_75, HiForest_V02_74, HiForest_V02_73, HiForest_V02_72, HiForest_V02_71, HiForest_V02_70, HiForest_V02_69, HiForest_V02_68, HiForest_V02_67, HiForest_V02_66, HiForest_V02_65, HiForest_V02_64, HiForest_V02_63, HiForest_V02_62, HiForest_V02_61, HiForest_V02_60, HiForest_V02_59, HiForest_V02_58, HiForest_V02_57, HiForest_V02_56, HiForest_V02_55, HiForest_V02_54, HiForest_V02_53, HiForest_V02_52, HiForest_V02_51, HiForest_V02_50, HiForest_V02_49, HiForest_V02_48, HiForest_V02_47, HiForest_V02_46, HiForest_V02_45, HiForest_V02_44, HiForest_V02_43, HiForest_V02_42, HiForest_V02_41, HiForest_V02_40, HiForest_V02_39, HiForest_V02_38, HiForest_V02_37, HiForest_V02_36, HiForest_V02_35, HiForest_V02_34, HiForest_V02_33, HiForest_V02_32, HiForest_V02_31, HiForest_V02_30, HiForest_V02_27, HiForest_V02_26, QM_2012, HiForest_V02_25, HiForest_V02_24, HiForest_V02_23, HiForest_V02_22, HiForest_V02_21, HiForest_V02_20, HiForest_V02_19, HiForest_V02_18, HiForest_V02_17, HiForest_V02_16, HiForest_V02_15, HiForest_V02_14, HiForest_V02_13, HiForest_V02_12, HiForest_V02_11, HiForest_V02_10, HiForest_V02_09, HiForest_V02_08, HiForest_V02_07, HiForest_V02_06, HiForest_V02_05, HiForest_V02_04, HiForest_V02_03, HiForest_V02_02, HiForest_V02_01, HiForest_V02_00, hi44X_02, hi413_03, hi441_1, hi441_0, hi413_11, hi413_10, hi413_09, hi413_08, hi413_07, hi413_06, hi413_05, hi413_04, hi413_02, hi39X_01, HEAD
Branch point for: branch_44x
Log Message:
muons

File Contents

# Content
1 #ifndef RecoParticleFlow_PFProducer_PFBlockAlgo_h
2 #define RecoParticleFlow_PFProducer_PFBlockAlgo_h
3
4 #include <set>
5 #include <vector>
6 #include <iostream>
7
8 // #include "FWCore/Framework/interface/Handle.h"
9 #include "DataFormats/Common/interface/Handle.h"
10 // #include "FWCore/Framework/interface/OrphanHandle.h"
11 #include "DataFormats/Common/interface/OrphanHandle.h"
12
13
14 #include "DataFormats/ParticleFlowReco/interface/PFClusterFwd.h"
15 #include "DataFormats/ParticleFlowReco/interface/PFRecTrackFwd.h"
16 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
17 #include "DataFormats/ParticleFlowReco/interface/PFRecTrack.h"
18 #include "DataFormats/ParticleFlowReco/interface/PFDisplacedTrackerVertex.h" // gouzevitch
19 #include "DataFormats/ParticleFlowReco/interface/PFConversionFwd.h"
20 #include "DataFormats/ParticleFlowReco/interface/PFConversion.h"
21 #include "DataFormats/ParticleFlowReco/interface/PFV0Fwd.h"
22 #include "DataFormats/ParticleFlowReco/interface/PFV0.h"
23 #include "DataFormats/ParticleFlowReco/interface/GsfPFRecTrack.h"
24 #include "DataFormats/ParticleFlowReco/interface/GsfPFRecTrackFwd.h"
25 #include "DataFormats/ParticleFlowReco/interface/PFBrem.h"
26 #include "DataFormats/ParticleFlowReco/interface/PFTrajectoryPoint.h"
27
28 #include "DataFormats/ParticleFlowReco/interface/PFBlockElement.h"
29 #include "DataFormats/ParticleFlowReco/interface/PFBlock.h"
30 #include "DataFormats/ParticleFlowReco/interface/PFBlockFwd.h"
31
32 // #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
33
34 #include "RecoParticleFlow/PFProducer/interface/PFMuonAlgo.h"
35 #include "RecoParticleFlow/PFClusterTools/interface/PFResolutionMap.h"
36
37 #include "DataFormats/MuonReco/interface/MuonFwd.h"
38 #include "DataFormats/MuonReco/interface/Muon.h"
39 #include "DataFormats/EgammaReco/interface/SuperClusterFwd.h"
40 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
41 #include "DataFormats/GsfTrackReco/interface/GsfTrackFwd.h"
42 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
43 #include "DataFormats/EgammaReco/interface/ElectronSeedFwd.h"
44 #include "DataFormats/EgammaReco/interface/ElectronSeed.h"
45 #include "RecoParticleFlow/PFClusterTools/interface/ClusterClusterMapping.h"
46
47 #include "RecoParticleFlow/PFProducer/interface/PFBlockLink.h"
48
49 #include <map>
50
51
52 /// \brief Particle Flow Algorithm
53 /*!
54 \author Colin Bernet
55 \date January 2006
56 */
57
58 class PFBlockAlgo {
59
60 public:
61
62 PFBlockAlgo();
63
64 ~PFBlockAlgo();
65
66
67 void setParameters( std::vector<double>& DPtovPtCut,
68 std::vector<unsigned>& NHitCut,
69 bool useConvBremPFRecTracks,
70 bool useIterTracking,
71 int nuclearInteractionsPurity);
72
73 typedef std::vector<bool> Mask;
74
75 /// set input collections of tracks and clusters
76 template< template<typename> class T>
77 void setInput(const T<reco::PFRecTrackCollection>& trackh,
78 const T<reco::GsfPFRecTrackCollection>& gsftrackh,
79 const T<reco::GsfPFRecTrackCollection>& convbremgsftrackh,
80 const T<reco::MuonCollection>& muonh,
81 const T<reco::PFDisplacedTrackerVertexCollection>& displacedh,
82 const T<reco::PFConversionCollection>& conv,
83 const T<reco::PFV0Collection>& v0,
84 const T<reco::PFClusterCollection>& ecalh,
85 const T<reco::PFClusterCollection>& hcalh,
86 const T<reco::PFClusterCollection>& hfemh,
87 const T<reco::PFClusterCollection>& hfhadh,
88 const T<reco::PFClusterCollection>& psh,
89 const Mask& trackMask = dummyMask_,
90 const Mask& gsftrackMask = dummyMask_,
91 const Mask& ecalMask = dummyMask_,
92 const Mask& hcalMask = dummyMask_,
93 const Mask& hfemMask = dummyMask_,
94 const Mask& hfhadMask = dummyMask_,
95 const Mask& psMask = dummyMask_ );
96
97 ///COLIN: I think this is for particle flow at HLT...
98 template< template<typename> class T >
99 void setInput(const T<reco::PFRecTrackCollection>& trackh,
100 const T<reco::MuonCollection>& muonh,
101 const T<reco::PFClusterCollection>& ecalh,
102 const T<reco::PFClusterCollection>& hcalh,
103 const T<reco::PFClusterCollection>& hfemh,
104 const T<reco::PFClusterCollection>& hfhadh,
105 const T<reco::PFClusterCollection>& psh,
106 const Mask& trackMask = dummyMask_,
107 const Mask& ecalMask = dummyMask_,
108 const Mask& hcalMask = dummyMask_,
109 const Mask& psMask = dummyMask_ ) {
110 T<reco::GsfPFRecTrackCollection> gsftrackh;
111 T<reco::GsfPFRecTrackCollection> convbremgsftrackh;
112 //T<reco::MuonCollection> muonh;
113 T<reco::PFDisplacedTrackerVertexCollection> displacedh;
114 T<reco::PFConversionCollection> convh;
115 T<reco::PFV0Collection> v0;
116 setInput<T>( trackh, gsftrackh, convbremgsftrackh, muonh, displacedh, convh, v0,
117 ecalh, hcalh, hfemh, hfhadh, psh,
118 trackMask, ecalMask, hcalMask, psMask);
119 }
120
121 ///COLIN: what is this setinput function for? can it be removed?
122 template< template<typename> class T >
123 void setInput(const T<reco::PFRecTrackCollection>& trackh,
124 const T<reco::GsfPFRecTrackCollection>& gsftrackh,
125 const T<reco::PFClusterCollection>& ecalh,
126 const T<reco::PFClusterCollection>& hcalh,
127 const T<reco::PFClusterCollection>& psh,
128 const Mask& trackMask = dummyMask_,
129 const Mask& gsftrackMask = dummyMask_,
130 const Mask& ecalMask = dummyMask_,
131 const Mask& hcalMask = dummyMask_,
132 const Mask& psMask = dummyMask_ ) {
133 T<reco::GsfPFRecTrackCollection> convbremgsftrackh;
134 T<reco::MuonCollection> muonh;
135 T<reco::PFDisplacedTrackerVertexCollection>& displacedh;
136 T<reco::PFConversionCollection> convh;
137 T<reco::PFV0Collection> v0;
138 setInput<T>( trackh, gsftrackh, convbremgsftrackh, muonh, displacedh, convh, v0, ecalh, hcalh, psh,
139 trackMask, gsftrackMask,ecalMask, hcalMask, psMask);
140 }
141
142
143 /// sets debug printout flag
144 void setDebug( bool debug ) {debug_ = debug;}
145
146 /// build blocks
147 void findBlocks();
148
149
150 /// \return collection of blocks
151 /* const reco::PFBlockCollection& blocks() const {return *blocks_;} */
152 const std::auto_ptr< reco::PFBlockCollection >& blocks() const
153 {return blocks_;}
154
155 /// \return auto_ptr to collection of blocks
156 std::auto_ptr< reco::PFBlockCollection > transferBlocks() {return blocks_;}
157
158 /// define these in *Fwd files in DataFormats/ParticleFlowReco?
159 typedef std::list< reco::PFBlockElement* >::iterator IE;
160 typedef std::list< reco::PFBlockElement* >::const_iterator IEC;
161 typedef reco::PFBlockCollection::const_iterator IBC;
162
163 private:
164
165 /// recursive procedure which adds elements from
166 /// elements_ to the current block, ie blocks_->back().
167 /// the resulting links between elements are stored in links,
168 /// not in the block. afterwards,
169 /// packLinks( reco::PFBlock& block, const vector<PFBlockLink>& links)
170 /// has to be called in order to pack the link information in the block.
171 IE associate(IE next, IE last, std::vector<PFBlockLink>& links);
172
173 /// compute missing links in the blocks
174 /// (the recursive procedure does not build all links)
175 void packLinks(reco::PFBlock& block,
176 const std::vector<PFBlockLink>& links) const;
177
178 /// remove extra links between primary track and clusters
179 void checkDisplacedVertexLinks( reco::PFBlock& block ) const;
180
181 ///COLIN: not used. Could be removed.
182 /// Could also be implemented, to produce a graph of a block,
183 /// Showing the connections between the elements
184 void buildGraph();
185
186 /// check whether 2 elements are linked. Returns distance and linktype
187 void link( const reco::PFBlockElement* el1,
188 const reco::PFBlockElement* el2,
189 PFBlockLink::Type& linktype,
190 reco::PFBlock::LinkTest& linktest,
191 double& dist) const;
192
193 /// tests association between a track and a PS cluster
194 /// returns distance
195 double testTrackAndPS(const reco::PFRecTrack& track,
196 const reco::PFCluster& ps) const;
197
198 /// tests association between an ECAL and an HCAL cluster
199 /// \returns distance
200 double testECALAndHCAL(const reco::PFCluster& ecal,
201 const reco::PFCluster& hcal) const;
202
203 /// tests association between a PS1 v cluster and a PS2 h cluster
204 /// returns distance
205 double testPS1AndPS2(const reco::PFCluster& ps1,
206 const reco::PFCluster& ps2) const;
207
208 //tests association between a track and a cluster by rechit
209 double testTrackAndClusterByRecHit( const reco::PFRecTrack& track,
210 const reco::PFCluster& cluster,
211 bool isBrem = false) const;
212
213 //tests association between ECAL and PS clusters by rechit
214 double testECALAndPSByRecHit( const reco::PFCluster& clusterECAL,
215 const reco::PFCluster& clusterPS) const;
216
217 /// test association between HFEM and HFHAD, by rechit
218 double testHFEMAndHFHADByRecHit( const reco::PFCluster& clusterHFEM,
219 const reco::PFCluster& clusterHFHAD) const;
220
221 /// test association by Supercluster between two ECAL
222 double testLinkBySuperCluster(const reco::PFClusterRef & elt1,
223 const reco::PFClusterRef & elt2) const;
224
225 /// computes a chisquare
226 double computeDist( double eta1, double phi1,
227 double eta2, double phi2 ) const;
228
229 /// checks size of the masks with respect to the vectors
230 /// they refer to. throws std::length_error if one of the
231 /// masks has the wrong size
232 void checkMaskSize( const reco::PFRecTrackCollection& tracks,
233 const reco::GsfPFRecTrackCollection& gsftracks,
234 const reco::PFClusterCollection& ecals,
235 const reco::PFClusterCollection& hcals,
236 const reco::PFClusterCollection& hfems,
237 const reco::PFClusterCollection& hfhads,
238 const reco::PFClusterCollection& pss,
239 const Mask& trackMask,
240 const Mask& gsftrackMask,
241 const Mask& ecalMask,
242 const Mask& hcalMask,
243 const Mask& hfemMask,
244 const Mask& hfhadMask,
245 const Mask& psMask ) const;
246
247 /// open a resolution map
248 // PFResolutionMap* openResolutionMap(const char* resMapName);
249
250 /// check the Pt resolution
251 bool goodPtResolution( const reco::TrackRef& trackref);
252
253 double testLinkByVertex(const reco::PFBlockElement* elt1,
254 const reco::PFBlockElement* elt2) const;
255
256 int muAssocToTrack( const reco::TrackRef& trackref,
257 const edm::Handle<reco::MuonCollection>& muonh) const;
258 int muAssocToTrack( const reco::TrackRef& trackref,
259 const edm::OrphanHandle<reco::MuonCollection>& muonh) const;
260
261
262 std::auto_ptr< reco::PFBlockCollection > blocks_;
263
264 /// actually, particles will be created by a separate producer
265 // std::vector< reco::PFCandidate > particles_;
266
267 // the test elements will be transferred to the blocks
268 std::list< reco::PFBlockElement* > elements_;
269
270 static const Mask dummyMask_;
271
272 /// DPt/Pt cut for creating atrack element
273 std::vector<double> DPtovPtCut_;
274
275 /// Number of layers crossed cut for creating atrack element
276 std::vector<unsigned> NHitCut_;
277
278 /// Flag to turn off quality cuts which require iterative tracking (for heavy-ions)
279 bool useIterTracking_;
280
281
282 // This parameters defines the level of purity of
283 // nuclear interactions choosen.
284 // Level 1 is only high Purity sample labeled as isNucl
285 // Level 2 isNucl + isNucl_Loose (2 secondary tracks vertices)
286 // Level 3 isNucl + isNucl_Loose + isNucl_Kink
287 // (low purity sample made of 1 primary and 1 secondary track)
288 // By default the level 1 is teh safest one.
289 int nuclearInteractionsPurity_;
290
291 /// switch on/off Conversions Brem Recovery with KF Tracks
292 bool useConvBremPFRecTracks_;
293
294 /// PS strip resolution
295 double resPSpitch_;
296
297 /// PS resolution along strip
298 double resPSlength_;
299
300 /// list of superclusters
301 std::vector<const reco::SuperCluster *> superClusters_;
302
303 /// SC corresponding to the PF cluster
304 // std::map<reco::PFClusterRef,int> pfcRefSCMap_;
305 std::vector<int> pfcSCVec_;
306
307 /// PF clusters corresponding to a given SC
308 std::vector<std::vector<reco::PFClusterRef> > scpfcRefs_;
309 /// if true, debug printouts activated
310 bool debug_;
311
312 friend std::ostream& operator<<(std::ostream&, const PFBlockAlgo&);
313
314 };
315
316 #include "DataFormats/ParticleFlowReco/interface/PFBlockElementGsfTrack.h"
317 #include "DataFormats/ParticleFlowReco/interface/PFBlockElementBrem.h"
318 #include "DataFormats/ParticleFlowReco/interface/PFBlockElementTrack.h"
319 #include "DataFormats/ParticleFlowReco/interface/PFBlockElementCluster.h"
320 #include "DataFormats/ParticleFlowReco/interface/PFLayer.h"
321
322 template< template<typename> class T >
323 void
324 PFBlockAlgo::setInput(const T<reco::PFRecTrackCollection>& trackh,
325 const T<reco::GsfPFRecTrackCollection>& gsftrackh,
326 const T<reco::GsfPFRecTrackCollection>& convbremgsftrackh,
327 const T<reco::MuonCollection>& muonh,
328 const T<reco::PFDisplacedTrackerVertexCollection>& displacedh,
329 const T<reco::PFConversionCollection>& convh,
330 const T<reco::PFV0Collection>& v0,
331 const T<reco::PFClusterCollection>& ecalh,
332 const T<reco::PFClusterCollection>& hcalh,
333 const T<reco::PFClusterCollection>& hfemh,
334 const T<reco::PFClusterCollection>& hfhadh,
335 const T<reco::PFClusterCollection>& psh,
336 const Mask& trackMask,
337 const Mask& gsftrackMask,
338 const Mask& ecalMask,
339 const Mask& hcalMask,
340 const Mask& hfemMask,
341 const Mask& hfhadMask,
342 const Mask& psMask ) {
343
344
345 checkMaskSize( *trackh,
346 *gsftrackh,
347 *ecalh,
348 *hcalh,
349 *hfemh,
350 *hfhadh,
351 *psh,
352 trackMask,
353 gsftrackMask,
354 ecalMask,
355 hcalMask,
356 hfemMask,
357 hfhadMask,
358 psMask );
359
360
361
362
363 /// -------------- GSF Primary tracks and brems ---------------------
364 std::vector<reco::PFRecTrackRef> convBremPFRecTracks;
365 convBremPFRecTracks.clear();
366 std::vector<reco::PFRecTrackRef> primaryKF_GSF;
367 primaryKF_GSF.clear();
368 // Super cluster mapping
369 superClusters_.clear();
370 scpfcRefs_.clear();
371 pfcSCVec_.clear();
372
373 if(gsftrackh.isValid() ) {
374 const reco::GsfPFRecTrackCollection PFGsfProd = *(gsftrackh.product());
375 for(unsigned i=0;i<gsftrackh->size(); i++) {
376 if( !gsftrackMask.empty() &&
377 !gsftrackMask[i] ) continue;
378 reco::GsfPFRecTrackRef refgsf(gsftrackh,i );
379
380 if((refgsf).isNull()) continue;
381 reco::GsfTrackRef gsf=refgsf->gsfTrackRef();
382
383 // retrieve and save the SC if ECAL-driven - Florian
384 if(gsf->extra().isAvailable() && gsf->extra()->seedRef().isAvailable()) {
385 reco::ElectronSeedRef seedRef= gsf->extra()->seedRef().castTo<reco::ElectronSeedRef>();
386 // check that the seed is valid
387 if(seedRef.isAvailable() && seedRef->isEcalDriven()) {
388 reco::SuperClusterRef scRef = seedRef->caloCluster().castTo<reco::SuperClusterRef>();
389 if(scRef.isNonnull()) {
390 superClusters_.push_back(&(*scRef));
391 }
392 }
393 }
394
395 reco::PFBlockElement* gsfEl;
396
397 const std::vector<reco::PFTrajectoryPoint>
398 PfGsfPoint = PFGsfProd[i].trajectoryPoints();
399
400 unsigned int c_gsf=0;
401 bool PassTracker = false;
402 bool GetPout = false;
403 unsigned int IndexPout = 0;
404
405 typedef std::vector<reco::PFTrajectoryPoint>::const_iterator IP;
406 for(IP itPfGsfPoint = PfGsfPoint.begin();
407 itPfGsfPoint!= PfGsfPoint.end();itPfGsfPoint++) {
408
409 if (itPfGsfPoint->isValid()){
410 int layGsfP = itPfGsfPoint->layer();
411 if (layGsfP == -1) PassTracker = true;
412 if (PassTracker && layGsfP > 0 && GetPout == false) {
413 IndexPout = c_gsf-1;
414 GetPout = true;
415 }
416 //const math::XYZTLorentzVector GsfMoment = itPfGsfPoint->momentum();
417 c_gsf++;
418 }
419 }
420 math::XYZTLorentzVector pin = PfGsfPoint[0].momentum();
421 math::XYZTLorentzVector pout = PfGsfPoint[IndexPout].momentum();
422
423 /// get tracks from converted brems
424 if(useConvBremPFRecTracks_) {
425 const std::vector<reco::PFRecTrackRef>& temp_convBremPFRecTracks(refgsf->convBremPFRecTrackRef());
426 if(temp_convBremPFRecTracks.size() > 0) {
427 for(unsigned int iconv = 0; iconv <temp_convBremPFRecTracks.size(); iconv++) {
428 convBremPFRecTracks.push_back(temp_convBremPFRecTracks[iconv]);
429 }
430 }
431 }
432
433 // get the kf track that seeded the gsf
434 if(refgsf->kfPFRecTrackRef().isNonnull())
435 primaryKF_GSF.push_back(refgsf->kfPFRecTrackRef());
436
437
438 gsfEl = new reco::PFBlockElementGsfTrack(refgsf, pin, pout);
439
440 elements_.push_back( gsfEl);
441
442 std::vector<reco::PFBrem> pfbrem = refgsf->PFRecBrem();
443
444 for (unsigned i2=0;i2<pfbrem.size(); i2++) {
445 const double DP = pfbrem[i2].DeltaP();
446 const double SigmaDP = pfbrem[i2].SigmaDeltaP();
447 const unsigned int TrajP = pfbrem[i2].indTrajPoint();
448 if(TrajP == 99) continue;
449
450 reco::PFBlockElement* bremEl;
451 bremEl = new reco::PFBlockElementBrem(refgsf,DP,SigmaDP,TrajP);
452 elements_.push_back(bremEl);
453
454 }
455 }
456 // set the vector to the right size so to allow random access
457 scpfcRefs_.resize(superClusters_.size());
458 }
459
460
461
462 /// -------------- conversions ---------------------
463
464 /// The tracks from conversions are filled into the elements collection
465
466 if(convh.isValid() ) {
467 reco::PFBlockElement* trkFromConversionElement;
468 for(unsigned i=0;i<convh->size(); i++) {
469 reco::PFConversionRef convRef(convh,i);
470
471 unsigned int trackSize=(convRef->pfTracks()).size();
472 if ( convRef->pfTracks().size() < 2) continue;
473 for(unsigned iTk=0;iTk<trackSize; iTk++) {
474
475 reco::PFRecTrackRef compPFTkRef = convRef->pfTracks()[iTk];
476 trkFromConversionElement = new reco::PFBlockElementTrack(convRef->pfTracks()[iTk]);
477 trkFromConversionElement->setConversionRef( convRef->originalConversion(), reco::PFBlockElement::T_FROM_GAMMACONV);
478
479 /// The primary KF tracks associated the GSF seed are not labeled secondary
480 /// This can be changed but PFElectronAlgo.cc needs to changed too. Contact Daniele
481 bool isPrimGSF = false;
482 for(unsigned ikfgsf =0; ikfgsf<primaryKF_GSF.size();ikfgsf++) {
483 if(compPFTkRef->trackRef() == primaryKF_GSF[ikfgsf]->trackRef()) {
484 isPrimGSF = true;
485 continue;
486 }
487 }
488 if(isPrimGSF) continue;
489
490 elements_.push_back( trkFromConversionElement );
491
492 trkFromConversionElement = new reco::PFBlockElementTrack(convRef->pfTracks()[iTk]);
493 trkFromConversionElement->setConversionRef( convRef->originalConversion(), reco::PFBlockElement::T_FROM_GAMMACONV);
494 elements_.push_back( trkFromConversionElement );
495
496 if (debug_){
497 std::cout << "PF Block Element from Conversion electron " <<
498 (*trkFromConversionElement).trackRef().key() << std::endl;
499 std::cout << *trkFromConversionElement << std::endl;
500 }
501
502 }
503 }
504 }
505
506
507 /// -------------- V0 ---------------------
508
509 /// The tracks from V0 are filled into the elements collection
510
511 if(v0.isValid() ) {
512 reco::PFBlockElement* trkFromV0Element = 0;
513 for(unsigned i=0;i<v0->size(); i++) {
514 reco::PFV0Ref v0Ref( v0, i );
515 unsigned int trackSize=(v0Ref->pfTracks()).size();
516 for(unsigned iTk=0;iTk<trackSize; iTk++) {
517
518 reco::PFRecTrackRef newPFRecTrackRef = (v0Ref->pfTracks())[iTk];
519 reco::TrackBaseRef newTrackBaseRef(newPFRecTrackRef->trackRef());
520 bool bNew = true;
521
522 /// One need to cross check if those tracks was not already filled
523 /// from the conversion collection
524 for(IE iel = elements_.begin(); iel != elements_.end(); iel++){
525 reco::TrackBaseRef elemTrackBaseRef((*iel)->trackRef());
526 if (newTrackBaseRef == elemTrackBaseRef){
527 trkFromV0Element = *iel;
528 bNew = false;
529 continue;
530 }
531 }
532
533 /// The primary KF tracks associated the GSF seed are not labeled secondary
534 /// This can be changed but PFElectronAlgo.cc needs to changed too. Contact Daniele
535 bool isPrimGSF = false;
536 for(unsigned ikfgsf =0; ikfgsf<primaryKF_GSF.size();ikfgsf++) {
537 reco::TrackBaseRef elemTrackBaseRef(primaryKF_GSF[ikfgsf]->trackRef());
538 if(newTrackBaseRef == elemTrackBaseRef){
539 isPrimGSF = true;
540 continue;
541 }
542 }
543 if(isPrimGSF) continue;
544
545 /// This is a new track not yet included into the elements collection
546 if (bNew) {
547 trkFromV0Element = new reco::PFBlockElementTrack(v0Ref->pfTracks()[iTk]);
548 elements_.push_back( trkFromV0Element );
549 }
550
551 trkFromV0Element->setV0Ref( v0Ref->originalV0(),
552 reco::PFBlockElement::T_FROM_V0 );
553
554 if (debug_){
555 std::cout << "PF Block Element from V0 track New = " << bNew
556 << (*trkFromV0Element).trackRef().key() << std::endl;
557 std::cout << *trkFromV0Element << std::endl;
558 }
559
560
561 }
562 }
563 }
564
565 /// -------------- Displaced Vertices ---------------------
566
567 /// The tracks from Displaced Vertices are filled into the elements collection
568
569 if(displacedh.isValid()) {
570 reco::PFBlockElement* trkFromDisplacedVertexElement = 0;
571 for(unsigned i=0;i<displacedh->size(); i++) {
572
573 const reco::PFDisplacedTrackerVertexRef dispacedVertexRef( displacedh, i );
574
575 // std::cout << "In PFBlockAlgo" << std::endl;
576
577 // dispacedVertexRef->displacedVertexRef()->Dump();
578 //bool bIncludeVertices = true;
579
580
581 bool bIncludeVertices = false;
582 bool bNucl = dispacedVertexRef->displacedVertexRef()->isNucl();
583 bool bNucl_Loose = dispacedVertexRef->displacedVertexRef()->isNucl_Loose();
584 bool bNucl_Kink = dispacedVertexRef->displacedVertexRef()->isNucl_Kink();
585
586 if (nuclearInteractionsPurity_ >= 1) bIncludeVertices = bNucl;
587 if (nuclearInteractionsPurity_ >= 2) bIncludeVertices = bIncludeVertices || bNucl_Loose;
588 if (nuclearInteractionsPurity_ >= 3) bIncludeVertices = bIncludeVertices || bNucl_Kink;
589
590 if (bIncludeVertices){
591
592 unsigned int trackSize= dispacedVertexRef->pfRecTracks().size();
593 if (debug_){
594 std::cout << "" << std::endl;
595 std::cout << "Displaced Vertex " << i << std::endl;
596 dispacedVertexRef->displacedVertexRef()->Dump();
597 }
598 for(unsigned iTk=0;iTk < trackSize; iTk++) {
599
600 reco::PFRecTrackRef newPFRecTrackRef = dispacedVertexRef->pfRecTracks()[iTk];
601 reco::TrackBaseRef newTrackBaseRef(newPFRecTrackRef->trackRef());
602 bool bNew = true;
603 reco::PFBlockElement::TrackType blockType;
604
605 /// One need to cross check if those tracks was not already filled
606 /// from the conversion or V0 collections
607 for(IE iel = elements_.begin(); iel != elements_.end(); iel++){
608 reco::TrackBaseRef elemTrackBaseRef((*iel)->trackRef());
609 if (newTrackBaseRef == elemTrackBaseRef){
610 trkFromDisplacedVertexElement = *iel;
611 bNew = false;
612 continue;
613 }
614 }
615
616 /// The primary KF tracks associated the GSF seed are not labeled secondary
617 /// This can be changed but PFElectronAlgo.cc needs to changed too. Contact Daniele
618 bool isPrimGSF = false;
619 for(unsigned ikfgsf =0; ikfgsf<primaryKF_GSF.size();ikfgsf++) {
620 reco::TrackBaseRef elemTrackBaseRef(primaryKF_GSF[ikfgsf]->trackRef());
621 if(newTrackBaseRef == elemTrackBaseRef){
622 isPrimGSF = true;
623 continue;
624 }
625 }
626 if(isPrimGSF) continue;
627
628 /// This is a new track not yet included into the elements collection
629 if (bNew) {
630 trkFromDisplacedVertexElement = new reco::PFBlockElementTrack(newPFRecTrackRef);
631 elements_.push_back( trkFromDisplacedVertexElement );
632 }
633
634 if (dispacedVertexRef->isIncomingTrack(newPFRecTrackRef))
635 blockType = reco::PFBlockElement::T_TO_DISP;
636 else if (dispacedVertexRef->isOutgoingTrack(newPFRecTrackRef))
637 blockType = reco::PFBlockElement::T_FROM_DISP;
638 else
639 blockType = reco::PFBlockElement::DEFAULT;
640
641 /// Fill the displaced vertex ref
642 trkFromDisplacedVertexElement->setDisplacedVertexRef( dispacedVertexRef, blockType );
643
644
645 if (debug_){
646 std::cout << "PF Block Element from DisplacedTrackingVertex track New = " << bNew
647 << (*trkFromDisplacedVertexElement).trackRef().key() << std::endl;
648 std::cout << *trkFromDisplacedVertexElement << std::endl;
649 }
650
651
652 }
653 }
654 }
655
656 if (debug_) std::cout << "" << std::endl;
657
658 }
659
660 /// -------------- Tracks ---------------------
661
662 /// Mask the tracks in trackh collection already included from Conversions
663 /// V0 and Displaced Vertices. Take care that all those collections come
664 /// from the same "generalTracks" collection.
665
666 if(trackh.isValid() ) {
667
668 if (debug_) std::cout << "Tracks already in from Displaced Vertices " << std::endl;
669
670 Mask trackMaskVertex;
671
672 for(unsigned i=0;i<trackh->size(); i++) {
673 reco::PFRecTrackRef pfRefTrack( trackh,i );
674 reco::TrackRef trackRef = pfRefTrack->trackRef();
675
676 bool bMask = true;
677 for(IE iel = elements_.begin(); iel != elements_.end(); iel++){
678 reco::TrackRef elemTrackRef = (*iel)->trackRef();
679 if( trackRef == elemTrackRef ) {
680 if (debug_) std::cout << " " << trackRef.key();
681 bMask = false; continue;
682 }
683 }
684
685 trackMaskVertex.push_back(bMask);
686
687
688
689 }
690
691 if (debug_) std::cout << "" << std::endl;
692
693 if (debug_) std::cout << "Additionnal tracks from main collection " << std::endl;
694
695 for(unsigned i=0;i<trackh->size(); i++) {
696
697
698 // this track has been disabled
699 if( (!trackMask.empty() && !trackMask[i])
700 ||
701 (!trackMaskVertex.empty() && !trackMaskVertex[i]) ) continue;
702
703
704
705 reco::PFRecTrackRef ref( trackh,i );
706
707 if (debug_) std::cout << " " << ref->trackRef().key();
708
709 // Get the eventual muon associated to this track
710 int muId_ = muAssocToTrack( ref->trackRef(), muonh );
711 bool thisIsAPotentialMuon = false;
712 if( muId_ != -1 ) {
713 reco::MuonRef muonref( muonh, muId_ );
714 thisIsAPotentialMuon =
715 PFMuonAlgo::isLooseMuon(muonref) ||
716 PFMuonAlgo::isMuon(muonref);
717 }
718 // Reject bad tracks (except if identified as muon
719 if( !thisIsAPotentialMuon && !goodPtResolution( ref->trackRef() ) ) continue;
720
721 if (thisIsAPotentialMuon && debug_) std::cout << "Potential Muon P " << ref->trackRef()->p()
722 << " pt " << ref->trackRef()->p() << std::endl;
723
724 reco::PFBlockElement* primaryElement = new reco::PFBlockElementTrack( ref );
725
726 if( muId_ != -1 ) {
727 // if a muon has been found
728 reco::MuonRef muonref( muonh, muId_ );
729 primaryElement->setMuonRef( muonref );
730 }
731
732 // set track type T_FROM_GAMMA for pfrectracks associated to conv brems
733 if(useConvBremPFRecTracks_) {
734 if(convBremPFRecTracks.size() > 0.) {
735 for(unsigned int iconv = 0; iconv < convBremPFRecTracks.size(); iconv++) {
736 if((*ref).trackRef() == (*convBremPFRecTracks[iconv]).trackRef()) {
737 bool value = true;
738 primaryElement->setTrackType(reco::PFBlockElement::T_FROM_GAMMACONV, value);
739 }
740 }
741 }
742 }
743 elements_.push_back( primaryElement );
744 }
745
746 if (debug_) std::cout << " " << std::endl;
747
748 }
749
750
751 // -------------- GSF tracks and brems for Conversion Recovery ----------
752
753 if(convbremgsftrackh.isValid() ) {
754
755
756 const reco::GsfPFRecTrackCollection ConvPFGsfProd = *(convbremgsftrackh.product());
757 for(unsigned i=0;i<convbremgsftrackh->size(); i++) {
758
759 reco::GsfPFRecTrackRef refgsf(convbremgsftrackh,i );
760
761 if((refgsf).isNull()) continue;
762
763 reco::PFBlockElement* gsfEl;
764
765 const std::vector<reco::PFTrajectoryPoint>
766 PfGsfPoint = ConvPFGsfProd[i].trajectoryPoints();
767
768 unsigned int c_gsf=0;
769 bool PassTracker = false;
770 bool GetPout = false;
771 unsigned int IndexPout = -1;
772
773 typedef std::vector<reco::PFTrajectoryPoint>::const_iterator IP;
774 for(IP itPfGsfPoint = PfGsfPoint.begin();
775 itPfGsfPoint!= PfGsfPoint.end();itPfGsfPoint++) {
776
777 if (itPfGsfPoint->isValid()){
778 int layGsfP = itPfGsfPoint->layer();
779 if (layGsfP == -1) PassTracker = true;
780 if (PassTracker && layGsfP > 0 && GetPout == false) {
781 IndexPout = c_gsf-1;
782 GetPout = true;
783 }
784 //const math::XYZTLorentzVector GsfMoment = itPfGsfPoint->momentum();
785 c_gsf++;
786 }
787 }
788 math::XYZTLorentzVector pin = PfGsfPoint[0].momentum();
789 math::XYZTLorentzVector pout = PfGsfPoint[IndexPout].momentum();
790
791
792
793 gsfEl = new reco::PFBlockElementGsfTrack(refgsf, pin, pout);
794
795 bool valuegsf = true;
796 // IMPORTANT SET T_FROM_GAMMACONV trackType() FOR CONVERSIONS
797 gsfEl->setTrackType(reco::PFBlockElement::T_FROM_GAMMACONV, valuegsf);
798
799
800
801 elements_.push_back( gsfEl);
802 std::vector<reco::PFBrem> pfbrem = refgsf->PFRecBrem();
803
804 for (unsigned i2=0;i2<pfbrem.size(); i2++) {
805 const double DP = pfbrem[i2].DeltaP();
806 const double SigmaDP = pfbrem[i2].SigmaDeltaP();
807 const unsigned int TrajP = pfbrem[i2].indTrajPoint();
808 if(TrajP == 99) continue;
809
810 reco::PFBlockElement* bremEl;
811 bremEl = new reco::PFBlockElementBrem(refgsf,DP,SigmaDP,TrajP);
812 elements_.push_back(bremEl);
813
814 }
815 }
816 }
817
818
819 // -------------- ECAL clusters ---------------------
820
821
822 if(ecalh.isValid() ) {
823 pfcSCVec_.resize(ecalh->size(),-1);
824 for(unsigned i=0;i<ecalh->size(); i++) {
825
826 // this ecal cluster has been disabled
827 if( !ecalMask.empty() &&
828 !ecalMask[i] ) continue;
829
830 reco::PFClusterRef ref( ecalh,i );
831 reco::PFBlockElement* te
832 = new reco::PFBlockElementCluster( ref,
833 reco::PFBlockElement::ECAL);
834 elements_.push_back( te );
835 // Now mapping with Superclusters
836 int scindex= ClusterClusterMapping::checkOverlap(*ref,superClusters_);
837
838 if(scindex>=0) {
839 pfcSCVec_[ref.key()]=scindex;
840 scpfcRefs_[scindex].push_back(ref);
841 }
842 }
843 }
844
845 // -------------- HCAL clusters ---------------------
846
847 if(hcalh.isValid() ) {
848
849 for(unsigned i=0;i<hcalh->size(); i++) {
850
851 // this hcal cluster has been disabled
852 if( !hcalMask.empty() &&
853 !hcalMask[i] ) continue;
854
855 reco::PFClusterRef ref( hcalh,i );
856 reco::PFBlockElement* th
857 = new reco::PFBlockElementCluster( ref,
858 reco::PFBlockElement::HCAL );
859 elements_.push_back( th );
860 }
861 }
862
863
864 // -------------- HFEM clusters ---------------------
865
866 if(hfemh.isValid() ) {
867
868 for(unsigned i=0;i<hfemh->size(); i++) {
869
870 // this hfem cluster has been disabled
871 if( !hfemMask.empty() &&
872 !hfemMask[i] ) continue;
873
874 reco::PFClusterRef ref( hfemh,i );
875 reco::PFBlockElement* th
876 = new reco::PFBlockElementCluster( ref,
877 reco::PFBlockElement::HFEM );
878 elements_.push_back( th );
879 }
880 }
881
882
883 // -------------- HFHAD clusters ---------------------
884
885 if(hfhadh.isValid() ) {
886
887 for(unsigned i=0;i<hfhadh->size(); i++) {
888
889 // this hfhad cluster has been disabled
890 if( !hfhadMask.empty() &&
891 !hfhadMask[i] ) continue;
892
893 reco::PFClusterRef ref( hfhadh,i );
894 reco::PFBlockElement* th
895 = new reco::PFBlockElementCluster( ref,
896 reco::PFBlockElement::HFHAD );
897 elements_.push_back( th );
898 }
899 }
900
901
902
903
904 // -------------- PS clusters ---------------------
905
906 if(psh.isValid() ) {
907 for(unsigned i=0;i<psh->size(); i++) {
908
909 // this ps cluster has been disabled
910 if( !psMask.empty() &&
911 !psMask[i] ) continue;
912 reco::PFBlockElement::Type type = reco::PFBlockElement::NONE;
913 reco::PFClusterRef ref( psh,i );
914 // two types of elements: PS1 (V) and PS2 (H)
915 // depending on layer: PS1 or PS2
916 switch(ref->layer()){
917 case PFLayer::PS1:
918 type = reco::PFBlockElement::PS1;
919 break;
920 case PFLayer::PS2:
921 type = reco::PFBlockElement::PS2;
922 break;
923 default:
924 break;
925 }
926 reco::PFBlockElement* tp
927 = new reco::PFBlockElementCluster( ref,
928 type );
929 elements_.push_back( tp );
930
931 }
932 }
933 }
934
935
936
937 #endif
938
939