ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/Selection/src/SelectionFuncs.cc
Revision: 1.12
Committed: Tue Oct 23 11:39:21 2012 UTC (12 years, 6 months ago) by dkralph
Content type: text/plain
Branch: MAIN
Changes since 1.11: +567 -86 lines
Log Message: 
*** empty log message ***

File Contents

# User Rev Content
1 dkralph 1.1 #include "SelectionFuncs.h"
2    
3     extern vector<SimpleLepton> failingLeptons;
4     extern vector<SimpleLepton> passingLeptons;
5    
6     extern vector<unsigned> cutvec;
7     extern vector<vector<unsigned> > zcutvec;
8     extern vector<vector<unsigned> > zzcutvec;
9     extern map<unsigned,float> evtrhoMap;
10     extern bool passes_HLT_MC;
11 dkralph 1.12 extern RunLumiRangeMap rlrm;
12 khahn 1.5
13 dkralph 1.1 //----------------------------------------------------------------------------
14 dkralph 1.12 bool setPV(ControlFlags &ctrl,
15 dkralph 1.1 const mithep::Array<mithep::Vertex> * vtxArr,
16     const mithep::Vertex* &vtx)
17     //----------------------------------------------------------------------------
18 dkralph 1.6 {
19 dkralph 1.1
20     const mithep::Vertex *bestPV = 0;
21     float best_sumpt=-1;
22    
23     // good PV requirements
24     const UInt_t fMinNTracksFit = 0;
25     const Double_t fMinNdof = 4;
26     const Double_t fMaxAbsZ = 24;
27     const Double_t fMaxRho = 2;
28    
29     for(int i=0; i<vtxArr->GetEntries(); ++i) {
30     const mithep::Vertex *pv = (mithep::Vertex*)(vtxArr->At(i));
31     if( ctrl.debug ) cout << "vertex :: " << i << "\tntrks: " << pv->NTracks() << endl;
32    
33     // Select best PV for corrected d0; if no PV passing cuts, the first PV in the collection will be used
34     if(!pv->IsValid()) continue;
35     if(pv->NTracksFit() < fMinNTracksFit) continue;
36     if(pv->Ndof() < fMinNdof) continue;
37     if(fabs(pv->Z()) > fMaxAbsZ) continue;
38     if(pv->Position().Rho() > fMaxRho) continue;
39    
40     // take the first ...
41     bestPV = pv;
42     break;
43    
44     // this never reached ...
45     float tmp_sumpt=0;
46     for( int t=0; t<pv->NTracks(); t++ )
47     tmp_sumpt += pv->Trk(t)->Pt();
48    
49     if( tmp_sumpt > best_sumpt ) {
50     bestPV = pv;
51     best_sumpt = tmp_sumpt;
52     if( ctrl.debug) cout << "new PV set, pt : " << best_sumpt << endl;
53     }
54     }
55    
56     // sync
57     if(!bestPV)
58     return false;
59     else {
60     vtx = bestPV;
61     return true;
62     }
63     }
64     //----------------------------------------------------------------------------------------
65     void writeEntries(ControlFlags ctrl, unsigned total_unskimmed)
66     {
67     ofstream ofs;
68     TString ofname(ctrl.outputfile);
69     ofs.open(ofname.ReplaceAll(".root",".nevents"));
70     assert(ofs.is_open());
71     ofs << total_unskimmed << endl;
72     ofs.close();
73     }
74     //----------------------------------------------------------------------------------------
75     void getEATargets(ControlFlags &ctrl, mithep::MuonTools::EMuonEffectiveAreaTarget &eraMu, mithep::ElectronTools::EElectronEffectiveAreaTarget &eraEle)
76     {
77     if( !ctrl.mc && ctrl.era == 2011 ) {
78     eraMu = mithep::MuonTools::kMuEAData2011;
79     eraEle = mithep::ElectronTools::kEleEAData2011;
80     } else if( !ctrl.mc && ctrl.era == 2012 ) {
81     eraMu = mithep::MuonTools::kMuEAData2012;
82     eraEle = mithep::ElectronTools::kEleEAData2012;
83     } else if( ctrl.mc && ctrl.era == 2011 ) {
84     eraMu = mithep::MuonTools::kMuEAFall11MC;
85     eraEle = mithep::ElectronTools::kEleEAFall11MC;
86     } else if( ctrl.mc && ctrl.era == 2012 ) {
87     eraMu = mithep::MuonTools::kMuEAData2012;
88     eraEle = mithep::ElectronTools::kEleEAData2012;
89     } else {
90     cerr << "unknown era for effective areas ... quitting." << endl;
91     exit(1);
92     }
93     }
94     //----------------------------------------------------------------------------
95 anlevin 1.10 unsigned makePFnoPUArray(const mithep::Array<PFCandidate> * fPFCandidates,
96 dkralph 1.1 vector<bool> &pfNoPileUpflag,
97     const mithep::Array<mithep::Vertex> * vtxArr )
98     //----------------------------------------------------------------------------
99     {
100     for(UInt_t i = 0; i < fPFCandidates->GetEntries(); i++) {
101     const PFCandidate *pf = fPFCandidates->At(i);
102     assert(pf);
103    
104     if(pf->PFType() == PFCandidate::eHadron) {
105     if(pf->HasTrackerTrk() &&
106     vtxArr->At(0)->HasTrack(pf->TrackerTrk()) &&
107     vtxArr->At(0)->TrackWeight(pf->TrackerTrk()) > 0) {
108    
109     pfNoPileUpflag.push_back(1);
110    
111     } else {
112    
113     Bool_t vertexFound = kFALSE;
114     const Vertex *closestVtx = 0;
115     Double_t dzmin = 10000;
116    
117     // loop over vertices
118     for(UInt_t j = 0; j < vtxArr->GetEntries(); j++) {
119     const Vertex *vtx = vtxArr->At(j);
120     assert(vtx);
121    
122     if(pf->HasTrackerTrk() &&
123     vtx->HasTrack(pf->TrackerTrk()) &&
124     vtx->TrackWeight(pf->TrackerTrk()) > 0) {
125     vertexFound = kTRUE;
126     closestVtx = vtx;
127     break;
128     }
129     Double_t dz = fabs(pf->SourceVertex().Z() - vtx->Z());
130     if(dz < dzmin) {
131     closestVtx = vtx;
132     dzmin = dz;
133     }
134     }
135    
136     // if(fCheckClosestZVertex) {
137     if(1) {
138     // Fallback: if track is not associated with any vertex,
139     // associate it with the vertex closest in z
140     if(vertexFound || closestVtx != vtxArr->At(0)) {
141     // pfPileUp->Add(pf);
142 dkralph 1.12 pfNoPileUpflag.push_back(0); // i.e., don't use it in the isolation sum
143 dkralph 1.1 } else {
144 dkralph 1.12 pfNoPileUpflag.push_back(1); // use it
145 dkralph 1.1 }
146     } else {
147     if(vertexFound && closestVtx != vtxArr->At(0)) {
148     // pfPileUp->Add(pf);
149     pfNoPileUpflag.push_back(0);
150     } else {
151     // PFCandidate * pfNoPileUp->AddNew(); // Ridiculous but that's how it is
152     pfNoPileUpflag.push_back(1);
153     }
154     }
155     } //hadron & trk stuff
156     } else { // hadron
157     // PFCandidate * ptr = pfNoPileUp->AddNew();
158     pfNoPileUpflag.push_back(1);
159     }
160     } // Loop over PF candidates
161    
162     return pfNoPileUpflag.size();
163     }
164 khahn 1.4 //--------------------------------------------------------------------------------
165 dkralph 1.12 void setEfficiencyWeights(ControlFlags &ctrl,
166     unsigned era,
167 khahn 1.8 EventData & evtdat,
168     std::bitset<1024> triggerBits,
169     mithep::TriggerTable *hltTable,
170     mithep::Array<mithep::TriggerObject> *hltObjArr,
171     mithep::TriggerObjectsTable *fTrigObjs,
172     WeightStruct & weights )
173 khahn 1.4 //--------------------------------------------------------------------------------
174 dkralph 1.1 {
175 dkralph 1.12 //
176     // offline weights
177     //
178 dkralph 1.1 vector<SimpleLepton> lepvec = evtdat.Z1leptons;
179     lepvec.insert(lepvec.end(), evtdat.Z2leptons.begin(), evtdat.Z2leptons.end());
180 khahn 1.8
181 dkralph 1.12 // get eff for each leg, for each lepton
182     vector< pair <double,double> > wlegs; // pair here is eff & err
183     for( int k=0; k<lepvec.size(); k++ ) {
184     double absEta = fabs(lepvec[k].vec.Eta());
185     double pt = lepvec[k].vec.Pt();
186     if( abs(lepvec[k].type) == 13 )
187     wlegs.push_back( muonPerLegOfflineEfficiencyWeight( ctrl, era, absEta, pt ) );
188     else
189     wlegs.push_back( elePerLegOfflineEfficiencyWeight( ctrl, era, absEta, pt, lepvec[k].isTight ) );
190     }
191     // combine weights and propagate uncertainties
192     pair<double,double> offpair = getOfflineEfficiencyWeight( wlegs );
193     weights.woff = offpair.first;
194     weights.werroff = offpair.second;
195    
196     //
197     // online weights
198     //
199 khahn 1.8 vector<SimpleLepton> muvec, elvec;
200     if(abs(evtdat.Z1leptons[0].type) == 11 ) {
201     elvec.push_back(evtdat.Z1leptons[0]);
202     elvec.push_back(evtdat.Z1leptons[1]);
203     } else {
204     muvec.push_back(evtdat.Z1leptons[0]);
205     muvec.push_back(evtdat.Z1leptons[1]);
206     }
207     if(abs(evtdat.Z2leptons[0].type) == 11 ) {
208     elvec.push_back(evtdat.Z2leptons[0]);
209     elvec.push_back(evtdat.Z2leptons[1]);
210     } else {
211     muvec.push_back(evtdat.Z2leptons[0]);
212     muvec.push_back(evtdat.Z2leptons[1]);
213     }
214 dkralph 1.12 pair<double,double> onpair = getOnlineEfficiencyWeight( ctrl, triggerBits, hltTable, hltObjArr, fTrigObjs, muvec, elvec );
215 dkralph 1.1 weights.won = onpair.first;
216     weights.werron = onpair.second;
217     }
218     //----------------------------------------------------------------------------
219     void initEvtRhoMap( map<unsigned, float> & evtrhoMap )
220     //----------------------------------------------------------------------------
221     {
222     unsigned evt;
223     float rho;
224    
225     cout << "initialzing EvtRhoMap ";
226     //FILE * f = fopen("evtrho.dat", "r");
227     // FILE * f = fopen("allrho.cali.uniq.dat", "r");
228     // FILE * f = fopen("/data/blue/khahn/CMSSW_4_4_2/src/MitHzz4l/allrhoAA.cali.uniq.dat", "r");
229     assert(fopen("/scratch/dkralph/khahn/CMSSW_4_4_2/src/MitHzz4l/allrhoAA.cali.uniq.dat","r")); FILE * f = fopen("/scratch/dkralph/khahn/CMSSW_4_4_2/src/MitHzz4l/allrhoAA.cali.uniq.dat", "r");
230     int lcount=0;
231     while( fscanf( f, "%u:%f", &evt, &rho ) ) {
232     if( feof(f) ) break;
233     evtrhoMap[evt] = rho;
234     lcount++;
235     if( !(lcount%1000) ) cout << "."; flush(cout);
236     }
237     cout << " done" << endl;
238     }
239    
240 khahn 1.5 //----------------------------------------------------------------------------------------
241 dkralph 1.1 void setEra(string fname, ControlFlags &ctrl)
242     {
243 khahn 1.5 //----------------------------------------------------------------------------------------
244 dkralph 1.1 if( ctrl.mc && (fname.find("f11-") != string::npos) ) ctrl.era=2011 ;
245 dkralph 1.9 else if( ctrl.mc && (fname.find("s12-") != string::npos) ) ctrl.era=2012 ;
246     else if( !ctrl.mc && (fname.find("r11") != string::npos) ) ctrl.era=2011 ;
247     else if( !ctrl.mc && (fname.find("r12") != string::npos) ) ctrl.era=2012 ;
248     else assert(0);
249 dkralph 1.1 }
250 khahn 1.5
251 dkralph 1.1 //----------------------------------------------------------------------------------------
252 dkralph 1.6 unsigned getNPU(mithep::Array<mithep::PileupInfo> *puArr, int bx)
253 khahn 1.5 //----------------------------------------------------------------------------------------
254 dkralph 1.1 {
255    
256     // Identify bunch crossing 0 in PU info
257     int ibx0=-1, ibxm=-1, ibxp=-1;
258     assert(puArr);
259     for(unsigned i=0; i<puArr->GetEntries(); ++i) {
260     if(puArr->At(i)->GetBunchCrossing() == 0) ibx0=i;
261     if(puArr->At(i)->GetBunchCrossing() ==-1) ibxm=i;
262     if(puArr->At(i)->GetBunchCrossing() == 1) ibxp=i;
263     }
264    
265     if(bx==0) return puArr->At(ibx0)->GetPU_NumInteractions();
266     else if(bx==-1) return puArr->At(ibxm)->GetPU_NumInteractions();
267     else if(bx== 1) return puArr->At(ibxp)->GetPU_NumInteractions();
268     else assert(0);
269    
270 dkralph 1.6 return 0;
271 dkralph 1.1 }
272 khahn 1.5
273    
274     //----------------------------------------------------------------------------
275     void initPUWeights()
276     //----------------------------------------------------------------------------
277     {
278 dkralph 1.12 // TFile * puf;
279 khahn 1.5
280 dkralph 1.12 // puf = new TFile("data/pileup/PUWeights_F11To2011.root");
281     // hpu_2011 = (TH1D*)(puf->Get("puWeights"));
282     // hpu_2011->SetDirectory(0);
283     // puf->Close();
284    
285     // puf = new TFile("data/pileup/PUWeights_S12To2012_190456-199011.root");
286     // hpu_2012 = (TH1D*)(puf->Get("puWeights"));
287     // hpu_2012->SetDirectory(0);
288     // puf->Close();
289 khahn 1.5 }
290    
291     //----------------------------------------------------------------------------
292 dkralph 1.12 double getPUWeight(int era, TString fname, unsigned npu)
293 khahn 1.5 //----------------------------------------------------------------------------
294     {
295 dkralph 1.12 if(npu > 59) npu = 59; // hardcoded crap in the pu functions
296    
297     if(era==2011)
298     return weightTrue2011(npu);
299     else assert(era==2012);
300    
301     // assume that if it's in /029/ that we want to reweigh to 53x data
302     TString vers,dataVers,mcVers;
303     if(fname.Contains("/028/")) dataVers = mcVers = vers = "52x";
304     else if(fname.Contains("/029/")) dataVers = mcVers = vers = "53x";
305     else {
306     // if(ctrl.debug)
307     // cout << "WARNING in getPUWeight(): couldn't figure out if this is 52x or 53x, so just using 53x" << endl;
308     dataVers = mcVers = vers = "53x";
309 khahn 1.7 }
310 dkralph 1.12
311     if(dataVers=="52x" && mcVers=="52x")
312     return weightTruePileupV07toIchep52X(npu);
313     else if(dataVers=="53x" && mcVers=="53x")
314     return weightTruePileupV10toHcp53X(npu);
315     // else if(dataVers=="53x" && mcVers=="52x") // umm... can't really do this a.t.m..., but this would be the proper function
316     // return weightTruePileupV07toHcp53X(npu);
317     else { cout << "Couldn't work out which npu to use!" << endl; assert(0); };
318    
319     return 1;
320 khahn 1.5 }
321 dkralph 1.12 //----------------------------------------------------------------------------------------
322     void initJets(ControlFlags ctrl,
323     string cmsswpath,
324     vector<TString> &fJetCorrParsv,
325     vector<JetCorrectorParameters> &correctionParameters,
326     FactorizedJetCorrector *&fJetCorrector,
327     JetIDMVA *&fJetIDMVA)
328     //----------------------------------------------------------------------------------------
329     {
330     //
331     // jet energy corrections
332     //
333     vector<TString> corrTypes;
334     corrTypes.push_back("L1FastJet");
335     corrTypes.push_back("L2Relative");
336     corrTypes.push_back("L3Absolute");
337     if(!ctrl.mc)
338     corrTypes.push_back("L2L3Residual");
339    
340     TString tag;
341     if(ctrl.era==2011) {
342     cout << "\n\nNOTE: dont have the right JEC tags yet!\n\n" << endl;
343     // // from h4l hcp twiki:
344     if(ctrl.mc) tag = "START42_V17"; //START44_V13
345     else tag = "GR_R_42_V23"; //GR_R_44_V15
346     } else if(ctrl.era==2012) {
347     cout << "\n\nNOTE: dont have the right JEC tags yet!\n\n" << endl;
348     // // from h4l hcp twiki:
349     if(ctrl.mc) tag = "START53_V7F"; //START53_V10
350     else tag = "GR_P_V41_AN2"; //GR_P_V41_AN1
351     assert(TString(cmsswpath).Contains("CMSSW_5_3"));
352     } else assert(0);
353    
354     TString dataPath = TString(cmsswpath) + "/MitPhysics/data/";
355     for(unsigned icorr=0; icorr<corrTypes.size(); icorr++) {
356     if(ctrl.era==2011) corrTypes[icorr] = "AK5PF_" + corrTypes[icorr];
357     else if(ctrl.era==2012) corrTypes[icorr] = corrTypes[icorr] + "_AK5PF";
358     else assert(0);
359     TString fname(dataPath + tag + "_" + corrTypes[icorr] + ".txt");
360     if( !fopen(fname, "r") ) { cout << " no such file! " << fname << endl; assert(0); }
361     addJetCorr(fname, fJetCorrParsv);
362     }
363 khahn 1.5
364 dkralph 1.12 for(unsigned icorr=0; icorr<fJetCorrParsv.size(); icorr++)
365     correctionParameters.push_back(JetCorrectorParameters(fJetCorrParsv[icorr].Data()));
366     fJetCorrector = new FactorizedJetCorrector(correctionParameters);
367    
368     //
369     // jet ID
370     //
371     fJetIDMVA = new JetIDMVA();
372     JetIDMVA::MVAType mvaType;
373     TString loPtXml(dataPath+"mva_JetID_lowpt.weights.xml");
374     TString cfi(cmsswpath+"/MitPhysics/Utils/python/JetIdParams_cfi.py");
375     TString Xml;
376     if(ctrl.era==2011) {
377     mvaType = JetIDMVA::k42;
378     Xml = dataPath+"TMVAClassification_PuJetIdOptMVA.weights.xml";
379     } else if(ctrl.era==2012) {
380     mvaType = JetIDMVA::k52;
381     Xml = dataPath+"TMVAClassification_5x_BDT_fullPlusRMS.weights.xml";
382     } else assert(0);
383    
384     fJetIDMVA->Initialize(JetIDMVA::kLoose, loPtXml, Xml, mvaType, cfi);
385     }
386     //----------------------------------------------------------------------------------------
387     int getGenSampleType(TString ofname)
388     {
389     if(ofname.Contains(TRegexp("-zz4[em]-"))) return ESampleType::kVV; // powheg ZZ
390     else if(ofname.Contains(TRegexp("-zz2e2m-"))) return ESampleType::kVV; // ggZZ
391     else if(ofname.Contains(TRegexp("-zz4l-gf"))) return -1; // damn ggZZ has differences so MC parsing code fails ESampleType::kVV; // ggZZ
392     else if(ofname.Contains(TRegexp("-zz2l2l-gf"))) return -1; // damn ggZZ has differences so MC parsing code fails ESampleType::kVV; // ggZZ
393     else if(ofname.Contains(TRegexp("-ggzz4l-"))) return -1; // damn ggZZ has differences so MC parsing code fails ESampleType::kVV; // ggZZ
394     else if(ofname.Contains(TRegexp("-ggzz2l2l-"))) return -1; // damn ggZZ has differences so MC parsing code fails ESampleType::kVV; // ggZZ
395     else if(ofname.Contains(TRegexp("-h[0-9][0-9]*zz4l-"))) return ESampleType::kHZZ; // powheg HZZ
396     else assert(0);
397     }
398     //----------------------------------------------------------------------------------------
399 anlevin 1.11 TString getChannel(int z1type, int z2type)
400     {
401     if( abs(z1type) == 11 && abs(z2type) == 11 ) return "4e";//0;
402     else if( abs(z1type) == 13 && abs(z2type) == 13 ) return "4m";//1;
403     else return "2e2m";//2;
404     }
405 dkralph 1.12 //----------------------------------------------------------------------------------------
406     // START: functions from Mangano
407     //----------------------------------------------------------------------------------------
408     float weightTruePileupV07toIchep52X(float input){
409     float w[60] = {
410     2.59113e-07,
411     1.77107e-05,
412     0.0263017,
413     1.66396,
414     0.16421,
415     0.754166,
416     2.84622,
417     5.57103,
418     8.66558,
419     11.5716,
420     11.8712,
421     12.8102,
422     10.3421,
423     8.91019,
424     7.614,
425     6.10397,
426     4.52745,
427     3.26031,
428     2.39558,
429     1.83297,
430     1.47821,
431     1.26728,
432     1.14716,
433     1.06707,
434     0.98066,
435     0.860877,
436     0.708281,
437     0.539789,
438     0.37652,
439     0.237298,
440     0.1338,
441     0.0671236,
442     0.0299236,
443     0.0118723,
444     0.00420968,
445     0.00134235,
446     0.000389563,
447     0.000104892,
448     2.69214e-05,
449     6.79674e-06,
450     1.73307e-06,
451     4.52553e-07,
452     1.21124e-07,
453     3.29924e-08,
454     9.10616e-09,
455     2.53998e-09,
456     7.16146e-10,
457     2.03786e-10,
458     5.84308e-11,
459     1.68192e-11,
460     4.8434e-12,
461     1.38959e-12,
462     3.96112e-13,
463     1.11358e-13,
464     3.17245e-14,
465     5.34916e-15,
466     0,
467     0,
468     0,
469     0};
470    
471     //return w[(int)floor(input+0.5)];
472     return w[(int)floor(input)];
473    
474     /*
475     TH1F h("boh","boh",60,0.,60.);
476     for(int k=0;k<60;k++){
477     h.SetBinContent(k+1,w[k]);
478     }
479     return h.GetBinContent(h.FindBin(input));
480     */
481     }
482    
483    
484     float weightTruePileupV07toHcp53X(float input){
485     float w[60] = {
486     0.0447136,
487     0.11785,
488     0.23825,
489     1.08447,
490     0.102575,
491     0.454605,
492     1.79761,
493     4.00271,
494     6.83281,
495     9.83701,
496     10.7966,
497     12.2356,
498     10.0247,
499     8.49395,
500     7.1125,
501     5.69527,
502     4.31256,
503     3.19305,
504     2.42035,
505     1.91666,
506     1.58485,
507     1.36297,
508     1.21166,
509     1.09466,
510     0.978941,
511     0.84653,
512     0.699235,
513     0.548996,
514     0.408673,
515     0.288194,
516     0.193367,
517     0.124653,
518     0.0781124,
519     0.0479268,
520     0.0287763,
521     0.0167744,
522     0.00941834,
523     0.00507877,
524     0.00264364,
525     0.00134612,
526     0.000682678,
527     0.000351412,
528     0.0001864,
529     0.00010259,
530     5.87818e-05,
531     3.5033e-05,
532     2.17116e-05,
533     1.39777e-05,
534     9.36123e-06,
535     6.53328e-06,
536     4.76598e-06,
537     3.64139e-06,
538     2.92018e-06,
539     2.4602e-06,
540     2.17291e-06,
541     2.01107e-06,
542     1.94392e-06,
543     1.9598e-06,
544     2.0583e-06,
545     2.24895e-06};
546    
547     return w[(int)floor(input)];
548     }
549    
550    
551     float weightTruePileupV10toIchep53X(float input){
552     float w[60] = {
553     2.35693e-06,
554     7.51928e-05,
555     0.0263529,
556     0.609947,
557     0.737917,
558     1.29365,
559     0.994503,
560     0.85454,
561     1.01559,
562     1.33243,
563     1.72454,
564     2.01264,
565     2.00573,
566     1.80333,
567     1.56328,
568     1.37452,
569     1.24753,
570     1.16481,
571     1.11738,
572     1.09701,
573     1.08843,
574     1.08796,
575     1.09768,
576     1.10763,
577     1.09328,
578     1.0339,
579     0.92408,
580     0.771537,
581     0.59283,
582     0.41266,
583     0.256892,
584     0.14188,
585     0.0692543,
586     0.029902,
587     0.0114564,
588     0.00391383,
589     0.00120625,
590     0.000341485,
591     9.09127e-05,
592     2.34008e-05,
593     5.95438e-06,
594     1.5122e-06,
595     3.82094e-07,
596     9.51794e-08,
597     2.32205e-08,
598     5.51698e-09,
599     1.27267e-09,
600     2.84346e-10,
601     6.12799e-11,
602     1.26731e-11,
603     2.50309e-12,
604     4.69797e-13,
605     8.35153e-14,
606     1.39452e-14,
607     2.24718e-15,
608     2.03841e-16,
609     0,
610     0,
611     0,
612     0};
613    
614     return w[(int)floor(input)];
615     }
616    
617    
618    
619     float weightTruePileupV10toHcp53X(float input){
620     float w[60] = {
621     0.409409,
622     0.527276,
623     0.39328,
624     0.507892,
625     0.48029,
626     0.787701,
627     0.632356,
628     0.618033,
629     0.806089,
630     1.14018,
631     1.5788,
632     1.93507,
633     1.957,
634     1.73004,
635     1.46737,
636     1.28278,
637     1.18189,
638     1.13388,
639     1.12578,
640     1.14415,
641     1.16048,
642     1.1618,
643     1.15318,
644     1.13405,
645     1.09239,
646     1.01915,
647     0.914837,
648     0.786744,
649     0.644879,
650     0.502039,
651     0.371688,
652     0.263586,
653     0.18067,
654     0.120472,
655     0.0780184,
656     0.0486113,
657     0.0289039,
658     0.0163367,
659     0.00879674,
660     0.00456046,
661     0.0023098,
662     0.00115977,
663     0.000583207,
664     0.000294815,
665     0.000149865,
666     7.62892e-05,
667     3.87537e-05,
668     1.96105e-05,
669     9.87744e-06,
670     4.95418e-06,
671     2.47913e-06,
672     1.23919e-06,
673     6.19751e-07,
674     3.10125e-07,
675     1.54934e-07,
676     7.71425e-08,
677     3.8182e-08,
678     1.87455e-08,
679     9.10765e-09,
680     9.19802e-09};
681     return w[(int)floor(input)];
682     }
683    
684    
685    
686    
687    
688    
689    
690    
691    
692     float weightTrue2011(float input){
693     if(input>50)
694     return 1;
695    
696    
697     float w[50];
698    
699    
700     w[0]= 0.212929;
701     w[1]= 0.0208114;
702     w[2]= 0.0584048;
703     w[3]= 0.538898;
704     w[4]= 1.357;
705     w[5]= 1.49913;
706     w[6]= 1.42247;
707     w[7]= 1.35904;
708     w[8]= 1.29946;
709     w[9]= 1.27925;
710     w[10]= 1.37845;
711     w[11]= 1.71246;
712     w[12]= 1.5291;
713     w[13]= 1.35234;
714     w[14]= 1.22215;
715     w[15]= 1.0155;
716     w[16]= 1.01137;
717     w[17]= 0.395465;
718     w[18]= 0.230984;
719     w[19]= 0.109883;
720     w[20]= 0.0433739;
721     w[21]= 0.0111497;
722     w[22]= 0.00408801;
723     w[23]= 0.00115678;
724     w[24]= 0.000365505;
725     w[25]= 0.000112391;
726     w[26]= 3.83894e-05;
727     w[27]= 1.60651e-05;
728     w[28]= 4.81412e-06;
729     w[29]= 1.39717e-06;
730     w[30]= 1.92368e-06;
731     w[31]= 4.10748e-06;
732     w[32]= 2.33157e-05;
733     w[33]= 4.0181e-05;
734     w[34]= 4.87786e-05;
735     w[35]= 0.00194128;
736     w[36]= 8.97414e-05;
737     w[37]= 1;
738     w[38]= 1;
739     w[39]= 0.000162709;
740     w[40]= 1;
741     w[41]= 1;
742     w[42]= 1;
743     w[43]= 1;
744     w[44]= 1;
745     w[45]= 1;
746     w[46]= 1;
747     w[47]= 1;
748     w[48]= 1;
749     w[49]= 1;
750    
751    
752     TH1F h("boh","boh",50,0.,50.);
753    
754     for(int k=0;k<50;k++){
755     h.SetBinContent(k+1,w[k]);
756     }
757    
758     return h.GetBinContent(h.FindBin(input));
759    
760     }
761     //----------------------------------------------------------------------------------------
762     float weightTrue2011to2012(float input){
763     if(input>50)
764     return 1;
765    
766     float w[50];
767    
768     w[0]= 0.000443112;
769     w[1]= 0.000248044;
770     w[2]= 0.000273111;
771     w[3]= 0.00109511;
772     w[4]= 0.00195699;
773     w[5]= 0.00480746;
774     w[6]= 0.027013;
775     w[7]= 0.074795;
776     w[8]= 0.166231;
777     w[9]= 0.309545;
778     w[10]= 0.577657;
779     w[11]= 1.12488;
780     w[12]= 1.36899;
781     w[13]= 1.56925;
782     w[14]= 1.89846;
783     w[15]= 2.20828;
784     w[16]= 3.14112;
785     w[17]= 1.87712;
786     w[18]= 1.97062;
787     w[19]= 2.07067;
788     w[20]= 2.17791;
789     w[21]= 1.7176;
790     w[22]= 2.10953;
791     w[23]= 2.0805;
792     w[24]= 2.29498;
793     w[25]= 2.42189;
794     w[26]= 2.80303;
795     w[27]= 3.94091;
796     w[28]= 3.67917;
797     w[29]= 2.26081;
798     w[30]= 2.99726;
799     w[31]= 3.76553;
800     w[32]= 11.285;
801     w[33]= 10.2781;
802     w[34]= 6.73407;
803     w[35]= 148.182;
804     w[36]= 3.88144;
805     w[37]= 1;
806     w[38]= 1;
807     w[39]= 1.48128;
808     w[40]= 1;
809     w[41]= 1;
810     w[42]= 1;
811     w[43]= 1;
812     w[44]= 1;
813     w[45]= 1;
814     w[46]= 1;
815     w[47]= 1;
816     w[48]= 1;
817     w[49]= 1;
818    
819    
820     TH1F h("boh","boh",50,0.,50.);
821    
822     for(int k=0;k<50;k++){
823     h.SetBinContent(k+1,w[k]);
824     }
825    
826     return h.GetBinContent(h.FindBin(input));
827    
828     }
829     //----------------------------------------------------------------------------------------
830     // END: functions from Mangano
831     //----------------------------------------------------------------------------------------