ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/Jeng/scripts/ikstest.cc

Comparing UserCode/Jeng/scripts/ikstest.cc (file contents):
Revision 1.3 by jengbou, Thu Sep 2 05:36:21 2010 UTC vs.
Revision 1.5 by jengbou, Tue Sep 7 03:59:59 2010 UTC

#	Line 1 | Line 1
1	<	#include <vector>
2	<	#include <TFile.h>
3	<	#include <TH1.h>
4	<	#include <TH2.h>
5	<	#include <TLatex.h>
6	<	#include <TCanvas.h>
7	<	#include <TMath.h>
8	<	#include <TLegend.h>
9	<	#include <TObject.h>
10	<	#include <iostream>
11	<	#include <sstream>
12	<	#include <fstream>
13	<	#include <iomanip>
14	<	#include <map>
1	>	#include "ikstest.h"
2	>
3	>	//=================================
4	>	// Program to run IKS test
5	>	// * Input directories:
6	>	//   * Data: skimmed_Data_x.xxpb-1
7	>	//      => default selection : Data_<suffix>.root
8	>	//      => n-1 selection     : Data_<suffix>_<invName>.root
9	>	//   * MC  : skimmed_MC {QCD,WJets,TTbar...}
10	>	//      => default selection : QCD_<suffix>.root
11	>	//      => n-1 selection     : Data_<suffix>_<invName>.root
12	>	// * Upmost output dir specified by <desDir>
13	>	//   subdirectory created according to <useInv> and <realData>
14	>	//      => <desDir'><xyz><suffix>.pdf (.txt)
15	>	// * e.g. to do KS test on data with MC QCD shape in signal region:
16	>	//   useInv=false; realData=true
17	>	//=================================
18
19		using namespace std;
20
21	<	// User defined parameters
22	<	bool useInv           = false; // whether to use n-1 QCD template
23	<	bool realData         = false;
24	<	const double stepsize = 0.001;
22	<
23	<	//define the constants: 2.78/pb
24	<	const double weight_[3]    = {0.0506107, //QCD
25	<	0.0088222, //WJets
26	<	0.000295   //TTbar
21	>	//define the constants: 2.88/pb
22	>	const double weight_[3]    = {0.0524313, //QCD
23	>	0.0091395, //WJets
24	>	0.000306   //TTbar
25		};
26		const Double_t procQCD     = 1.46;
27		const Double_t procWjets   = 1.03;
28		const Double_t procttjets  = 1.0;
29
30	<	struct testMC {
31	<	testMC(Double_t p = 0., Double_t sb = 0., Double_t ss = 0.){prob=p; scaleF_backg = sb; scaleF_sample = ss;}
32	<	Double_t prob;
33	<	Double_t scaleF_backg;
34	<	Double_t scaleF_sample;
35	<	};
36	<
37	<	//function for doing KS test
38	<	vector<testMC> doKStest(Double_t NsS, Double_t Ns1, Double_t Ns2, TH1F* mixS, TH1F* s1, TH1F* s2) {
41	<	vector<testMC> output;
42	<	//define the scale factors
43	<	Double_t sf1 = 0.0; // QCD
44	<	Double_t sf2 = 0.0; // Wjets
45	<	//KS test
46	<	do {
47	<	sf1 = (NsS - Ns2*sf2)/Ns1;
48	<	if (sf1 < 0) break;
49	<	//cout << "..........sf1 = " << sf1 << endl;
50	<	int nbins = mixS->GetNbinsX();
51	<	double xmin = mixS->GetXaxis()->GetBinLowEdge(1);
52	<	double xmax = mixS->GetXaxis()->GetBinUpEdge(nbins);
53	<	TH1F *test = new TH1F("test", "", nbins, xmin, xmax);
54	<	test -> Sumw2();
55	<	test -> Add(s1,s2,sf1,sf2);
56	<	//test->Scale(1./(1.*test->Integral()));
57	<	Double_t probability = mixS -> KolmogorovTest(test,"");
58	<	testMC temp = testMC(probability,sf1,sf2);
59	<	output.push_back(temp);
60	<	//     cout << "probability = " << setw(15) << temp.prob
61	<	//       << "; sfQCD = "     << setw(10) << temp.scaleF_backg
62	<	//       << "; sfWjets = "   << setw(6)  << temp.scaleF_sample << endl;
63	<	delete test;
64	<	sf2 = sf2 + stepsize;
65	<	} while(sf1 > 0 && sf2 <= 2.0);
66	<	return output;
67	<	}
68	<
69	<	//get the maximum KS test result
70	<	testMC getMax(vector<testMC> vec) {
71	<	testMC maxKSRes;
72	<	Double_t maximum = 0.0;
73	<	for (size_t i = 0; i < vec.size(); i++) {
74	<	if (maximum < vec.at(i).prob)  {
75	<	maximum = vec.at(i).prob;
76	<	maxKSRes = vec.at(i);
77	<	}
78	<	}
79	<	cout << "for maximum: " << setw(12) << maxKSRes.prob
80	<	<< "; sb = " << setw(10) << maxKSRes.scaleF_backg
81	<	<< "; ss = " << setw(5) << maxKSRes.scaleF_sample << endl;
82	<	return maxKSRes;
83	<	}
84	<
30	>	// Output directory
31	>	TString desDir = "Results_2.88pb-1/";
32	>	// User defined parameters
33	>	bool useInv    = false; // whether to use n-1 QCD template
34	>	bool realData  = false;
35	>	// Ntuples to use
36	>	TString suffix = "Sel0"; // Suffix of selection
37	>	TString invNames[2] = {"RelIsogt0p1","D0gt0p02"};
38	>	map<TString,TCanvas*> cvs; // map of usual histogram
39
40		//=================================
41		// Main program
42		//=================================
43		void ikstest() {
44	<	//Style();
44	>	gROOT->SetStyle("CMS");
45		TLatex *latex = new TLatex();
46		latex->SetNDC();
47	<
48	<	ofstream outprint( "ikstest_results_20100901.txt" );
49	<	//open the files with histograms
50	<	map<string,TFile*> mfile;
51	<	mfile["Data"] = TFile::Open("skimmed_Data_20100901/Data_RefSel_v3.root");
52	<	// n-1 cuts
53	<	if (useInv) {
54	<	if (realData)
55	<	//       mfile["InvSel"] = TFile::Open("skimmed_Data_20100825/Data_D0ge0p03.root");
56	<	mfile["InvSel"] = TFile::Open("skimmed_Data_20100901/Data_RelIsoge0p1_v3.root");
103	<	else
104	<	//       mfile["InvSel"] = TFile::Open("skimmed_MC/QCD_D0ge0p03.root");
105	<	mfile["InvSel"] = TFile::Open("skimmed_MC/QCD_RelIsoge0p1_v3.root");
106	<	}
107	<
108	<	// RefSel MC
109	<	mfile["0"] = TFile::Open("skimmed_MC/QCD_RefSel_v3.root");
110	<	mfile["1"] = TFile::Open("skimmed_MC/WJets_RefSel_v3.root");
111	<	mfile["2"] = TFile::Open("skimmed_MC/TTbar_RefSel_v3.root");
112	<
113	<	//define histograms and related parameters
114	<	string histoName[3] = {"h_mu_pt_calo","h_met_calo","h_mt_calo"};
115	<	string histoLabelX[3] = {"p_{T}^{good Muons}", "E_{T}^{#nu}", "m_{T}^{W}"};
116	<	Int_t xbins[3] = {20,20,40};
117	<	Double_t xlow[3] = {0.,0.,0.};
118	<	Double_t xhigh[3] = {100.,100.,200.};
119	<	string sample[3] = {"QCD","Wjets","ttjets"};
120	<
121	<	TH1F* h_[9];
122	<	TH1F* mixh_[3];
123	<	TH1F* hQCD_NEW[3];
124	<	TH1F* hKSres_[3];
125	<	TH1F* hKSvalues_[3];
126	<
127	<	//load the histograms from the root files
128	<	for (int i = 0; i < 3; i++) {// 3 variables
129	<	//cout << "file[" << i << "] : " << endl;
130	<	string nameNewHisto = "mix_"+histoName[i];
131	<	string nameNewHistoSFKS = "finalSF_"+histoName[i];
132	<	string nameNewHistoKSvalues = "KSvalues_"+histoLabelX[i];
133	<
134	<	mixh_[i] = new TH1F(nameNewHisto.c_str(),"",xbins[i],xlow[i],xhigh[i]);
135	<	hKSres_[i] = new TH1F(nameNewHistoSFKS.c_str(),"",xbins[i],xlow[i],xhigh[i]);
136	<	hKSvalues_[i] = new TH1F(nameNewHistoKSvalues.c_str(),"",2./stepsize, stepsize, 2.+stepsize);
137	<
138	<	if (!useInv) {//use QCD MC sample
139	<	hQCD_NEW[i] = (TH1F*) mfile["0"]->Get(TString(histoName[i]))->Clone();
140	<	hQCD_NEW[i] -> Scale(weight_[0]);
141	<	hQCD_NEW[i] -> SetName((histoName[i]).c_str());
47	>	TString invName;
48	>	int size_ninv = (useInv ? 2 : 1);
49	>	for (int ninv = 0;ninv < size_ninv; ++ninv) {
50	>	invName = invNames[ninv];
51	>	if (!useInv) {
52	>	if (!realData) desDir += "MC/";
53	>	else desDir += "Data_MC/";
54	>	} else {
55	>	if (!realData) desDir += "MC_"+invName+"/";
56	>	else desDir += "Data_"+invName+"/";
57		}
58	<	else {
59	<	hQCD_NEW[i] = (TH1F*) mfile["InvSel"]->Get(TString(histoName[i]));
60	<	if (!realData) hQCD_NEW[i] -> Scale(weight_[0]);
61	<	hQCD_NEW[i] -> SetName((histoName[i]).c_str());
58	>	struct stat stDir;
59	>	if (!stat(desDir,&stDir)){
60	>	cout << "Output folder exists! Continues? (enter to continue; 'q' for quit)" << endl;
61	>	char incmd;
62	>	cin.get(incmd);
63	>	if (incmd == 'q') return;
64	>	} else {
65	>	cout << "Creating folder : " << desDir << endl;
66	>	if (mkdir(desDir,0755) == -1){
67	>	std::cerr << "Error creating folder" << endl;
68	>	return;
69	>	}
70		}
71
72	<	mixh_[i] -> Sumw2();
73	<	hKSres_[i] -> Sumw2();
74	<	hKSvalues_[i] -> Sumw2();
75	<	}
72	>	ofstream outprint(TString(desDir+"Results_"+suffix+".txt"));
73	>	//open the files with histograms
74	>	map<string,TFile*> mfile;
75	>	mfile["Data"] = TFile::Open(TString("skimmed_Data_2.88pb-1/Data_"+suffix+".root"));
76	>	// n-1 cuts
77	>	if (useInv) {
78	>	if (realData)
79	>	mfile["InvSel"] = TFile::Open(TString("skimmed_Data_2.88pb-1/Data_"+suffix+"_"+invName+".root"));
80	>	else
81	>	mfile["InvSel"] = TFile::Open(TString("skimmed_MC/v5/QCD_"+suffix+"_"+invName+".root"));
82	>	}
83	>	// RefSel MC
84	>	mfile["0"] = TFile::Open(TString("skimmed_MC/v5/QCD_"+suffix+".root"));
85	>	mfile["1"] = TFile::Open(TString("skimmed_MC/v5/WJets_"+suffix+".root"));
86	>	mfile["2"] = TFile::Open(TString("skimmed_MC/v5/TTbar_"+suffix+".root"));
87	>
88	>	//define histograms and related parameters
89	>	string histoName[3] = {"h_mu_pt_calo","h_met_calo","h_mt_calo"};
90	>	string histoLabelX[3] = {"p_{T}^{#mu} [GeV/c]", "#slash{E}_{T} [GeV/c]", "M_{T}^{W} [GeV/c^{2}]"};
91	>	Int_t xbins[3] = {20,20,40};
92	>	Double_t xlow[3] = {0.,0.,0.};
93	>	Double_t xhigh[3] = {100.,100.,200.};
94	>	string sample[3] = {"QCD","Wjets","ttjets"};
95	>
96	>	TH1F* h_[9];
97	>	TH1F* mixh_[3];
98	>	TH1F* hQCD_NEW[3];
99	>	TH1F* hKSres_[3];
100	>	TH1F* hKSvalues_[3];
101	>
102	>	//load the histograms from the root files
103	>	for (int i = 0; i < 3; i++) {// 3 variables
104	>	//cout << "file[" << i << "] : " << endl;
105	>	string nameNewHisto = "mix_"+histoName[i];
106	>	string nameNewHistoSFKS = "finalSF_"+histoName[i];
107	>	string nameNewHistoKSvalues = "KSvalues_"+histoLabelX[i];
108	>
109	>	mixh_[i] = new TH1F(nameNewHisto.c_str(),"",xbins[i],xlow[i],xhigh[i]);
110	>	hKSres_[i] = new TH1F(nameNewHistoSFKS.c_str(),"",xbins[i],xlow[i],xhigh[i]);
111	>	hKSvalues_[i] = new TH1F(nameNewHistoKSvalues.c_str(),"",2./stepsize, stepsize, 2.+stepsize);
112	>
113	>	if (!useInv) {//use QCD MC sample
114	>	hQCD_NEW[i] = (TH1F*) mfile["0"]->Get(TString(histoName[i]))->Clone();
115	>	hQCD_NEW[i] -> Scale(weight_[0]);
116	>	hQCD_NEW[i] -> SetName((histoName[i]).c_str());
117	>	}
118	>	else {
119	>	hQCD_NEW[i] = (TH1F*) mfile["InvSel"]->Get(TString(histoName[i]));
120	>	if (!realData) hQCD_NEW[i] -> Scale(weight_[0]);
121	>	hQCD_NEW[i] -> SetName((histoName[i]).c_str());
122	>	}
123	>
124	>	mixh_[i] -> Sumw2();
125	>	hKSres_[i] -> Sumw2();
126	>	hKSvalues_[i] -> Sumw2();
127	>	}
128
129	<	for (int n = 0; n < 3; ++n) {// 3 MC samples
130	<	for (int ihisto = 0; ihisto < 3; ihisto++) {// 3 variables
131	<	//cout << "Variable[" << ihisto << "]" << endl;
132	<	string histo_name = histoName[ihisto]+sample[n];
133	<	ostringstream ss;
134	<	ss << n;
135	<	h_[n*3+ihisto] = (TH1F*) mfile[ss.str()]->Get(TString(histoName[ihisto]))->Clone();
136	<	h_[n*3+ihisto] -> Scale(weight_[n]);
137	<	h_[n*3+ihisto] -> SetName(histo_name.c_str());
129	>	for (int n = 0; n < 3; ++n) {// 3 MC samples
130	>	for (int ihisto = 0; ihisto < 3; ihisto++) {// 3 variables
131	>	//cout << "Variable[" << ihisto << "]" << endl;
132	>	string histo_name = histoName[ihisto]+sample[n];
133	>	ostringstream ss;
134	>	ss << n;
135	>	h_[n*3+ihisto] = (TH1F*) mfile[ss.str()]->Get(TString(histoName[ihisto]))->Clone();
136	>	h_[n*3+ihisto] -> Scale(weight_[n]);
137	>	h_[n*3+ihisto] -> SetName(histo_name.c_str());
138	>	}
139		}
164	–	}
140
141	<	//create the mixed samples = "data"
142	<	TCanvas *canvas0 = new TCanvas("Data","Data distributions");
143	<	canvas0->Divide(3,1);
144	<	for (int i = 0; i < 3; i++) {
145	<	canvas0->cd(i+1);
171	<	if (!realData) {
172	<	mixh_[i] -> Add(h_[i],h_[i+3], procQCD,procWjets);
173	<	//mixh_[i] -> Add(mixh_[i],h_[i+6], 1,procttjets);
174	<	//cout << "histo_name: " << mixh_[0]->GetNbinsX() << endl;
141	>	//create the mixed samples = "data"
142	>	TString cvsName0 = "Data";
143	>	if (useInv) {
144	>	cvsName0 += "_";
145	>	cvsName0 += invName;
146		}
147	<	else {
148	<	TH1F *htmp = (TH1F*) mfile["Data"]->Get(TString(histoName[i]));
149	<	mixh_[i] -> Add(htmp,1.);
147	>	cvs[cvsName0] = new TCanvas(cvsName0,"Data distributions",600,700);
148	>	cvs[cvsName0]->Divide(3,1);
149	>	for (int i = 0; i < 3; i++) {
150	>	cvs[cvsName0]->cd(i+1);
151	>	if (!realData) {
152	>	mixh_[i] -> Add(h_[i],h_[i+3], procQCD,procWjets);
153	>	//mixh_[i] -> Add(mixh_[i],h_[i+6], 1,procttjets);
154	>	//cout << "histo_name: " << mixh_[0]->GetNbinsX() << endl;
155	>	}
156	>	else {
157	>	TH1F *htmp = (TH1F*) mfile["Data"]->Get(TString(histoName[i]));
158	>	mixh_[i] -> Add(htmp,1.);
159	>	}
160	>	mixh_[i]->GetXaxis()->SetTitle(histoLabelX[i].c_str());
161	>	mixh_[i]->GetYaxis()->SetTitle("Entries");
162	>	mixh_[i]->DrawClone();
163		}
164	<	mixh_[i]->DrawClone();
181	<	}
182	<	canvas0->SaveAs("Data_distributions.pdf");
164	>	cvs[cvsName0]->SaveAs(TString(desDir+"Data_distributions.pdf"));
165
166	<	//define the weight corrections for each sample
167	<	double NevData                = mixh_[2]->Integral();
168	<	double corr_NevQCD            = h_[2]->Integral();
169	<	double corr_NevQCD_NEW        = hQCD_NEW[2]->Integral();
170	<	double corr_NevWjets          = h_[5]->Integral();
171	<	double corr_Nevttjets         = h_[8]->Integral();
172	<	double corr_Nevmix            = procQCD*corr_NevQCD+procWjets*corr_NevWjets;
173	<	//double corr_Nevmix            = procQCD*corr_NevQCD+procWjets*corr_NevWjets+procttjets*corr_Nevttjets;
174	<	if (!realData)
175	<	outprint << "Events mix sample    = " << corr_Nevmix << endl;
176	<	else
177	<	outprint << "Events in Data       = " << NevData << endl;
178	<	outprint << "Events QCD sample    = " << corr_NevQCD << endl;
179	<	outprint << "Events Wjets sample  = " << corr_NevWjets << endl;
180	<	outprint << "Events InvSel sample = " << corr_NevQCD_NEW << endl;
181	<
182	<	//define the containers for chosen numbers (coressponding to the max KStest result)
183	<	testMC maxProb[3];
184	<
185	<	//define the scale factors calculated using information obtained from all parameters
186	<	Double_t SFbackg = 0.0;
187	<	Double_t sumSFbackg = 0.0;
188	<	Double_t SFsample = 0.0;
189	<	Double_t sumSFsample = 0.0;
190	<	Double_t allKS = 0.0;
191	<
192	<	//do the KS test by varying the scale factors
193	<	for (int i = 0; i < 3; i++) { // 3 variables
194	<	TH1F *data = (TH1F*)mixh_[i]->Clone();
195	<	data -> SetName("dataClone");
196	<	//data -> Scale(1./data->Integral());
197	<	vector<testMC> resultsKS = doKStest((realData ? NevData : corr_Nevmix),
198	<	(useInv ? corr_NevQCD_NEW : corr_NevQCD),
199	<	corr_NevWjets,
200	<	data, hQCD_NEW[i], h_[i+3]);
201	<	testMC tksmax = getMax(resultsKS);
202	<	maxProb[i] = tksmax;
203	<	outprint << "\nFor the plot " << histoLabelX[i] << " the results are:"<< endl;
204	<	outprint << "\tmax Probability = " << maxProb[i].prob << endl;
205	<	outprint << "\tproc_background = " << maxProb[i].scaleF_backg << endl;
206	<	outprint << "\tproc_sample     = " << maxProb[i].scaleF_sample << endl;
207	<
208	<	outprint << "\n\tpercent_B of Data = "
209	<	<< maxProb[i].scaleF_backg*corr_NevQCD_NEW*100/(realData ? NevData : corr_Nevmix) << endl;
210	<	outprint << "\tpercent_S of Data = "
211	<	<< maxProb[i].scaleF_sample*corr_NevWjets*100/(realData ? NevData : corr_Nevmix) << endl;
212	<	outprint << "---------------------------" << endl;
213	<
214	<	//create the mixed samples with KS test results
215	<	sumSFbackg += maxProb[i].prob*maxProb[i].scaleF_backg;
216	<	sumSFsample += maxProb[i].prob*maxProb[i].scaleF_sample;
217	<	allKS += maxProb[i].prob;
218	<
219	<	//fill a histogram with the results from the KS test for each variable
220	<	for (int jiter = 0; jiter < resultsKS.size(); jiter++) {
221	<	if (resultsKS.at(jiter).prob == 1.)
222	<	cout << "variable [" << i << "]: prob[" << jiter << "]= " << resultsKS.at(jiter).prob << endl;
223	<	hKSvalues_[i]->SetBinContent(jiter,resultsKS.at(jiter).prob);
166	>	//define the weight corrections for each sample
167	>	double NevData                = mixh_[2]->Integral();
168	>	double corr_NevQCD          = h_[2]->Integral();
169	>	double corr_NevQCD_NEW        = hQCD_NEW[2]->Integral();
170	>	double corr_NevWjets        = h_[5]->Integral();
171	>	double corr_Nevttjets       = h_[8]->Integral();
172	>	double corr_Nevmix          = procQCD*corr_NevQCD+procWjets*corr_NevWjets;
173	>	//double corr_Nevmix            = procQCD*corr_NevQCD+procWjets*corr_NevWjets+procttjets*corr_Nevttjets;
174	>	if (!realData)
175	>	outprint << "Events mix sample    = " << corr_Nevmix << endl;
176	>	else
177	>	outprint << "Events in Data       = " << NevData << endl;
178	>	outprint << "Events QCD sample    = " << corr_NevQCD << endl;
179	>	outprint << "Events Wjets sample  = " << corr_NevWjets << endl;
180	>	outprint << "Events InvSel sample = " << corr_NevQCD_NEW << endl;
181	>
182	>	//define the containers for chosen numbers (coressponding to the max KStest result)
183	>	testMC maxProb[3];
184	>
185	>	//define the scale factors calculated using information obtained from all parameters
186	>	Double_t SFbackg = 0.0;
187	>	Double_t sumSFbackg = 0.0;
188	>	Double_t SFsample = 0.0;
189	>	Double_t sumSFsample = 0.0;
190	>	Double_t allKS = 0.0;
191	>
192	>	//do the KS test by varying the scale factors
193	>	for (int i = 0; i < 3; i++) { // 3 variables
194	>	TH1F *data = (TH1F*)mixh_[i]->Clone();
195	>	data -> SetName("dataClone");
196	>	//data -> Scale(1./data->Integral());
197	>	vector<testMC> resultsKS = doKStest((realData ? NevData : corr_Nevmix),
198	>	(useInv ? corr_NevQCD_NEW : corr_NevQCD),
199	>	corr_NevWjets,
200	>	data, hQCD_NEW[i], h_[i+3]);
201	>	testMC tksmax = getMax(resultsKS);
202	>	maxProb[i] = tksmax;
203	>	outprint << "\nFor the plot " << histoLabelX[i] << " the results are:"<< endl;
204	>	outprint << "\tmax Probability = " << maxProb[i].prob << endl;
205	>	outprint << "\tproc_background = " << maxProb[i].scaleF_backg << endl;
206	>	outprint << "\tproc_sample     = " << maxProb[i].scaleF_sample << endl;
207	>
208	>	outprint << "\n\tpercent_B of Data = "
209	>	<< maxProb[i].scaleF_backg*corr_NevQCD_NEW*100/(realData ? NevData : corr_Nevmix) << endl;
210	>	outprint << "\tpercent_S of Data = "
211	>	<< maxProb[i].scaleF_sample*corr_NevWjets*100/(realData ? NevData : corr_Nevmix) << endl;
212	>	outprint << "---------------------------" << endl;
213	>
214	>	//create the mixed samples with KS test results
215	>	sumSFbackg += maxProb[i].prob*maxProb[i].scaleF_backg;
216	>	sumSFsample += maxProb[i].prob*maxProb[i].scaleF_sample;
217	>	allKS += maxProb[i].prob;
218	>
219	>	//fill a histogram with the results from the KS test for each variable
220	>	for (int jiter = 0; jiter < resultsKS.size(); jiter++) {
221	>	if (resultsKS.at(jiter).prob == 1.)
222	>	cout << "variable [" << i << "]: prob[" << jiter << "]= " << resultsKS.at(jiter).prob << endl;
223	>	hKSvalues_[i]->SetBinContent(jiter,resultsKS.at(jiter).prob);
224	>	}
225	>	delete data;
226		}
243	–	delete data;
244	–	}
245	–
246	–	//calculate the final scale factors
247	–	SFbackg = sumSFbackg/allKS;
248	–	SFsample = sumSFsample/allKS;
249	–	outprint << "allKS = " << allKS << "\tbackground = " << SFbackg << "\tsample = " << SFsample << endl;
250	–	outprint << "==> Scale Factor for QCD MC = " << SFbackg*corr_NevQCD_NEW/corr_NevQCD << endl;
251	–	outprint << "\tcombined percent_B of Data = "
252	–	<< SFbackg*corr_NevQCD_NEW*100/(realData ? NevData : corr_Nevmix) << endl;
253	–	outprint << "\tcombined percent_S of Data = "
254	–	<< SFsample*corr_NevWjets*100/(realData ? NevData : corr_Nevmix) << endl;
255	–	outprint << "\n" << endl;
256	–	outprint << "=================================" << endl;
257	–	outprint << "\n" << endl;
258	–
259	–
260	–	//=================================
261	–	// Plots
262	–	//=================================
263	–	for (int i = 0; i < 3; i++) {// 3 variables
264	–	hKSres_[i] -> Add(hQCD_NEW[i],h_[i+3],SFbackg,SFsample);
265	–	outprint << "hKSres->Integral() = " << hKSres_[i]->Integral() << endl;
266	–	outprint << "Data->Integral()   = " << mixh_[i]->Integral() << endl;
267	–
268	–	mixh_[i]->Rebin(2);
269	–	hQCD_NEW[i]->Rebin(2);
270	–	h_[i]->Rebin(2);
271	–	h_[i+3]->Rebin(2);
272	–	hKSres_[i]->Rebin(2);
273	–	//hKSvalues_[i]->Rebin(2);
274	–
275	–	mixh_[i]     ->SetLineColor(1);
276	–	hQCD_NEW[i]  ->SetLineColor(3);
277	–	h_[i]        ->SetLineColor(6);
278	–	h_[i+3]      ->SetLineColor(4);
279	–	hKSres_[i]   ->SetLineColor(2);
280	–	hKSvalues_[i]->SetLineColor(i+1);
281	–
282	–	mixh_[i]     ->SetMarkerColor(1);
283	–	hQCD_NEW[i]  ->SetMarkerColor(3);
284	–	h_[i]        ->SetMarkerColor(6);
285	–	h_[i+3]      ->SetMarkerColor(4);
286	–	hKSres_[i]   ->SetMarkerColor(2);
287	–	hKSvalues_[i]->SetMarkerColor(i+1);
288	–
289	–	mixh_[i]     ->SetMarkerStyle(24);
290	–	hQCD_NEW[i]  ->SetMarkerStyle(20);
291	–	h_[i]        ->SetMarkerStyle(20);
292	–	h_[i+3]      ->SetMarkerStyle(20);
293	–	hKSres_[i]   ->SetMarkerStyle(20);
294	–	hKSvalues_[i]->SetMarkerStyle(20);
295	–
296	–	mixh_[i]     ->SetMarkerSize(1.4);
297	–	hQCD_NEW[i]  ->SetMarkerSize(1.1);
298	–	h_[i]        ->SetMarkerSize(1.1);
299	–	h_[i+3]      ->SetMarkerSize(1.1);
300	–	hKSres_[i]   ->SetMarkerSize(0.9);
301	–	hKSvalues_[i]->SetMarkerSize(1.1);
302	–
303	–	mixh_[i]     ->SetStats(0);
304	–	hQCD_NEW[i]  ->SetStats(0);
305	–	h_[i]        ->SetStats(0);
306	–	h_[i+3]      ->SetStats(0);
307	–	hKSres_[i]   ->SetStats(0);
308	–	hKSvalues_[i]->SetStats(0);
309	–
310	–	mixh_[i]->GetXaxis()->SetTitle(histoLabelX[i].c_str());
311	–	mixh_[i]->GetYaxis()->SetTitle("Entries");
312	–	hKSres_[i]->GetXaxis()->SetTitle(histoLabelX[i].c_str());
313	–	hKSres_[i]->GetYaxis()->SetTitle("Entries");
314	–	hKSvalues_[i]->GetXaxis()->SetTitle("iteration #");
315	–	hKSvalues_[i]->GetYaxis()->SetTitle("KS test values");
316	–
317	–	string nameCanvas1 = histoName[i]+"_QCD.pdf";
318	–	TCanvas *canvas1 = new TCanvas(nameCanvas1.c_str(), "");
319	–	hQCD_NEW[i] -> Scale(1./hQCD_NEW[i]->Integral());
320	–	h_[i] -> Scale(1./h_[i]->Integral());
321	–	h_[i+3] -> Scale(1./h_[i+3]->Integral());
322	–	outprint << "For " << histoName[i] << " , the KStest result btw MC_QCD/InvSel is = "
323	–	<< h_[i] -> KolmogorovTest(hQCD_NEW[i],"") << endl;
324	–	hQCD_NEW[i]->Draw("P");
325	–	h_[i]->Draw("sameP");
326	–	h_[i+3]->Draw("sameP");
327	–	TLegend *legend1 = new TLegend(0.7, 0.70, 0.9, 0.85);
328	–	legend1->AddEntry(h_[i], "default");
329	–	legend1->AddEntry(h_[i+3], "W+jets");
330	–	legend1->AddEntry(hQCD_NEW[i], "new");
331	–	legend1->Draw();
332	–	legend1->SetFillColor(kWhite);
333	–	latex->DrawLatex(0.22,0.91,histoName[i].c_str());
334	–	canvas1->SetLogy();
335	–	canvas1->SaveAs(nameCanvas1.c_str());
336	–
337	–	string nameCanvas2 = histoName[i]+"_dataKS.pdf";
338	–	TCanvas *canvas2 = new TCanvas(nameCanvas2.c_str(), "");
339	–	hKSres_[i]->Draw("P");
340	–	mixh_[i]->Draw("sameP");
341	–	TLegend *legend2 = new TLegend(0.7, 0.70, 0.9, 0.85);
342	–	legend2->AddEntry(mixh_[i], "Data");
343	–	legend2->AddEntry(hKSres_[i], "KS result");
344	–	legend2->Draw();
345	–	legend2->SetFillColor(kWhite);
346	–	latex->DrawLatex(0.22,0.91,histoName[i].c_str());
347	–	canvas2->SetLogy();
348	–	canvas2->SaveAs(nameCanvas2.c_str());
227
228	<	}
228	>	//calculate the final scale factors
229	>	SFbackg = sumSFbackg/allKS;
230	>	SFsample = sumSFsample/allKS;
231	>	outprint << "allKS = " << allKS << "\tbackground = " << SFbackg << "\tsample = " << SFsample << endl;
232	>	outprint << "==> Scale Factor for QCD MC = " << SFbackg*corr_NevQCD_NEW/corr_NevQCD << endl;
233	>	outprint << "\tcombined percent_B of Data = "
234	>	<< SFbackg*corr_NevQCD_NEW*100/(realData ? NevData : corr_Nevmix) << endl;
235	>	outprint << "\tcombined percent_S of Data = "
236	>	<< SFsample*corr_NevWjets*100/(realData ? NevData : corr_Nevmix) << endl;
237	>	outprint << "\n" << endl;
238	>	outprint << "=================================" << endl;
239	>	outprint << "\n" << endl;
240	>
241	>
242	>	//=================================
243	>	// Plots
244	>	//=================================
245	>	for (int i = 0; i < 3; i++) {// 3 variables
246	>	hKSres_[i] -> Add(hQCD_NEW[i],h_[i+3],SFbackg,SFsample);
247	>	outprint << "hKSres->Integral() = " << hKSres_[i]->Integral() << endl;
248	>	outprint << "Data->Integral()   = " << mixh_[i]->Integral() << endl;
249	>
250	>	mixh_[i]->Rebin(2);
251	>	hQCD_NEW[i]->Rebin(2);
252	>	h_[i]->Rebin(2);
253	>	h_[i+3]->Rebin(2);
254	>	hKSres_[i]->Rebin(2);
255	>	//hKSvalues_[i]->Rebin(2);
256	>
257	>	mixh_[i]     ->SetLineColor(1);
258	>	hQCD_NEW[i]  ->SetLineColor(2);
259	>	h_[i]        ->SetLineColor(4);
260	>	h_[i+3]      ->SetLineColor(3);
261	>	hKSres_[i]   ->SetLineColor(2);
262	>	hKSvalues_[i]->SetLineColor(i+1);
263	>
264	>	mixh_[i]     ->SetMarkerColor(1);
265	>	hQCD_NEW[i]  ->SetMarkerColor(2);
266	>	h_[i]        ->SetMarkerColor(4);
267	>	h_[i+3]      ->SetMarkerColor(3);
268	>	hKSres_[i]   ->SetMarkerColor(2);
269	>	hKSvalues_[i]->SetMarkerColor(i+1);
270	>
271	>	mixh_[i]     ->SetMarkerStyle(24);
272	>	hQCD_NEW[i]  ->SetMarkerStyle(20);
273	>	h_[i]        ->SetMarkerStyle(20);
274	>	h_[i+3]      ->SetMarkerStyle(20);
275	>	hKSres_[i]   ->SetMarkerStyle(20);
276	>	hKSvalues_[i]->SetMarkerStyle(20);
277	>
278	>	mixh_[i]     ->SetMarkerSize(1.4);
279	>	hQCD_NEW[i]  ->SetMarkerSize(1.1);
280	>	h_[i]        ->SetMarkerSize(1.1);
281	>	h_[i+3]      ->SetMarkerSize(1.1);
282	>	hKSres_[i]   ->SetMarkerSize(0.9);
283	>	hKSvalues_[i]->SetMarkerSize(1.1);
284	>
285	>	mixh_[i]     ->SetStats(0);
286	>	hQCD_NEW[i]  ->SetStats(0);
287	>	h_[i]        ->SetStats(0);
288	>	h_[i+3]      ->SetStats(0);
289	>	hKSres_[i]   ->SetStats(0);
290	>	hKSvalues_[i]->SetStats(0);
291	>
292	>	hKSres_[i]->GetXaxis()->SetTitle(histoLabelX[i].c_str());
293	>	hKSres_[i]->GetYaxis()->SetTitle("Entries");
294	>	hKSvalues_[i]->GetXaxis()->SetTitle("iteration #");
295	>	hKSvalues_[i]->GetYaxis()->SetTitle("KS test values");
296	>	h_[i]->GetXaxis()->SetTitle(histoLabelX[i].c_str());
297	>	h_[i]->GetYaxis()->SetTitle("A.U.");
298	>
299	>	TString nameCanvas1 = desDir+histoName[i]+"_QCD_"+suffix+".pdf";
300	>	TString cvsName1 = histoName[i]+"_QCD";
301	>	if(useInv) cvsName1 = cvsName1 + "_" + invName;
302	>	cvs[cvsName1] = new TCanvas(cvsName1,"",600,700);
303	>	hQCD_NEW[i] -> Scale(1./hQCD_NEW[i]->Integral());
304	>	h_[i] -> Scale(1./h_[i]->Integral());
305	>	h_[i+3] -> Scale(1./h_[i+3]->Integral());
306	>	outprint << "For " << histoName[i] << " , the KStest result btw MC_QCD/InvSel is = "
307	>	<< h_[i] -> KolmogorovTest(hQCD_NEW[i],"") << endl;
308	>	h_[i]->Draw("P");
309	>	if (useInv)
310	>	hQCD_NEW[i]->Draw("sameP");
311	>	h_[i+3]->Draw("sameP");
312	>	TLegend *legend1 = new TLegend(0.7, 0.70, 0.9, 0.85);
313	>	legend1->AddEntry(h_[i], "QCD");
314	>	if (useInv)
315	>	legend1->AddEntry(hQCD_NEW[i], "QCD - InvSel");
316	>	legend1->AddEntry(h_[i+3], "W+jets");
317	>	legend1->Draw();
318	>	legend1->SetFillColor(kWhite);
319	>	//latex->DrawLatex(0.22,0.91,histoName[i].c_str());
320	>	//cvs[cvsName1]->SetLogy();
321	>	cvs[cvsName1]->SaveAs(nameCanvas1);
322	>
323	>	TString nameCanvas2 = desDir+histoName[i]+"_dataKS_"+suffix+".pdf";
324	>	TString cvsName2 = histoName[i]+"_dataKS";
325	>	if(useInv) cvsName2 = cvsName2 + "_" + invName;
326	>	cvs[cvsName2] = new TCanvas(cvsName2,"",600,700);
327	>	hKSres_[i]->Draw("P");
328	>	mixh_[i]->Draw("sameP");
329	>	TLegend *legend2 = new TLegend(0.7, 0.70, 0.9, 0.85);
330	>	legend2->AddEntry(mixh_[i], "Data");
331	>	legend2->AddEntry(hKSres_[i], "KS result");
332	>	legend2->Draw();
333	>	legend2->SetFillColor(kWhite);
334	>	//latex->DrawLatex(0.22,0.91,histoName[i].c_str());
335	>	//cvs[cvsName2]->SetLogy();
336	>	cvs[cvsName2]->SaveAs(nameCanvas2);
337
338	+	}
339
353	–	TCanvas *canvas3 = new TCanvas("KStestValues", "");
354	–	//hKSvalues_[0]->GetXaxis()->SetRangeUser(0.9,1.1);
355	–	//hKSvalues_[0]->GetYaxis()->SetRangeUser(1e-3,1.1);
356	–	hKSvalues_[0]->Draw();
357	–	hKSvalues_[1]->SetLineColor(2);
358	–	hKSvalues_[1]->Draw("same");
359	–	hKSvalues_[2]->SetLineColor(4);
360	–	hKSvalues_[2]->Draw("same");
361	–	TLegend *legend3 = new TLegend(0.7, 0.70, 0.9, 0.85);
362	–	legend3->AddEntry(hKSvalues_[0], "muon_pT");
363	–	legend3->AddEntry(hKSvalues_[1], "MET");
364	–	legend3->AddEntry(hKSvalues_[2], "W_mT");
365	–	legend3->Draw();
366	–	legend3->SetFillColor(kWhite);
367	–	latex->DrawLatex(0.22,0.91,"KS test values");
368	–	canvas3->SetLogy();
369	–	string nameCanvas3 = "KStestValues_newQCD.pdf";
370	–	canvas3->SaveAs(nameCanvas3.c_str());
340
341	+	TString cvsName3 = "KStestValues";
342	+	if(useInv) cvsName3 = cvsName3 + "_" + invName;
343	+	cvs[cvsName3] = new TCanvas(cvsName3,"",600,700);
344	+	//hKSvalues_[0]->GetXaxis()->SetRangeUser(0.9,1.2);
345	+	hKSvalues_[0]->GetYaxis()->SetRangeUser(1e-36,1.2);
346	+	hKSvalues_[0]->Draw();
347	+	hKSvalues_[1]->Draw("same");
348	+	hKSvalues_[2]->Draw("same");
349	+	TLegend *legend3 = new TLegend(0.7, 0.70, 0.9, 0.85);
350	+	legend3->AddEntry(hKSvalues_[0], "muon_pT");
351	+	legend3->AddEntry(hKSvalues_[1], "MET");
352	+	legend3->AddEntry(hKSvalues_[2], "W_mT");
353	+	legend3->Draw();
354	+	legend3->SetFillColor(kWhite);
355	+	//latex->DrawLatex(0.22,0.91,"KS test values");
356	+	cvs[cvsName3]->SetLogy();
357	+	TString nameCanvas3 = desDir+"KStestValues_newQCD"+suffix+".pdf";
358	+	cvs[cvsName3]->SaveAs(nameCanvas3);
359	+	}
360		}

Diff Legend

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