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.3 by bortigno, Thu Jan 10 14:21:00 2013 UTC vs.
Revision 1.7 by nmohr, Thu Mar 21 17:25:17 2013 UTC

# Line 1 | Line 1
1   #include "TLorentzVector.h"
2   #include "TVector3.h"
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*/
9  
10   namespace VHbb {
11 +  
12  
13    double deltaPhi(double phi1,double phi2)
14    {
# Line 38 | Line 44 | namespace VHbb {
44      TVector3 n1(H1);
45      TVector3 n2(H2);
46      
47 <    float det= n1.Px() * n2.Py() - n2.Px() * n1.Py();
47 >    double det= n1.Px() * n2.Py() - n2.Px() * n1.Py();
48      
49      H1.SetMag( (  - n2.Py() * V.Px() + n2.Px() * V.Py() )  / (sin(n1.Theta()) *det ) );
50      H2.SetMag( ( + n1.Py() * V.Px() - n1.Px() * V.Py() )  / (sin(n2.Theta())  *det ) );
51      
52 <    float mass=TMath::Sqrt( TMath::Power( (H1.Mag()+H2.Mag()),2 ) - TMath::Power(( ( H1+H2 ).Mag()),2) );
52 >    double mass=TMath::Sqrt( TMath::Power( (H1.Mag()+H2.Mag()),2 ) - TMath::Power(( ( H1+H2 ).Mag()),2) );
53      
54      return mass;
55      
# Line 71 | Line 77 | namespace VHbb {
77  
78  
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;
86 +
87 +  TVector3 bZ =  msum.BoostVector();
88 +
89 +  m1.Boost(-bZ);  
90 +  m2.Boost(-bZ);  
91 +
92 +  TVector3 b1;
93 +
94 +
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);
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;
110 +
111 + TVector3 bZ =  msum.BoostVector();
112 +
113 + m1.Boost(-bZ);
114 + m2.Boost(-bZ);  
115 +
116 + TVector3 b1;
117 +
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 +   }
126 +
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;
136 +    } else {
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 + }
145 +
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;
156 +    } else {
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;
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;
170 +    else
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 + {
176 +    double metx = met * cos(metphi);
177 +    double mety = met * sin(metphi);
178 +    TLorentzVector j1;
179 +    TLorentzVector j2;
180 +    j1.SetPtEtaPhiE(pt1,eta1,phi1, e1 );
181 +    j2.SetPtEtaPhiE(pt2,eta2,phi2, e2 );
182 +    metx += j1.Px()*(1-corr1);
183 +    metx += j2.Px()*(1-corr2);
184 +    mety += j1.Py()*(1-corr1);
185 +    mety += j2.Py()*(1-corr2);
186 +    TVector2 corrMET(metx, mety);
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){
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){
195 +    return metType1Reg(met, metphi, corr1, corr2, pt1, eta1, phi1, e1, pt2, eta2, phi2, e2).Mod();
196 +
197 + }
198 +
199 +
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 + }
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 + }
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;  
222 +  }
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();
231 +
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();
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 + }*/
250 +
251 + double ptWeightDY( double lheV_pt)
252 + {
253 +    double SF = 1.;
254 +    if (50. < lheV_pt && lheV_pt < 100.){
255 +        SF = 0.873885+0.00175853*lheV_pt;
256 +    }
257 +    else if (lheV_pt > 100){
258 +        SF = 1.10651-0.000705265*lheV_pt;
259 +    }
260 +    return SF;
261 + }
262 +
263  
264  
265   }

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines