AnaTools/src/SFWeight.cc

#include "OSUT3Analysis/AnaTools/interface/SFWeight.h"

MuonSFWeight::MuonSFWeight (const string &sfFile, const string &dataOverMC)
{

  TFile *fin = TFile::Open (sfFile.c_str ());
  TGraphAsymmErrors *data = (TGraphAsymmErrors *) fin->Get (dataOverMC.c_str ());

  double bins[data->GetN () + 1];
  for (int i = 0; i < data->GetN (); i++)
    {
      double x, y, exl, exh;

      data->GetPoint (i, x, y);
      exl = data->GetErrorXlow (i);
      exh = data->GetErrorXhigh (i);
      if (i == 0)
        bins[i] = x - exl;
      bins[i + 1] = x + exh;
    }

  muonSFWeight_ = new TH1D ("muonSFWeight", "", data->GetN (), bins);
  muonSFWeight_->SetDirectory (0);
  for (int i = 0; i < data->GetN (); i++)
    {
      double dataX, dataY;

      data->GetPoint (i, dataX, dataY);
      muonSFWeight_->SetBinContent (i + 1, dataY);
    }

  fin->Close ();
  delete data;
}

MuonSFWeight::MuonSFWeight (const string &sfFile, const string &dataEff, const string &mcEff)
{
  TFile *fin = TFile::Open (sfFile.c_str ());
  TGraphAsymmErrors *mc = (TGraphAsymmErrors *) fin->Get (mcEff.c_str ());
  TGraphAsymmErrors *data = (TGraphAsymmErrors *) fin->Get (dataEff.c_str ());

  double bins[mc->GetN () + 1];
  for (int i = 0; i < mc->GetN (); i++)
    {
      double x, y, exl, exh;

      mc->GetPoint (i, x, y);
      exl = mc->GetErrorXlow (i);
      exh = mc->GetErrorXhigh (i);
      if (i == 0)
        bins[i] = x - exl;
      bins[i + 1] = x + exh;
    }

  muonSFWeight_ = new TH1D ("muonSFWeight", "", mc->GetN (), bins);
  muonSFWeight_->SetDirectory (0);
  for (int i = 0; i < mc->GetN (); i++)
    {
      double mcX, mcY, dataX, dataY;

      mc->GetPoint (i, mcX, mcY);
      data->GetPoint (i, dataX, dataY);
      muonSFWeight_->SetBinContent (i + 1, dataY / mcY);
    }

  fin->Close ();
  delete mc;
  delete data;
}

MuonSFWeight::~MuonSFWeight ()
{
  delete muonSFWeight_;
}

ElectronSFWeight::ElectronSFWeight (const string &cmsswRelease, const string &id) :
  cmsswRelease_ (cmsswRelease),
  id_ (id)
{
}

double
ElectronSFWeight::at (const double &eta, const double &pt)
{
  //cout << eta << " " << pt << endl;
  //cout << "doing it\n";
  if (cmsswRelease_ == "53X")
    {
      if (id_ == "loose")
        {
          if (fabs (eta) < 0.8)
            {
              if (pt < 15.0 && pt >= 10.0)
                return 0.855;
              if (pt < 20.0 && pt >= 15.0)
                return 0.962;
              if (pt < 30.0 && pt >= 20.0)
                return 1.005;
              if (pt < 40.0 && pt >= 30.0)
                return 1.004;
              if (pt < 50.0 && pt >= 40.0)
                return 1.008;
              if (pt < 200.0 && pt >= 50.0)
                return 1.008;
            }
          if (fabs (eta) < 1.442 && fabs (eta) >= 0.8)
            {
              if (pt < 15.0 && pt >= 10.0)
                return 0.858;
              if (pt < 20.0 && pt >= 15.0)
                return 0.962;
              if (pt < 30.0 && pt >= 20.0)
                return 0.981;
              if (pt < 40.0 && pt >= 30.0)
                return 0.991;
              if (pt < 50.0 && pt >= 40.0)
                return 0.994;
              if (pt < 200.0 && pt >= 50.0)
                return 0.999;
            }
          if (fabs (eta) < 1.556 && fabs (eta) >= 1.442)
            {
              if (pt < 15.0 && pt >= 10.0)
                return 1.109;
              if (pt < 20.0 && pt >= 15.0)
                return 0.903;
              if (pt < 30.0 && pt >= 20.0)
                return 1.044;
              if (pt < 40.0 && pt >= 30.0)
                return 0.998;
              if (pt < 50.0 && pt >= 40.0)
                return 0.989;
              if (pt < 200.0 && pt >= 50.0)
                return 0.994;
            }
          if (fabs (eta) < 2.0 && fabs (eta) >= 1.556)
            {
              if (pt < 15.0 && pt >= 10.0)
                return 0.838;
              if (pt < 20.0 && pt >= 15.0)
                return 0.939;
              if (pt < 30.0 && pt >= 20.0)
                return 0.980;
              if (pt < 40.0 && pt >= 30.0)
                return 0.992;
              if (pt < 50.0 && pt >= 40.0)
                return 1.004;
              if (pt < 200.0 && pt >= 50.0)
                return 1.006;
            }
          if (fabs (eta) < 2.5 && fabs (eta) >= 2.0)
            {
              if (pt < 15.0 && pt >= 10.0)
                return 1.034;
              if (pt < 20.0 && pt >= 15.0)
                return 0.970;
              if (pt < 30.0 && pt >= 20.0)
                return 1.017;
              if (pt < 40.0 && pt >= 30.0)
                return 1.019;
              if (pt < 50.0 && pt >= 40.0)
                return 1.005;
              if (pt < 200.0 && pt >= 50.0)
                return 1.009;
            }
        }
      if (id_ == "tight")
        {
          if (fabs (eta) < 0.8)
            {
              if (pt < 15.0 && pt >= 10.0)
                return 0.818;
              if (pt < 20.0 && pt >= 15.0)
                return 0.928;
              if (pt < 30.0 && pt >= 20.0)
                return 0.973;
              if (pt < 40.0 && pt >= 30.0)
                return 0.979;
              if (pt < 50.0 && pt >= 40.0)
                return 0.984;
              if (pt < 200.0 && pt >= 50.0)
                return 0.983;
            }
          if (fabs (eta) < 1.442 && fabs (eta) >= 0.8)
            {
              if (pt < 15.0 && pt >= 10.0)
                return 0.840;
              if (pt < 20.0 && pt >= 15.0)
                return 0.914;
              if (pt < 30.0 && pt >= 20.0)
                return 0.948;
              if (pt < 40.0 && pt >= 30.0)
                return 0.961;
              if (pt < 50.0 && pt >= 40.0)
                return 0.972;
              if (pt < 200.0 && pt >= 50.0)
                return 0.977;
            }
          if (fabs (eta) < 1.556 && fabs (eta) >= 1.442)
            {
              if (pt < 15.0 && pt >= 10.0)
                return 1.008;
              if (pt < 20.0 && pt >= 15.0)
                return 0.877;
              if (pt < 30.0 && pt >= 20.0)
                return 0.983;
              if (pt < 40.0 && pt >= 30.0)
                return 0.983;
              if (pt < 50.0 && pt >= 40.0)
                return 0.957;
              if (pt < 200.0 && pt >= 50.0)
                return 0.978;
            }
          if (fabs (eta) < 2.0 && fabs (eta) >= 1.556)
            {
              if (pt < 15.0 && pt >= 10.0)
                return 0.906;
              if (pt < 20.0 && pt >= 15.0)
                return 0.907;
              if (pt < 30.0 && pt >= 20.0)
                return 0.957;
              if (pt < 40.0 && pt >= 30.0)
                return 0.962;
              if (pt < 50.0 && pt >= 40.0)
                return 0.985;
              if (pt < 200.0 && pt >= 50.0)
                return 0.986;
            }
          if (fabs (eta) < 2.5 && fabs (eta) >= 2.0)
            {
              if (pt < 15.0 && pt >= 10.0)
                return 0.991;
              if (pt < 20.0 && pt >= 15.0)
                return 0.939;
              if (pt < 30.0 && pt >= 20.0)
                return 1.001;
              if (pt < 40.0 && pt >= 30.0)
                return 1.002;
              if (pt < 50.0 && pt >= 40.0)
                return 0.999;
              if (pt < 200.0 && pt >= 50.0)
                return 0.995;
            }
        }
      if (id_ == "mvaNonTrig")
        {
          //cout << "mvaNonTrig\n";
          if (fabs (eta) < 0.8)
            {
              if (pt < 20.0)
                return 1.0;
              if (pt < 30.0 && pt >= 20.0)
                return 0.994;
              if (pt < 40.0 && pt >= 30.0)
                return 0.974;
              if (pt < 50.0 && pt >= 40.0)
                return 0.988;
              if (pt < 200.0 && pt >= 50.0)
                return 0.999;
            }
          if (fabs (eta) < 1.478 && fabs (eta) >= 0.8)
            {
              if (pt < 20.0)
                return 1.0;
              if (pt < 30.0 && pt >= 20.0)
                return 0.971;
              if (pt < 40.0 && pt >= 30.0)
                return 0.967;
              if (pt < 50.0 && pt >= 40.0)
                return 0.984;
              if (pt < 200.0 && pt >= 50.0)
                return 0.998;
            }
          if (fabs (eta) <= 2.5 && fabs (eta) >= 1.478)
            {
              if (pt < 20.0)
                return 1.0;
              if (pt < 30.0 && pt >= 20.0)
                return 0.921;
              if (pt < 40.0 && pt >= 30.0)
                return 0.932;
              if (pt < 50.0 && pt >= 40.0)
                return 0.968;
              if (pt < 200.0 && pt >= 50.0)
                return 0.983;
            }
        }
    }

  return 0.0;
}

ElectronSFWeight::~ElectronSFWeight ()
{
}
Revision:	1.5
Committed:	Tue Jul 9 17:22:09 2013 UTC (11 years, 10 months ago) by ahart
Content type:	text/plain
Branch:	MAIN
CVS Tags:	V02-03-02
Changes since 1.4:	+78 -0 lines
Log Message:	Added official electron scale factors for loose ID.
#	User	Rev	Content
1	ahart	1.1	#include "OSUT3Analysis/AnaTools/interface/SFWeight.h"
2
3			MuonSFWeight::MuonSFWeight (const string &sfFile, const string &dataOverMC)
4			{
5	lantonel	1.2
6	ahart	1.1	TFile *fin = TFile::Open (sfFile.c_str ());
7			TGraphAsymmErrors data = (TGraphAsymmErrors ) fin->Get (dataOverMC.c_str ());
8
9			double bins[data->GetN () + 1];
10			for (int i = 0; i < data->GetN (); i++)
11			{
12			double x, y, exl, exh;
13
14			data->GetPoint (i, x, y);
15			exl = data->GetErrorXlow (i);
16			exh = data->GetErrorXhigh (i);
17			if (i == 0)
18			bins[i] = x - exl;
19			bins[i + 1] = x + exh;
20			}
21
22			muonSFWeight_ = new TH1D ("muonSFWeight", "", data->GetN (), bins);
23			muonSFWeight_->SetDirectory (0);
24			for (int i = 0; i < data->GetN (); i++)
25			{
26			double dataX, dataY;
27
28			data->GetPoint (i, dataX, dataY);
29			muonSFWeight_->SetBinContent (i + 1, dataY);
30			}
31
32			fin->Close ();
33			delete data;
34			}
35
36			MuonSFWeight::MuonSFWeight (const string &sfFile, const string &dataEff, const string &mcEff)
37			{
38			TFile *fin = TFile::Open (sfFile.c_str ());
39			TGraphAsymmErrors mc = (TGraphAsymmErrors ) fin->Get (mcEff.c_str ());
40			TGraphAsymmErrors data = (TGraphAsymmErrors ) fin->Get (dataEff.c_str ());
41
42			double bins[mc->GetN () + 1];
43			for (int i = 0; i < mc->GetN (); i++)
44			{
45			double x, y, exl, exh;
46
47			mc->GetPoint (i, x, y);
48			exl = mc->GetErrorXlow (i);
49			exh = mc->GetErrorXhigh (i);
50			if (i == 0)
51			bins[i] = x - exl;
52			bins[i + 1] = x + exh;
53			}
54
55			muonSFWeight_ = new TH1D ("muonSFWeight", "", mc->GetN (), bins);
56			muonSFWeight_->SetDirectory (0);
57			for (int i = 0; i < mc->GetN (); i++)
58			{
59			double mcX, mcY, dataX, dataY;
60
61			mc->GetPoint (i, mcX, mcY);
62			data->GetPoint (i, dataX, dataY);
63			muonSFWeight_->SetBinContent (i + 1, dataY / mcY);
64			}
65
66			fin->Close ();
67			delete mc;
68			delete data;
69			}
70
71			MuonSFWeight::~MuonSFWeight ()
72			{
73			delete muonSFWeight_;
74			}
75
76			ElectronSFWeight::ElectronSFWeight (const string &cmsswRelease, const string &id) :
77			cmsswRelease_ (cmsswRelease),
78			id_ (id)
79			{
80			}
81
82			double
83			ElectronSFWeight::at (const double &eta, const double &pt)
84			{
85	ahart	1.3	//cout << eta << " " << pt << endl;
86			//cout << "doing it\n";
87	ahart	1.1	if (cmsswRelease_ == "53X")
88			{
89	ahart	1.5	if (id_ == "loose")
90			{
91			if (fabs (eta) < 0.8)
92			{
93			if (pt < 15.0 && pt >= 10.0)
94			return 0.855;
95			if (pt < 20.0 && pt >= 15.0)
96			return 0.962;
97			if (pt < 30.0 && pt >= 20.0)
98			return 1.005;
99			if (pt < 40.0 && pt >= 30.0)
100			return 1.004;
101			if (pt < 50.0 && pt >= 40.0)
102			return 1.008;
103			if (pt < 200.0 && pt >= 50.0)
104			return 1.008;
105			}
106			if (fabs (eta) < 1.442 && fabs (eta) >= 0.8)
107			{
108			if (pt < 15.0 && pt >= 10.0)
109			return 0.858;
110			if (pt < 20.0 && pt >= 15.0)
111			return 0.962;
112			if (pt < 30.0 && pt >= 20.0)
113			return 0.981;
114			if (pt < 40.0 && pt >= 30.0)
115			return 0.991;
116			if (pt < 50.0 && pt >= 40.0)
117			return 0.994;
118			if (pt < 200.0 && pt >= 50.0)
119			return 0.999;
120			}
121			if (fabs (eta) < 1.556 && fabs (eta) >= 1.442)
122			{
123			if (pt < 15.0 && pt >= 10.0)
124			return 1.109;
125			if (pt < 20.0 && pt >= 15.0)
126			return 0.903;
127			if (pt < 30.0 && pt >= 20.0)
128			return 1.044;
129			if (pt < 40.0 && pt >= 30.0)
130			return 0.998;
131			if (pt < 50.0 && pt >= 40.0)
132			return 0.989;
133			if (pt < 200.0 && pt >= 50.0)
134			return 0.994;
135			}
136			if (fabs (eta) < 2.0 && fabs (eta) >= 1.556)
137			{
138			if (pt < 15.0 && pt >= 10.0)
139			return 0.838;
140			if (pt < 20.0 && pt >= 15.0)
141			return 0.939;
142			if (pt < 30.0 && pt >= 20.0)
143			return 0.980;
144			if (pt < 40.0 && pt >= 30.0)
145			return 0.992;
146			if (pt < 50.0 && pt >= 40.0)
147			return 1.004;
148			if (pt < 200.0 && pt >= 50.0)
149			return 1.006;
150			}
151			if (fabs (eta) < 2.5 && fabs (eta) >= 2.0)
152			{
153			if (pt < 15.0 && pt >= 10.0)
154			return 1.034;
155			if (pt < 20.0 && pt >= 15.0)
156			return 0.970;
157			if (pt < 30.0 && pt >= 20.0)
158			return 1.017;
159			if (pt < 40.0 && pt >= 30.0)
160			return 1.019;
161			if (pt < 50.0 && pt >= 40.0)
162			return 1.005;
163			if (pt < 200.0 && pt >= 50.0)
164			return 1.009;
165			}
166			}
167	ahart	1.1	if (id_ == "tight")
168			{
169			if (fabs (eta) < 0.8)
170			{
171			if (pt < 15.0 && pt >= 10.0)
172			return 0.818;
173			if (pt < 20.0 && pt >= 15.0)
174			return 0.928;
175			if (pt < 30.0 && pt >= 20.0)
176			return 0.973;
177			if (pt < 40.0 && pt >= 30.0)
178			return 0.979;
179			if (pt < 50.0 && pt >= 40.0)
180			return 0.984;
181			if (pt < 200.0 && pt >= 50.0)
182			return 0.983;
183			}
184			if (fabs (eta) < 1.442 && fabs (eta) >= 0.8)
185			{
186			if (pt < 15.0 && pt >= 10.0)
187			return 0.840;
188			if (pt < 20.0 && pt >= 15.0)
189			return 0.914;
190			if (pt < 30.0 && pt >= 20.0)
191			return 0.948;
192			if (pt < 40.0 && pt >= 30.0)
193			return 0.961;
194			if (pt < 50.0 && pt >= 40.0)
195			return 0.972;
196			if (pt < 200.0 && pt >= 50.0)
197			return 0.977;
198			}
199			if (fabs (eta) < 1.556 && fabs (eta) >= 1.442)
200			{
201			if (pt < 15.0 && pt >= 10.0)
202			return 1.008;
203			if (pt < 20.0 && pt >= 15.0)
204			return 0.877;
205			if (pt < 30.0 && pt >= 20.0)
206			return 0.983;
207			if (pt < 40.0 && pt >= 30.0)
208			return 0.983;
209			if (pt < 50.0 && pt >= 40.0)
210			return 0.957;
211			if (pt < 200.0 && pt >= 50.0)
212			return 0.978;
213			}
214			if (fabs (eta) < 2.0 && fabs (eta) >= 1.556)
215			{
216			if (pt < 15.0 && pt >= 10.0)
217			return 0.906;
218			if (pt < 20.0 && pt >= 15.0)
219			return 0.907;
220			if (pt < 30.0 && pt >= 20.0)
221			return 0.957;
222			if (pt < 40.0 && pt >= 30.0)
223			return 0.962;
224			if (pt < 50.0 && pt >= 40.0)
225			return 0.985;
226			if (pt < 200.0 && pt >= 50.0)
227			return 0.986;
228			}
229			if (fabs (eta) < 2.5 && fabs (eta) >= 2.0)
230			{
231			if (pt < 15.0 && pt >= 10.0)
232			return 0.991;
233			if (pt < 20.0 && pt >= 15.0)
234			return 0.939;
235			if (pt < 30.0 && pt >= 20.0)
236			return 1.001;
237			if (pt < 40.0 && pt >= 30.0)
238			return 1.002;
239			if (pt < 50.0 && pt >= 40.0)
240			return 0.999;
241			if (pt < 200.0 && pt >= 50.0)
242			return 0.995;
243			}
244			}
245	lantonel	1.2	if (id_ == "mvaNonTrig")
246			{
247	ahart	1.4	//cout << "mvaNonTrig\n";
248	lantonel	1.2	if (fabs (eta) < 0.8)
249			{
250			if (pt < 20.0)
251			return 1.0;
252			if (pt < 30.0 && pt >= 20.0)
253			return 0.994;
254			if (pt < 40.0 && pt >= 30.0)
255			return 0.974;
256			if (pt < 50.0 && pt >= 40.0)
257			return 0.988;
258			if (pt < 200.0 && pt >= 50.0)
259			return 0.999;
260			}
261			if (fabs (eta) < 1.478 && fabs (eta) >= 0.8)
262			{
263			if (pt < 20.0)
264			return 1.0;
265			if (pt < 30.0 && pt >= 20.0)
266			return 0.971;
267			if (pt < 40.0 && pt >= 30.0)
268			return 0.967;
269			if (pt < 50.0 && pt >= 40.0)
270			return 0.984;
271			if (pt < 200.0 && pt >= 50.0)
272			return 0.998;
273			}
274			if (fabs (eta) <= 2.5 && fabs (eta) >= 1.478)
275			{
276			if (pt < 20.0)
277			return 1.0;
278			if (pt < 30.0 && pt >= 20.0)
279			return 0.921;
280			if (pt < 40.0 && pt >= 30.0)
281			return 0.932;
282			if (pt < 50.0 && pt >= 40.0)
283			return 0.968;
284			if (pt < 200.0 && pt >= 50.0)
285			return 0.983;
286			}
287			}
288	ahart	1.1	}
289
290			return 0.0;
291			}
292
293			ElectronSFWeight::~ElectronSFWeight ()
294			{
295			}