Go to the documentation of this file.00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012 #include "RecoMET/METAlgorithms/interface/PFSpecificAlgo.h"
00013 #include "RecoMET/METAlgorithms/interface/SignAlgoResolutions.h"
00014 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
00015
00016
00017 using namespace reco;
00018
00019
00020 reco::PFMET PFSpecificAlgo::addInfo(edm::Handle<edm::View<Candidate> > PFCandidates, CommonMETData met)
00021 {
00022 SpecificPFMETData specific = mkSpecificPFMETData(PFCandidates);
00023
00024 const LorentzVector p4(met.mex , met.mey, 0.0, met.met);
00025 const Point vtx(0.0,0.0,0.0);
00026 PFMET pfMET(specific, met.sumet, p4, vtx );
00027
00028 if(doSignificance) pfMET.setSignificanceMatrix(mkSignifMatrix(PFCandidates));
00029
00030 return pfMET;
00031 }
00032
00033
00034 void PFSpecificAlgo::runSignificance(metsig::SignAlgoResolutions &resolutions, edm::Handle<edm::View<reco::PFJet> > jets)
00035 {
00036 doSignificance = true;
00037 pfsignalgo_.setResolutions( &resolutions );
00038 pfsignalgo_.addPFJets(jets);
00039 }
00040
00041
00042 void PFSpecificAlgo::initializeSpecificPFMETData(SpecificPFMETData &specific)
00043 {
00044 specific.NeutralEMFraction = 0.0;
00045 specific.NeutralHadFraction = 0.0;
00046 specific.ChargedEMFraction = 0.0;
00047 specific.ChargedHadFraction = 0.0;
00048 specific.MuonFraction = 0.0;
00049 specific.Type6Fraction = 0.0;
00050 specific.Type7Fraction = 0.0;
00051 }
00052
00053
00054 SpecificPFMETData PFSpecificAlgo::mkSpecificPFMETData(edm::Handle<edm::View<reco::Candidate> > &PFCandidates)
00055 {
00056 if(!PFCandidates->size())
00057 {
00058 SpecificPFMETData specific;
00059 initializeSpecificPFMETData(specific);
00060 return specific;
00061 }
00062
00063 double NeutralEMEt = 0.0;
00064 double NeutralHadEt = 0.0;
00065 double ChargedEMEt = 0.0;
00066 double ChargedHadEt = 0.0;
00067 double MuonEt = 0.0;
00068 double type6Et = 0.0;
00069 double type7Et = 0.0;
00070
00071 for( edm::View<reco::Candidate>::const_iterator iParticle = (PFCandidates.product())->begin(); iParticle != (PFCandidates.product())->end(); ++iParticle )
00072 {
00073 const PFCandidate* pfCandidate = dynamic_cast<const PFCandidate*> (&(*iParticle));
00074 if (!pfCandidate) continue;
00075 const double theta = pfCandidate->theta();
00076 const double e = pfCandidate->energy();
00077 const double et = e*sin(theta);
00078
00079 if (pfCandidate->particleId() == 1) ChargedHadEt += et;
00080 if (pfCandidate->particleId() == 2) ChargedEMEt += et;
00081 if (pfCandidate->particleId() == 3) MuonEt += et;
00082 if (pfCandidate->particleId() == 4) NeutralEMEt += et;
00083 if (pfCandidate->particleId() == 5) NeutralHadEt += et;
00084 if (pfCandidate->particleId() == 6) type6Et += et;
00085 if (pfCandidate->particleId() == 7) type7Et += et;
00086
00087
00088 }
00089
00090 const double Et_total = NeutralEMEt + NeutralHadEt + ChargedEMEt + ChargedHadEt + MuonEt + type6Et + type7Et;
00091 SpecificPFMETData specific;
00092 initializeSpecificPFMETData(specific);
00093 if (Et_total!=0.0)
00094 {
00095 specific.NeutralEMFraction = NeutralEMEt/Et_total;
00096 specific.NeutralHadFraction = NeutralHadEt/Et_total;
00097 specific.ChargedEMFraction = ChargedEMEt/Et_total;
00098 specific.ChargedHadFraction = ChargedHadEt/Et_total;
00099 specific.MuonFraction = MuonEt/Et_total;
00100 specific.Type6Fraction = type6Et/Et_total;
00101 specific.Type7Fraction = type7Et/Et_total;
00102 }
00103 return specific;
00104 }
00105
00106
00107 TMatrixD PFSpecificAlgo::mkSignifMatrix(edm::Handle<edm::View<reco::Candidate> > &PFCandidates)
00108 {
00109 pfsignalgo_.useOriginalPtrs(PFCandidates.id());
00110 for(edm::View<reco::Candidate>::const_iterator iParticle = (PFCandidates.product())->begin(); iParticle != (PFCandidates.product())->end(); ++iParticle )
00111 {
00112 const PFCandidate* pfCandidate = dynamic_cast<const PFCandidate*> (&(*iParticle));
00113 if (!pfCandidate) continue;
00114 reco::CandidatePtr dau(PFCandidates, iParticle - PFCandidates->begin());
00115 if(dau.isNull()) continue;
00116 if(!dau.isAvailable()) continue;
00117 reco::PFCandidatePtr pf(dau.id(), pfCandidate, dau.key());
00118 pfsignalgo_.addPFCandidate(pf);
00119 }
00120 return pfsignalgo_.getSignifMatrix();
00121 }
00122
00123