ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitProd/TreeFiller/src/FillerElectrons.cc
Revision: 1.65
Committed: Sat May 5 16:49:59 2012 UTC (13 years ago) by paus
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_027a, Mit_027
Changes since 1.64: +90 -26 lines
Log Message:
Version 027 - complete version for ICHEP 2012.

File Contents

# User Rev Content
1 paus 1.65 // $Id: $
2    
3 pharris 1.59 #include "FWCore/MessageLogger/interface/MessageLogger.h"
4 loizides 1.1 #include "MitProd/TreeFiller/interface/FillerElectrons.h"
5     #include "DataFormats/TrackReco/interface/Track.h"
6     #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
7     #include "DataFormats/TrackReco/interface/TrackFwd.h"
8 sixie 1.11 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
9     #include "DataFormats/EgammaCandidates/interface/GsfElectronFwd.h"
10     #include "DataFormats/EgammaReco/interface/ClusterShape.h"
11     #include "DataFormats/EgammaReco/interface/BasicClusterShapeAssociation.h"
12 bendavid 1.18 #include "DataFormats/Common/interface/RefToPtr.h"
13 loizides 1.35 #include "DataFormats/Common/interface/ValueMap.h"
14     #include "AnalysisDataFormats/Egamma/interface/ElectronID.h"
15     #include "AnalysisDataFormats/Egamma/interface/ElectronIDAssociation.h"
16 sixie 1.11 #include "RecoEgamma/EgammaIsolationAlgos/interface/ElectronTkIsolation.h"
17     #include "RecoEgamma/EgammaIsolationAlgos/interface/EgammaEcalIsolation.h"
18 sixie 1.14 #include "RecoEgamma/EgammaIsolationAlgos/interface/EgammaTowerIsolation.h"
19 loizides 1.35 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterLazyTools.h"
20 bendavid 1.43 #include "DataFormats/VertexReco/interface/VertexFwd.h"
21     #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
22     #include "TrackingTools/TransientTrack/plugins/TransientTrackBuilderESProducer.h"
23 bendavid 1.44 #include "RecoVertex/GaussianSumVertexFit/interface/GsfVertexTrackCompatibilityEstimator.h"
24 bendavid 1.43 #include "TrackingTools/IPTools/interface/IPTools.h"
25     #include "RecoEgamma/EgammaTools/interface/ConversionFinder.h"
26     #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
27     #include "MagneticField/Engine/interface/MagneticField.h"
28 loizides 1.35 #include "MitAna/DataTree/interface/ElectronCol.h"
29     #include "MitAna/DataTree/interface/Names.h"
30     #include "MitAna/DataTree/interface/Track.h"
31 loizides 1.33 #include "MitEdm/DataFormats/interface/RefToBaseToPtr.h"
32 loizides 1.35 #include "MitProd/ObjectService/interface/ObjectService.h"
33 bendavid 1.45 #include "MitEdm/DataFormats/interface/DecayPart.h"
34     #include "MitEdm/ConversionRejection/interface/ConversionMatcher.h"
35 bendavid 1.51 #include "RecoVertex/VertexTools/interface/LinearizedTrackStateFactory.h"
36     #include "RecoVertex/VertexTools/interface/VertexTrackFactory.h"
37     #include "RecoVertex/VertexPrimitives/interface/VertexTrack.h"
38     #include "RecoVertex/VertexPrimitives/interface/CachingVertex.h"
39     #include "RecoVertex/KalmanVertexFit/interface/KalmanVertexUpdator.h"
40 ksung 1.57 #include "MitEdm/Tools/interface/VertexReProducer.h"
41 paus 1.65 #include "RecoEgamma/EgammaTools/interface/ggPFClusters.h"
42     #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
43     #include "Geometry/CaloTopology/interface/CaloSubdetectorTopology.h"
44     #include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
45     #include "Geometry/Records/interface/CaloGeometryRecord.h"
46     #include "FWCore/Framework/interface/ESHandle.h"
47     #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
48     #include "DataFormats/Math/interface/deltaR.h"
49 sixie 1.13
50 loizides 1.1 using namespace std;
51     using namespace edm;
52     using namespace mithep;
53    
54 loizides 1.5 //--------------------------------------------------------------------------------------------------
55 loizides 1.28 FillerElectrons::FillerElectrons(const edm::ParameterSet &cfg, const char *name, bool active) :
56 paus 1.65 BaseFiller (cfg,name,active),
57     edmName_ (Conf().getUntrackedParameter<string>("edmName","pixelMatchGsfElectrons")),
58     expectedHitsName_ (Conf().getUntrackedParameter<string>("expectedHitsName","")),
59     mitName_ (Conf().getUntrackedParameter<string>("mitName",Names::gkElectronBrn)),
60     gsfTrackMapName_ (Conf().getUntrackedParameter<string>("gsfTrackMapName","")),
61     trackerTrackMapName_ (Conf().getUntrackedParameter<string>("trackerTrackMapName","")),
62 sixie 1.12 barrelSuperClusterMapName_(Conf().getUntrackedParameter<string>("barrelSuperClusterMapName","")),
63     endcapSuperClusterMapName_(Conf().getUntrackedParameter<string>("endcapSuperClusterMapName","")),
64 paus 1.65 checkClusterActive_ (Conf().getUntrackedParameter<bool>("requireClusterAndGsfMap",true)),
65     pfSuperClusterMapName_ (Conf().getUntrackedParameter<string>("pfSuperClusterMapName","")),
66     pfClusterMapName_ (Conf().getUntrackedParameter<string>("pfClusterMapName","")),
67     electronMapName_ (Conf().getUntrackedParameter<string>("electronMapName","")),
68     eIDCutBasedTightName_ (Conf().getUntrackedParameter<string>("eIDCutBasedTightName","eidTight")),
69     eIDCutBasedLooseName_ (Conf().getUntrackedParameter<string>("eIDCutBasedLooseName","eidLoose")),
70     eIDLikelihoodName_ (Conf().getUntrackedParameter<string>("eIDLikelihoodName","")),
71     pvEdmName_ (Conf().getUntrackedParameter<string>("pvEdmName","offlinePrimaryVertices")),
72     pvBSEdmName_ (Conf().getUntrackedParameter<string>("pvBSEdmName","offlinePrimaryVerticesWithBS")),
73     recomputeConversionInfo_ (Conf().getUntrackedParameter<bool>("recomputeConversionInfo",false)),
74     fitUnbiasedVertex_ (Conf().getUntrackedParameter<bool>("fitUnbiasedVertex",true)),
75     electronMap_ (new mithep::ElectronMap),
76     electrons_ (new mithep::ElectronArr(16)),
77     gsfTrackMap_ (0),
78     trackerTrackMap_ (0),
79     barrelSuperClusterMap_ (0),
80     endcapSuperClusterMap_ (0),
81     pfSuperClusterMap_ (0),
82     pfClusterMap_ (0)
83 loizides 1.1 {
84     // Constructor.
85     }
86    
87 loizides 1.5 //--------------------------------------------------------------------------------------------------
88 loizides 1.1 FillerElectrons::~FillerElectrons()
89     {
90     // Destructor.
91 loizides 1.6
92     delete electrons_;
93 pharris 1.62 delete electronMap_;
94 loizides 1.1 }
95    
96 loizides 1.5 //--------------------------------------------------------------------------------------------------
97 bendavid 1.42 void FillerElectrons::BookDataBlock(TreeWriter &tws)
98 loizides 1.1 {
99 loizides 1.10 // Add electron branch to our tree and get our maps.
100 loizides 1.1
101 loizides 1.29 tws.AddBranch(mitName_,&electrons_);
102     OS()->add<mithep::ElectronArr>(electrons_,mitName_);
103     if (!gsfTrackMapName_.empty()) {
104     gsfTrackMap_ = OS()->get<TrackMap>(gsfTrackMapName_);
105     if (gsfTrackMap_)
106     AddBranchDep(mitName_,gsfTrackMap_->GetBrName());
107     }
108     if (!trackerTrackMapName_.empty()) {
109     trackerTrackMap_ = OS()->get<TrackMap>(trackerTrackMapName_);
110     if (trackerTrackMap_)
111     AddBranchDep(mitName_,trackerTrackMap_->GetBrName());
112     }
113     if (!barrelSuperClusterMapName_.empty()) {
114     barrelSuperClusterMap_ = OS()->get<SuperClusterMap>(barrelSuperClusterMapName_);
115     if (barrelSuperClusterMap_)
116     AddBranchDep(mitName_,barrelSuperClusterMap_->GetBrName());
117     }
118     if (!endcapSuperClusterMapName_.empty()) {
119     endcapSuperClusterMap_ = OS()->get<SuperClusterMap>(endcapSuperClusterMapName_);
120     if (endcapSuperClusterMap_)
121     AddBranchDep(mitName_,endcapSuperClusterMap_->GetBrName());
122     }
123 bendavid 1.36 if (!pfSuperClusterMapName_.empty()) {
124     pfSuperClusterMap_ = OS()->get<SuperClusterMap>(pfSuperClusterMapName_);
125     if (pfSuperClusterMap_)
126     AddBranchDep(mitName_,pfSuperClusterMap_->GetBrName());
127     }
128 paus 1.65 if (!pfClusterMapName_.empty()) {
129     pfClusterMap_ = OS()->get<BasicClusterMap>(pfClusterMapName_);
130     if (pfClusterMap_)
131     AddBranchDep(mitName_,pfClusterMap_->GetBrName());
132     }
133 pharris 1.62 if (!electronMapName_.empty()) {
134     electronMap_->SetBrName(mitName_);
135     OS()->add<ElectronMap>(electronMap_,electronMapName_);
136     }
137 loizides 1.1 }
138    
139 loizides 1.5 //--------------------------------------------------------------------------------------------------
140 loizides 1.29 void FillerElectrons::FillDataBlock(const edm::Event &event, const edm::EventSetup &setup)
141 loizides 1.1 {
142     // Fill electrons from edm collection into our collection.
143    
144 pharris 1.62 electrons_ ->Delete();
145     electronMap_->Reset();
146    
147 sixie 1.11 Handle<reco::GsfElectronCollection> hElectronProduct;
148 loizides 1.7 GetProduct(edmName_, hElectronProduct, event);
149 loizides 1.1
150 loizides 1.29 // handles to get the electron ID information
151 peveraer 1.47 Handle<edm::ValueMap<float> > eidLooseMap;
152     GetProduct(eIDCutBasedLooseName_, eidLooseMap, event);
153     Handle<edm::ValueMap<float> > eidTightMap;
154 peveraer 1.48 GetProduct(eIDCutBasedTightName_, eidTightMap, event);
155 bendavid 1.49 edm::Handle<edm::ValueMap<float> > eidLikelihoodMap;
156     if (!eIDLikelihoodName_.empty()) {
157     GetProduct(eIDLikelihoodName_, eidLikelihoodMap, event);
158     }
159    
160 bendavid 1.43 edm::Handle<reco::VertexCollection> hVertex;
161     event.getByLabel(pvEdmName_, hVertex);
162     const reco::VertexCollection *pvCol = hVertex.product();
163     edm::Handle<reco::VertexCollection> hVertexBS;
164     event.getByLabel(pvBSEdmName_, hVertexBS);
165     const reco::VertexCollection *pvBSCol = hVertexBS.product();
166     edm::Handle<reco::TrackCollection> hGeneralTracks;
167     event.getByLabel("generalTracks", hGeneralTracks);
168     //const reco::VertexCollection *trackCol = hGeneralTracks.product();
169 bendavid 1.45
170 bendavid 1.54 edm::Handle<reco::GsfTrackCollection> hGsfTracks;
171     event.getByLabel("electronGsfTracks", hGsfTracks);
172    
173 bendavid 1.45 edm::Handle<std::vector<mitedm::DecayPart> > hConversions;
174     event.getByLabel("mvfConversionRemoval", hConversions);
175    
176     mitedm::ConversionMatcher convMatcher;
177 bendavid 1.50
178 bendavid 1.43 edm::ESHandle<TransientTrackBuilder> hTransientTrackBuilder;
179     setup.get<TransientTrackRecord>().get("TransientTrackBuilder",hTransientTrackBuilder);
180     const TransientTrackBuilder *transientTrackBuilder = hTransientTrackBuilder.product();
181 bendavid 1.19
182 bendavid 1.44 GsfVertexTrackCompatibilityEstimator gsfEstimator;
183    
184 bendavid 1.51 LinearizedTrackStateFactory lTrackFactory;
185     VertexTrackFactory<5> vTrackFactory;
186     KalmanVertexUpdator<5> updator;
187    
188 paus 1.65 //pf photon stuff
189     edm::Handle< EcalRecHitCollection > pEBRecHits;
190     event.getByLabel("reducedEcalRecHitsEB", pEBRecHits );
191     edm::Handle< EcalRecHitCollection > pEERecHits;
192     event.getByLabel( "reducedEcalRecHitsEE", pEERecHits );
193    
194     edm::ESHandle<CaloGeometry> pGeometry;
195     setup.get<CaloGeometryRecord>().get(pGeometry);
196    
197     const CaloSubdetectorGeometry *geometryEB = pGeometry->getSubdetectorGeometry(DetId::Ecal, EcalBarrel);
198     const CaloSubdetectorGeometry *geometryEE = pGeometry->getSubdetectorGeometry(DetId::Ecal, EcalEndcap);
199    
200     ggPFClusters pfclusters(pEBRecHits, pEERecHits, geometryEB, geometryEE);
201    
202 pharris 1.61 //Get Magnetic Field from event setup, taking value at (0,0,0)
203 bendavid 1.54 edm::ESHandle<MagneticField> magneticField;
204     setup.get<IdealMagneticFieldRecord>().get(magneticField);
205     const double bfield = magneticField->inTesla(GlobalPoint(0.,0.,0.)).z();
206 bendavid 1.43
207 sixie 1.11 const reco::GsfElectronCollection inElectrons = *(hElectronProduct.product());
208 loizides 1.29 // loop over electrons
209 sixie 1.11 for (reco::GsfElectronCollection::const_iterator iM = inElectrons.begin();
210 loizides 1.1 iM != inElectrons.end(); ++iM) {
211 sixie 1.13
212 loizides 1.29 // the index and Ref are needed for the eID association Map
213 sixie 1.13 unsigned int iElectron = iM - inElectrons.begin();
214     reco::GsfElectronRef eRef(hElectronProduct, iElectron);
215    
216 loizides 1.8 mithep::Electron *outElectron = electrons_->AddNew();
217 bendavid 1.31
218     outElectron->SetPtEtaPhi(iM->pt(),iM->eta(),iM->phi());
219 sixie 1.12
220 bendavid 1.40 outElectron->SetCharge(iM->charge());
221     outElectron->SetScPixCharge(iM->scPixCharge());
222 pharris 1.62 outElectron->SetESuperClusterOverP(iM->eSuperClusterOverP());
223 loizides 1.30 outElectron->SetESeedClusterOverPout(iM->eSeedClusterOverPout());
224 sixie 1.63 outElectron->SetEEleClusterOverPout(iM->eEleClusterOverPout());
225 loizides 1.30 outElectron->SetPIn(iM->trackMomentumAtVtx().R());
226 loizides 1.29 outElectron->SetPOut(iM->trackMomentumOut().R());
227 loizides 1.30 outElectron->SetDeltaEtaSuperClusterTrackAtVtx(iM->deltaEtaSuperClusterTrackAtVtx());
228     outElectron->SetDeltaEtaSeedClusterTrackAtCalo(iM->deltaEtaSeedClusterTrackAtCalo());
229     outElectron->SetDeltaPhiSuperClusterTrackAtVtx(iM->deltaPhiSuperClusterTrackAtVtx());
230     outElectron->SetDeltaPhiSeedClusterTrackAtCalo(iM->deltaPhiSeedClusterTrackAtCalo());
231     outElectron->SetIsEnergyScaleCorrected(iM->isEnergyScaleCorrected());
232 bendavid 1.55 //outElectron->SetIsMomentumCorrected(iM->isMomentumCorrected());
233 bendavid 1.36 outElectron->SetNumberOfClusters(iM->basicClustersSize());
234 loizides 1.30 outElectron->SetClassification(iM->classification());
235 bendavid 1.37 outElectron->SetFBrem(iM->fbrem());
236    
237 loizides 1.38 // pflow electron stuff
238 bendavid 1.41 outElectron->SetIsEcalDriven(iM->ecalDrivenSeed());
239     outElectron->SetIsTrackerDriven(iM->trackerDrivenSeed());
240 bendavid 1.37 outElectron->SetMva(iM->mva());
241    
242 loizides 1.29 // shower shape variables
243 bendavid 1.37 outElectron->SetE15(iM->e1x5());
244     outElectron->SetE25Max(iM->e2x5Max());
245     outElectron->SetE55(iM->e5x5());
246     outElectron->SetCovEtaEta(iM->sigmaEtaEta());
247     outElectron->SetCoviEtaiEta(iM->sigmaIetaIeta());
248     outElectron->SetHadronicOverEm(iM->hcalOverEcal());
249     outElectron->SetHcalDepth1OverEcal(iM->hcalDepth1OverEcal());
250     outElectron->SetHcalDepth2OverEcal(iM->hcalDepth2OverEcal());
251 bendavid 1.64 outElectron->SetHadOverEmTow(iM->hcalOverEcalBc());
252    
253    
254 loizides 1.38 // fill isolation variables for both cone sizes
255 bendavid 1.37 outElectron->SetEcalRecHitIsoDr04(iM->dr04EcalRecHitSumEt());
256     outElectron->SetHcalDepth1TowerSumEtDr04(iM->dr04HcalDepth1TowerSumEt());
257     outElectron->SetHcalDepth2TowerSumEtDr04(iM->dr04HcalDepth2TowerSumEt());
258     outElectron->SetTrackIsolationDr04(iM->dr04TkSumPt());
259 bendavid 1.64 outElectron->SetHCalIsoTowDr04( iM->dr04HcalTowerSumEtBc() );
260 bendavid 1.37 outElectron->SetEcalRecHitIsoDr03(iM->dr03EcalRecHitSumEt());
261     outElectron->SetHcalTowerSumEtDr03(iM->dr03HcalTowerSumEt());
262     outElectron->SetHcalDepth1TowerSumEtDr03(iM->dr03HcalDepth1TowerSumEt());
263     outElectron->SetHcalDepth2TowerSumEtDr03(iM->dr03HcalDepth2TowerSumEt());
264     outElectron->SetTrackIsolationDr03(iM->dr03TkSumPt());
265 bendavid 1.64 outElectron->SetHCalIsoTowDr03( iM->dr03HcalTowerSumEtBc() );
266    
267    
268 bendavid 1.56 //pflow isolation
269     outElectron->SetPFChargedHadronIso(iM->pfIsolationVariables().chargedHadronIso);
270     outElectron->SetPFChargedHadronIso(iM->pfIsolationVariables().neutralHadronIso);
271     outElectron->SetPFChargedHadronIso(iM->pfIsolationVariables().photonIso);
272    
273 loizides 1.38 // fiducial flags
274 bendavid 1.37 outElectron->SetIsEB(iM->isEB());
275     outElectron->SetIsEE(iM->isEE());
276     outElectron->SetIsEBEEGap(iM->isEBEEGap());
277     outElectron->SetIsEBEtaGap(iM->isEBEtaGap());
278     outElectron->SetIsEBPhiGap(iM->isEBPhiGap());
279     outElectron->SetIsEEDeeGap(iM->isEEDeeGap());
280     outElectron->SetIsEERingGap(iM->isEERingGap());
281 sixie 1.12
282 loizides 1.38 // gsf-tracker match quality
283 bendavid 1.37 outElectron->SetFracSharedHits(iM->shFracInnerHits());
284 sixie 1.15
285 loizides 1.29 // make proper links to Tracks and Super Clusters
286 bendavid 1.49 if (gsfTrackMap_ && iM->gsfTrack().isNonnull()) {
287 pharris 1.59 try{outElectron->SetGsfTrk(gsfTrackMap_->GetMit(refToPtr(iM->gsfTrack())));}
288     catch(...) {
289     if(checkClusterActive_) throw edm::Exception(edm::errors::Configuration, "FillerElectrons:FillDataBlock()\n")
290     << "Error! GSF track unmapped collection " << edmName_ << endl;
291     }
292 bendavid 1.49 }
293     // make links to ambigous gsf tracks
294     if (gsfTrackMap_) {
295     for (reco::GsfTrackRefVector::const_iterator agsfi = iM->ambiguousGsfTracksBegin(); agsfi != iM->ambiguousGsfTracksEnd(); ++agsfi) {
296 pharris 1.59 try{outElectron->AddAmbiguousGsfTrack(gsfTrackMap_->GetMit(refToPtr(*agsfi)));}
297     catch(...) {
298     if(checkClusterActive_) throw edm::Exception(edm::errors::Configuration, "FillerElectrons:FillDataBlock()\n")
299     << "Error! GSF track unmapped collection " << edmName_ << endl;
300     }
301 bendavid 1.49 }
302     }
303    
304     // make tracker track links,
305 bendavid 1.37 if (trackerTrackMap_ && iM->closestCtfTrackRef().isNonnull()) {
306 bendavid 1.36 outElectron->SetTrackerTrk(trackerTrackMap_->GetMit(refToPtr(iM->closestCtfTrackRef())));
307 bendavid 1.19 }
308 pharris 1.61
309 loizides 1.38 if (barrelSuperClusterMap_ && endcapSuperClusterMap_ &&
310     pfSuperClusterMap_ && iM->superCluster().isNonnull()) {
311 bendavid 1.37 if(barrelSuperClusterMap_->HasMit(iM->superCluster())) {
312 sixie 1.12 outElectron->SetSuperCluster(barrelSuperClusterMap_->GetMit(iM->superCluster()));
313 bendavid 1.36 }
314     else if (endcapSuperClusterMap_->HasMit(iM->superCluster())) {
315 sixie 1.12 outElectron->SetSuperCluster(endcapSuperClusterMap_->GetMit(iM->superCluster()));
316     }
317 bendavid 1.36 else if (pfSuperClusterMap_->HasMit(iM->superCluster())) {
318     outElectron->SetSuperCluster(pfSuperClusterMap_->GetMit(iM->superCluster()));
319     }
320 pharris 1.59 else if(checkClusterActive_) {
321     throw edm::Exception(edm::errors::Configuration, "FillerElectrons:FillDataBlock()\n")
322 bendavid 1.36 << "Error! SuperCluster reference in unmapped collection " << edmName_ << endl;
323 pharris 1.59 }
324 bendavid 1.36 }
325 paus 1.65
326     if (pfSuperClusterMap_ && iM->pflowSuperCluster().isNonnull()) {
327     outElectron->SetPFSuperCluster(pfSuperClusterMap_->GetMit(iM->pflowSuperCluster()));
328     }
329    
330     //find matching egamma supercluster first by ref, or by geometric matching if only pfsupercluster is linked
331     const reco::SuperCluster *egammasc = 0;
332     double mindr = 999.;
333     if (iM->superCluster().isNonnull() && iM->superCluster() != iM->pflowSuperCluster()) {
334     egammasc = iM->superCluster().get();
335     }
336     else {
337     for (SuperClusterMap::fwdMapType::const_iterator scit = barrelSuperClusterMap_->FwdMap().begin(); scit != barrelSuperClusterMap_->FwdMap().end(); ++scit) {
338     double dr = reco::deltaR(*iM->pflowSuperCluster(),*scit->first);
339     if (dr<0.1 && dr<mindr) {
340     egammasc = scit->first.get();
341     mindr = dr;
342     }
343     }
344     for (SuperClusterMap::fwdMapType::const_iterator scit = endcapSuperClusterMap_->FwdMap().begin(); scit != endcapSuperClusterMap_->FwdMap().end(); ++scit) {
345     double dr = reco::deltaR(*iM->pflowSuperCluster(),*scit->first);
346     if (dr<0.1 && dr<mindr) {
347     egammasc = scit->first.get();
348     mindr = dr;
349     }
350     }
351     }
352    
353     //tag overlapping energy of pflow clusters
354     if (pfClusterMap_ && iM->pflowSuperCluster().isNonnull() && egammasc) {
355     for (reco::CaloCluster_iterator pfcit = iM->pflowSuperCluster()->clustersBegin();
356     pfcit!=iM->pflowSuperCluster()->clustersEnd(); ++pfcit) {
357     float eoverlap = pfclusters.getPFSuperclusterOverlap( **pfcit, *egammasc );
358     const_cast<mithep::BasicCluster*>(pfClusterMap_->GetMit(*pfcit))->SetMatchedEnergy(eoverlap);
359     }
360     }
361 pharris 1.61
362 sixie 1.63 //compute NLayersWithMeasurement for associated ctf track
363     if (iM->closestCtfTrackRef().isNonnull()) {
364     outElectron->SetCTFTrkNLayersWithMeasurement(iM->closestCtfTrackRef()->hitPattern().trackerLayersWithMeasurement());
365     } else {
366     outElectron->SetCTFTrkNLayersWithMeasurement(-1);
367     }
368    
369 bendavid 1.43 //compute impact parameter with respect to PV
370     if (iM->gsfTrack().isNonnull()) {
371     const reco::TransientTrack &tt = transientTrackBuilder->build(iM->gsfTrack());
372 bendavid 1.51
373     reco::TransientTrack ttckf;
374    
375 bendavid 1.58 reco::Vertex thevtx = pvCol->at(0);
376     reco::Vertex thevtxbs = pvBSCol->at(0);
377 bendavid 1.52
378 bendavid 1.58 reco::Vertex thevtxub = pvCol->at(0);
379     reco::Vertex thevtxubbs = pvBSCol->at(0);
380 pharris 1.61
381 bendavid 1.51 //check if closest ctf track is included in PV and if so, remove it before computing impact parameters and uncertainties
382     if (iM->closestCtfTrackRef().isNonnull()) {
383     ttckf = transientTrackBuilder->build(iM->closestCtfTrackRef());
384 ksung 1.57
385     //check if closest ctf track is included in PV and if so, remove it from the collection of
386     //tracks associated with the PV and perform a refit before computing impact parameters and uncertainties
387     reco::TrackCollection newTkCollection;
388     bool foundMatch = false;
389     for(reco::Vertex::trackRef_iterator itk = thevtx.tracks_begin(); itk!=thevtx.tracks_end(); itk++) {
390     if(iM->closestCtfTrack().ctfTrack.isNonnull()) {
391     bool refMatching = (itk->get() == &*(iM->closestCtfTrack().ctfTrack));
392     float shFraction = iM->closestCtfTrack().shFracInnerHits;
393     if(refMatching && shFraction > 0.5) {
394 pharris 1.61 foundMatch = true;
395     continue;
396 ksung 1.57 }
397     }
398     newTkCollection.push_back(*itk->get());
399     }
400 bendavid 1.51
401 pharris 1.61 if(foundMatch && fitUnbiasedVertex_) {
402 ksung 1.57 edm::Handle<reco::BeamSpot> bs;
403     event.getByLabel(edm::InputTag("offlineBeamSpot"),bs);
404 pharris 1.61
405     VertexReProducer revertex(hVertex,event); //Needs to be fixed
406    
407     edm::Handle<reco::BeamSpot> pvbeamspot;
408 ksung 1.57 event.getByLabel(revertex.inputBeamSpot(),pvbeamspot);
409 pharris 1.61 vector<TransientVertex> pvs = revertex.makeVertices(newTkCollection,*pvbeamspot,setup);
410     if(pvs.size()>0) {
411 bendavid 1.58 thevtxub = pvs.front(); // take the first in the list
412     }
413 ksung 1.57 VertexReProducer revertexbs(hVertexBS,event);
414     edm::Handle<reco::BeamSpot> pvbsbeamspot;
415     event.getByLabel(revertexbs.inputBeamSpot(),pvbsbeamspot);
416     vector<TransientVertex> pvbss = revertexbs.makeVertices(newTkCollection,*pvbsbeamspot,setup);
417 pharris 1.61 if(pvbss.size()>0) {
418 bendavid 1.58 thevtxubbs = pvbss.front(); // take the first in the list
419     }
420 ksung 1.57 }
421 bendavid 1.51 }
422     //preserve sign of transverse impact parameter (cross-product definition from track, not lifetime-signing)
423     const double gsfsign = ( (-iM->gsfTrack()->dxy(thevtx.position())) >=0 ) ? 1. : -1.;
424     const double gsfsignbs = ( (-iM->gsfTrack()->dxy(thevtxbs.position())) >=0 ) ? 1. : -1.;
425     const std::pair<bool,Measurement1D> &d0pv = IPTools::absoluteTransverseImpactParameter(tt,thevtx);
426 bendavid 1.43 if (d0pv.first) {
427 bendavid 1.51 outElectron->SetD0PV(gsfsign*d0pv.second.value());
428 bendavid 1.43 outElectron->SetD0PVErr(d0pv.second.error());
429     }
430 bendavid 1.44 else {
431 bendavid 1.51 outElectron->SetD0PV(-999.0);
432 bendavid 1.44 }
433    
434 bendavid 1.51 const std::pair<bool,Measurement1D> &ip3dpv = IPTools::absoluteImpactParameter3D(tt,thevtx);
435 bendavid 1.43 if (ip3dpv.first) {
436 bendavid 1.51 outElectron->SetIp3dPV(gsfsign*ip3dpv.second.value());
437 bendavid 1.43 outElectron->SetIp3dPVErr(ip3dpv.second.error());
438     }
439 bendavid 1.44 else {
440 bendavid 1.51 outElectron->SetIp3dPV(-999.0);
441 bendavid 1.44 }
442 bendavid 1.43
443 bendavid 1.51 const std::pair<bool,Measurement1D> &d0pvbs = IPTools::absoluteTransverseImpactParameter(tt,thevtxbs);
444 bendavid 1.43 if (d0pvbs.first) {
445 bendavid 1.51 outElectron->SetD0PVBS(gsfsignbs*d0pvbs.second.value());
446 bendavid 1.43 outElectron->SetD0PVBSErr(d0pvbs.second.error());
447     }
448 bendavid 1.44 else {
449 bendavid 1.51 outElectron->SetD0PVBS(-999.0);
450 bendavid 1.44 }
451 bendavid 1.43
452 bendavid 1.51 const std::pair<bool,Measurement1D> &ip3dpvbs = IPTools::absoluteImpactParameter3D(tt,thevtxbs);
453 bendavid 1.43 if (ip3dpvbs.first) {
454 bendavid 1.51 outElectron->SetIp3dPVBS(gsfsignbs*ip3dpvbs.second.value());
455 bendavid 1.43 outElectron->SetIp3dPVBSErr(ip3dpvbs.second.error());
456     }
457 bendavid 1.44 else {
458 bendavid 1.51 outElectron->SetIp3dPVBS(-999.0);
459 bendavid 1.44 }
460    
461 bendavid 1.58 const std::pair<bool,Measurement1D> &d0pvub = IPTools::absoluteTransverseImpactParameter(tt,thevtxub);
462     if (d0pvub.first) {
463     outElectron->SetD0PVUB(gsfsign*d0pvub.second.value());
464     outElectron->SetD0PVUBErr(d0pvub.second.error());
465     }
466     else {
467     outElectron->SetD0PVUB(-999.0);
468     }
469    
470    
471     const std::pair<bool,Measurement1D> &ip3dpvub = IPTools::absoluteImpactParameter3D(tt,thevtxub);
472     if (ip3dpvub.first) {
473     outElectron->SetIp3dPVUB(gsfsign*ip3dpvub.second.value());
474     outElectron->SetIp3dPVUBErr(ip3dpvub.second.error());
475     }
476     else {
477     outElectron->SetIp3dPVUB(-999.0);
478     }
479    
480     const std::pair<bool,Measurement1D> &d0pvubbs = IPTools::absoluteTransverseImpactParameter(tt,thevtxubbs);
481     if (d0pvubbs.first) {
482     outElectron->SetD0PVUBBS(gsfsignbs*d0pvubbs.second.value());
483     outElectron->SetD0PVUBBSErr(d0pvubbs.second.error());
484     }
485     else {
486     outElectron->SetD0PVUBBS(-999.0);
487     }
488    
489     const std::pair<bool,Measurement1D> &ip3dpvubbs = IPTools::absoluteImpactParameter3D(tt,thevtxubbs);
490     if (ip3dpvubbs.first) {
491     outElectron->SetIp3dPVUBBS(gsfsignbs*ip3dpvubbs.second.value());
492     outElectron->SetIp3dPVUBBSErr(ip3dpvubbs.second.error());
493     }
494     else {
495     outElectron->SetIp3dPVUBBS(-999.0);
496     }
497    
498 bendavid 1.51 if (iM->closestCtfTrackRef().isNonnull()) {
499    
500     const double ckfsign = ( (-iM->closestCtfTrackRef()->dxy(thevtx.position())) >=0 ) ? 1. : -1.;
501     const double ckfsignbs = ( (-iM->closestCtfTrackRef()->dxy(thevtxbs.position())) >=0 ) ? 1. : -1.;
502    
503     const std::pair<bool,Measurement1D> &d0pvckf = IPTools::absoluteTransverseImpactParameter(ttckf,thevtx);
504     if (d0pvckf.first) {
505     outElectron->SetD0PVCkf(ckfsign*d0pvckf.second.value());
506     outElectron->SetD0PVCkfErr(d0pvckf.second.error());
507     }
508     else {
509     outElectron->SetD0PVCkf(-999.0);
510     }
511    
512    
513     const std::pair<bool,Measurement1D> &ip3dpvckf = IPTools::absoluteImpactParameter3D(ttckf,thevtx);
514     if (ip3dpvckf.first) {
515     outElectron->SetIp3dPVCkf(ckfsign*ip3dpvckf.second.value());
516     outElectron->SetIp3dPVCkfErr(ip3dpvckf.second.error());
517     }
518     else {
519     outElectron->SetIp3dPVCkf(-999.0);
520     }
521    
522     const std::pair<bool,Measurement1D> &d0pvbsckf = IPTools::absoluteTransverseImpactParameter(ttckf,thevtxbs);
523     if (d0pvbsckf.first) {
524     outElectron->SetD0PVBSCkf(ckfsignbs*d0pvbsckf.second.value());
525     outElectron->SetD0PVBSCkfErr(d0pvbsckf.second.error());
526     }
527     else {
528     outElectron->SetD0PVBSCkf(-999.0);
529     }
530    
531     const std::pair<bool,Measurement1D> &ip3dpvbsckf = IPTools::absoluteImpactParameter3D(ttckf,thevtxbs);
532     if (ip3dpvbsckf.first) {
533     outElectron->SetIp3dPVBSCkf(ckfsignbs*ip3dpvbsckf.second.value());
534     outElectron->SetIp3dPVBSCkfErr(ip3dpvbsckf.second.error());
535     }
536     else {
537     outElectron->SetIp3dPVBSCkf(-999.0);
538     }
539 bendavid 1.58
540     const std::pair<bool,Measurement1D> &d0pvubckf = IPTools::absoluteTransverseImpactParameter(ttckf,thevtxub);
541     if (d0pvubckf.first) {
542     outElectron->SetD0PVUBCkf(ckfsign*d0pvubckf.second.value());
543     outElectron->SetD0PVUBCkfErr(d0pvubckf.second.error());
544     }
545     else {
546     outElectron->SetD0PVUBCkf(-999.0);
547     }
548    
549    
550     const std::pair<bool,Measurement1D> &ip3dpvubckf = IPTools::absoluteImpactParameter3D(ttckf,thevtxub);
551     if (ip3dpvubckf.first) {
552     outElectron->SetIp3dPVUBCkf(ckfsign*ip3dpvubckf.second.value());
553     outElectron->SetIp3dPVUBCkfErr(ip3dpvubckf.second.error());
554     }
555     else {
556     outElectron->SetIp3dPVUBCkf(-999.0);
557     }
558    
559     const std::pair<bool,Measurement1D> &d0pvubbsckf = IPTools::absoluteTransverseImpactParameter(ttckf,thevtxubbs);
560     if (d0pvubbsckf.first) {
561     outElectron->SetD0PVUBBSCkf(ckfsignbs*d0pvubbsckf.second.value());
562     outElectron->SetD0PVUBBSCkfErr(d0pvubbsckf.second.error());
563     }
564     else {
565     outElectron->SetD0PVUBBSCkf(-999.0);
566     }
567    
568     const std::pair<bool,Measurement1D> &ip3dpvubbsckf = IPTools::absoluteImpactParameter3D(ttckf,thevtxubbs);
569     if (ip3dpvubbsckf.first) {
570     outElectron->SetIp3dPVUBBSCkf(ckfsignbs*ip3dpvubbsckf.second.value());
571     outElectron->SetIp3dPVUBBSCkfErr(ip3dpvubbsckf.second.error());
572     }
573     else {
574     outElectron->SetIp3dPVUBBSCkf(-999.0);
575     }
576    
577 bendavid 1.51 }
578     else {
579     outElectron->SetD0PVCkf(-999.0);
580     outElectron->SetIp3dPVCkf(-999.0);
581     outElectron->SetD0PVBSCkf(-999.0);
582     outElectron->SetIp3dPVBSCkf(-999.0);
583 bendavid 1.58
584     outElectron->SetD0PVUBCkf(-999.0);
585     outElectron->SetIp3dPVUBCkf(-999.0);
586     outElectron->SetD0PVUBBSCkf(-999.0);
587     outElectron->SetIp3dPVUBBSCkf(-999.0);
588 bendavid 1.51 }
589    
590 bendavid 1.43 if (verbose_>1) {
591     printf("gsf track pt = %5f\n",iM->gsfTrack()->pt());
592     printf("gsf track mode pt = %5f\n",iM->gsfTrack()->ptMode());
593     printf("ttrack pt = %5f\n",tt.initialFreeState().momentum().perp());
594     //printf("ttrackgsf pt = %5f\n",ttgsf.innermostMeasurementState().globalMomentum().perp());
595 bendavid 1.44 printf("ip3dpv reduced chisquared = %5f, probability = %5f\n", ip3dpv.second.value()/ip3dpv.second.error(), TMath::Prob(ip3dpv.second.value()/ip3dpv.second.error(),1));
596     //printf("gsf reduced chisquared = %5f, probability = %5f\n", pvGsfCompat.second/2, TMath::Prob(pvGsfCompat.second,2));
597 bendavid 1.43 }
598     }
599 bendavid 1.53
600 bendavid 1.43 //fill conversion partner track info
601 bendavid 1.54 if (recomputeConversionInfo_) {
602     ConversionFinder convFinder; outElectron->SetConvPartnerDCotTheta(iM->convDcot());
603     ConversionInfo convInfo = convFinder.getConversionInfo(*iM, hGeneralTracks, hGsfTracks, bfield);
604    
605     outElectron->SetConvFlag(convInfo.flag());
606     outElectron->SetConvPartnerDCotTheta(convInfo.dcot());
607     outElectron->SetConvPartnerDist(convInfo.dist());
608     outElectron->SetConvPartnerRadius(convInfo.radiusOfConversion());
609     reco::TrackRef ckfconvTrackRef = convInfo.conversionPartnerCtfTk();
610     reco::GsfTrackRef gsfconvTrackRef = convInfo.conversionPartnerGsfTk();
611    
612    
613     if ( gsfconvTrackRef.isNonnull() && gsfTrackMap_ ) {
614 pharris 1.59 try {outElectron->SetConvPartnerTrk(gsfTrackMap_->GetMit(edm::refToPtr(gsfconvTrackRef)));}
615     catch(...) {
616     if(checkClusterActive_) throw edm::Exception(edm::errors::Configuration, "FillerElectrons:FillDataBlock()\n")
617     << "Error! GSF track unmapped collection " << edmName_ << endl;
618     }
619    
620 bendavid 1.54 }
621     else if (ckfconvTrackRef.isNonnull() && trackerTrackMap_) {
622     outElectron->SetConvPartnerTrk(trackerTrackMap_->GetMit(edm::refToPtr(ckfconvTrackRef)));
623 bendavid 1.50 }
624 bendavid 1.54 }
625     else {
626     outElectron->SetConvFlag(iM->convFlags());
627     outElectron->SetConvPartnerDCotTheta(iM->convDcot());
628     outElectron->SetConvPartnerDist(iM->convDist());
629     outElectron->SetConvPartnerRadius(iM->convRadius());
630     reco::TrackBaseRef convTrackRef = iM->convPartner();
631     if (convTrackRef.isNonnull()) {
632     if ( dynamic_cast<const reco::GsfTrack*>(convTrackRef.get()) && gsfTrackMap_ ) {
633 pharris 1.59 try{outElectron->SetConvPartnerTrk(gsfTrackMap_->GetMit(mitedm::refToBaseToPtr(convTrackRef)));}
634     catch(...) {
635     if(checkClusterActive_) throw edm::Exception(edm::errors::Configuration, "FillerElectrons:FillDataBlock()\n")
636     << "Error! GSF track unmapped collection " << edmName_ << endl;
637     }
638 bendavid 1.54 }
639     else if (trackerTrackMap_) {
640     outElectron->SetConvPartnerTrk(trackerTrackMap_->GetMit(mitedm::refToBaseToPtr(convTrackRef)));
641     }
642 bendavid 1.50 }
643 bendavid 1.45 }
644 bendavid 1.50
645 loizides 1.29 // fill Electron ID information
646 peveraer 1.47 outElectron->SetPassLooseID((*eidLooseMap)[eRef]);
647     outElectron->SetPassTightID((*eidTightMap)[eRef]);
648 bendavid 1.49 if (!eIDLikelihoodName_.empty()) {
649     outElectron->SetIDLikelihood((*eidLikelihoodMap)[eRef]);
650     }
651     // fill corrected expected inner hits
652 bendavid 1.50 if(iM->gsfTrack().isNonnull()) {
653     outElectron->SetCorrectedNExpectedHitsInner(iM->gsfTrack()->trackerExpectedHitsInner().numberOfHits());
654 bendavid 1.49 }
655 bendavid 1.45 //fill additional conversion flag
656     outElectron->SetMatchesVertexConversion(convMatcher.matchesGoodConversion(*iM,hConversions));
657    
658 pharris 1.62 // add electron to map
659     edm::Ptr<reco::GsfElectron> thePtr(hElectronProduct, iM - inElectrons.begin());
660     electronMap_->Add(thePtr, outElectron);
661    
662 bendavid 1.31 if (verbose_>1) {
663     double recomass = sqrt(iM->energy()*iM->energy() - iM->p()*iM->p());
664 loizides 1.32 printf(" mithep::Electron, pt=%5f, eta=%5f, phi=%5f, energy=%5f, p=%5f, mass=%5f\n",
665     outElectron->Pt(), outElectron->Eta(), outElectron->Phi(),
666     outElectron->E(), outElectron->P(), outElectron->Mass());
667     printf("reco::GsfElectron , pt=%5f, eta=%5f, phi=%5f, energy=%5f, p=%5f, mass=%5f\n",
668     iM->pt(), iM->eta(), iM->phi(), iM->energy(), iM->p(), recomass);
669 bendavid 1.31 }
670 pharris 1.59 }
671 loizides 1.1 electrons_->Trim();
672     }