ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/cbrown/Development/Plotting/Modules/ZbTools.C

Comparing UserCode/cbrown/Development/Plotting/Modules/ZbTools.C (file contents):
Revision 1.22 by buchmann, Mon Jan 21 14:43:36 2013 UTC vs.
Revision 1.26 by buchmann, Mon Mar 18 16:05:04 2013 UTC

#	Line 17 | Line 17
17		#include <TProfile.h>
18		#include <TLegendEntry.h>
19
20	+	#include <RooRealVar.h>
21	+	#include <RooArgSet.h>
22	+	#include <RooDataSet.h>
23	+	#include <RooDataHist.h>
24	+	#include <RooMCStudy.h>
25	+	#include <RooCategory.h>
26	+
27	+	#include <RooPlot.h>
28	+	#include <RooSimultaneous.h>
29	+	#include <RooAddPdf.h>
30	+	#include <RooFitResult.h>
31	+	#include <RooVoigtian.h>
32	+	#include <RooChebychev.h>
33	+	#include <RooMsgService.h>
34	+
35	+	#include <RooStats/ModelConfig.h>
36	+	#include "RooStats/ProfileLikelihoodCalculator.h"
37	+	#include "RooStats/LikelihoodInterval.h"
38	+	#include "RooStats/HypoTestResult.h"
39	+	#include "RooStats/SimpleLikelihoodRatioTestStat.h"
40	+	#include "RooStats/ProfileLikelihoodTestStat.h"
41	+	#include "RooStats/HybridCalculatorOriginal.h"
42	+	#include "RooStats/HypoTestInverterOriginal.h"
43	+
44	+
45	+
46		using namespace std;
47
48		using namespace PlottingSetup;
#	Line 26 | Line 52 | namespace ZbData {
52		TCut LeadingB("ZbCHS3010_bTagProbCSVBP[0]>0.898");
53		TCut EtaB("abs(ZbCHS3010_pfJetGoodEta[0])<1.3");
54		TCut PhiZcut("abs(ZbCHS3010_pfJetDphiZ[0])>2.7");
55	<
55	>	float LowerHistoBoundary=0.5; // where to cut on MPF (or Rabs)
56	>	float UpperHistoBoundary=1.5; // where to cut on MPF (or Rabs)
57	>
58		}
59
60		TCut STDWEIGHT;
#	Line 64 | Line 92 | public:
92		Value Rabs_Result_FaceValue;
93		Value MPF_Result_Extrapolated;
94		Value Rabs_Result_Extrapolated;
95	+	void Print();
96
97		ZbResultContainer(string name);
69	–	ZbResultContainer(const ZbResultContainer &old);
98		};
99
100
101		ZbResultContainer::ZbResultContainer(string _name) {
102		name=_name;
103	<	cout<< "Result Container for \"" << name << "\" has been initialized" << endl;
103	>	dout<< "Result Container for \"" << name << "\" has been initialized" << endl;
104
105		}
106
107	<	std::ostream &operator<<(std::ostream &ostr, ZbResultContainer &b)
107	>	void ZbResultContainer::Print()
108		{
109	<	ostr << "Results for " << b.name << endl;
110	<	ostr << " MPF and Rabs locations: " << endl;
111	<	for(int i=0;i<b.MPF_Rabs_location.size();i++) ostr << "   " << b.MPF_Rabs_location[i] << endl;
112	<	ostr << " Results: " << endl;
113	<	ostr << "   Face Value: " << endl;
114	<	ostr << "     MPF " << b.MPF_Result_FaceValue << endl;
115	<	ostr << "     Rabs " << b.Rabs_Result_FaceValue << endl;
116	<	ostr << "   Extrapolated: " << endl;
117	<	ostr << "     MPF " << b.MPF_Result_Extrapolated << endl;
118	<	ostr << "             Rabs " << b.Rabs_Result_Extrapolated << endl;
91	<
92	<	return ostr;
109	>	dout << "Results for " << this->name << endl;
110	>	dout << " MPF and Rabs locations: " << endl;
111	>	for(int i=0;i<this->MPF_Rabs_location.size();i++) dout << "   " << this->MPF_Rabs_location[i] << endl;
112	>	dout << " Results: " << endl;
113	>	dout << "   Face Value: " << endl;
114	>	dout << "     MPF " << this->MPF_Result_FaceValue << endl;
115	>	dout << "     Rabs " << this->Rabs_Result_FaceValue << endl;
116	>	dout << "   Extrapolated: " << endl;
117	>	dout << "     MPF " << this->MPF_Result_Extrapolated << endl;
118	>	dout << "             Rabs " << this->Rabs_Result_Extrapolated << endl;
119		}
120
121
96	–
122		using namespace ZbData;
123
124
#	Line 161 | Line 186 | CutYields print_yield(TCut cut) {
186		if(Contains(name,"Z_b_")) Yield.Zb+=hist->Integral();
187		}
188
189	<	cout << Yield << endl;
190	<	cout << endl << endl;
189	>	dout << Yield << endl;
190	>	dout << endl << endl;
191
192		delete data;
193		CleanLegends();
#	Line 225 | Line 250 | void draw_kin_variable(string variable,
250		speciall->Draw();
251		}
252
253	<	save_with_ratio(datahisto,mcstack,kinpad->cd(),saveas);
253	>	Save_With_Ratio(datahisto,mcstack,kinpad->cd(),saveas);
254
255		datahisto->Delete();
256		delete kinleg;
#	Line 236 | Line 261 | void draw_kin_variable(string variable,
261
262		void print_all_b_yields() {
263		vector<CutYields> yields;
264	<	cout << "Basic selection with a b jet" << endl;
264	>	dout << "Basic selection with a b jet" << endl;
265		yields.push_back(print_yield(ZplusBsel&&"ZbCHS3010_pfJetGoodNumBtag>0"));
266	<	cout << "Leading jet is a b-jet " << endl;
266	>	dout << "Leading jet is a b-jet " << endl;
267		yields.push_back(print_yield(ZplusBsel&&LeadingB));
268	<	cout << "|#eta_{b}|<1.3 " << endl;
268	>	dout << "|#eta_{b}|<1.3 " << endl;
269		yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB));
270	<	cout << "|#delta#phi(b<Z)|>2.7" << endl;
270	>	dout << "|#delta#phi(b<Z)|>2.7" << endl;
271		yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut));
272	<	cout << "30<ptZ<1000 GeV" << endl;
272	>	dout << "30<ptZ<1000 GeV" << endl;
273		yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000"));
274	<	cout << "#alpha < 0.35" << endl;
274	>	dout << "#alpha < 0.35" << endl;
275		yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000"&&"(ZbCHS3010_alpha)<0.35"));
276	<	cout << "#alpha < 0.3" << endl;
276	>	dout << "#alpha < 0.3" << endl;
277		yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000"&&"(ZbCHS3010_alpha)<0.3"));
278	<	/*  cout << "#alpha < 0.2" << endl;
278	>	/*  dout << "#alpha < 0.2" << endl;
279		yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000"&&"(ZbCHS3010_alpha)<0.2"));
280	<	cout << "#alpha < 0.15" << endl;
280	>	dout << "#alpha < 0.15" << endl;
281		yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000"&&"(ZbCHS3010_alpha)<0.15"));
282	<	cout << "#alpha < 0.1" << endl;
282	>	dout << "#alpha < 0.1" << endl;
283		yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000"&&"(ZbCHS3010_alpha)<0.1"));*/
284
285		for(int iy=0;iy<yields.size();iy++) {
286	<	cout << yields[iy] << endl;
286	>	dout << yields[iy] << endl;
287		}
288		}
289
#	Line 344 | Line 369 | void DrawEvilCutFlow() {
369
370		void draw_Zb_kin_vars() {
371
372	<	draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000"),25,0,200,"Z p_{T}","Official/Zpt",1);
372	>	draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000"),25,0,200,"Z p_{T}","Official/Zpt",0);
373	>	draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&id1==0"),25,0,200,"Z p_{T}","Official/Zpt_ee",1);
374	>	draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&id1==1"),25,0,200,"Z p_{T}","Official/Zpt_mm",1);
375	>	draw_kin_variable("pt",ZplusBsel&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000"),25,0,200,"Z p_{T}","Official/Zpt__nonBtagged",1);
376	>	draw_kin_variable("pt",ZplusBsel&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&id1==0"),25,0,200,"Z p_{T}","Official/Zpt__nonBtagged_ee",1);
377	>	draw_kin_variable("pt",ZplusBsel&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&id1==1"),25,0,200,"Z p_{T}","Official/Zpt__nonBtagged_mm",1);
378	>
379		draw_kin_variable("ZbCHS3010_pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000"),25,0,200,"Sub-Leading Jet Pt","Official/Jet2Pt",1);
380	+	draw_kin_variable("ZbCHS3010_pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&id1==0"),25,0,200,"Sub-Leading Jet Pt","Official/Jet2Pt_ee",1);
381	+	draw_kin_variable("ZbCHS3010_pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&id1==1"),25,0,200,"Sub-Leading Jet Pt","Official/Jet2Pt_mm",1);
382	+	draw_kin_variable("ZbCHS3010_pfJetGoodPt[1]",ZplusBsel&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000"),25,0,200,"Sub-Leading Jet Pt","Official/Jet2Pt__nonBtagged",1);
383	+	draw_kin_variable("ZbCHS3010_pfJetGoodPt[1]",ZplusBsel&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&id1==0"),25,0,200,"Sub-Leading Jet Pt","Official/Jet2Pt__nonBtagged_ee",1);
384	+	draw_kin_variable("ZbCHS3010_pfJetGoodPt[1]",ZplusBsel&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&id1==1"),25,0,200,"Sub-Leading Jet Pt","Official/Jet2Pt__nonBtagged_mm",1);
385	+
386		draw_kin_variable("ZbCHS3010_pfJetGoodPt[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000"),25,0,200,"Leading Jet Pt (B)","Official/LeadingJetPt",1);
387	+	draw_kin_variable("ZbCHS3010_pfJetGoodPt[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&id1==0"),25,0,200,"Leading Jet Pt (B)","Official/LeadingJetPt_ee",1);
388	+	draw_kin_variable("ZbCHS3010_pfJetGoodPt[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&id1==1"),25,0,200,"Leading Jet Pt (B)","Official/LeadingJetPt_mm",1);
389	+	draw_kin_variable("ZbCHS3010_pfJetGoodPt[0]",ZplusBsel&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000"),25,0,200,"Leading Jet Pt (B)","Official/LeadingJetPt__nonBtagged",1);
390	+	draw_kin_variable("ZbCHS3010_pfJetGoodPt[0]",ZplusBsel&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&id1==0"),25,0,200,"Leading Jet Pt (B)","Official/LeadingJetPt__nonBtagged_ee",1);
391	+	draw_kin_variable("ZbCHS3010_pfJetGoodPt[0]",ZplusBsel&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&id1==1"),25,0,200,"Leading Jet Pt (B)","Official/LeadingJetPt__nonBtagged_mm",1);
392	+
393	+
394		draw_kin_variable("ZbCHS3010_pfJetGoodPt[1]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000"),20,0,2,"p_{T}^{J2} / p_{T}^{Z}","Official/SubLeadingJetPt_Over_ZPt",1);
395		draw_kin_variable("ZbCHS3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<50"),   40,0,2,"#alpha","Official/alpha_pt_30_to_50",1);
396		draw_kin_variable("ZbCHS3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>50&&pt<75"),  40,0,2,"#alpha","Official/alpha_50_to_75",1);
#	Line 371 | Line 415 | void draw_mpf_vars() {
415		draw_kin_variable("ZbCHS3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.05"),20,0.0,2.0,"p_{T} b-jet / p_{T} Z","ptb_over_ptz__mpf_alpha_smaller_0p05",0,"#alpha<0.05, 30 GeV < p_{T}^{Z} < 1000 GeV");
416		}
417
418	+	void make_RooFit(string ContainerName, TH1F *h, float &peak, float &error_peak, string saveas) {
419	+	TCanvas *fcan = new TCanvas("fcan","fcan");
420	+	RooRealVar mpf("mpf","mpf", -5, 5, "");
421	+	RooDataHist AllData("AllData", "CMS Real Dataset", mpf , h);
422	+
423	+	RooRealVar fraction("fraction", "fraction", 0.8, 0, 1);
424	+	RooRealVar mean("mean", "mean", 1, 0, 3);
425	+	RooRealVar sigma("sigma", "sigma", 1, 0, 5);
426	+	RooRealVar width("width", "width", 2.4, 0, 5);
427	+	RooRealVar chevvar1("chevvar1", "chevvar1", 0.2,-1.0, 1.0);
428	+	RooRealVar chevvar2("chevvar2", "chevvar2", 0.2,-1.0, 1.0);
429	+
430	+	RooVoigtian signal("signal", "signal", mpf, mean, width, sigma);
431	+	RooChebychev background("background_", "background", mpf, RooArgList(chevvar1, chevvar2));
432	+
433	+	RooAddPdf pdf("pdf",  "pdf", RooArgList(signal, background), fraction);
434	+
435	+	RooFitResult *result = pdf.fitTo(AllData, RooFit::Save());
436	+	RooPlot* frame = mpf.frame(RooFit::Bins(30), RooFit::Title("MPF")) ;
437	+	AllData.plotOn(frame);
438	+	pdf.plotOn(frame, RooFit::LineColor(kBlue)) ;
439	+	pdf.plotOn(frame, RooFit::Components(RooArgSet(background)), RooFit::LineStyle(kDashed), RooFit::LineColor(kRed)) ;
440	+	pdf.plotOn(frame, RooFit::Components(RooArgSet(signal)), RooFit::LineStyle(kDashed), RooFit::LineColor(kGreen)) ;
441	+
442	+	RooRealVar* thePeak = (RooRealVar*) result->floatParsFinal().find("mean");
443	+	peak = thePeak->getVal();
444	+	error_peak = thePeak->getError();
445	+
446	+	///******************* just making it beautiful ....
447	+	RooPlot* frame2 = mpf.frame(RooFit::Bins(30), RooFit::Title("MPF")) ;
448	+	TH1F *hcopy = (TH1F*)h->Clone("hcopy");
449	+	hcopy->Rebin(40);
450	+	RooDataHist AllRData("AllRData", "CMS Dataset", mpf , hcopy);
451	+
452	+	AllRData.plotOn(frame2);
453	+	pdf.plotOn(frame2, RooFit::LineColor(kBlue)) ;
454	+	pdf.plotOn(frame2, RooFit::Components(RooArgSet(background)), RooFit::LineStyle(kDashed), RooFit::LineColor(kRed)) ;
455	+	pdf.plotOn(frame2, RooFit::Components(RooArgSet(signal)), RooFit::LineStyle(kDashed), RooFit::LineColor(kGreen)) ;
456	+
457	+	string saveasFinal = saveas + "_Fit";
458	+	fcan->cd() ; gPad->SetLeftMargin(0.15) ; frame->GetYaxis()->SetTitleOffset(1.4) ; frame2->Draw();
459	+
460	+	if(h->GetName()=="hdata") DrawPrelim();
461	+	else DrawMCPrelim();
462	+
463	+	CompleteSave(fcan, ContainerName+"/Fit/"+saveasFinal);
464	+
465	+	delete hcopy;
466	+	delete fcan;
467	+	delete frame;
468	+	delete thePeak;
469	+	}
470	+
471		Value get_Zb_data_over_mc(string ContainerName, ZbResultContainer &res, string variable, string variable2, TCut cut, string saveas, float pt1, float pt2, float a1, float a2) {
472		//  write_warning(__FUNCTION__,"Debugging this function, therefore always returning 3.1415+/-1.010101"); return Value(3.1415,1.010101);
473		TCanvas *cn = new TCanvas("cn","cn");
#	Line 378 | Line 475 | Value get_Zb_data_over_mc(string Contain
475
476		bool DoFit=false;
477		bool UseFit=false;
478	+	bool DoRooFit=false;
479
480
481		if(Contains(variable,"pfJetGood")) varname="R_{abs}";
482		//  TH1F *hdata  = allsamples.Draw("hdata",variable,1200,0,30, varname, "events",  cut,data,luminosity);
483		//  TH1F *hmc    = allsamples.Draw("hmc" , variable,1200,0,30, varname, "events", cut, mc,  luminosity);
484	<	TH1F *hdata  = allsamples.Draw("hdata",variable,1200,0.5,1.5, varname, "events",  cut,data,luminosity); // making it more precise
485	<	TH1F *hmc    = allsamples.Draw("hmc" , variable,1200,0.5,1.5, varname, "events", cut, mc,  luminosity); // making it more precise
484	>	TH1F *hdata, hmc;
485	>	if(!DoRooFit) {
486	>	hdata = allsamples.Draw("hdata",variable,1200,ZbData::LowerHistoBoundary,ZbData::UpperHistoBoundary, varname, "events",  cut,data,luminosity); // making it more precise
487	>	hmc    = allsamples.Draw("hmc" , variable,1200,ZbData::LowerHistoBoundary,ZbData::UpperHistoBoundary, varname, "events", cut, mc,  luminosity); // making it more precise
488	>	}
489	>
490		hmc->SetLineColor(kBlue);
491
492		float a=hdata->GetMean();
#	Line 407 | Line 509 | Value get_Zb_data_over_mc(string Contain
509		gaus_mc->SetParameter(2,hmc->GetRMS());
510		hmc->Fit(gaus_mc);
511		}
512	<
512	>	if(DoRooFit) {
513	>	hdata  = allsamples.Draw("hdata",variable,7200,0.0,6.0, varname, "events",  cut,data,luminosity); // making it more precise
514	>	hmc    = allsamples.Draw("hmc" , variable,7200,0.0,6.0, varname, "events", cut, mc,  luminosity); // making it more precise
515	>	if(hdata->Integral()>10) {
516	>	make_RooFit(ContainerName, hdata, a, da, saveas+"_data");
517	>	make_RooFit(ContainerName, hmc, b, db, saveas+"_mc");
518	>	factor = a / b;
519	>	error = TMath::Sqrt(  (1/(b*b)) * (da*da) + ((a*a)/(b*b))*db*db);
520	>	} else {
521	>	factor=NAN;//set to nan
522	>	error=NAN;//set to nan
523	>	}
524	>	}
525	>
526		cn->cd();
527	<	hdata->Rebin(0.05*1200/(1.5-0.5));
528	<	hmc->Rebin(0.05*1200/(1.5-0.5));
527	>	hdata->Rebin(int(0.05*1200/(ZbData::UpperHistoBoundary-ZbData::LowerHistoBoundary)));
528	>	hmc->Rebin(int(0.05*1200/(ZbData::UpperHistoBoundary-ZbData::LowerHistoBoundary)));
529		hdata->GetYaxis()->SetTitle("events / 0.05");
530		hdata->SetMarkerColor(kRed);
531		if(DoFit) {
#	Line 463 | Line 578 | Value get_Zb_data_over_mc(string Contain
578		CompleteSave(cn,ContainerName+"/ResponsePlots/"+saveas);
579
580		res.MPF_Rabs_location.push_back(ZbPointResult(variable,pt1,pt2,a1,a2,Value(a,da),Value(b,db)));
581	<	cout << "Data;" << a << ";"<<da<<";MC;"<<b<<";"<<db<<endl;
581	>	dout << "Data;" << a << ";"<<da<<";MC;"<<b<<";"<<db<<endl;
582		delete cn;
583		delete hdata;
584		delete hmc;
#	Line 479 | Line 594 | float PurityMedium=0.54;
594		float PurityTight=0.82;
595
596
597	+	void GetPtPointPosition(float AvPos[], const int nptcuts,string AlphaVariable, TCut basecut, float FaceValueToBeTakenAt, float ptcuts[]) {
598	+	for(int i=0;i<nptcuts-1;i++) {
599	+	stringstream specialcut;
600	+	specialcut << "((pt>30&& pt< 1000) && (" << AlphaVariable << "<" << FaceValueToBeTakenAt << "))";
601	+	TH1F *hdata  = allsamples.Draw("hada", "pt",100,ptcuts[i],ptcuts[i+1], "MPF", "events", TCut(basecut && specialcut.str().c_str()), data,  luminosity);
602	+	AvPos[i]=hdata->GetMean();
603	+	dout << AvPos[i] << " : ";
604	+	delete hdata;
605	+	}
606	+	dout << endl;
607	+	}
608	+
609	+
610		void DeterminePurity(TCut basecut, int nptcuts, float ptcuts[], string AlphaVariable, string ContainerName, float FaceValueToBeTakenAt) {
611
612		bool ComputePurity=false;
#	Line 515 | Line 643 | void DeterminePurity(TCut basecut, int n
643		}
644
645		float purity=Zb/total;
646	<	dout << "Determined a purity of " << 100*purity << " % for the container " << ContainerName << endl;
646	>	dout << "Determined a purity of " << 100*purity << " % for the container \"" << ContainerName << "\"" << endl;
647
648		if(ContainerName=="inclusive") PurityInclusive=purity;
649		if(ContainerName=="loose") PurityLoose=purity;
#	Line 534 | Line 662 | string NiceAlphaPrint(float alpha) {
662
663		ZbResultContainer new_data_mc_agreement_2d(bool gofast, string AlphaVariable, string ContainerName, TCut BTagCut, TCut BTagWgt, float FaceValueToBeTakenAt) {
664
537	–
665		if(gofast) write_info(__FUNCTION__,"Fast mode activated for "+ContainerName);
666		ZbResultContainer results(ContainerName);
667
#	Line 550 | Line 677 | ZbResultContainer new_data_mc_agreement_
677
678		const int nptcuts=5;
679		float ptcuts[nptcuts]={30,50,75,125,1000};
680	+	float ptPointPosition[nptcuts];
681	+	GetPtPointPosition(ptPointPosition,nptcuts,AlphaVariable,basecut,FaceValueToBeTakenAt,ptcuts);
682	+
683
684		DeterminePurity(basecut, nptcuts,ptcuts, AlphaVariable, ContainerName, FaceValueToBeTakenAt);
685
#	Line 572 | Line 702 | ZbResultContainer new_data_mc_agreement_
702		alow=0;
703		ahigh=alphacuts[ia];
704		}
705	<	cout << specialcut.str() << endl;
705	>	dout << specialcut.str() << endl;
706		Value MPF_data_over_mc = get_Zb_data_over_mc(ContainerName,results,"mpf","",TCut(basecut && specialcut.str().c_str()),"MPF___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__alpha_"+any2string(NiceAlphaPrint(alphacuts[ia])), ptcuts[ipt], ptcuts[ipt+1],alow,ahigh);
707		Value RABS_data_over_mc = get_Zb_data_over_mc(ContainerName,results,"ZbCHS3010_pfJetGoodPt[0]/pt","",TCut(basecut && specialcut.str().c_str()),"Rabs___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__alpha_"+any2string(NiceAlphaPrint(alphacuts[ia])), ptcuts[ipt], ptcuts[ipt+1],alow,ahigh);
708
#	Line 582 | Line 712 | ZbResultContainer new_data_mc_agreement_
712		RABS_Results[ia][ipt]=RABS_data_over_mc.getValue();
713		RABS_Errors[ia][ipt]=RABS_data_over_mc.getError();
714
715	<	cout << "alpha cut at " << alphacuts[ia+1] << " and pt range " << ptcuts[ipt] << " , " << ptcuts[ipt+1] << " \t : \t" << specialcut.str() << " \t \t --> " << MPF_data_over_mc << " (MPF) " << RABS_data_over_mc << " (RABS)" << endl;
715	>	dout << "alpha cut at " << alphacuts[ia+1] << " and pt range " << ptcuts[ipt] << " , " << ptcuts[ipt+1] << " \t : \t" << specialcut.str() << " \t \t --> " << MPF_data_over_mc << " (MPF) " << RABS_data_over_mc << " (RABS)" << endl;
716		}
717		}
718
#	Line 601 | Line 731 | ZbResultContainer new_data_mc_agreement_
731
732		stringstream filename;
733		filename << "Extrapolation/Zplusb_data_over_mc___" << ptcuts[ipt] << "_to_" << ptcuts[ipt+1];
734	<	cout << "Going to create histo for " << filename.str() << endl;
734	>	dout << "Going to create histo for " << filename.str() << endl;
735		TGraphErrors *mpf_gr =new TGraphErrors();
736		TGraphErrors *rabs_gr =new TGraphErrors();
737
#	Line 609 | Line 739 | ZbResultContainer new_data_mc_agreement_
739		int mpf_counter=0;
740
741		for(int ia=0;ia<nalphacuts && !gofast;ia++) {
742	<	cout << "Setting a point for an alpha cut at " << alphacuts[ia] << " with value " << MPF_Results[ia][ipt] << " (MPF) and " << RABS_Results[ia][ipt] << " (RABS)" << endl;
742	>	dout << "Setting a point for an alpha cut at " << alphacuts[ia] << " with value " << MPF_Results[ia][ipt] << " (MPF) and " << RABS_Results[ia][ipt] << " (RABS)" << endl;
743		if(MPF_Results[ia][ipt]==MPF_Results[ia][ipt] && MPF_Errors[ia][ipt]==MPF_Errors[ia][ipt]) { // remove nan's
744		mpf_gr->SetPoint(mpf_counter,alphacuts[ia],MPF_Results[ia][ipt]);
745		mpf_gr->SetPointError(mpf_counter,0,MPF_Errors[ia][ipt]);
#	Line 617 | Line 747 | ZbResultContainer new_data_mc_agreement_
747		} else {
748		dout << "Detected BS for MPF method (bin ignored) : " << endl;
749		dout << "       alpha: " << alphacuts[ia] << " , pt range: " << ptcuts[ipt] << " < pt < " << ptcuts[ipt+1] << endl;
750	<	cout << "       value (MPF) : " << MPF_Results[ia][ipt] << " +/- " << MPF_Errors[ia][ipt] << endl;
750	>	dout << "       value (MPF) : " << MPF_Results[ia][ipt] << " +/- " << MPF_Errors[ia][ipt] << endl;
751		}
752		if(RABS_Results[ia][ipt]==RABS_Results[ia][ipt] && RABS_Errors[ia][ipt]==RABS_Errors[ia][ipt]) {
753		rabs_gr->SetPoint(rabs_counter,alphacuts[ia],RABS_Results[ia][ipt]);
#	Line 626 | Line 756 | ZbResultContainer new_data_mc_agreement_
756		} else {
757		dout << "Detected BS for RABS method (bin ignored) : " << endl;
758		dout << "       alpha: " << alphacuts[ia] << " , pt range: " << ptcuts[ipt] << " < pt < " << ptcuts[ipt+1] << endl;
759	<	cout << "       value (RABS) : " << RABS_Results[ia][ipt] << " +/- " << RABS_Errors[ia][ipt] << endl;
759	>	dout << "       value (RABS) : " << RABS_Results[ia][ipt] << " +/- " << RABS_Errors[ia][ipt] << endl;
760		}
761		}
762
#	Line 634 | Line 764 | ZbResultContainer new_data_mc_agreement_
764		specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (" << AlphaVariable << "<" << FaceValueToBeTakenAt << "))";
765		Value MPF_data_over_mc = get_Zb_data_over_mc(ContainerName,results,"mpf","",TCut(basecut && specialcut.str().c_str()),"MPF___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__INCLUSIVEalpha_"+any2string(FaceValueToBeTakenAt),ptcuts[ipt],ptcuts[ipt+1],0,FaceValueToBeTakenAt);
766		Value RABS_data_over_mc = get_Zb_data_over_mc(ContainerName,results,"ZbCHS3010_pfJetGoodPt[0]/pt","",TCut(basecut && specialcut.str().c_str()),"Rabs___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__INCLUSIVEalpha_"+any2string(FaceValueToBeTakenAt),ptcuts[ipt],ptcuts[ipt+1],0,FaceValueToBeTakenAt);
767	<	cout << "About to add a point in pt (" << (ptcuts[ipt]+ptcuts[ipt+1])/2 << " ) to face value : " << MPF_data_over_mc << " (MPF) , " <<  RABS_data_over_mc << "(RABS)" << endl;
767	>	dout << "About to add a point in pt (" << (ptcuts[ipt]+ptcuts[ipt+1])/2 << " ) to face value : " << MPF_data_over_mc << " (MPF) , " <<  RABS_data_over_mc << "(RABS)" << endl;
768		MPF_FaceValueAtPoint3->SetPoint(ipt,(ptcuts[ipt]+ptcuts[ipt+1])/2,MPF_data_over_mc.getValue());
769		MPF_FaceValueAtPoint3->SetPointError(ipt,0,MPF_data_over_mc.getError());
770		RABS_FaceValueAtPoint3->SetPoint(ipt,(ptcuts[ipt]+ptcuts[ipt+1])/2,RABS_data_over_mc.getValue());
771		RABS_FaceValueAtPoint3->SetPointError(ipt,0,RABS_data_over_mc.getError());
772
773
774	+	TH1D *mpf_band, *rabs_band;
775
776		if(!gofast) {
777		can->cd();
778		mpf_gr->SetMarkerStyle(21);
779		mpf_gr->SetMarkerColor(kBlue);
780		mpf_gr->Draw("AP");
781	<	mpf_gr->Fit("pol1","");
781	>
782	>	TF1 *mpf_pol = new TF1("mpf_pol","[0]+[1]*x",ZbData::LowerHistoBoundary,ZbData::UpperHistoBoundary);
783	>
784	>	mpf_gr->Fit(mpf_pol,"E"); // better error estimation using MINOS (E), quiet (Q) and using loglikelihood method
785		mpf_gr->GetXaxis()->SetTitle("#alpha");
786		mpf_gr->GetXaxis()->CenterTitle();
787		mpf_gr->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
788		mpf_gr->GetYaxis()->CenterTitle();
789	<	TF1 *mpf_pol=(TF1*)mpf_gr->GetFunction("pol1");
789	>	//TF1 *mpf_pol=(TF1*)mpf_gr->GetFunction("pol1");
790	>	mpf_band = getBand(mpf_pol,"MPF_band");
791		mpf_pol->SetLineWidth(0);
792		TF1 *mpf_pol_complete = new TF1("mpf_pol_complete","[0]+[1]*x");
793		mpf_pol_complete->SetParameters(mpf_pol->GetParameters());
#	Line 672 | Line 807 | ZbResultContainer new_data_mc_agreement_
807		histo->GetYaxis()->SetRangeUser(min,max);
808		mpf_gr->GetYaxis()->SetRangeUser(min,max);
809
810	<	TPolyLine *mpf_fit_uncert = GetFitUncertaintyShape(mpf_gr, "pol1", min, max,0.0,0.4);
810	>	//TPolyLine *mpf_fit_uncert = GetFitUncertaintyShape(mpf_gr, "pol1", min, max,0.0,0.4);
811	>
812		mpf_pol->SetLineColor(TColor::GetColor("#0101DF"));
813	<	mpf_fit_uncert->SetFillColor(TColor::GetColor("#5353FC"));
813	>	//mpf_fit_uncert->SetFillColor(TColor::GetColor("#A2A2FA"));
814	>	mpf_band->SetFillColor(TColor::GetColor("#A2A2FA"));
815		histo->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
816		histo->Draw();
817	<	mpf_fit_uncert->Draw("F");
818	<	mpf_fit_uncert->Draw("");
817	>	//mpf_fit_uncert->Draw("F");
818	>	//mpf_fit_uncert->Draw("");
819	>	mpf_gr->Draw("P");
820	>	mpf_band->Draw("C E3 SAME");
821		mpf_gr->Draw("P");
822		mpf_pol_complete->Draw("same");
823
#	Line 702 | Line 841 | ZbResultContainer new_data_mc_agreement_
841		filename << "__MPF";
842		DrawPrelim();
843		CompleteSave(can,ContainerName+"/"+filename.str());
844	<	cout << "MPF : " << mpf_a << " + " << mpf_b << " * alpha " << endl;
844	>	dout << "MPF : " << mpf_a << " + " << mpf_b << " * alpha " << endl;
845
846		filename.str("");
847
848	+	TF1 *rabs_pol = new TF1("rabs_pol","[0]+[1]*x",ZbData::LowerHistoBoundary,ZbData::UpperHistoBoundary);
849		rabs_gr->SetMarkerStyle(21);
850		rabs_gr->SetMarkerColor(kRed);
851		rabs_gr->Draw("AP");
852	<	rabs_gr->Fit("pol1","");
852	>	//rabs_gr->Fit("pol1","E");
853	>	rabs_gr->Fit(rabs_pol,"E");
854	>	rabs_band = getBand(rabs_pol,"Rabs_band");
855	>
856	>
857		rabs_gr->GetXaxis()->SetTitle("#alpha");
858		rabs_gr->GetXaxis()->CenterTitle();
859		rabs_gr->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc}");
860		rabs_gr->GetYaxis()->CenterTitle();
861	<	TF1 *rabs_pol=(TF1*)rabs_gr->GetFunction("pol1");
861	>	//TF1 *rabs_pol=(TF1*)rabs_gr->GetFunction("pol1");
862		rabs_pol->SetLineWidth(0);
863		TF1 *rabs_pol_complete = new TF1("rabs_pol_complete","[0]+[1]*x");
864		rabs_pol_complete->SetParameters(rabs_pol->GetParameters());
#	Line 723 | Line 867 | ZbResultContainer new_data_mc_agreement_
867		FindMinMax(rabs_gr,min,max);
868		rabs_gr->GetYaxis()->SetRangeUser(min,max);
869		histo->GetYaxis()->SetRangeUser(min,max);
870	<	TPolyLine *rabs_fit_uncert = GetFitUncertaintyShape(rabs_gr, "pol1", min, max,0.0,0.4);
871	<	rabs_pol->SetLineColor(TColor::GetColor("#FA5858"));
872	<	rabs_pol->SetFillColor(TColor::GetColor("#FA5858"));
873	<	rabs_fit_uncert->SetLineColor(TColor::GetColor("#FA5858"));
874	<	rabs_fit_uncert->SetFillColor(TColor::GetColor("#FA5858"));
870	>	//TPolyLine *rabs_fit_uncert = GetFitUncertaintyShape(rabs_gr, "pol1", min, max,0.0,0.4);
871	>	rabs_pol->SetLineColor(TColor::GetColor("#F99C9C"));
872	>	rabs_pol->SetFillColor(TColor::GetColor("#F99C9C"));
873	>	rabs_band->SetLineColor(TColor::GetColor("#F99C9C"));
874	>	rabs_band->SetFillColor(TColor::GetColor("#F99C9C"));
875	>	//rabs_fit_uncert->SetLineColor(TColor::GetColor("#F99C9C"));
876	>	//rabs_fit_uncert->SetFillColor(TColor::GetColor("#F99C9C"));
877		histo->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc}");
878		histo->Draw();
879	<	rabs_fit_uncert->Draw("F");
880	<	rabs_fit_uncert->Draw("");
879	>	rabs_band->Draw("C E3 SAME");
880	>	//rabs_fit_uncert->Draw("F");
881	>	//rabs_fit_uncert->Draw("");
882		rabs_gr->Draw("P");
883		rabs_pol_complete->Draw("same");
884
#	Line 740 | Line 887 | ZbResultContainer new_data_mc_agreement_
887		RABS_FinalGraph->SetPoint(ipt,0.5*(ptcuts[ipt]+ptcuts[ipt+1]),rabs_a);
888		RABS_FinalGraph->SetPointError(ipt,0,rabs_pol->GetParError(0));
889
890	<	cout << "!+!+!+!+!+!+!+!+!+ Just added a point to the final plot : " << rabs_a << " +/- " << rabs_pol->GetParError(0) << endl;
890	>	dout << "!+!+!+!+!+!+!+!+!+ Just added a point to the final plot : " << rabs_a << " +/- " << rabs_pol->GetParError(0) << endl;
891
892		stringstream rabs_info;
893		rabs_info << "#splitline{#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (R_{abs}) }}{Extrapolated value: " << std::setprecision(4) << rabs_a << "}}{#chi^{2}/NDF : " << rabs_pol->GetChisquare() << " / " << rabs_pol->GetNDF() << "}";
#	Line 754 | Line 901 | ZbResultContainer new_data_mc_agreement_
901		filename << filenamebkp << "__RABS";
902		DrawPrelim();
903		CompleteSave(can,ContainerName+"/"+filename.str());
904	<	cout << "RABS : " << rabs_a << " + " << rabs_b << " * alpha " << endl;
904	>	dout << "RABS : " << rabs_a << " + " << rabs_b << " * alpha " << endl;
905		}
906		}
907
#	Line 890 | Line 1037 | ZbResultContainer new_data_mc_agreement_
1037		filename2 << "Extrapolation/FINAL_RESULT_RABS_FaceValp3";
1038		CompleteSave(can,ContainerName+"/"+filename2.str());
1039
1040	<	cout << "FINAL RESULTS: " << endl;
1041	<	if(!gofast) cout << "MPF : " << MPFResult  << " +/- " << MPFResultError  << endl;
1042	<	if(!gofast) cout << "Rabs: " << RabsResult << " +/- " << RabsResultError << endl;
1043	<	cout << "Face value at " << FaceValueToBeTakenAt << " :" << FaceValResult << " +/- " << FaceValResultError << "  (MPF) " << endl;
1044	<	cout << "Face value at " << FaceValueToBeTakenAt << " :" << RABSFaceValResult << " +/- " << RABSFaceValResultError << "  (RABS) " << endl;
1040	>	dout << "FINAL RESULTS: " << endl;
1041	>	if(!gofast) dout << "MPF : " << MPFResult  << " +/- " << MPFResultError  << endl;
1042	>	if(!gofast) dout << "Rabs: " << RabsResult << " +/- " << RabsResultError << endl;
1043	>	dout << "Face value at " << FaceValueToBeTakenAt << " :" << FaceValResult << " +/- " << FaceValResultError << "  (MPF) " << endl;
1044	>	dout << "Face value at " << FaceValueToBeTakenAt << " :" << RABSFaceValResult << " +/- " << RABSFaceValResultError << "  (RABS) " << endl;
1045
1046		delete can;
1047
#	Line 902 | Line 1049 | ZbResultContainer new_data_mc_agreement_
1049		results.Rabs_Result_FaceValue=Value(RABSFaceValResult,RABSFaceValResultError);
1050		if(!gofast) results.MPF_Result_Extrapolated=Value(MPFResult,MPFResultError);
1051		if(!gofast) results.Rabs_Result_Extrapolated=Value(RabsResult,RabsResultError);
1052	<	cout << results << endl;
1052	>	results.Print();
1053		return results;
1054		}
1055
#	Line 914 | Line 1061 | void DoPUStudy(string identifier) {
1061		TCanvas *can = new TCanvas("can","can");
1062		TH1F *malpha[8];
1063		TH1F *dalpha[8];
1064	<	cout << identifier << ";";
1065	<	cout << endl;
1064	>	dout << identifier << ";";
1065	>	dout << endl;
1066		for(int i=1;i<5;i++) {
1067		stringstream sSpecialCut;
1068		if(i<4) {
1069		sSpecialCut << "numVtx>" << numVtxcuts[i-1] << " && numVtx<" << numVtxcuts[i];
1070	<	cout << numVtxcuts[i-1] << "  < nVtx < " << numVtxcuts[i] << ";";
1070	>	dout << numVtxcuts[i-1] << "  < nVtx < " << numVtxcuts[i] << ";";
1071		} else {
1072		sSpecialCut << "numVtx>" << numVtxcuts[i-1];
1073	<	cout << numVtxcuts[i-1] << ";";
1073	>	dout << numVtxcuts[i-1] << ";";
1074		}
1075
1076		sSpecialCut << " && " << identifier << "_alpha<0.3";
#	Line 941 | Line 1088 | void DoPUStudy(string identifier) {
1088		float factor = a / b;
1089		float error = TMath::Sqrt(  (1/(b*b)) * (da*da) + ((a*a)/(b*b))*db*db);
1090
1091	<	cout << malpha[i]->GetMean() << ";" << malpha[i]->GetMeanError() << ";" << dalpha[i]->GetMean() << ";" << dalpha[i]->GetMeanError() << ";" << malpha[i]->Integral()<<";"<<dalpha[i]->Integral() << ";"<<factor << ";" << error <<  endl;
1091	>	dout << malpha[i]->GetMean() << ";" << malpha[i]->GetMeanError() << ";" << dalpha[i]->GetMean() << ";" << dalpha[i]->GetMeanError() << ";" << malpha[i]->Integral()<<";"<<dalpha[i]->Integral() << ";"<<factor << ";" << error <<  endl;
1092		delete malpha[i];
1093		delete dalpha[i];
1094		}
#	Line 1022 | Line 1169 | void ScenarioComparison() {
1169		}
1170
1171
1172	+	gStyle->SetOptFit(0);
1173		gr[i]->SetTitle();
1174		gr[i]->GetXaxis()->SetTitle("N(Vtx)");
1175		gr[i]->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
#	Line 1154 | Line 1302 | void ScenarioComparisonInclusive() {
1302		void compare_selection(string identifier) {
1303		bool recognized_scenario=false;
1304
1305	<	cout << "Running with identifier " << identifier << endl;
1305	>	dout << "Running with identifier " << identifier << endl;
1306		if(identifier=="Zb30") {
1307		LeadingB=TCut ("Zb30_bTagProbCSVBP[0]>0.898");
1308		EtaB=TCut("abs(Zb30_pfJetGoodEta[0])<1.3");
#	Line 1266 | Line 1414 | void compare_selection(string identifier
1414
1415		assert(recognized_scenario);
1416
1417	<	cout << "Cuts: " << endl;
1418	<	cout << "   " << (const char*) LeadingB << endl;
1419	<	cout << "   " << (const char*) EtaB << endl;
1420	<	cout << "   " << (const char*) PhiZcut << endl;
1417	>	dout << "Cuts: " << endl;
1418	>	dout << "   " << (const char*) LeadingB << endl;
1419	>	dout << "   " << (const char*) EtaB << endl;
1420	>	dout << "   " << (const char*) PhiZcut << endl;
1421
1422	<	cout << endl << endl;
1422	>	dout << endl << endl;
1423
1424		//  ScenarioComparison();
1425		//  ScenarioComparisonInclusive();
#	Line 1326 | Line 1474 | void TEST_GetNumberEventsInsideOutsideAl
1474		//   GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2&&id1==1"), "mm");
1475		//   GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2"), "ee/mm");
1476
1477	<	cout << "BASE SELECTION" << endl;
1477	>	dout << "BASE SELECTION" << endl;
1478		GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2&&id1==0"), "ee");
1479		GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2&&id1==1"), "mm");
1480		GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2"), "ee/mm");
1481
1482	<	cout << "BASE SELECTION: Z+b" << endl;
1482	>	dout << "BASE SELECTION: Z+b" << endl;
1483		GetNumberEventsInsideOutsideAlphaWindow(ZplusBsel&&TCut("id1==id2&&id1==0"), "ee");
1484		GetNumberEventsInsideOutsideAlphaWindow(ZplusBsel&&TCut("id1==id2&&id1==1"), "mm");
1485		GetNumberEventsInsideOutsideAlphaWindow(ZplusBsel&&TCut("id1==id2"), "ee/mm");
1486
1487	<	cout << "Leading B" << endl;
1487	>	dout << "Leading B" << endl;
1488		GetNumberEventsInsideOutsideAlphaWindow(LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==0"), "ee");
1489		GetNumberEventsInsideOutsideAlphaWindow(LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==1"), "mm");
1490		GetNumberEventsInsideOutsideAlphaWindow(LeadingB&&ZplusBsel&&TCut("id1==id2"), "ee/mm");
1491
1492	<	cout << "PhiZcut" << endl;
1492	>	dout << "PhiZcut" << endl;
1493		GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==0"), "ee");
1494		GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==1"), "mm");
1495		GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2"), "ee/mm");
1496
1497	<	cout << "alpha cut" << endl;
1497	>	dout << "alpha cut" << endl;
1498		GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==0")&&TCut("ZbCHS3010_alpha<0.3"), "ee");
1499		GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==1")&&TCut("ZbCHS3010_alpha<0.3"), "mm");
1500		GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2")&&TCut("ZbCHS3010_alpha<0.3"), "ee/mm");
#	Line 1357 | Line 1505 | void TEST_GetNumberEventsInsideOutsideAl
1505		}
1506
1507
1508	<
1508	>	void MetSpectrumDepletionIllustration(int isdata) {
1509	>
1510	>	TCanvas *can = new TCanvas("can","can");
1511	>	TCut basecut(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
1512	>
1513	>	cutWeight=STDWEIGHT*TCut("ZbCHS3010_BTagWgtTUp*(1.0)*(ZbCHS3010_bTagProbCSVBP[0]>0.898)");
1514	>	TH1F *standard  = allsamples.Draw("standard","met[1]",80,0,400, "MET [GeV]", "events", basecut,isdata,luminosity);
1515	>
1516	>	cutWeight=STDWEIGHT*TCut("ZbCHS3010_BTagWgtTUp*(1.0*softMuon)*(ZbCHS3010_bTagProbCSVBP[0]>0.898)");
1517	>	TH1F *Neutrinos  = allsamples.Draw("Neutrinos","met[1]",80,0,400, "MET [GeV]", "events", basecut,isdata,luminosity);
1518	>
1519	>	cutWeight=STDWEIGHT*TCut("ZbCHS3010_BTagWgtTUp*(1.0*(1.0-softMuon))*(ZbCHS3010_bTagProbCSVBP[0]>0.898)");
1520	>	TH1F *NoNeutrinos  = allsamples.Draw("NoNeutrinos","met[1]",80,0,400, "MET [GeV]", "events", basecut,isdata,luminosity);
1521	>
1522	>	standard->Scale(1.0/standard->Integral());
1523	>	Neutrinos->Scale(1.0/Neutrinos->Integral());
1524	>	NoNeutrinos->Scale(1.0/NoNeutrinos->Integral());
1525	>
1526	>	Neutrinos->SetLineColor(kGreen);
1527	>	NoNeutrinos->SetLineColor(kRed);
1528	>
1529	>	can->SetLogy(1);
1530	>	TLegend *leg = make_legend();
1531	>	leg->AddEntry(standard,"Standard sample","ep");
1532	>	leg->AddEntry(Neutrinos,"#nu enriched sample","l");
1533	>	leg->AddEntry(NoNeutrinos,"#nu depleted sample","l");
1534	>	leg->SetX1(0.9*leg->GetX1());
1535	>	standard->Draw("e1");
1536	>	Neutrinos->Draw("histo,same");
1537	>	NoNeutrinos->Draw("histo,same");
1538	>	if(isdata==data) leg->SetHeader("Data:");
1539	>	else leg->SetHeader("MC:");
1540	>	leg->Draw();
1541	>	DrawPrelim();
1542	>
1543	>	if(isdata==data) CompleteSave(can,"NeutrinoMetSpectrum__Data");
1544	>	else CompleteSave(can,"NeutrinoMetSpectrum__MC");
1545	>	delete standard;
1546	>	delete Neutrinos;
1547	>	delete NoNeutrinos;
1548	>	delete can;
1549	>	CleanLegends();
1550	>	}
1551
1552		void do_basic_ZB_analysis(bool DoCompleteAnalysis) {
1553		STDWEIGHT=TCut(cutWeight);
1554		cutWeight=TCut(STDWEIGHT*TCut("ZbCHS3010_BTagWgtT"));
1555	<	cout << "The lepton requirement is " << (const char*) leptoncut << endl;
1555	>	dout << "The lepton requirement is " << (const char*) leptoncut << endl;
1556		bool doquick=false;
1557
1558		TCanvas *zbcanvas = new TCanvas("zbcanvas","zbcanvas");
1559
1560	<	//   ScenarioComparison();
1561	<	//   ScenarioComparisonInclusive();
1562	<
1560	>
1561	>	//   MetSpectrumDepletionIllustration(data);
1562	>	//   MetSpectrumDepletionIllustration(mc);
1563	>
1564		//   compare_selection("ZbCHS3010_p5Clean");
1565		//   compare_selection("ZbCHS3010");
1566		//   compare_selection("Zb30");
#	Line 1388 | Line 1579 | void do_basic_ZB_analysis(bool DoComplet
1579		//   compare_selection("Zb40");
1580
1581		if(!doquick) {
1582	<	//draw_kin_variable("mll","",28,60.0,200.0,"m_{ll}","mll",1,"");//nice and quick test to see if the normalization is ok (i.e. all data is processed etc.)
1583	<	//print_all_b_yields();
1582	>	//    draw_kin_variable("mll","",28,60.0,200.0,"m_{ll}","mll",1,"");//nice and quick test to see if the normalization is ok (i.e. all data is processed etc.)
1583	>	//    print_all_b_yields();
1584		//     draw_mpf_vars();
1585	<	//DrawEvilCutFlow();
1585	>	//    DrawEvilCutFlow();
1586
1587		//     draw_Zb_kin_vars();
1588
#	Line 1410 | Line 1601 | void do_basic_ZB_analysis(bool DoComplet
1601
1602		bool fast=true;bool slow=false;//don't change these
1603
1604	<	ZbResultContainer inclusive = new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","inclusive",TCut("mll>0"),TCut("1.0"),0.3);
1604	>	ZbResultContainer inclusive = new_data_mc_agreement_2d(slow,"ZbCHS3010_alpha","inclusive",TCut("mll>0"),TCut("1.0"),0.3);
1605		ZbResultContainer Reference = new_data_mc_agreement_2d(slow,"ZbCHS3010_alpha","reference",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),TCut("ZbCHS3010_BTagWgtT"),0.3);
1606
1607	+
1608		ZbResultContainer NeutrinoQ = new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","NeutrinoQ",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),TCut("ZbCHS3010_BTagWgtTUp*(1.0*softMuon)"),0.3);
1609		ZbResultContainer ANeutrinoQ = new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","ANeutrinoQ",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),TCut("ZbCHS3010_BTagWgtTUp*(1.0*(1.0-softMuon))"),0.3);
1610

Diff Legend

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