ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/HbbAnalysis/src/HistosMuons.cc
Revision: 1.4
Committed: Mon Jun 1 13:57:48 2009 UTC (15 years, 11 months ago) by amagnan
Content type: text/plain
Branch: MAIN
CVS Tags: v00-03-00, v00-02-01, v00-02-00, v00-01-00
Changes since 1.3: +1 -1 lines
Log Message:
add histos for elecs+MET selection

File Contents

# User Rev Content
1 amagnan 1.1 #include <iostream>
2     #include <fstream>
3    
4     #include "FWCore/MessageLogger/interface/MessageLogger.h"
5    
6     #include "DataFormats/MuonReco/interface/Muon.h"
7     #include "DataFormats/MuonReco/interface/MuonIsolation.h"
8    
9    
10     #include "UserCode/HbbAnalysis/interface/HistosMuons.hh"
11    
12     namespace HbbAnalysis {//namespace
13    
14     void HistosMuons::Initialise(TFileDirectory & aDir, std::string aName, bool aDoGenMatched){
15    
16     doGenMatched_ = aDoGenMatched;
17     CreateHistos(aName,aDir);
18    
19     p_nMuons = aDir.make<TH1F>("p_nMuons",";N_{muons};N_{entries}",20,0,20);
20    
21     p_caloCompat = aDir.make<TH1F>("p_caloCompat",";caloCompatibility;N_{entries}/0.01",100,0,1);
22     p_segCompat = aDir.make<TH1F>("p_segCompat",";segmentCompatibility;N_{entries}/0.01",100,0,1);
23     p_nChambers = aDir.make<TH1F>("p_nChambers",";N_{chambers};N_{entries}",50,0,50);
24     p_nMatchesLoose = aDir.make<TH1F>("p_nMatchesLoose",";N_{matches} (NoArbitration);N_{entries}",50,0,50);
25     p_nMatchesMedium = aDir.make<TH1F>("p_nMatchesMedium",";N_{matches} (SegmentArbitration);N_{entries}",50,0,50);
26     p_nMatchesTight = aDir.make<TH1F>("p_nMatchesTight",";N_{matches} (SegmentAndTrackArbitration);N_{entries}",50,0,50);
27     p_type = aDir.make<TH1F>("p_type",";type;N_{good}/type",15,0,15);
28     p_muonType = aDir.make<TH1F>("p_muonType",";muon type;N_{muons}/type",5,0,5);
29     p_muonID = aDir.make<TH1F>("p_muonID",";muon ID;N_{good}/ID",15,0,15);
30    
31     p_caloCompatvsPt = aDir.make<TH2F>("p_caloCompatvsPt",";p_{T} (GeV);caloCompatibility",20,0,100,100,0,1);
32     p_segCompatvsPt = aDir.make<TH2F>("p_segCompatvsPt",";p_{T} (GeV);segmentCompatibility",20,0,100,100,0,1);
33     p_nChambersvsPt = aDir.make<TH2F>("p_nChambersvsPt",";p_{T} (GeV);N_{chambers}",20,0,100,50,0,50);
34     p_nMatchesvsPt = aDir.make<TH2F>("p_nMatchesvsPt",";p_{T} (GeV);N_{matches} (SegmentArbitration)",20,0,100,50,0,50);
35     p_muonTypevsPt = aDir.make<TH2F>("p_muonTypevsPt",";p_{T} (GeV);muon type",20,0,100,5,0,5);
36     p_muonIDvsPt = aDir.make<TH2F>("p_muonIDvsPt",";p_{T} (GeV);muon ID",20,0,100,15,0,15);
37    
38     p_caloCompatvsEta = aDir.make<TH2F>("p_caloCompatvsEta",";#eta;caloCompatibility",50,-2.5,2.5,100,0,1);
39     p_segCompatvsEta = aDir.make<TH2F>("p_segCompatvsEta",";#eta;segmentCompatibility",50,-2.5,2.5,100,0,1);
40     p_nChambersvsEta = aDir.make<TH2F>("p_nChambersvsEta",";#eta;N_{chambers}",50,-2.5,2.5,50,0,50);
41     p_nMatchesvsEta = aDir.make<TH2F>("p_nMatchesvsEta",";#eta;N_{matches} (SegmentArbitration)",50,-2.5,2.5,50,0,50);
42     p_muonTypevsEta = aDir.make<TH2F>("p_muonTypevsEta",";#eta;muon type",50,-2.5,2.5,5,0,5);
43     p_muonIDvsEta = aDir.make<TH2F>("p_muonIDvsEta",";#eta;muon ID",50,-2.5,2.5,15,0,15);
44    
45     //efficiency of ID cut
46     peff_muID = aDir.make<TH1F>("peff_muID",";muon ID;N_{ID}/N_{tot}",15,0,15);
47     peff_muIDvsEta = aDir.make<TH2F>("peff_muIDvsEta",";#eta;muon ID;N_{ID}/N_{tot}",50,-2.5,2.5,15,0,15);
48     idEff_.initialise(peff_muID->GetNbinsX(),peff_muID->GetXaxis()->GetXmin(),peff_muID->GetXaxis()->GetXmax());
49     for (unsigned id(0); id<15;id++){
50     idEffEta_[id].initialise(peff_muIDvsEta->GetNbinsX(),peff_muIDvsEta->GetXaxis()->GetXmin(),peff_muIDvsEta->GetXaxis()->GetXmax());
51     }
52    
53    
54     p_isoR03_emEt = aDir.make<TH1F>("p_isoR03_emEt",";emEt (#DeltaR=0.3, GeV);N_{entries}/0.05 GeV",500,0,25);
55     p_isoR05_emEt = aDir.make<TH1F>("p_isoR05_emEt",";emEt (#DeltaR=0.5, GeV);N_{entries}/0.05 GeV",500,0,25);
56     p_isoR03_hadEt = aDir.make<TH1F>("p_isoR03_hadEt",";hadEt (#DeltaR=0.3, GeV);N_{entries}/0.05 GeV",500,0,25);
57     p_isoR05_hadEt = aDir.make<TH1F>("p_isoR05_hadEt",";hadEt (#DeltaR=0.5, GeV);N_{entries}/0.05 GeV",500,0,25);
58     p_isoR03_nTracks = aDir.make<TH1F>("p_isoR03_nTracks",";N_{tracks} (#DeltaR=0.3);N_{entries}",20,0,20);
59     p_isoR05_nTracks = aDir.make<TH1F>("p_isoR05_nTracks",";N_{tracks} (#DeltaR=0.5);N_{entries}",20,0,20);
60     p_isoR03_nJets = aDir.make<TH1F>("p_isoR03_nJets",";N_{Jets} (#DeltaR=0.3);N_{entries}",10,0,10);
61     p_isoR05_nJets = aDir.make<TH1F>("p_isoR05_nJets",";N_{Jets} (#DeltaR=0.5);N_{entries}",10,0,10);
62    
63    
64     //TString lIsoVar[4] = {"sumPt","sumPtOverPt","combIso","combIsoOverPt"};
65     assert (getNVar() == 4);
66     p_isoR03[0] = aDir.make<TH1F>("p_isoR03_sumPt",";sumPt (#DeltaR=0.3, GeV);N_{entries}/0.05 GeV",300,0,30);
67     p_isoR03[1] = aDir.make<TH1F>("p_isoR03_sumPtOverPt",";sumPt/Pt (#DeltaR=0.3);N_{entries}",100,0,1);
68     p_isoR03[2] = aDir.make<TH1F>("p_isoR03_combIso",";combIso (#DeltaR=0.3, GeV);N_{entries}/0.05 GeV",300,0,30);
69     p_isoR03[3] = aDir.make<TH1F>("p_isoR03_combIsoOverPt",";combIso/Pt (#DeltaR=0.3);N_{entries}",100,0,1);
70    
71     p_isoR05[0] = aDir.make<TH1F>("p_isoR05_sumPt",";sumPt (#DeltaR=0.5, GeV);N_{entries}/0.05 GeV",300,0,30);
72     p_isoR05[1] = aDir.make<TH1F>("p_isoR05_sumPtOverPt",";sumPt/Pt (#DeltaR=0.5);N_{entries}",100,0,1);
73     p_isoR05[2] = aDir.make<TH1F>("p_isoR05_combIso",";combIso (#DeltaR=0.5, GeV);N_{entries}/0.05 GeV",300,0,30);
74     p_isoR05[3] = aDir.make<TH1F>("p_isoR05_combIsoOverPt",";combIso/Pt (#DeltaR=0.5);N_{entries}",100,0,1);
75    
76     //efficiency curves vs cut value
77     peff_muIsoR03[0] = aDir.make<TH1F>("peff_muIsoR03_sumPt",";sumPt (#DeltaR=0.3, GeV);N_{sel}/N_{tot}",50,0,5);
78     peff_muIsoR03[1] = aDir.make<TH1F>("peff_muIsoR03_sumPtOverPt",";sumPt/Pt (#DeltaR=0.3);N_{sel}/N_{tot}",50,0,1);
79     peff_muIsoR03[2] = aDir.make<TH1F>("peff_muIsoR03_combIso",";combIso (#DeltaR=0.3, GeV);N_{sel}/N_{tot}",50,0,10);
80     peff_muIsoR03[3] = aDir.make<TH1F>("peff_muIsoR03_combIsoOverPt",";combIso/Pt (#DeltaR=0.3);N_{sel}/N_{tot}",50,0,1);
81    
82     peff_muIsoR05[0] = aDir.make<TH1F>("peff_muIsoR05_sumPt",";sumPt (#DeltaR=0.5, GeV);N_{sel}/N_{tot}",50,0,5);
83     peff_muIsoR05[1] = aDir.make<TH1F>("peff_muIsoR05_sumPtOverPt",";sumPt/Pt (#DeltaR=0.5);N_{sel}/N_{tot}",50,0,1);
84     peff_muIsoR05[2] = aDir.make<TH1F>("peff_muIsoR05_combIso",";combIso (#DeltaR=0.5, GeV);N_{sel}/N_{tot}",50,0,10);
85     peff_muIsoR05[3] = aDir.make<TH1F>("peff_muIsoR05_combIsoOverPt",";combIso/Pt (#DeltaR=0.5);N_{sel}/N_{tot}",50,0,1);
86    
87     //efficiency curves vs pT muon
88     peff_muPt[0] = aDir.make<TH1F>("peff_muPt_sumPt",";p_{T} (GeV);N_{sel}/N_{tot}",18,10,100);
89     const unsigned int lNBinsPt = peff_muPt[0]->GetNbinsX();
90     const double lMinPt = peff_muPt[0]->GetXaxis()->GetXmin();
91     const double lMaxPt = peff_muPt[0]->GetXaxis()->GetXmax();
92     peff_muPt[1] = aDir.make<TH1F>("peff_muPt_sumPtOverPt",";p_{T} (GeV);N_{sel}/N_{tot}",lNBinsPt,lMinPt,lMaxPt);
93     peff_muPt[2] = aDir.make<TH1F>("peff_muPt_combIso",";p_{T} (GeV);N_{sel}/N_{tot}",lNBinsPt,lMinPt,lMaxPt);
94     peff_muPt[3] = aDir.make<TH1F>("peff_muPt_combIsoOverPt",";p_{T} (GeV);N_{sel}/N_{tot}",lNBinsPt,lMinPt,lMaxPt);
95    
96    
97     //efficiency curves vs eta muon
98     peff_muEta[0] = aDir.make<TH1F>("peff_muEta_sumPt",";#eta;N_{sel}/N_{tot}",42,-2.1,2.1);
99     const unsigned int lNBinsEta = peff_muEta[0]->GetNbinsX();
100     const double lMinEta = peff_muEta[0]->GetXaxis()->GetXmin();
101     const double lMaxEta = peff_muEta[0]->GetXaxis()->GetXmax();
102     peff_muEta[1] = aDir.make<TH1F>("peff_muEta_sumPtOverPt",";#eta;N_{sel}/N_{tot}",lNBinsEta,lMinEta,lMaxEta);
103     peff_muEta[2] = aDir.make<TH1F>("peff_muEta_combIso",";#eta;N_{sel}/N_{tot}",lNBinsEta,lMinEta,lMaxEta);
104     peff_muEta[3] = aDir.make<TH1F>("peff_muEta_combIsoOverPt",";#eta;N_{sel}/N_{tot}",lNBinsEta,lMinEta,lMaxEta);
105    
106     //efficiency calculation tool
107     for (unsigned int i(0); i<getNVar(); i++){
108     const unsigned int lNBinsR03 = peff_muIsoR03[i]->GetNbinsX();
109     const double lMinR03 = peff_muIsoR03[i]->GetXaxis()->GetXmin();
110     const double lMaxR03 = peff_muIsoR03[i]->GetXaxis()->GetXmax();
111     const unsigned int lNBinsR05 = peff_muIsoR05[i]->GetNbinsX();
112     const double lMinR05 = peff_muIsoR05[i]->GetXaxis()->GetXmin();
113     const double lMaxR05 = peff_muIsoR05[i]->GetXaxis()->GetXmax();
114     isoR03Eff_[i].initialise(lNBinsR03,lMinR03,lMaxR03);
115     isoR05Eff_[i].initialise(lNBinsR05,lMinR05,lMaxR05);
116     isoEffPt_[i].initialise(lNBinsPt,lMinPt,lMaxPt);
117     isoEffEta_[i].initialise(lNBinsEta,lMinEta,lMaxEta);
118     }
119    
120    
121     }//Initialise
122    
123 amagnan 1.3 void HistosMuons::FillEventHistograms(const edm::Handle<std::vector<pat::Muon> > & aCol){
124 amagnan 1.1
125     if (doGenMatched_) {
126 amagnan 1.2 unsigned int nMatched = 0;
127     for (std::vector<pat::Muon>::const_iterator iMuon = aCol->begin();
128     iMuon != aCol->end();
129 amagnan 1.1 iMuon++)
130     {
131 amagnan 1.2 if (MatchesGenMuon(*iMuon)) nMatched++;
132 amagnan 1.1 }
133 amagnan 1.2 p_nMuons->Fill(nMatched);
134 amagnan 1.1 }
135 amagnan 1.2 else p_nMuons->Fill(aCol->size());
136 amagnan 1.1
137     }
138    
139     void HistosMuons::FillHistograms(const pat::Muon & aMuon, bool isLead){//FillHistograms
140    
141     const reco::Muon* recoMuon = dynamic_cast<const reco::Muon*>(aMuon.originalObject());
142     bool lIsGenMatched = !doGenMatched_ || (doGenMatched_ && MatchesGenMuon(aMuon));
143     if (lIsGenMatched) {
144     if (isLead) {
145 amagnan 1.4 FillBaseHistograms(aMuon.pt(),aMuon.eta(),aMuon.phi(),aMuon.charge());
146 amagnan 1.1 p_caloCompat->Fill(recoMuon->caloCompatibility());
147     p_segCompat->Fill(recoMuon->segmentCompatibility());
148     p_nChambers->Fill(recoMuon->numberOfChambers());
149     p_nMatchesLoose->Fill(recoMuon->numberOfMatches(reco::Muon::NoArbitration));
150     p_nMatchesMedium->Fill(recoMuon->numberOfMatches(reco::Muon::SegmentArbitration));
151     p_nMatchesTight->Fill(recoMuon->numberOfMatches(reco::Muon::SegmentAndTrackArbitration));
152     p_type->Fill(recoMuon->type());
153    
154     if (recoMuon->isGlobalMuon()) {
155     p_muonType->Fill(1);
156     p_muonTypevsPt->Fill(aMuon.pt(),1);
157     p_muonTypevsEta->Fill(aMuon.eta(),1);
158     }
159     else if (recoMuon->isTrackerMuon()) {
160     p_muonType->Fill(2);
161     p_muonTypevsPt->Fill(aMuon.pt(),2);
162     p_muonTypevsEta->Fill(aMuon.eta(),2);
163     }
164     else if (recoMuon->isStandAloneMuon()) {
165     p_muonType->Fill(3);
166     p_muonTypevsPt->Fill(aMuon.pt(),3);
167     p_muonTypevsEta->Fill(aMuon.eta(),3);
168     }
169     else if (recoMuon->isCaloMuon()) {
170     p_muonType->Fill(4);
171     p_muonTypevsPt->Fill(aMuon.pt(),4);
172     p_muonTypevsEta->Fill(aMuon.eta(),4);
173     }
174     else {
175     p_muonType->Fill(0);
176     p_muonTypevsPt->Fill(aMuon.pt(),0);
177     p_muonTypevsEta->Fill(aMuon.eta(),0);
178     }
179    
180     for (unsigned int id(0); id<idEff_.numberOfBins(); id++){
181     idEff_.incrementTotal(id);
182     for (unsigned ieta(0); ieta<idEffEta_[id].numberOfBins(); ieta++){
183     double lCutMin = idEffEta_[id].xMin()+idEffEta_[id].stepSize()*ieta;
184     double lCutMax = idEffEta_[id].xMin()+idEffEta_[id].stepSize()*(ieta+1);
185    
186     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
187     idEffEta_[id].incrementTotal(ieta);
188     }
189     }
190     }
191    
192     for (unsigned ieta(0); ieta<idEffEta_[0].numberOfBins(); ieta++){
193     double lCutMin = idEffEta_[0].xMin()+idEffEta_[0].stepSize()*ieta;
194     double lCutMax = idEffEta_[0].xMin()+idEffEta_[0].stepSize()*(ieta+1);
195    
196    
197     if (recoMuon->isGood(reco::Muon::AllGlobalMuons)){
198     if (ieta == 0) {
199     p_muonID->Fill(1);
200     p_muonIDvsPt->Fill(aMuon.pt(),1);
201     p_muonIDvsEta->Fill(aMuon.eta(),1);
202     idEff_.incrementPass(1);
203     }
204     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
205     idEffEta_[1].incrementPass(ieta);
206     }
207     }
208     if (recoMuon->isGood(reco::Muon::AllStandAloneMuons)){
209     if (ieta == 0) {
210     p_muonID->Fill(2);
211     p_muonIDvsPt->Fill(aMuon.pt(),2);
212     p_muonIDvsEta->Fill(aMuon.eta(),2);
213     idEff_.incrementPass(2);
214     }
215     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
216     idEffEta_[2].incrementPass(ieta);
217     }
218     }
219     if (recoMuon->isGood(reco::Muon::AllTrackerMuons)){
220     if (ieta == 0) {
221     p_muonID->Fill(3);
222     p_muonIDvsPt->Fill(aMuon.pt(),3);
223     p_muonIDvsEta->Fill(aMuon.eta(),3);
224     idEff_.incrementPass(3);
225     }
226     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
227     idEffEta_[3].incrementPass(ieta);
228     }
229     }
230     if (recoMuon->isGood(reco::Muon::TrackerMuonArbitrated)){
231     if (ieta == 0) {
232     p_muonID->Fill(4);
233     p_muonIDvsPt->Fill(aMuon.pt(),4);
234     p_muonIDvsEta->Fill(aMuon.eta(),4);
235     idEff_.incrementPass(4);
236     }
237     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
238     idEffEta_[4].incrementPass(ieta);
239     }
240     }
241     if (recoMuon->isGood(reco::Muon::AllArbitrated)){
242     if (ieta == 0) {
243     p_muonID->Fill(5);
244     p_muonIDvsPt->Fill(aMuon.pt(),5);
245     p_muonIDvsEta->Fill(aMuon.eta(),5);
246     idEff_.incrementPass(5);
247     }
248     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
249     idEffEta_[5].incrementPass(ieta);
250     }
251     }
252     if (recoMuon->isGood(reco::Muon::GlobalMuonPromptTight)){
253     if (ieta == 0) {
254     p_muonID->Fill(6);
255     p_muonIDvsPt->Fill(aMuon.pt(),6);
256     p_muonIDvsEta->Fill(aMuon.eta(),6);
257     idEff_.incrementPass(6);
258     }
259     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
260     idEffEta_[6].incrementPass(ieta);
261     }
262     }
263     if (recoMuon->isGood(reco::Muon::TMLastStationLoose)){
264     if (ieta == 0) {
265     p_muonID->Fill(7);
266     p_muonIDvsPt->Fill(aMuon.pt(),7);
267     p_muonIDvsEta->Fill(aMuon.eta(),7);
268     idEff_.incrementPass(7);
269     }
270     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
271     idEffEta_[7].incrementPass(ieta);
272     }
273     }
274     if (recoMuon->isGood(reco::Muon::TMLastStationTight)){
275     if (ieta == 0) {
276     p_muonID->Fill(8);
277     p_muonIDvsPt->Fill(aMuon.pt(),8);
278     p_muonIDvsEta->Fill(aMuon.eta(),8);
279     idEff_.incrementPass(8);
280     }
281     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
282     idEffEta_[8].incrementPass(ieta);
283     }
284     }
285     if (recoMuon->isGood(reco::Muon::TM2DCompatibilityLoose)){
286     if (ieta == 0) {
287     p_muonID->Fill(9);
288     p_muonIDvsPt->Fill(aMuon.pt(),9);
289     p_muonIDvsEta->Fill(aMuon.eta(),9);
290     idEff_.incrementPass(9);
291     }
292     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
293     idEffEta_[9].incrementPass(ieta);
294     }
295     }
296     if (recoMuon->isGood(reco::Muon::TM2DCompatibilityTight)){
297     if (ieta == 0) {
298     p_muonID->Fill(10);
299     p_muonIDvsPt->Fill(aMuon.pt(),10);
300     p_muonIDvsEta->Fill(aMuon.eta(),10);
301     idEff_.incrementPass(10);
302     }
303     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
304     idEffEta_[10].incrementPass(ieta);
305     }
306     }
307     if (recoMuon->isGood(reco::Muon::TMOneStationLoose)){
308     if (ieta == 0) {
309     p_muonID->Fill(11);
310     p_muonIDvsPt->Fill(aMuon.pt(),11);
311     p_muonIDvsEta->Fill(aMuon.eta(),11);
312     idEff_.incrementPass(11);
313     }
314     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
315     idEffEta_[11].incrementPass(ieta);
316     }
317     }
318     if (recoMuon->isGood(reco::Muon::TMOneStationTight)){
319     if (ieta == 0) {
320     p_muonID->Fill(12);
321     p_muonIDvsPt->Fill(aMuon.pt(),12);
322     p_muonIDvsEta->Fill(aMuon.eta(),12);
323     idEff_.incrementPass(12);
324     }
325     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
326     idEffEta_[12].incrementPass(ieta);
327     }
328     }
329     if (recoMuon->isGood(reco::Muon::TMLastStationOptimizedLowPtLoose)){
330     if (ieta == 0) {
331     p_muonID->Fill(13);
332     p_muonIDvsPt->Fill(aMuon.pt(),13);
333     p_muonIDvsEta->Fill(aMuon.eta(),13);
334     idEff_.incrementPass(13);
335     }
336     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
337     idEffEta_[13].incrementPass(ieta);
338     }
339     }
340     if (recoMuon->isGood(reco::Muon::TMLastStationOptimizedLowPtTight)){
341     if (ieta == 0) {
342     p_muonID->Fill(14);
343     p_muonIDvsPt->Fill(aMuon.pt(),14);
344     p_muonIDvsEta->Fill(aMuon.eta(),14);
345     idEff_.incrementPass(14);
346     }
347     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
348     idEffEta_[14].incrementPass(ieta);
349     }
350     }
351    
352     }
353    
354     p_caloCompatvsPt->Fill(aMuon.pt(),recoMuon->caloCompatibility());
355     p_segCompatvsPt->Fill(aMuon.pt(),recoMuon->segmentCompatibility());
356     p_nChambersvsPt->Fill(aMuon.pt(),recoMuon->numberOfChambers());
357     p_nMatchesvsPt->Fill(aMuon.pt(),recoMuon->numberOfMatches());
358    
359     p_caloCompatvsEta->Fill(aMuon.eta(),recoMuon->caloCompatibility());
360     p_segCompatvsEta->Fill(aMuon.eta(),recoMuon->segmentCompatibility());
361     p_nChambersvsEta->Fill(aMuon.eta(),recoMuon->numberOfChambers());
362     p_nMatchesvsEta->Fill(aMuon.eta(),recoMuon->numberOfMatches());
363    
364     }
365    
366     p_isoR03_emEt->Fill(recoMuon->isolationR03().emEt);
367     p_isoR05_emEt->Fill(recoMuon->isolationR05().emEt);
368     p_isoR03_hadEt->Fill(recoMuon->isolationR03().hadEt);
369     p_isoR05_hadEt->Fill(recoMuon->isolationR05().hadEt);
370     p_isoR03_nTracks->Fill(recoMuon->isolationR03().nTracks);
371     p_isoR05_nTracks->Fill(recoMuon->isolationR05().nTracks);
372     p_isoR03_nJets->Fill(recoMuon->isolationR03().nJets);
373     p_isoR05_nJets->Fill(recoMuon->isolationR05().nJets);
374    
375     double lIsoR03Var[4] = {
376     recoMuon->isolationR03().sumPt,
377     recoMuon->isolationR03().sumPt/aMuon.pt(),
378     recoMuon->isolationR03().sumPt+recoMuon->isolationR03().emEt,
379     (recoMuon->isolationR03().sumPt+recoMuon->isolationR03().emEt)/aMuon.pt()
380     };
381     double lIsoR05Var[4] = {
382     recoMuon->isolationR05().sumPt,
383     recoMuon->isolationR05().sumPt/aMuon.pt(),
384     recoMuon->isolationR05().sumPt+recoMuon->isolationR05().emEt,
385     (recoMuon->isolationR05().sumPt+recoMuon->isolationR05().emEt)/aMuon.pt()
386     };
387    
388     double lIsoR03Cut[4] = {3.,0.09,6.,0.18};
389    
390     for (unsigned int i(0); i<getNVar(); i++){//loop on iso variables
391     p_isoR03[i]->Fill(lIsoR03Var[i]);
392     p_isoR05[i]->Fill(lIsoR05Var[i]);
393    
394     //R03
395     for (unsigned int lBin(0); lBin<isoR03Eff_[i].numberOfBins(); lBin++){
396     double lCut = isoR03Eff_[i].xMin()+isoR03Eff_[i].stepSize()*lBin;
397    
398     //std::cout << "******** var R03 #" << i << ", bin " << lBin << " isolation = " << lIsoR03Var[i] << ", cut = " << lCut << std::endl;
399    
400     isoR03Eff_[i].incrementTotal(lBin);
401     if (lIsoR03Var[i] <= lCut){
402     isoR03Eff_[i].incrementPass(lBin);
403     }
404     }
405    
406     //R05
407     for (unsigned int lBin(0); lBin<isoR05Eff_[i].numberOfBins(); lBin++){
408     double lCut = isoR05Eff_[i].xMin()+isoR05Eff_[i].stepSize()*lBin;
409    
410     //std::cout << "******** var R05 #" << i << ", bin " << lBin << " isolation = " << lIsoR03Var[i] << ", cut = " << lCut << std::endl;
411    
412     isoR05Eff_[i].incrementTotal(lBin);
413     if (lIsoR05Var[i] <= lCut){
414     isoR05Eff_[i].incrementPass(lBin);
415     }
416     }
417    
418     //eff vs pT
419     for (unsigned int lBin(0); lBin<isoEffPt_[i].numberOfBins(); lBin++){
420     double lCutMin = isoEffPt_[i].xMin()+isoEffPt_[i].stepSize()*lBin;
421     double lCutMax = isoEffPt_[i].xMin()+isoEffPt_[i].stepSize()*(lBin+1);
422    
423     if (aMuon.pt() >= lCutMin && aMuon.pt() < lCutMax){
424     isoEffPt_[i].incrementTotal(lBin);
425    
426     if (lIsoR03Var[i] <= lIsoR03Cut[i]){
427     isoEffPt_[i].incrementPass(lBin);
428     }
429     }
430     }
431    
432     //eff vs eta
433     for (unsigned int lBin(0); lBin<isoEffEta_[i].numberOfBins(); lBin++){
434     double lCutMin = isoEffEta_[i].xMin()+isoEffEta_[i].stepSize()*lBin;
435     double lCutMax = isoEffEta_[i].xMin()+isoEffEta_[i].stepSize()*(lBin+1);
436    
437     if (aMuon.eta() >= lCutMin && aMuon.eta() < lCutMax){
438     isoEffEta_[i].incrementTotal(lBin);
439    
440     if (lIsoR03Var[i] <= lIsoR03Cut[i]){
441     isoEffEta_[i].incrementPass(lBin);
442     }
443     }
444     }
445    
446     }//loop on iso variables
447     }//lIsGenMatched
448    
449     }//FillHistograms
450    
451    
452     void HistosMuons::FillEffHistograms(){//FillEffHistograms
453    
454    
455     for (unsigned int id(0); id<idEff_.numberOfBins(); id++){
456     peff_muID->SetBinContent(id+1,idEff_.getRatio(id));
457     peff_muID->SetBinError(id+1,idEff_.getRatioError(id));
458    
459     for (unsigned int lBin(0); lBin<idEffEta_[id].numberOfBins(); lBin++){
460     peff_muIDvsEta->SetBinContent(lBin+1,id+1,idEffEta_[id].getRatio(lBin));
461     peff_muIDvsEta->SetBinError(lBin+1,id+1,idEffEta_[id].getRatioError(lBin));
462     }
463    
464     }
465    
466     for (unsigned int i(0); i<getNVar(); i++){//loop on iso variables
467    
468     for (unsigned int lBin(0); lBin<isoR03Eff_[i].numberOfBins(); lBin++){
469     peff_muIsoR03[i]->SetBinContent(lBin+1,isoR03Eff_[i].getRatio(lBin));
470     peff_muIsoR03[i]->SetBinError(lBin+1,isoR03Eff_[i].getRatioError(lBin));
471     }
472     for (unsigned int lBin(0); lBin<isoR05Eff_[i].numberOfBins(); lBin++){
473     peff_muIsoR05[i]->SetBinContent(lBin+1,isoR05Eff_[i].getRatio(lBin));
474     peff_muIsoR05[i]->SetBinError(lBin+1,isoR05Eff_[i].getRatioError(lBin));
475     }
476    
477     //eff vs pT
478     for (unsigned int lBin(0); lBin<isoEffPt_[i].numberOfBins(); lBin++){
479     peff_muPt[i]->SetBinContent(lBin+1,isoEffPt_[i].getRatio(lBin));
480     peff_muPt[i]->SetBinError(lBin+1,isoEffPt_[i].getRatioError(lBin));
481     }
482    
483     //eff vs eta
484     for (unsigned int lBin(0); lBin<isoEffEta_[i].numberOfBins(); lBin++){
485     peff_muEta[i]->SetBinContent(lBin+1,isoEffEta_[i].getRatio(lBin));
486     peff_muEta[i]->SetBinError(lBin+1,isoEffEta_[i].getRatioError(lBin));
487     }
488    
489     }//loop on iso variables
490    
491     }//FillEffHistograms
492    
493    
494     bool HistosMuons::MatchesGenMuon(const pat::Muon& aPatMuon)
495     {
496    
497 amagnan 1.2 bool isGenMatched = false;
498 amagnan 1.1
499     //if (aPatMuon.genParticlesSize() == 0) {
500     // edm::LogWarning("MatchesGenMuon") << "No ref to gen particles for pat::Muon ! ";
501     // return false;
502     // }
503    
504     const std::vector<reco::GenParticleRef> & lVec = aPatMuon.genParticleRefs();
505     //for ( std::vector<reco::GenParticleRef>::const_iterator it = aPatMuon.genParticleRefs().begin();
506     // it != aPatMuon.genParticleRefs().end(); ++it ) {
507     for ( std::vector<reco::GenParticleRef>::const_iterator it = lVec.begin();
508     it != lVec.end(); ++it ) {
509     if ( it->ref().isNonnull() && it->ref().isValid() ) {
510     const reco::GenParticleRef & genParticle = (*it);
511     //assert (genParticle->isNonnull() && genParticle->isValid());
512 amagnan 1.2 if ( genParticle->pdgId() == -13 || genParticle->pdgId() == +13 ) isGenMatched = true;
513 amagnan 1.1 } else {
514     edm::LogWarning("MatchesGenMuon") << " edm::Ref of genParticle associated to pat::Muon is invalid !!";
515     }
516     }
517 amagnan 1.2 return isGenMatched;
518 amagnan 1.1 }
519    
520     }//namespace
521    
522