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.9 by buchmann, Wed Nov 14 10:46:45 2012 UTC vs.
Revision 1.14 by buchmann, Mon Nov 26 12:06:26 2012 UTC

#	Line 23 | Line 23 | using namespace PlottingSetup;
23		namespace ZbData {
24		vector<float> data_over_mc;
25
26	<	TCut ZplusBsel("pt1>20&&pt2>20&&mll>60&&mll<120");
26	>	TCut ZplusBsel("pt1>20&&pt2>20&&mll>60&&mll<120&&id1==id2");
27		TCut LeadingB("Zb3010_bTagProbCSVBP[0]>0.898");
28		TCut EtaB("abs(Zb3010_pfJetGoodEta[0])<1.3");
29	<	TCut PhiZcut("abs(pfJetDphiZ[0])>2.7");
29	>	TCut PhiZcut("abs(Zb3010_pfJetDphiZ[0])>2.7");
30
31		}
32
#	Line 120 | Line 120 | void print_all_b_yields() {
120		print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
121		cout << "10<ptZ<1000 GeV" << endl;
122		print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000");
123	+	cout << "#alpha < 0.35" << endl;
124	+	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.35");
125		cout << "#alpha < 0.3" << endl;
126		print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.3");
127		cout << "#alpha < 0.2" << endl;
#	Line 220 | Line 222 | void draw_Zb_kin_vars() {
222		draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>125&&pt<1000",40,0,2,"#alpha","Official/alpha_pt_125_to_1000",1);
223		draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000",40,0,2,"#alpha","Official/alpha_pt_30_to_1000",1);
224
225	<	draw_kin_variable("pfBJetDphiZ[0]",ZplusBsel&&LeadingB&&"pt>10&&pt<1000&&pfBJetDphiZ[0]>0",30,0,3.2,"#delta#phi (Z,b lead)","DeltaPhiZBlead",0);
225	>	draw_kin_variable("Zb3010_pfBJetDphiZ[0]",ZplusBsel&&LeadingB&&"pt>10&&pt<1000&&pfBJetDphiZ[0]>0",30,0,3.2,"#delta#phi (Z,b lead)","DeltaPhiZBlead",0);
226
227		/*
228		draw_kin_variable("pt1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,200,"p_{T}^{l1}","Lep/pt1",1);
#	Line 254 | Line 256 | void draw_mpf_vars() {
256		draw_kin_variable("Zb3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.05",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz__mpf_alpha_smaller_0p05",0,"#alpha<0.05, 10 GeV < p_{T}^{Z} < 1000 GeV");
257		}
258
259	<	Value get_Zb_data_over_mc(string variable, TCut cut) {
259	>	/*double GetMedian(TH1F *histo) {
260	>	int numBins = histo->GetXaxis()->GetNbins();
261	>	Double_t x[numBins];
262	>	Double_t y[numBins];
263	>	for (int i = 1; i <= numBins; i++) {
264	>	x[i] = histo->GetBinCenter(i);
265	>	y[i] = histo->GetBinContent(i);
266	>	}
267	>	return TMath::Median(numBins, &x, &y);
268	>	}*/
269	>
270	>	Value get_Zb_data_over_mc(string variable, string variable2, TCut cut, string saveas) {
271		//  write_warning(__FUNCTION__,"Debugging this function, therefore always returning 3.1415+/-1.010101"); return Value(3.1415,1.010101);
272	<	TH1F *hdata  = allsamples.Draw("hdata",variable,1000,0,30, "#alpha [GeV]", "events", cut,data,luminosity);
273	<	TH1F *hmc    = allsamples.Draw("hmc" , variable,1000,0,30, "#alpha [GeV]", "events", cut,mc,  luminosity);
272	>	TCanvas *cn = new TCanvas("cn","cn");
273	>	string varname="MPF";
274	>	if(Contains(variable,"pfJetGood")) varname="R_{abs}";
275	>	TH1F *hdata  = allsamples.Draw("hdata",variable,1200,0,30, varname, "events",  cut,data,luminosity);
276	>	TH1F *hmc    = allsamples.Draw("hmc" , variable,1200,0,30, varname, "events", cut, mc,  luminosity);
277	>	hmc->SetLineColor(kBlue);
278
279	<	hdata->Fit("gaus","QE");
280	<	hmc->Fit("gaus","QE");
281	<	TF1 *fdata = (TF1*) hdata->GetFunction("gaus");
282	<	TF1 *fmc   = (TF1*) hmc->GetFunction("gaus");
266	<
267	<	/*  float a=hdata->GetMean();
279	>
280	>
281	>
282	>	float a=hdata->GetMean();
283		float b=hmc->GetMean();
284		float da=hdata->GetMeanError();
285	<	float db=hmc->GetMeanError();*/
271	<	float a = fdata->GetParameter(1);
272	<	float da = fdata->GetParError(1);
273	<	float b = fmc->GetParameter(1);
274	<	float db = fmc->GetParError(1);
285	>	float db=hmc->GetMeanError();
286		float factor = a / b;
287		float error = TMath::Sqrt(  (1/(b*b)) * (da*da) + ((a*a)/(b*b))*db*db);
288	+
289	+	TF1 *gaus_data = new TF1("gaus_data","gaus",0.7,1.3);
290	+	gaus_data->SetParameter(0,hdata->GetMaximum());
291	+	gaus_data->SetParameter(2,hdata->GetMean());
292	+	gaus_data->SetParameter(1,hdata->GetRMS());
293	+	hdata->Fit(gaus_data);
294	+	TF1 *gaus_mc = new TF1("gaus_mc","gaus",0.7,1.3);
295	+	gaus_mc->SetParameter(0,hmc->GetMaximum());
296	+	gaus_mc->SetParameter(1,hmc->GetMean());
297	+	gaus_mc->SetParameter(2,hmc->GetRMS());
298	+	hmc->Fit(gaus_mc);
299	+
300	+
301	+	cn->cd();
302	+	hdata->GetYaxis()->SetTitle("events / 0.1");
303	+	hdata->SetMarkerColor(kRed);
304	+	/*  hdata->Rebin(4);
305	+	hmc->Rebin(4);
306	+
307	+	gaus_mc->SetParameter(0,gaus_mc->GetParameter(0)*4);//rebinning
308	+	gaus_data->SetParameter(0,gaus_data->GetParameter(0)*4);//rebinning
309	+	*/
310	+	gaus_mc->SetLineColor(kBlue);
311	+	gaus_data->SetLineColor(kRed);
312	+
313	+	hdata->GetXaxis()->SetRangeUser(0,2);
314	+	hdata->Draw("e1");
315	+	hmc->Draw("histo,same");
316	+	hdata->Draw("e1,same");
317	+	gaus_data->Draw("same");
318	+	gaus_mc->Draw("same");
319	+	stringstream summary;
320	+	summary << "#splitline{" << varname << "}{#splitline{data: " << std::setprecision(4) << a << " #pm " << std::setprecision(4) << da << "}{";
321	+	summary << "#splitline{mc:   " << std::setprecision(4) << b << " #pm " << std::setprecision(4) << db << "}{#splitline{ratio: " << factor;
322	+	summary << " #pm " << error << "}{#splitline{}{#splitline{data fit: " << std::setprecision(4) << gaus_data->GetParameter(1) << "+/-";
323	+	summary << std::setprecision(4) << gaus_data->GetParError(1) << "}{#splitline{";
324	+	summary << "mc fit:" << std::setprecision(4) << gaus_mc->GetParameter(1) << "+/-" << std::setprecision(4) << gaus_mc->GetParError(1) << "}{ratio: ";
325	+	summary << gaus_data->GetParameter(1)/gaus_mc->GetParameter(1) << "+/-" << gaus_data->GetParameter(1)/gaus_mc->GetParameter(1) * sqrt((gaus_data->GetParError(1) * gaus_data->GetParError(1))/(gaus_data->GetParameter(1)*gaus_data->GetParameter(1)) + (gaus_mc->GetParError(1) * gaus_mc->GetParError(1))/(gaus_mc->GetParameter(1)*gaus_mc->GetParameter(1)) ) << "}}}}}}}";
326	+	TText *infobox = write_title(summary.str());
327	+	infobox->SetX(0.75);
328	+	infobox->SetTextSize(0.03);
329	+	infobox->SetY(0.75);
330	+	infobox->Draw();
331	+
332	+
333	+	DrawPrelim();
334	+	CompleteSave(cn,"ResponsePlots/"+saveas);
335	+
336	+	cout << "Data;" << a << ";"<<da<<";MC"<<b<<";"<<db<<endl;
337	+	delete cn;
338		delete hdata;
339		delete hmc;
340
#	Line 285 | Line 346 | Value get_Zb_data_over_mc(string variabl
346
347
348		void new_data_mc_agreement_2d() {
349	+
350		gStyle->SetOptFit(0);
351		TCut basecut(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
352	<	const int nalphacuts=5;
353	<	float alphacuts[nalphacuts] = {0.1,0.15,0.2,0.25,0.3};
292	<	int PointThree=-1;
293	<	for(int i=0;i<nalphacuts;i++) {
294	<	if((alphacuts[i]-0.3)<0.01) PointThree=i-1;
295	<	}
296	<	if(PointThree<0) {
297	<	write_warning(__FUNCTION__,"alpha=0.3 MUST BE ONE OF THE ALPHA VALUES!");
298	<	assert(PointThree<0);
299	<	}
352	>	const int nalphacuts=6;
353	>	float alphacuts[nalphacuts] = {0.1,0.15,0.2,0.25,0.3,0.35};
354
355		const int nptcuts=5;
356		float ptcuts[nptcuts]={30,50,75,125,1000};
#	Line 310 | Line 364 | void new_data_mc_agreement_2d() {
364		for(int ia=0;ia<nalphacuts;ia++) {
365		for(int ipt=0;ipt<nptcuts-1;ipt++) {
366		stringstream specialcut;
367	<	specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (Zb3010_alpha<" << alphacuts[ia] << "))";
367	>	if(ia>0) specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (Zb3010_alpha<" << alphacuts[ia] << " && Zb3010_alpha>" << alphacuts[ia-1] << "))";
368	>	else     specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (Zb3010_alpha<" << alphacuts[ia] << "))";
369		cout << specialcut.str() << endl;
370	<	Value MPF_data_over_mc = get_Zb_data_over_mc("mpf",TCut(basecut && specialcut.str().c_str()));
371	<	Value RABS_data_over_mc = get_Zb_data_over_mc("Zb3010_pfJetGoodPt[0]/pt",TCut(basecut && specialcut.str().c_str()));
370	>	Value MPF_data_over_mc = get_Zb_data_over_mc("mpf","",TCut(basecut && specialcut.str().c_str()),"MPF___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__alpha_"+any2string(alphacuts[ia]));
371	>	Value RABS_data_over_mc = get_Zb_data_over_mc("Zb3010_pfJetGoodPt[0]/pt","",TCut(basecut && specialcut.str().c_str()),"Rabs___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__alpha_"+any2string(alphacuts[ia]));
372
373		MPF_Results[ia][ipt]=MPF_data_over_mc.getValue();
374		MPF_Errors[ia][ipt]=MPF_data_over_mc.getError();
#	Line 329 | Line 384 | void new_data_mc_agreement_2d() {
384		float RABS_ExtraPolatedResults[nptcuts];
385
386		TCanvas *can = new TCanvas("can");
387	+	can->cd(1)->SetLeftMargin(0.15);
388
389		TGraphErrors *MPF_FinalGraph = new TGraphErrors();
390		TGraphErrors *RABS_FinalGraph = new TGraphErrors();
391	<	TGraphErrors *FaceValueAtPoint3 = new TGraphErrors();
391	>	TGraphErrors *MPF_FaceValueAtPoint3 = new TGraphErrors();
392	>	TGraphErrors *RABS_FaceValueAtPoint3 = new TGraphErrors();
393
394		for(int ipt=0;ipt<nptcuts-1;ipt++) {
395
#	Line 341 | Line 398 | void new_data_mc_agreement_2d() {
398		cout << "Going to create histo for " << filename.str() << endl;
399		TGraphErrors *mpf_gr =new TGraphErrors();
400		TGraphErrors *rabs_gr =new TGraphErrors();
401	+
402	+	int rabs_counter=0;
403	+	int mpf_counter=0;
404	+
405		for(int ia=0;ia<nalphacuts;ia++) {
406		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;
407	<	mpf_gr->SetPoint(ia,alphacuts[ia],MPF_Results[ia][ipt]);     // avoiding empty bins (nan's)
408	<	mpf_gr->SetPointError(ia,0,MPF_Errors[ia][ipt]);             // avoiding empty bins (nan's)
409	<	rabs_gr->SetPoint(ia,alphacuts[ia],RABS_Results[ia][ipt]); // avoiding empty bins (nan's)
410	<	rabs_gr->SetPointError(ia,0,RABS_Errors[ia][ipt]);         // avoiding empty bins (nan's)
411	<	}
412	<	cout << "About to add a point in pt (" << ptcuts[ipt]+ptcuts[ipt+1] << " ) to face value : " << MPF_Results[PointThree][ipt] << endl;
413	<	FaceValueAtPoint3->SetPoint(ipt,(ptcuts[ipt]+ptcuts[ipt+1]),MPF_Results[PointThree][ipt]);
414	<	FaceValueAtPoint3->SetPointError(ipt,0,MPF_Errors[PointThree][ipt]);
407	>	if(MPF_Results[ia][ipt]==MPF_Results[ia][ipt] && MPF_Errors[ia][ipt]==MPF_Errors[ia][ipt]) { // remove nan's
408	>	mpf_gr->SetPoint(mpf_counter,alphacuts[ia],MPF_Results[ia][ipt]);
409	>	mpf_gr->SetPointError(mpf_counter,0,MPF_Errors[ia][ipt]);
410	>	mpf_counter++;
411	>	} else {
412	>	dout << "Detected BS for MPF method (bin ignored) : " << endl;
413	>	dout << "       alpha: " << alphacuts[ia] << " , pt range: " << ptcuts[ipt] << " < pt < " << ptcuts[ipt+1] << endl;
414	>	cout << "       value (MPF) : " << MPF_Results[ia][ipt] << " +/- " << MPF_Errors[ia][ipt] << endl;
415	>	}
416	>	if(RABS_Results[ia][ipt]==RABS_Results[ia][ipt] && RABS_Errors[ia][ipt]==RABS_Errors[ia][ipt]) {
417	>	rabs_gr->SetPoint(rabs_counter,alphacuts[ia],RABS_Results[ia][ipt]);
418	>	rabs_gr->SetPointError(rabs_counter,0,RABS_Errors[ia][ipt]);
419	>	rabs_counter++;
420	>	} else {
421	>	dout << "Detected BS for RABS method (bin ignored) : " << endl;
422	>	dout << "       alpha: " << alphacuts[ia] << " , pt range: " << ptcuts[ipt] << " < pt < " << ptcuts[ipt+1] << endl;
423	>	cout << "       value (RABS) : " << RABS_Results[ia][ipt] << " +/- " << RABS_Errors[ia][ipt] << endl;
424	>	}
425	>	}
426
427	+	stringstream specialcut;
428	+	specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (Zb3010_alpha<0.3))";
429	+	Value MPF_data_over_mc = get_Zb_data_over_mc("mpf","",TCut(basecut && specialcut.str().c_str()),"MPF___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__INCLUSIVEalpha_0.3");
430	+	Value RABS_data_over_mc = get_Zb_data_over_mc("Zb3010_pfJetGoodPt[0]/pt","",TCut(basecut && specialcut.str().c_str()),"Rabs___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__INCLUSIVEalpha_0.3");
431	+	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;
432	+	MPF_FaceValueAtPoint3->SetPoint(ipt,(ptcuts[ipt]+ptcuts[ipt+1])/2,MPF_data_over_mc.getValue());
433	+	MPF_FaceValueAtPoint3->SetPointError(ipt,0,MPF_data_over_mc.getError());
434	+	RABS_FaceValueAtPoint3->SetPoint(ipt,(ptcuts[ipt]+ptcuts[ipt+1])/2,RABS_data_over_mc.getValue());
435	+	RABS_FaceValueAtPoint3->SetPointError(ipt,0,RABS_data_over_mc.getError());
436	+
437	+	can->cd();
438		mpf_gr->SetMarkerStyle(21);
439		mpf_gr->SetMarkerColor(kBlue);
440		mpf_gr->Draw("AP");
441	<	mpf_gr->Fit("pol1","QE"); // quiet, use minos
441	>	mpf_gr->Fit("pol1","");
442		mpf_gr->GetXaxis()->SetTitle("#alpha");
443		mpf_gr->GetXaxis()->CenterTitle();
444		mpf_gr->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
445		mpf_gr->GetYaxis()->CenterTitle();
446		TF1 *mpf_pol=(TF1*)mpf_gr->GetFunction("pol1");
447	+	mpf_pol->SetLineWidth(0);
448	+	TF1 *mpf_pol_complete = new TF1("mpf_pol_complete","[0]+[1]*x");
449	+	mpf_pol_complete->SetParameters(mpf_pol->GetParameters());
450	+	mpf_pol_complete->SetLineWidth(1);
451	+	mpf_pol_complete->SetLineColor(kBlue);
452	+	float min,max;
453	+	FindMinMax(mpf_gr,min,max);
454	+	TH1F *histo = new TH1F("histo","histo",1,0,0.4);
455	+	histo->GetXaxis()->SetTitle("#alpha");
456	+	histo->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
457	+	histo->GetXaxis()->CenterTitle();
458	+	histo->GetYaxis()->CenterTitle();
459	+	histo->SetStats(0);
460	+	histo->GetXaxis()->SetRangeUser(0,0.4);
461	+	histo->GetYaxis()->SetRangeUser(min,max);
462	+	mpf_gr->GetYaxis()->SetRangeUser(min,max);
463	+
464	+	TPolyLine *mpf_fit_uncert = GetFitUncertaintyShape(mpf_gr, "pol1", min, max,0.0,0.4);
465	+	mpf_pol->SetLineColor(TColor::GetColor("#0101DF"));
466	+	mpf_fit_uncert->SetFillColor(TColor::GetColor("#5353FC"));
467	+	histo->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
468	+	histo->Draw();
469	+	mpf_fit_uncert->Draw("F");
470	+	mpf_fit_uncert->Draw("");
471	+	mpf_gr->Draw("P");
472	+	mpf_pol_complete->Draw("same");
473	+
474		float mpf_a=mpf_pol->GetParameter(0);
475		float mpf_b=mpf_pol->GetParameter(1);
476		MPF_FinalGraph->SetPoint(ipt,0.5*(ptcuts[ipt]+ptcuts[ipt+1]),mpf_a);
#	Line 369 | Line 479 | void new_data_mc_agreement_2d() {
479		MPF_ExtraPolatedResults[ipt]=mpf_a;
480
481		stringstream mpf_info;
482	<	mpf_info <<  "#splitline{#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (MPF) }}{Extrapolated value: " << std::setprecision(3) << mpf_a << "}}{#chi^{2}/NDF : " << mpf_pol->GetChisquare() << " / " << mpf_pol->GetNDF() << "}";
482	>	mpf_info <<  "#splitline{#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (MPF) }}{Extrapolated value: " << std::setprecision(4) << mpf_a << "}}{#chi^{2}/NDF : " << mpf_pol->GetChisquare() << " / " << mpf_pol->GetNDF() << "}";
483
484		TText *mpf_mark = write_title(mpf_info.str());
485		mpf_mark->SetX(0.75);
#	Line 386 | Line 496 | void new_data_mc_agreement_2d() {
496		filename.str("");
497
498		rabs_gr->SetMarkerStyle(21);
499	<	rabs_gr->SetMarkerColor(kBlue);
390	<	rabs_gr->Fit("pol1","QE"); // quiet, use minos
499	>	rabs_gr->SetMarkerColor(kRed);
500		rabs_gr->Draw("AP");
501	+	rabs_gr->Fit("pol1","");
502		rabs_gr->GetXaxis()->SetTitle("#alpha");
503		rabs_gr->GetXaxis()->CenterTitle();
504		rabs_gr->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc}");
505		rabs_gr->GetYaxis()->CenterTitle();
396	–
506		TF1 *rabs_pol=(TF1*)rabs_gr->GetFunction("pol1");
507	+	rabs_pol->SetLineWidth(0);
508	+	TF1 *rabs_pol_complete = new TF1("rabs_pol_complete","[0]+[1]*x");
509	+	rabs_pol_complete->SetParameters(rabs_pol->GetParameters());
510	+	rabs_pol_complete->SetLineColor(kRed);
511	+	rabs_pol_complete->SetLineWidth(1);
512	+	FindMinMax(rabs_gr,min,max);
513	+	rabs_gr->GetYaxis()->SetRangeUser(min,max);
514	+	histo->GetYaxis()->SetRangeUser(min,max);
515	+	TPolyLine *rabs_fit_uncert = GetFitUncertaintyShape(rabs_gr, "pol1", min, max,0.0,0.4);
516	+	rabs_pol->SetLineColor(TColor::GetColor("#FA5858"));
517	+	rabs_pol->SetFillColor(TColor::GetColor("#FA5858"));
518	+	rabs_fit_uncert->SetLineColor(TColor::GetColor("#FA5858"));
519	+	rabs_fit_uncert->SetFillColor(TColor::GetColor("#FA5858"));
520	+	histo->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc}");
521	+	histo->Draw();
522	+	rabs_fit_uncert->Draw("F");
523	+	rabs_fit_uncert->Draw("");
524	+	rabs_gr->Draw("P");
525	+	rabs_pol_complete->Draw("same");
526	+
527	+
528		float rabs_a=rabs_pol->GetParameter(0);
529		float rabs_b=rabs_pol->GetParameter(1);
530		RABS_FinalGraph->SetPoint(ipt,0.5*(ptcuts[ipt]+ptcuts[ipt+1]),rabs_a);
#	Line 402 | Line 532 | void new_data_mc_agreement_2d() {
532
533		cout << "!+!+!+!+!+!+!+!+!+ Just added a point to the final plot : " << rabs_a << " +/- " << rabs_pol->GetParError(0) << endl;
534
405	–
535		stringstream rabs_info;
536	<	rabs_info << "#splitline{#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (R_{abs}) }}{Extrapolated value: " << std::setprecision(3) << rabs_a << "}}{#chi^{2}/NDF : " << rabs_pol->GetChisquare() << " / " << rabs_pol->GetNDF() << "}";
536	>	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() << "}";
537		TText *rabs_mark = write_title(rabs_info.str());
538		rabs_mark->SetX(0.75);
539		rabs_mark->SetTextSize(0.03);
#	Line 412 | Line 541 | void new_data_mc_agreement_2d() {
541		rabs_mark->Draw();
542
543
415	–
544		RABS_ExtraPolatedResults[ipt]=rabs_a;
545		filename << filenamebkp << "__RABS";
546		DrawPrelim();
#	Line 421 | Line 549 | void new_data_mc_agreement_2d() {
549		}
550
551		MPF_FinalGraph->SetMarkerStyle(21);
552	<	MPF_FinalGraph->SetMarkerColor(kRed);
553	<	MPF_FinalGraph->Fit("pol0","QE"); // quiet, use minos
552	>	MPF_FinalGraph->SetMarkerColor(kBlue);
553	>	MPF_FinalGraph->Fit("pol0",""); // quiet, use minos
554
555		TF1 *mpf_pol0=(TF1*)MPF_FinalGraph->GetFunction("pol0");
556	+	TF1 *mpf_pol0_complete = new TF1("mpf_pol0_complete","[0]+[1]*x");
557	+	mpf_pol0_complete->SetLineWidth(1);
558	+	mpf_pol0_complete->SetLineColor(kBlue);
559	+	mpf_pol0->SetLineColor(kBlue);
560		MPF_FinalGraph->Draw("AP");
561		MPF_FinalGraph->GetXaxis()->SetTitle("p_{T}^{Z}");
562		MPF_FinalGraph->GetXaxis()->CenterTitle();
563		MPF_FinalGraph->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
564		MPF_FinalGraph->GetYaxis()->CenterTitle();
565		MPF_FinalGraph->Draw("AP");
566	+	mpf_pol0_complete->Draw("same");
567
568		stringstream mpf_result;
569		mpf_result << "#splitline{C_{abs}= " << std::setprecision(5) << 1/mpf_pol0->GetParameter(0) << " #pm " << std::setprecision(2) << mpf_pol0->GetParError(0)/(mpf_pol0->GetParameter(0)*mpf_pol0->GetParameter(0)) << "}{ (MPF)}";
#	Line 476 | Line 609 | void new_data_mc_agreement_2d() {
609		CompleteSave(can,filename2.str());
610
611		can->cd();
612	<	TCanvas *FaceValueCanvas = new TCanvas("FaceValueCanvas","FaceValueCanvas");
613	<	FaceValueAtPoint3->SetMarkerStyle(21);
614	<	FaceValueAtPoint3->SetMarkerStyle(21);
615	<	FaceValueAtPoint3->Fit("pol0","QE"); // quiet, use minos
616	<	FaceValueAtPoint3->Draw("AP");
617	<	FaceValueAtPoint3->GetXaxis()->SetTitle("p_{T}^{Z}");
618	<	FaceValueAtPoint3->GetXaxis()->CenterTitle();
619	<	FaceValueAtPoint3->GetYaxis()->SetTitle("data/mc for #alpha<0.3");
620	<	FaceValueAtPoint3->GetYaxis()->CenterTitle();
621	<	FaceValueAtPoint3->Draw("AP");
622	<	TF1 *faceval_pol0=(TF1*)FaceValueAtPoint3->GetFunction("pol0");
612	>	MPF_FaceValueAtPoint3->SetMarkerStyle(21);
613	>	MPF_FaceValueAtPoint3->SetMarkerColor(kBlue);
614	>	MPF_FaceValueAtPoint3->Fit("pol0","QE"); // quiet, use minos
615	>	MPF_FaceValueAtPoint3->Draw("AP");
616	>	MPF_FaceValueAtPoint3->GetXaxis()->SetTitle("p_{T}^{Z}");
617	>	MPF_FaceValueAtPoint3->GetXaxis()->CenterTitle();
618	>	MPF_FaceValueAtPoint3->GetYaxis()->SetTitle("<MPF>__{data}/<MPF>_{mc} for #alpha<0.3");
619	>	MPF_FaceValueAtPoint3->GetYaxis()->CenterTitle();
620	>	MPF_FaceValueAtPoint3->Draw("AP");
621	>	TF1 *faceval_pol0=(TF1*)MPF_FaceValueAtPoint3->GetFunction("pol0");
622	>	faceval_pol0->SetLineColor(kBlue);
623	>	faceval_pol0->SetLineWidth(1);
624	>	faceval_pol0->Draw("same");
625	>	TF1 *faceval_pol0_complete = new TF1("faceval_pol0_complete","[0]+[1]*x");
626	>	faceval_pol0_complete->SetParameters(faceval_pol0->GetParameters());
627		stringstream faceval_result;
628		faceval_result << "#splitline{C_{abs}= " << std::setprecision(5) << 1/faceval_pol0->GetParameter(0) << " #pm " << std::setprecision(2) << faceval_pol0->GetParError(0)/(faceval_pol0->GetParameter(0)*faceval_pol0->GetParameter(0)) << "}{ (MPF at #alpha<0.3)}";;
629		TText *faceval_final_mark = write_title(faceval_result.str());
#	Line 503 | Line 640 | void new_data_mc_agreement_2d() {
640		filename2 << "Extrapolation/FINAL_RESULT_MPF_FaceValp3";
641		CompleteSave(can,filename2.str());
642
643	+	can->cd();
644	+	RABS_FaceValueAtPoint3->SetMarkerStyle(21);
645	+	RABS_FaceValueAtPoint3->SetMarkerColor(kRed);
646	+	RABS_FaceValueAtPoint3->Fit("pol0","QE"); // quiet, use minos
647	+	RABS_FaceValueAtPoint3->Draw("AP");
648	+	RABS_FaceValueAtPoint3->GetXaxis()->SetTitle("p_{T}^{Z}");
649	+	RABS_FaceValueAtPoint3->GetXaxis()->CenterTitle();
650	+	RABS_FaceValueAtPoint3->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc} for #alpha<0.3");
651	+	RABS_FaceValueAtPoint3->GetYaxis()->CenterTitle();
652	+	RABS_FaceValueAtPoint3->Draw("AP");
653	+	TF1 *faceval_pol0RABS=(TF1*)RABS_FaceValueAtPoint3->GetFunction("pol0");
654	+	TF1 *faceval_pol0RABS_complete = new TF1("faceval_pol0_complete","[0]+[1]*x");
655	+	faceval_pol0RABS_complete->SetParameters(faceval_pol0RABS->GetParameters());
656	+	faceval_pol0RABS->SetLineWidth(1);
657	+	faceval_pol0RABS->SetLineColor(kRed);
658	+	faceval_pol0RABS_complete->Draw("same");
659	+	stringstream RABSfaceval_result;
660	+	RABSfaceval_result << "#splitline{C_{abs}= " << std::setprecision(5) << 1/faceval_pol0RABS->GetParameter(0) << " #pm " << std::setprecision(2) << faceval_pol0RABS->GetParError(0)/(faceval_pol0RABS->GetParameter(0)*faceval_pol0RABS->GetParameter(0)) << "}{ (R_{abs} at #alpha<0.3)}";;
661	+	TText *RABSfaceval_final_mark = write_title(RABSfaceval_result.str());
662	+	RABSfaceval_final_mark->SetX(0.75);
663	+	RABSfaceval_final_mark->SetY(0.75);
664	+	RABSfaceval_final_mark->SetTextSize(0.03);
665	+	RABSfaceval_final_mark->Draw();
666	+	DrawPrelim();
667	+
668	+	float RABSFaceValResult=1/faceval_pol0RABS->GetParameter(0);
669	+	float RABSFaceValResultError=faceval_pol0RABS->GetParError(0)/(faceval_pol0RABS->GetParameter(0)*faceval_pol0RABS->GetParameter(0));
670	+
671	+	filename2.str("");
672	+	filename2 << "Extrapolation/FINAL_RESULT_RABS_FaceValp3";
673	+	CompleteSave(can,filename2.str());
674	+
675		cout << "FINAL RESULTS: " << endl;
676		cout << "MPF : " << MPFResult  << " +/- " << MPFResultError  << endl;
677		cout << "Rabs: " << RabsResult << " +/- " << RabsResultError << endl;
678	<	cout << "Face value at 0.3 :" << FaceValResult << " +/- " << FaceValResultError << endl;
678	>	cout << "Face value at 0.3 :" << FaceValResult << " +/- " << FaceValResultError << "  (MPF) " << endl;
679	>	cout << "Face value at 0.3 :" << RABSFaceValResult << " +/- " << RABSFaceValResultError << "  (RABS) " << endl;
680
681		delete can;
682
#	Line 565 | Line 735 | void DoPUStudy(string identifier) {
735		}
736
737		void ScenarioComparison() {
738	<
738	>	//very dumb way to do this.
739		TGraphAsymmErrors *gr[5];
740		TF1 *fit[5];
741		bool dofit=false;
#	Line 598 | Line 768 | void ScenarioComparison() {
768		double dyl[4] = {0.011,0.009,0.013,0.020};
769		double dyh[4] = {0.011,0.009,0.013,0.020};
770		gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
771	<	name="30/10";
771	>	name="30/15";
772		color="#2464bc";
773		}
774		if(i==2) {
#	Line 655 | Line 825 | void ScenarioComparison() {
825		}
826
827		void ScenarioComparisonInclusive() {
828	<
828	>	//dumbest way ever to do this. but ok we only need to do it once.
829		TGraphAsymmErrors *gr[5];
830		TF1 *fit[5];
831		bool dofit=true;
#	Line 688 | Line 858 | void ScenarioComparisonInclusive() {
858		double dyl[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
859		double dyh[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
860		gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
861	<	name="30/10";
861	>	name="30/15";
862		color="#2464bc";
863		}
864		if(i==2) {
#	Line 874 | Line 1044 | void compare_selection(string identifier
1044
1045		}
1046
1047	<	void do_basic_ZB_analysis() {
1047	>	void GetNumberEventsInsideOutsideAlphaWindow(TCut specialcut, string title) {
1048	>
1049	>	//  TCut cut =  ZplusBsel&&LeadingB&&PhiZcut&&specialcut&&TCut("Zb3010_alpha<0.3");
1050	>	TCut cut = specialcut;
1051	>	TCut in  =  EtaB;
1052	>	TCut out =  TCut("abs(Zb3010_pfJetGoodEta[0])>1.3");
1053	>
1054	>	TH1F *data_IN  = allsamples.Draw("data_IN",    "mll",1,0,150, "nothing [GeV]", "events", cut&&in, data,luminosity);
1055	>	TH1F *data_OUT = allsamples.Draw("data_OUT",   "mll",1,0,150, "nothing [GeV]", "events", cut&&out,data,luminosity);
1056	>	TH1F *mc_IN    = allsamples.Draw("mc_IN",      "mll",1,0,150, "nothing [GeV]", "events", cut&&in, mc  ,luminosity);
1057	>	TH1F *mc_OUT   = allsamples.Draw("mc_OUT",     "mll",1,0,150, "nothing [GeV]", "events", cut&&out,mc  ,luminosity);
1058	>
1059	>	dout << title << " : " << endl;
1060	>	dout << "  Data: In " << data_IN->Integral() << " vs out " << data_OUT->Integral() << endl;
1061	>	dout << "  MC  : In " << mc_IN->Integral()   << " vs out " << mc_OUT->Integral()   << endl;
1062	>	dout << "       ratio in  : " << data_IN->Integral()/mc_IN->Integral() << " +/- " << (data_IN->Integral()/mc_IN->Integral()) * sqrt(1.0/mc_IN->Integral()+ 1.0/data_IN->Integral()) << endl;
1063	>	dout << "       ratio out : " << data_OUT->Integral()/mc_OUT->Integral() << " +/- " << (data_OUT->Integral()/mc_OUT->Integral()) * sqrt(1.0/mc_OUT->Integral()+ 1.0/data_OUT->Integral()) << endl;
1064	>
1065	>	delete data_IN;
1066	>	delete data_OUT;
1067	>	delete mc_IN;
1068	>	delete mc_OUT;
1069	>	}
1070	>
1071	>	void GetNumberEventsInsideOutsideAlphaWindow() {
1072	>	TCanvas *ca = new TCanvas("ca","ca");
1073	>	//   GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2&&id1==0"), "ee");
1074	>	//   GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2&&id1==1"), "mm");
1075	>	//   GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2"), "ee/mm");
1076	>
1077	>	cout << "BASE SELECTION" << endl;
1078	>	GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2&&id1==0"), "ee");
1079	>	GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2&&id1==1"), "mm");
1080	>	GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2"), "ee/mm");
1081	>
1082	>	cout << "BASE SELECTION: Z+b" << endl;
1083	>	GetNumberEventsInsideOutsideAlphaWindow(ZplusBsel&&TCut("id1==id2&&id1==0"), "ee");
1084	>	GetNumberEventsInsideOutsideAlphaWindow(ZplusBsel&&TCut("id1==id2&&id1==1"), "mm");
1085	>	GetNumberEventsInsideOutsideAlphaWindow(ZplusBsel&&TCut("id1==id2"), "ee/mm");
1086	>
1087	>	cout << "Leading B" << endl;
1088	>	GetNumberEventsInsideOutsideAlphaWindow(LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==0"), "ee");
1089	>	GetNumberEventsInsideOutsideAlphaWindow(LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==1"), "mm");
1090	>	GetNumberEventsInsideOutsideAlphaWindow(LeadingB&&ZplusBsel&&TCut("id1==id2"), "ee/mm");
1091	>
1092	>	cout << "PhiZcut" << endl;
1093	>	GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==0"), "ee");
1094	>	GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==1"), "mm");
1095	>	GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2"), "ee/mm");
1096	>
1097	>	cout << "alpha cut" << endl;
1098	>	GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==0")&&TCut("Zb3010_alpha<0.3"), "ee");
1099	>	GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==1")&&TCut("Zb3010_alpha<0.3"), "mm");
1100	>	GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2")&&TCut("Zb3010_alpha<0.3"), "ee/mm");
1101	>
1102	>	delete ca;
1103	>
1104	>
1105	>	}
1106	>
1107	>
1108	>
1109	>
1110	>	void do_basic_ZB_analysis(int do_inclusive) {
1111
1112	+	//https://twiki.cern.ch/twiki/bin/viewauth/CMS/BtagPOG#2011_Data_and_MC
1113	+	//https://twiki.cern.ch/twiki/pub/CMS/BtagPOG/SFb-ttbar_payload.txt
1114	+
1115	+	cout << "The lepton requirement is " << (const char*) leptoncut << endl;
1116	+	if(do_inclusive) {
1117	+	LeadingB=TCut("Zb3010_bTagProbCSVBP[0]>=-10000");
1118	+	cutWeight=TCut("(weight*(weight<1000)*(is_data+(!is_data)))");
1119	+	write_info(__FUNCTION__,"Inclusive selection - no btag requirement, no btag efficiency/mistag correction (weight)");
1120	+	} else {
1121	+	LeadingB=TCut("Zb3010_bTagProbCSVBP[0]>0.898");
1122	+	cutWeight=TCut("(weight*(weight<1000)*(is_data+(!is_data)*((id1==id2&&id1==0)*0.95+(id1==id2&&id1==1)*0.88+(id1!=id2)*0.92)*Zb3010_BTagWgt))");
1123	+	write_info(__FUNCTION__,"Exclusive selection - switched on btag requirement and btag efficiency/mistag correction (weight)");
1124	+	}
1125	+
1126	+	bool doquick=true;
1127
1128		TCanvas *zbcanvas = new TCanvas("zbcanvas","zbcanvas");
1129
#	Line 893 | Line 1141 | void do_basic_ZB_analysis() {
1141		//     compare_selection("Zb3010");
1142		//     compare_selection("Zb40");
1143
1144	+	if(!doquick) {
1145	+	print_all_b_yields();
1146	+	draw_mpf_vars();
1147	+	//  DrawEvilCutFlow();
1148	+
1149	+	draw_Zb_kin_vars();
1150	+
1151	+	}
1152
1153	<	//  print_all_b_yields();
898	<	//  draw_mpf_vars();
899	<	//  DrawEvilCutFlow();
900	<
901	<	//  draw_Zb_kin_vars();
902	<
1153	>	//  GetNumberEventsInsideOutsideAlphaWindow();
1154		new_data_mc_agreement_2d();
1155
905	–
1156		delete zbcanvas;
1157		}
1158

Diff Legend

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