ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/VHbb/interface/VHbbNameSpace.h
(Generate patch)

Comparing UserCode/VHbb/interface/VHbbNameSpace.h (file contents):
Revision 1.7 by nmohr, Thu Mar 21 17:25:17 2013 UTC vs.
Revision 1.13 by peller, Mon Apr 15 11:45:39 2013 UTC

# Line 3 | Line 3
3   #include "TVector2.h"
4   #include "TMath.h"
5   /*#if !defined(__CINT__) && !defined(__MAKECINT__)
6 < #include "CondFormats/JetMETObjects/interface/JetCorrectorParameters.h"
7 < #include "CondFormats/JetMETObjects/interface/JetCorrectionUncertainty.h"
8 < #endif*/
6 >  #include "CondFormats/JetMETObjects/interface/JetCorrectorParameters.h"
7 >  #include "CondFormats/JetMETObjects/interface/JetCorrectionUncertainty.h"
8 >  #endif*/
9  
10   namespace VHbb {
11    
# Line 79 | Line 79 | namespace VHbb {
79    }
80    
81    double ANGLELZ(double pt, double eta, double phi, double mass, double pt2, double eta2, double phi2, double mass2){
82 <  TLorentzVector m1, m2, msum;
83 <  m1.SetPtEtaPhiM(pt, eta, phi, mass);
84 <  m2.SetPtEtaPhiM(pt2, eta2, phi2, mass2);
85 <  msum = m1 + m2;
82 >    TLorentzVector m1, m2, msum;
83 >    m1.SetPtEtaPhiM(pt, eta, phi, mass);
84 >    m2.SetPtEtaPhiM(pt2, eta2, phi2, mass2);
85 >    msum = m1 + m2;
86  
87 <  TVector3 bZ =  msum.BoostVector();
87 >    TVector3 bZ =  msum.BoostVector();
88  
89 <  m1.Boost(-bZ);  
90 <  m2.Boost(-bZ);  
89 >    m1.Boost(-bZ);  
90 >    m2.Boost(-bZ);  
91  
92 <  TVector3 b1;
92 >    TVector3 b1;
93  
94  
95 <  if((int) (pt) % 2 == 0)
96 <    b1 =  m1.BoostVector();
97 <  else
98 <    b1 =  m2.BoostVector();
95 >    if((int) (pt) % 2 == 0)
96 >      b1 =  m1.BoostVector();
97 >    else
98 >      b1 =  m2.BoostVector();
99  
100 < double cosTheta = b1.Dot(msum.BoostVector()) / (b1.Mag()*msum.BoostVector().Mag());
101 < return(cosTheta);
100 >    double cosTheta = b1.Dot(msum.BoostVector()) / (b1.Mag()*msum.BoostVector().Mag());
101 >    return(cosTheta);
102    }
103  
104  
105 <   double ANGLEHB(double pt, double eta, double phi, double e, double pt2, double eta2, double phi2, double e2){
106 < TLorentzVector m1, m2, msum;
107 < m1.SetPtEtaPhiE(pt, eta, phi, e);
108 < m2.SetPtEtaPhiE(pt2, eta2, phi2, e2);
109 < msum = m1 + m2;
105 >  double ANGLEHB(double pt, double eta, double phi, double e, double pt2, double eta2, double phi2, double e2){
106 >    TLorentzVector m1, m2, msum;
107 >    m1.SetPtEtaPhiE(pt, eta, phi, e);
108 >    m2.SetPtEtaPhiE(pt2, eta2, phi2, e2);
109 >    msum = m1 + m2;
110  
111 < TVector3 bZ =  msum.BoostVector();
111 >    TVector3 bZ =  msum.BoostVector();
112  
113 < m1.Boost(-bZ);
114 < m2.Boost(-bZ);  
113 >    m1.Boost(-bZ);
114 >    m2.Boost(-bZ);  
115  
116 < TVector3 b1;
116 >    TVector3 b1;
117  
118 <  if((int) (pt) % 2 == 0)
119 <    b1 =  m1.BoostVector();
120 <  else
121 <    b1 =  m2.BoostVector();
118 >    if((int) (pt) % 2 == 0)
119 >      b1 =  m1.BoostVector();
120 >    else
121 >      b1 =  m2.BoostVector();
122  
123 < double cosTheta = b1.Dot(msum.BoostVector()) / (b1.Mag()*msum.BoostVector().Mag());
124 < return(cosTheta);
125 <   }
123 >    double cosTheta = b1.Dot(msum.BoostVector()) / (b1.Mag()*msum.BoostVector().Mag());
124 >    return(cosTheta);
125 >  }
126  
127 <   double metCorSysShift(double met, double metphi, int Nvtx, int EVENT_run)
128 < {
127 >  double metCorSysShift(double met, double metphi, int Nvtx, int EVENT_run)
128 >  {
129      double metx = met * cos(metphi);
130      double mety = met * sin(metphi);
131      double px = 0.0, py = 0.0;
132      if (EVENT_run!=1) {
133 <        //pfMEtSysShiftCorrParameters_2012runAplusBvsNvtx_data
134 <        px = +1.68804e-01 + 3.37139e-01*Nvtx;
135 <        py = -1.72555e-01 - 1.79594e-01*Nvtx;
133 >      //pfMEtSysShiftCorrParameters_2012runAplusBvsNvtx_data
134 >      px = +1.68804e-01 + 3.37139e-01*Nvtx;
135 >      py = -1.72555e-01 - 1.79594e-01*Nvtx;
136      } else {
137 <        //pfMEtSysShiftCorrParameters_2012runAplusBvsNvtx_mc
138 <        px = +2.22335e-02 - 6.59183e-02*Nvtx;
139 <        py = +1.52720e-01 - 1.28052e-01*Nvtx;
137 >      //pfMEtSysShiftCorrParameters_2012runAplusBvsNvtx_mc
138 >      px = +2.22335e-02 - 6.59183e-02*Nvtx;
139 >      py = +1.52720e-01 - 1.28052e-01*Nvtx;
140      }
141      metx -= px;
142      mety -= py;
143      return std::sqrt(metx*metx + mety*mety);
144 < }
144 >  }
145  
146 <    double metphiCorSysShift(double met, double metphi, int Nvtx, int EVENT_run)
147 < {
146 >  double metphiCorSysShift(double met, double metphi, int Nvtx, int EVENT_run)
147 >  {
148      double metx = met * cos(metphi);
149      double mety = met * sin(metphi);
150      double px = 0.0, py = 0.0;
151      if (EVENT_run!=1) {
152  
153 <        //pfMEtSysShiftCorrParameters_2012runAplusBvsNvtx_data
154 <        px = +1.68804e-01 + 3.37139e-01*Nvtx;
155 <        py = -1.72555e-01 - 1.79594e-01*Nvtx;
153 >      //pfMEtSysShiftCorrParameters_2012runAplusBvsNvtx_data
154 >      px = +1.68804e-01 + 3.37139e-01*Nvtx;
155 >      py = -1.72555e-01 - 1.79594e-01*Nvtx;
156      } else {
157 <        //pfMEtSysShiftCorrParameters_2012runAplusBvsNvtx_mc
158 <        px = +2.22335e-02 - 6.59183e-02*Nvtx;
159 <        py = +1.52720e-01 - 1.28052e-01*Nvtx;
157 >      //pfMEtSysShiftCorrParameters_2012runAplusBvsNvtx_mc
158 >      px = +2.22335e-02 - 6.59183e-02*Nvtx;
159 >      py = +1.52720e-01 - 1.28052e-01*Nvtx;
160      }
161      metx -= px;
162      mety -= py;
163      if (metx == 0.0 && mety == 0.0)
164 <        return 0.0;
164 >      return 0.0;
165  
166      double phi1 = std::atan2(mety,metx);
167      double phi2 = std::atan2(mety,metx)-2.0*M_PI;
168      if (std::abs(phi1-metphi) < std::abs(phi2-metphi)+0.5*M_PI)
169 <        return phi1;
169 >      return phi1;
170      else
171 <        return phi2;
172 < }
171 >      return phi2;
172 >  }
173  
174 < TVector2 metType1Reg(double met, double metphi, double corr1, double corr2, double pt1, double eta1, double phi1, double e1, double pt2, double eta2, double phi2, double e2)
175 < {
174 >  TVector2 metType1Reg(double met, double metphi, double corr1, double corr2, double pt1, double eta1, double phi1, double e1, double pt2, double eta2, double phi2, double e2)
175 >  {
176      double metx = met * cos(metphi);
177      double mety = met * sin(metphi);
178      TLorentzVector j1;
# Line 184 | Line 184 | TVector2 metType1Reg(double met, double
184      mety += j1.Py()*(1-corr1);
185      mety += j2.Py()*(1-corr2);
186      TVector2 corrMET(metx, mety);
187 <     return corrMET;
188 < }
187 >    return corrMET;
188 >  }
189  
190 < double metType1Phi(double met, double metphi, double corr1, double corr2, double pt1, double eta1, double phi1, double e1, double pt2, double eta2, double phi2, double e2){
190 >  double metType1Phi(double met, double metphi, double corr1, double corr2, double pt1, double eta1, double phi1, double e1, double pt2, double eta2, double phi2, double e2){
191      return metType1Reg(met, metphi, corr1, corr2, pt1, eta1, phi1, e1, pt2, eta2, phi2, e2).Phi();
192  
193 < }
194 < double metType1Et(double met, double metphi, double corr1, double corr2, double pt1, double eta1, double phi1, double e1, double pt2, double eta2, double phi2, double e2){
193 >  }
194 >  double metType1Et(double met, double metphi, double corr1, double corr2, double pt1, double eta1, double phi1, double e1, double pt2, double eta2, double phi2, double e2){
195      return metType1Reg(met, metphi, corr1, corr2, pt1, eta1, phi1, e1, pt2, eta2, phi2, e2).Mod();
196  
197 < }
197 >  }
198  
199  
200 <   double met_MPF(double met, double metphi, double pt, double phi)
201 < {
200 >  double met_MPF(double met, double metphi, double pt, double phi)
201 >  {
202      return 1.+met*pt*std::cos( deltaPhi(metphi,phi) ) / (pt*pt);
203  
204 < }
204 >  }
205  
206 < double resolutionBias(double eta)
207 < {
208 < // return 0;//Nominal!
209 <  if(eta< 0.5) return 0.052;
210 <  if(eta< 1.1) return 0.057;
211 <  if(eta< 1.7) return 0.096;
212 <  if(eta< 2.3) return 0.134;
213 <  if(eta< 5) return 0.28;
214 <  return 0;
215 < }
206 >  double resolutionBias(double eta)
207 >  {
208 >    // return 0;//Nominal!
209 >    if(eta< 0.5) return 0.052;
210 >    if(eta< 1.1) return 0.057;
211 >    if(eta< 1.7) return 0.096;
212 >    if(eta< 2.3) return 0.134;
213 >    if(eta< 5) return 0.28;
214 >    return 0;
215 >  }
216  
217 < double evalJERBias( double ptreco, double ptgen, double eta1){
218 <  double eta = fabs(eta1);
219 <  double cor =1;  
220 <  if ((fabs(ptreco - ptgen)/ ptreco)<0.5) { //Limit the effect to the core
221 <     cor = (ptreco +resolutionBias(eta) *(ptreco-ptgen))/ptreco;  
217 >  double evalJERBias( double ptreco, double ptgen, double eta1){
218 >    double eta = fabs(eta1);
219 >    double cor =1;  
220 >    if ((fabs(ptreco - ptgen)/ ptreco)<0.5) { //Limit the effect to the core
221 >      cor = (ptreco +resolutionBias(eta) *(ptreco-ptgen))/ptreco;  
222 >    }
223 >    if (ptgen > 0.) return ptreco*cor;
224 >    else return ptreco;
225    }
223  if (ptgen > 0.) return ptreco*cor;
224  else return ptreco;
225 }
226  
227 < double evalEt( double pt, double eta, double phi, double e){
228 <  TLorentzVector j;
229 <  j.SetPtEtaPhiE(pt,eta,phi, e );
230 <  return j.Et();
227 >  double evalEt( double pt, double eta, double phi, double e){
228 >    TLorentzVector j;
229 >    j.SetPtEtaPhiE(pt,eta,phi, e );
230 >    return j.Et();
231  
232 < }
232 >  }
233  
234 < double evalMt( double pt, double eta, double phi, double e){
235 <  TLorentzVector j;
236 <  j.SetPtEtaPhiE(pt,eta,phi, e );
237 <  return j.Mt();
234 >  double evalMt( double pt, double eta, double phi, double e){
235 >    TLorentzVector j;
236 >    j.SetPtEtaPhiE(pt,eta,phi, e );
237 >    return j.Mt();
238  
239 < }
240 < /*double evalJECUnc( double pt, double eta){
241 < // Total uncertainty for reference
242 < JetCorrectionUncertainty *total = new JetCorrectionUncertainty("/shome/nmohr/CMSSW_5_2_6_patch1/src/UserCode/VHbb/data/START53_V15MC_Uncertainty_AK5PFchs.txt");
243 <
244 < total->setJetPt(pt);
245 < total->setJetEta(eta);
246 < double uncert =  total->getUncertainty(true);
247 < delete total;
248 < return uncert;
249 < }*/
239 >  }
240 >  /*double evalJECUnc( double pt, double eta){
241 >  // Total uncertainty for reference
242 >  JetCorrectionUncertainty *total = new JetCorrectionUncertainty("/shome/nmohr/CMSSW_5_2_6_patch1/src/UserCode/VHbb/data/START53_V15MC_Uncertainty_AK5PFchs.txt");
243 >
244 >  total->setJetPt(pt);
245 >  total->setJetEta(eta);
246 >  double uncert =  total->getUncertainty(true);
247 >  delete total;
248 >  return uncert;
249 >  }*/
250  
251 < double ptWeightDY( double lheV_pt)
252 < {
251 >  double ptWeightDY( double lheV_pt, double sign = 1.)
252 >  {
253      double SF = 1.;
254      if (50. < lheV_pt && lheV_pt < 100.){
255 <        SF = 0.873885+0.00175853*lheV_pt;
255 >      SF = 0.873885+0.00175853*sign*lheV_pt;
256      }
257      else if (lheV_pt > 100){
258 <        SF = 1.10651-0.000705265*lheV_pt;
258 >      SF = 1.10651-0.000705265*sign*lheV_pt;
259      }
260      return SF;
261 < }
261 >  }
262 >
263 >  // weights correction for EWK NLO correction
264 >  double ptWeightZllH( double genHPt){
265 >    double SF = 1.;
266 >    if ( genHPt > 1.) SF = 0.999757-0.000174527*genHPt;
267 >    return SF>0?SF:0;
268 >  }
269 >
270  
271 +  double minCSVold(int EVENT_run, double hJet_csv0, double hJet_csv1, double hJet_csvOld0, double hJet_csvOld1){
272  
273 +    if( EVENT_run < 2 )
274 +      return std::min(hJet_csvOld0,hJet_csvOld1);
275 +    else
276 +      return std::min(hJet_csv0,hJet_csv1);
277 +  }
278 +
279 +  double maxCSVold(int EVENT_run, double hJet_csv0, double hJet_csv1, double hJet_csvOld0, double hJet_csvOld1){
280 +
281 +    if( EVENT_run < 2 )
282 +      return std::max(hJet_csvOld0,hJet_csvOld1);
283 +    else
284 +      return std::max(hJet_csv0,hJet_csv1);
285 +  }
286 +
287 +
288 +  // weights correction for EWK NLO correction from Atlas
289 +  double ewkAtlas8TeVZllH( double genHPt, double genZPt){
290 +    if (genHPt < 1.) return 1.;
291 +    double hll8_contents[95] = {0.000664024, -0.00357095, -0.00767076, -0.00967366, -0.0134844, -0.0157148, -0.0181885, -0.0209647, -0.0232788, -0.0252373, -0.0265634, -0.0275069, -0.0285776, -0.0281683, -0.0294206, -0.0299975, -0.0308047, -0.0311716, -0.030913, -0.0324821, -0.0323192, -0.0324639, -0.0319356, -0.0322621, -0.0331146, -0.0338905, -0.0345189, -0.0358591, -0.0358407, -0.040018, -0.0396389, -0.0407177, -0.0445103, -0.0441406, -0.0471215, -0.0463301, -0.0513777, -0.0536773, -0.0546446, -0.0568508, -0.0590333, -0.0612157, -0.0633981, -0.0655805, -0.067763, -0.0699454, -0.0721278, -0.0743103, -0.0764927, -0.0786751, -0.0808575, -0.08304, -0.0852224, -0.0874048, -0.0895872, -0.0917697, -0.0939521, -0.0961345, -0.098317, -0.100499, -0.102682, -0.104864, -0.107047, -0.109229, -0.111412, -0.113594, -0.115776, -0.117959, -0.120141, -0.122324, -0.124506, -0.126689, -0.128871, -0.131053, -0.133236, -0.135418, -0.137601, -0.139783, -0.141965, -0.144148, -0.14633, -0.148513, -0.150695, -0.152878, -0.15506, -0.157242, -0.159425, -0.161607, -0.16379, -0.165972, -0.168155, -0.170337, -0.172519, -0.174702, -0.176884};
292 +    int corrBin = (int) ( (genZPt - 25.) / 5. );
293 +    if (corrBin < 0) corrBin = 0;
294 +    if (corrBin > 94) corrBin = 94;
295 +    double SF = 1.+hll8_contents[corrBin];
296 +    return SF;
297 +  }
298 +  
299 +  double mueEff(int Vtype, double eta0, double eta1, double pt0, double pt1){
300 +      if (Vtype == 0) return 1.;
301 +      if (Vtype == 1) {
302 +            double corr = 1.0;
303 +            if (abs(eta0) < 1.3) corr*=0.95;
304 +            if (abs(eta1) < 1.3) corr*=0.95;
305 +            if (min(pt0,pt1) < 25.) corr*=0.6;
306 +            else if (min(pt0,pt1) < 30.) corr*=0.8;
307 +            // scale to lumi ABC/D
308 +            corr = 1.-(1.-corr)*6.9/12.1;
309 +            return corr;
310 +      }
311 +      return 1.;
312 +  }
313  
314   }
315  

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines