ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/cbrown/AnalysisFramework/Plotting/Modules/Systematics.C

Comparing UserCode/cbrown/AnalysisFramework/Plotting/Modules/Systematics.C (file contents):
Revision 1.29 by buchmann, Wed Sep 7 06:42:21 2011 UTC vs.
Revision 1.36 by buchmann, Wed Sep 28 15:58:09 2011 UTC

#	Line 1 | Line 1
1		#include <iostream>
2		#include <vector>
3		#include <sys/stat.h>
4	+	#include <algorithm>
5	+	#include <cmath>
6
7		#include <TMath.h>
8		#include <TColor.h>
#	Line 276 | Line 278 | void PeakError(TTree *events,float &resu
278
279		//____________________________________________________________________________________
280		// Total selection efficiency (MC)
281	<	void MCefficiency(TTree *events,float &result, float &resulterr,string mcjzb,bool requireZ,int Neventsinfile, string addcut="", int k = 0) {
281	>	//returns the efficiency WITHOUT signal contamination, and the result and resulterr contain the result and the corresponding error
282	>	Value MCefficiency(TTree *events,float &result, float &resulterr,string mcjzb,bool requireZ,int Neventsinfile, string addcut="", int k = 0) {
283	>
284	>	if(!events) {
285	>	write_error(__FUNCTION__,"Tree passed for efficiency calculation is invalid!");
286	>	return Value(0,0);
287	>	}
288
289		char jzbSelStr[256]; sprintf(jzbSelStr,"%f",jzbSel);
290		// All acceptance cuts at gen. level
#	Line 298 | Line 306 | void MCefficiency(TTree *events,float &r
306		string newNegSide = "(" + snegSide + ")*" + svar;
307
308		TH1F *effh= new TH1F("effh","effh",1,-14000,14000);
309	<	if(k>=0)events->Draw((mcjzbexpression+">>effh").c_str(), newPosSide.c_str(),"goff");
310	<	else events->Draw((mcjzbexpression+">>effh").c_str(), sposSide.c_str(),"goff");
309	>	cout << "right before drawing" << endl;
310	>	if(k>=0)events->Draw((mcjzbexpression+">>effh").c_str(), newPosSide.c_str(),"");
311	>	else events->Draw((mcjzbexpression+">>effh").c_str(), sposSide.c_str(),"");
312	>	cout << "right after drawing" << endl;
313	>	cout << "consider sigcont?" << ConsiderSignalContaminationForLimits << endl;
314		Float_t sel = effh->Integral();
315		Float_t nsel=0;
316		if(ConsiderSignalContaminationForLimits) {
317	<	if(k>=0)events->Draw((mcjzbexpression+">>effh").c_str(), newNegSide.c_str(),"goff");
318	<	else events->Draw((mcjzbexpression+">>effh").c_str(), snegSide.c_str(),"goff");
317	>	if(k>=0)events->Draw((mcjzbexpression+">>effh").c_str(), newNegSide.c_str(),"");
318	>	else events->Draw((mcjzbexpression+">>effh").c_str(), snegSide.c_str(),"");
319		nsel = effh->Integral();
320		}
321		//Corrections due to normalization in the PDF. This has to be applied as well to the number of events in a file if the definition changes at some point.
#	Line 317 | Line 328 | void MCefficiency(TTree *events,float &r
328		//      Float_t tot = events->GetSelectedRows();
329		Float_t tot = Neventsinfile;
330
331	+	Value result_wo_signalcont;
332	+
333		if(ConsiderSignalContaminationForLimits) {
334		result=(sel-nsel)/tot;
335		resulterr=(1.0/tot)*TMath::Sqrt(sel+nsel+(sel-nsel)*(sel-nsel)/tot);
336	+	result_wo_signalcont=Value(sel/tot,TMath::Sqrt(sel/tot*(1+sel/tot)/tot));
337		} else {//no signal contamination considered:
338		result=(sel)/tot;
339		resulterr=TMath::Sqrt(sel/tot*(1+sel/tot)/tot);
340	+	result_wo_signalcont=Value(result,resulterr);
341		}
342	+	if(!automatized && k>0 ) dout << "PDF assessment: ";
343		if(!automatized) dout << "  MC efficiency: " << result << "+-" << resulterr << "  ( JZB>" << jzbSel << " : " << sel << " , JZB<-" << jzbSel << " : " << nsel << " and nevents=" << tot << ") with normFact=" << normFactor << std::endl;
344		delete effh;
345	+	return result_wo_signalcont;
346		}
347
348
349		//____________________________________________________________________________________
350	<	// Total selection efficiency (MC)
350	>	// Selection efficiency for one process (MC)
351		vector<float> processMCefficiency(TTree *events,string mcjzb,bool requireZ,int Neventsinfile, string addcut) {
352		vector<float> process_efficiencies;
353		for(int iprocess=0;iprocess<=10;iprocess++) {
#	Line 412 | Line 429 | void doJZBscale(TTree *events, float &do
429		Float_t eff  = Interpolate(jzbSel,hist);
430		Float_t effp = Interpolate(jzbSel*(1.+systematic),hist);
431		Float_t effm = Interpolate(jzbSel*(1.-systematic),hist);
432	<	if(!automatized) dout << "  efficiency at JZB==" << jzbSel*(1.+systematic)  << "(-"<<syst*100<<"%) : " << effp << " (" << ((effp-eff)/eff)*100. << "%)"  << std::endl;
432	>	if(!automatized) dout << "  efficiency at JZB==" << jzbSel*(1.+systematic)  << "(-"<<systematic*100<<"%) : " << effp << " (" << ((effp-eff)/eff)*100. << "%)"  << std::endl;
433		if(!automatized) dout << "  efficiency at JZB==" << jzbSel  << ": " << eff << std::endl;
434	<	if(!automatized) dout << "  efficiency at JZB==" << jzbSel*(1.-systematic)  << "(-"<<syst*100<<"%) : " << effm << " (" << ((effm-eff)/eff)*100. << "%)"  << std::endl;
434	>	if(!automatized) dout << "  efficiency at JZB==" << jzbSel*(1.-systematic)  << "(-"<<systematic*100<<"%) : " << effm << " (" << ((effm-eff)/eff)*100. << "%)"  << std::endl;
435		up=((effp-eff)/eff);
436		down=((effm-eff)/eff);
437		}
#	Line 438 | Line 455 | void JZBresponse(TTree *events, bool req
455		hJzbResp->SetMinimum(0.2);
456		hJzbResp->Fit("pol0","Q");
457		TF1 *fittedfunction = hJzbResp->GetFunction("pol0");
458	<	resp=fittedfunction->GetParameter(0);
459	<	resperr=fittedfunction->GetParError(0);
460	<	if(!automatized) dout << "  Response: " << resp << " +/- " << resperr << endl;
458	>	if(!fittedfunction) {
459	>	// in case there are not enough points passing our selection
460	>	cout << "OOPS response function invalid, assuming 100% error !!!!" << endl;
461	>	resp=1;
462	>	resperr=1;
463	>	} else {
464	>	resp=fittedfunction->GetParameter(0);
465	>	resperr=fittedfunction->GetParError(0);
466	>	if(!automatized) dout << "  Response: " << resp << " +/- " << resperr << endl;
467	>	}
468		delete hJzbResp;
469		}
470
#	Line 478 | Line 502 | void do_systematics_for_one_file(TTree *
502		float triggereff=5.0/100;// in range [0,1]
503		dout << "Trigger efficiency not implemented in this script  yet, still using external one" << endl;
504		float leptonseleff=2.0/100;// in range [0,1]
505	+	leptonseleff=TMath::Sqrt(leptonseleff*leptonseleff+leptonseleff*leptonseleff); // because the 2% is per lepton
506		dout << "Lepton selection efficiency not implemented in this script  yet, still using external one" << endl;
507
508		int NPdfs=0;
#	Line 485 | Line 510 | void do_systematics_for_one_file(TTree *
510
511		float mceff,mcefferr,jzbeff,jzbefferr;
512		if(!automatized) dout << "MC efficiencies:" << endl;
513	<	MCefficiency(events,mceff,mcefferr,mcjzb,requireZ,Neventsinfile,addcut,-1);
514	<	JZBefficiency(events,informalname,jzbeff,jzbefferr,requireZ,addcut);
513	>	Value mceff_nosigcont = MCefficiency(events,mceff,mcefferr,mcjzb,requireZ,Neventsinfile,addcut,-1);
514	>	if(!automatized) cout << "   Without signal contamination, we find an efficiency of " << mceff_nosigcont << endl;
515	>
516	>	if(PlottingSetup::computeJZBefficiency) JZBefficiency(events,informalname,jzbeff,jzbefferr,requireZ,addcut);
517		if(!automatized) dout << "JZB efficiency: " << jzbeff << "+/-" << jzbefferr << endl;
518
519		if(!automatized) dout << "Error from Peak position:" << endl;
#	Line 503 | Line 530 | void do_systematics_for_one_file(TTree *
530
531		if(!automatized) dout << "JZB response: " << std::endl;
532		float resp,resperr;
533	<	JZBresponse(events,requireZ,resp,resperr,addcut);
533	>	if(PlottingSetup::computeJZBresponse) {
534	>	if(!automatized) dout << "JZB response: " << std::endl;
535	>	JZBresponse(events,requireZ,resp,resperr,addcut);
536	>	}
537
538		if(!automatized) dout << "Pileup: " << std::endl;
539	<	float resolution=pileup(events,requireZ,informalname,addcut);
539	>	float resolution;
540	>	resolution=pileup(events,requireZ,informalname,addcut);
541
542		float PDFuncert=0;
543	+	if(!automatized) dout << "Assessing PDF uncertainty: " << std::endl;
544		if(ismSUGRA) PDFuncert = get_pdf_uncertainty(events, mcjzb, requireZ, Neventsinfile, NPdfs, addcut);
545
546		dout << "_______________________________________________" << endl;
#	Line 523 | Line 555 | void do_systematics_for_one_file(TTree *
555		dout << "Resolution : " << resolution << endl; // in range [0,1]
556		dout << "From peak : " << sysfrompeak << endl; // in range [0,1]
557		if(ismSUGRA) dout << "PDF uncertainty  : " << PDFuncert << endl; // in range [0,1]
558	<	dout << "JZB efficiency: " << jzbeff << "+/-" << jzbefferr << " (not yet included below) " << endl; // in range [0,1]
559	<	dout << "JZB response  : " << resp << " +/-" << resperr << " (not yet included below) " << endl; // in range [0,1]
558	>	if(PlottingSetup::computeJZBefficiency) dout << "JZB efficiency: " << jzbeff << "+/-" << jzbefferr << " (not yet included below) " << endl; // in range [0,1]
559	>	if(PlottingSetup::computeJZBresponse)dout << "JZB response  : " << resp << " +/-" << resperr << " (not yet included below) " << endl; // in range [0,1]
560
561		float toterr=0;
562		toterr+=(triggereff)*(triggereff);
#	Line 540 | Line 572 | void do_systematics_for_one_file(TTree *
572
573		dout << "FINAL RESULT : " << 100*mceff << " +/- "<< 100*mcefferr << " (stat) +/- " << 100*systerr << " (syst)   %" << endl;
574		dout << "     we thus use the sqrt of the sum of the squares of the stat & syst err, which is : " << 100*toterr << endl;
575	+	dout << "_______________________________________________" << endl;
576
577		//Do not modify the lines below or mess with the order; this order is expected by all limit calculating functions!
578		vector<float> res;
#	Line 551 | Line 584 | void do_systematics_for_one_file(TTree *
584		if(fabs(jesup)>fabs(jesdown)) res.push_back(fabs(jesup)); else res.push_back(fabs(jesdown));
585		if(fabs(scaleup)>fabs(scaledown)) res.push_back(fabs(scaleup)); else res.push_back(fabs(scaledown));
586		res.push_back(fabs(resolution));
587	+	res.push_back(mceff_nosigcont.getValue());
588	+	res.push_back(mceff_nosigcont.getError());
589		if(ismSUGRA) res.push_back(PDFuncert);
590		results.push_back(res);
591		}

Diff Legend

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