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.48 by buchmann, Wed Nov 9 14:18:17 2011 UTC vs.
Revision 1.54 by pablom, Wed Nov 23 16:39:56 2011 UTC

#	Line 149 | Line 149 | float allcontributionsplot(TTree* events
149		TH1F* plotEff(TTree* events, TCut kbase, TString informalname, int flipped) {
150		iplot++;
151		int count=iplot;
152	+	iplot++;
153	+	int count2=iplot;
154		// Define new histogram
155		char hname[30]; sprintf(hname,"hJzbEff%d",count);
156	<	TH1F* hJzbEff = new TH1F(hname,"JZB selection efficiency ; JZB (GeV/c); Efficiency",
157	<	nBins,jzbMin,jzbMax);
156	>	char hname2[30]; sprintf(hname2,"hJzbEff%d",count2);
157	>	TH1F* hJzbEff = new TH1F(hname,"JZB selection efficiency ; JZB [GeV]; Efficiency",nBins,jzbMin,jzbMax);
158	>	TH1F* hJzbEff2= new TH1F(hname2,"JZB selection efficiency ; JZB [GeV]; Efficiency",1,-14000,14000);
159		Float_t step = (jzbMax-jzbMin)/static_cast<Float_t>(nBins);
160	<
161	<	if(flipped==0) events->Draw(mcjzbexpression.c_str(),"genJZB>-400"&&kbase,"goff");
162	<	else events->Draw(("-"+mcjzbexpression).c_str(),"genJZB>-14000"&&kbase,"goff");
163	<	Float_t maxEff = events->GetSelectedRows();
160	>
161	>	if(flipped==0) events->Draw((mcjzbexpression+">>"+(string)hname2).c_str(),("genJZB>-400"&&kbase),"goff");
162	>	else events->Draw(("(-"+mcjzbexpression+")>>"+(string)hname2).c_str(),("genJZB>-400"&&kbase),"goff");
163	>	Float_t maxEff = hJzbEff2->Integral();
164		if(verbose>0) dout << hname << " (" << informalname <<") " << maxEff <<  std::endl;
165
166		if(verbose>0) dout <<  "JZB max = " << jzbMax << std::endl;
#	Line 165 | Line 168 | TH1F* plotEff(TTree* events, TCut kbase,
168		char cut[256];
169		for ( Int_t iBin = 0; iBin<nBins; ++iBin ) {
170		sprintf(cut,"genJZB>%3f",jzbMin+iBin*step);
171	<	events->Draw(mcjzbexpression.c_str(),TCut(cut)&&kbase,"goff");
172	<	Float_t eff = static_cast<Float_t>(events->GetSelectedRows())/maxEff;
173	<	//     dout << "COUCOU " << __LINE__ << std::endl;
171	>	if(flipped==0) events->Draw((mcjzbexpression+">>"+(string)hname2).c_str(),(TCut(cut)&&kbase),"goff");
172	>	if(flipped>0) events->Draw(("(-"+mcjzbexpression+")>>"+(string)hname2).c_str(),(TCut(cut)&&kbase),"goff");
173	>	Float_t eff = static_cast<Float_t>(hJzbEff2->Integral())/maxEff;
174		hJzbEff->SetBinContent(iBin+1,eff);
175		hJzbEff->SetBinError(iBin+1,TMath::Sqrt(eff*(1-eff)/maxEff));
176		}
177	+	delete hJzbEff2;
178		return hJzbEff;
175	–
176	–
179		}
180
181
182	+
183		//________________________________________________________________________________________
184		// Master Formula
185		void master_formula(std::vector<float> eff, float &errHi, float &errLo) {
#	Line 294 | Line 297 | void PeakError(TTree *events,float &resu
297		}
298		if(TMath::Abs(rescent-resup)>TMath::Abs(rescent-resdown)) result=(TMath::Abs(rescent-resup)/(float)rescent);
299		else result=(TMath::Abs(rescent-resdown)/(float)rescent);
300	+	cout << "    " << result << endl;
301		}
302
303
#	Line 305 | Line 309 | void  MCPartialefficiency(TTree *events,
309		}
310
311		char jzbSelStr[256]; sprintf(jzbSelStr,"%f",jzbSel);
312	+	char metSel[256]; sprintf(metSel, "met[4] > %f", jzbSel);
313	+	string metSelection(metSel);
314		// All acceptance cuts at gen. level
315		//TCut kbase("abs(genMll-91.2)<20&&genNjets>2&&genZPt>0&&genJZB"+geq_or_leq()+TString(jzbSelStr)+"&&genId1==-genId2");
316		TCut kbase("");
#	Line 317 | Line 323 | void  MCPartialefficiency(TTree *events,
323		TCut njets(cutnJets);
324		TCut jzbp;
325		TCut jzbn;
326	+	TCut met(("pfJetGoodNum > 1 && abs(mll-91.2) < 10.0 && id1 == id2 &&" + metSelection).c_str());
327		if(flipped==0) {
328		jzbp=TCut((TString(mcjzb)+geq_or_leq()+TString(jzbSelStr)));
329		jzbn=TCut((TString(mcjzb)+ngeq_or_leq()+TString("-")+TString(jzbSelStr)));
#	Line 324 | Line 331 | void  MCPartialefficiency(TTree *events,
331		jzbp=TCut(TString(mcjzb)+ngeq_or_leq()+TString("-")+TString(jzbSelStr));
332		jzbn=TCut(TString(mcjzb)+geq_or_leq()+TString(jzbSelStr));
333		}
334	<	float ntotal = events->Draw("pt1", kbase, "goff");
334	>	float ntotal = events->Draw("pt1", addcut.c_str(), "goff");
335		TCut theCut;
336		switch(type) {
337		case 1:
#	Line 339 | Line 346 | void  MCPartialefficiency(TTree *events,
346		case 4:
347		theCut = kbase+massId+njets+jzbn;
348		break;
349	+	case 5:
350	+	theCut = kbase + met;
351	+	break;
352		default:
353		theCut = kbase+massId+njets+jzbn;
354		break;
355		}
356	<
356	>
357		string stheCut(theCut);
358		char var[20];
359		sprintf(var, "pdfW[%d]", k);
#	Line 366 | Line 376 | void  MCPartialefficiency(TTree *events,
376		result=(sel)/ntotal;
377		resulterr=TMath::Sqrt(sel/ntotal*(1+sel/ntotal)/ntotal);
378
379	+	delete effh;
380		}
381
382		//____________________________________________________________________________________
#	Line 424 | Line 435 | Value MCefficiency(TTree *events,float &
435		string emnewNegSide = "((id1!=id2)&&(" + snegSide + "))*" + svar; // only used for off peak analysis
436
437		TH1F *effh= new TH1F("effh","effh",1,-14000,14000);
438	<	if(k>=0)events->Draw((mcjzbexpression+">>effh").c_str(), newPosSide.c_str(),"goff");
439	<	else events->Draw((mcjzbexpression+">>effh").c_str(), (sposSide+"&&(id1==id2)").c_str(),"goff");//the OSSF condition is added for the offpeak analysis, in onpeak case it's there already but doesn't change anything.
438	>	if(k>=0)events->Draw((mcjzbexpression+">>effh").c_str(), TCut(newPosSide.c_str())*PlottingSetup::Weight,"goff");
439	>	else events->Draw((mcjzbexpression+">>effh").c_str(), TCut((sposSide+"&&(id1==id2)").c_str())*PlottingSetup::Weight,"goff");//the OSSF condition is added for the offpeak analysis, in onpeak case it's there already but doesn't change anything.
440	>
441		Float_t sel = effh->Integral();
442		Float_t nsel=0;
443
#	Line 434 | Line 446 | Value MCefficiency(TTree *events,float &
446		if(ConsiderSignalContaminationForLimits) {
447		flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
448		if(PlottingSetup::RestrictToMassPeak) {
449	<	events->Draw((mcjzbexpression+">>effh").c_str(), newNegSide.c_str(),"goff");
449	>	events->Draw((mcjzbexpression+">>effh").c_str(), (TCut(newNegSide.c_str()))*PlottingSetup::Weight,"goff");
450		nsel += effh->Integral();
451		} else {
452	<	events->Draw((mcjzbexpression+">>effh").c_str(), newNegSide.c_str(),"goff");
452	>	events->Draw((mcjzbexpression+">>effh").c_str(), (TCut(newNegSide.c_str()))*PlottingSetup::Weight,"goff");
453		nsel += effh->Integral();
454	<	events->Draw((mcjzbexpression+">>effh").c_str(), emnewPosSide.c_str(),"goff");
454	>	events->Draw((mcjzbexpression+">>effh").c_str(), (TCut(emnewPosSide.c_str()))*PlottingSetup::Weight,"goff");
455		nsel += effh->Integral();
456	<	events->Draw((mcjzbexpression+">>effh").c_str(), emnewNegSide.c_str(),"goff");
456	>	events->Draw((mcjzbexpression+">>effh").c_str(), (TCut(emnewNegSide.c_str()))*PlottingSetup::Weight,"goff");
457		nsel -= effh->Integral();
458		}
459		}
#	Line 452 | Line 464 | Value MCefficiency(TTree *events,float &
464		sel = sel/normFactor;
465		nsel = nsel/normFactor;
466
455	–	//      events->Draw(mcjzbexpression.c_str(),kbase,"goff");
456	–	//      Float_t tot = events->GetSelectedRows();
467		Float_t tot = Neventsinfile;
468
469		Value result_wo_signalcont;
#	Line 524 | Line 534 | void JZBjetScale(TTree *events, float &j
534		if ( !(plotName.Length()>1) ) plotName = informalname;
535
536		nBins = 1; jzbMin = jzbSel*0.95; jzbMax = jzbSel*1.05;
537	<	TH1F* hist = plotEff(events,(kbase&&ksel&&nJets),informalname,flipped);
528	<
537	>	TH1F* hist  = plotEff(events,(kbase&&ksel&&nJets),informalname,flipped);
538		TH1F* histp = plotEff(events,(kbase&&ksel&&nJetsP),informalname,flipped);
530	–
539		TH1F* histm = plotEff(events,(kbase&&ksel&&nJetsM),informalname,flipped);
540
541		// Dump some information
#	Line 571 | Line 579 | void doJZBscale(TTree *events, float &do
579
580		//________________________________________________________________________
581		// JZB response (true/reco. vs. true)
582	<	void JZBresponse(TTree *events, bool requireZ, float &resp, float &resperr, int flipped, string addcut="", bool isMET = kFALSE, Float_t myJzbMax = 200., Int_t nPeriods = 9 ) {
582	>	void JZBresponse(TTree *events, string name, bool requireZ, float &resp, float &resperr, int flipped, string addcut="", bool isMET = kFALSE, Float_t myJzbMax = 200., Int_t nPeriods = 9 ) {
583
584		jzbMin = 20;
585		flag_this_change(__FUNCTION__,__LINE__,false);//PlottingSetup::RestrictToMassPeak
#	Line 657 | Line 665 | void do_systematics_for_one_file(TTree *
665		if(PlottingSetup::computeJZBefficiency) JZBefficiency(events,informalname,jzbeff,jzbefferr,flipped,requireZ,addcut);
666		if(!automatized) dout << "JZB efficiency: " << jzbeff << "+/-" << jzbefferr << endl;
667
668	<	if(!automatized) dout << "Error from Peak position:" << endl;
668	>	if(!automatized) dout << "Error from Peak position:";
669		float sysfrompeak=0;
670		PeakError(events,sysfrompeak,mcjzb,peakerror,flipped,addcut);
671
672	<	if(!automatized) dout << "Jet energy scale: " << std::endl;
672	>	if(!automatized) dout << "Jet energy scale (JES): " << std::endl;
673		float jesup,jesdown;
674		JZBjetScale(events,jesdown,jesup,informalname,flipped,requireZ,addcut);
675	<
675	>
676		if(!automatized) dout << "JZB scale: " << std::endl;
677		float scaleup,scaledown,scalesyst;
678		doJZBscale(events,scaledown,scaleup,scalesyst,JZBScaleUncert,informalname,flipped,requireZ,addcut);
#	Line 673 | Line 681 | void do_systematics_for_one_file(TTree *
681		float resp,resperr;
682		if(PlottingSetup::computeJZBresponse) {
683		if(!automatized) dout << "JZB response: " << std::endl;
684	<	if(!ismSUGRA) JZBresponse(events,requireZ,resp,resperr,flipped,addcut);
684	>	if(!ismSUGRA) JZBresponse(events,informalname,requireZ,resp,resperr,flipped,addcut);
685		}
686
687		if(!automatized) dout << "Pileup: " << std::endl;
#	Line 725 | Line 733 | void do_systematics_for_one_file(TTree *
733		if(fabs(jesup)>fabs(jesdown)) res.push_back(fabs(jesup)); else res.push_back(fabs(jesdown));
734		if(fabs(scaleup)>fabs(scaledown)) res.push_back(fabs(scaleup)); else res.push_back(fabs(scaledown));
735		//  res.push_back(fabs(resolution));
736	+	res.push_back(0.0);
737		res.push_back(mceff_nosigcont.getValue());
738		res.push_back(mceff_nosigcont.getError());
739		if(ismSUGRA) res.push_back(PDFuncert);
#	Line 735 | Line 744 | vector<vector<float> > compute_systemati
744		automatized=true;
745		vector< vector<float> > systematics;
746		for (int isignal=0; isignal<signalsamples.collection.size();isignal++) {
747	<	dout << "Looking at signal " << (signalsamples.collection)[isignal].filename << endl;
748	<	for(int ibin=0;ibin<bins.size();ibin++) {
749	<	jzbSel=bins[ibin];
750	<	geqleq="geq";
751	<	do_systematics_for_one_file((signalsamples.collection)[isignal].events,(signalsamples.collection)[isignal].Nentries,(signalsamples.collection)[isignal].samplename,systematics,flipped,mcjzb,datajzb,mcpeakerror,requireZ);
752	<	}//end of bin loop
747	>	dout << "Looking at signal " << (signalsamples.collection)[isignal].filename << endl;
748	>	for(int ibin=0;ibin<bins.size();ibin++) {
749	>	jzbSel=bins[ibin];
750	>	geqleq="geq";
751	>	do_systematics_for_one_file((signalsamples.collection)[isignal].events,(signalsamples.collection)[isignal].Nentries,(signalsamples.collection)[isignal].samplename,systematics,flipped,mcjzb,datajzb,mcpeakerror,requireZ);
752	>	}//end of bin loop
753		}//end of signal loop
754		return systematics;
755		}

Diff Legend

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