ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

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

Comparing UserCode/cbrown/Development/Plotting/Modules/ShapeLimit.C (file contents):
Revision 1.7 by buchmann, Fri Apr 27 14:07:39 2012 UTC vs.
Revision 1.13 by buchmann, Thu May 17 19:52:55 2012 UTC

#	Line 18 | Line 18
18		#include <TSystem.h>
19		#include <TRandom3.h>
20
21	–	//#include "TTbar_stuff.C"
21		using namespace std;
22
23		using namespace PlottingSetup;
#	Line 44 | Line 43 | void EliminateNegativeEntries(TH1F *hist
43		}
44		}
45
46	<	void prepare_limit_datacard(string RunDirectory, TFile *f, float uJSU, float uPDF) {
47	<	TH1F *dataob = (TH1F*)f->Get("data_obs");
48	<	TH1F *signal = (TH1F*)f->Get("signal");
49	<	TH1F *Tbackground = (TH1F*)f->Get("TTbarBackground");
50	<	TH1F *Zbackground  = (TH1F*)f->Get("ZJetsBackground");
51	<
52	<	assert(dataob);
53	<	assert(signal);
54	<	assert(Tbackground);
55	<	assert(Zbackground);
46	>	vector<float> get_expected_points(float observed,int npoints, string RunDirectory, bool doPlot=false) {
47	>	TF1 *gaus = new TF1("gaus","gaus(0)",0,10*observed);
48	>	gaus->SetParameters(1,observed,2*observed);
49	>	gaus->SetParameter(0,1/(gaus->Integral(0,10*observed)));
50	>	float currentpoint=0.01;
51	>	float stepsize=observed/100;
52	>	vector<float> points;
53	>	while(currentpoint<25*observed && points.size()<npoints) {
54	>
55	>	if(gaus->Integral(0,currentpoint)>((points.size()+1.0)/(npoints+2.0))) {
56	>	if(Verbosity>0) dout << "Added points for calculation at " << currentpoint << " (integral is " << gaus->Integral(0,currentpoint) << ")" << endl;
57	>	points.push_back(currentpoint);
58	>	if(points.size()==npoints) break;
59	>	}
60	>	currentpoint+=stepsize;
61	>	}
62	>	if(doPlot) {
63	>	gaus->SetName("Expected limit computation points");
64	>	gaus->SetTitle("Expected limit computation points");
65	>	TCanvas *can = new TCanvas("can","can");
66	>	gaus->Draw();
67	>	TLine *lines[npoints];
68	>	for(int i=0;i<npoints;i++) {
69	>	lines[i] = new TLine(points[i],0,points[i],gaus->GetMaximum());
70	>	lines[i]->SetLineColor(kBlue);
71	>	lines[i]->Draw();
72	>	}
73	>	can->SaveAs((RunDirectory+"/DistributionOfComputedPoints.png").c_str());
74	>	delete can;
75	>	for(int i=0;i<npoints;i++) delete lines[i];
76	>	}
77	>	delete gaus;
78	>
79	>	return points;
80	>	}
81	>
82	>	void prepare_MC_datacard(string RunDirectory, TFile *f, float uJSU, float uPDF) {
83	>	TH1F *dataob = (TH1F*)f->Get("data_obs");
84	>	TH1F *signal = (TH1F*)f->Get("signal");
85	>	assert(dataob);
86	>	assert(signal);
87	>
88	>	write_warning(__FUNCTION__,"The following two defs for th1's are fake");TH1F *Tbackground; TH1F *Zbackground;
89	>	write_warning(__FUNCTION__,"Not correctly implemented yet for MC");
90	>	ofstream datacard;
91	>	write_warning(__FUNCTION__,"Need to rethink systematics we want to use !");
92	>	datacard.open ((RunDirectory+"/susydatacard.txt").c_str());
93	>	dout << "Writing this to a file.\n";
94	>	datacard << "imax 1\n"; // number of channels
95	>	datacard << "jmax 2\n"; // number of backgrounds (Z+Jets prediction and TTbar prediction)
96	>	datacard << "kmax *\n"; // number of nuisance parameters (sources of systematic uncertainties)
97	>	datacard << "---------------\n";
98	>	datacard << "shapes * * limitfile.root $PROCESS $PROCESS_$SYSTEMATIC\n";
99	>	datacard << "---------------\n";
100	>	datacard << "bin 1\n";
101	>	datacard << "observation "<<dataob->Integral()<<"\n";
102	>	datacard << "------------------------------\n";
103	>	datacard << "bin             1          1          1\n";
104	>	datacard << "process         signal     TTbarBackground     ZJetsBackground\n";
105	>	datacard << "process         0          1          2\n";
106	>	datacard << "rate            "<<signal->Integral()<<"         "<<Tbackground->Integral()<<"         "<<Zbackground->Integral()<<"\n";
107	>	datacard << "--------------------------------\n";
108	>	datacard << "lumi     lnN    " << 1+ PlottingSetup::lumiuncert << "       -       -    luminosity uncertainty\n"; // only affects MC -> signal!
109	>	datacard << "Stat   shape    -          1          -    statistical uncertainty (ttbar)\n";
110	>	datacard << "Stat   shape    -          -          1    statistical uncertainty (zjets)\n";
111	>	datacard << "Sys    shape    -          1          -    systematic uncertainty on ttbar\n";
112	>	datacard << "Sys    shape    -          -          1    systematic uncertainty on zjets\n";
113	>	datacard << "Stat   shape    1          -          -    statistical uncertainty (signal)\n";
114	>	datacard << "JES    shape    1          -          -    uncertainty on Jet Energy Scale\n";
115	>	datacard << "JSU      lnN    " << 1+uJSU << "          -          -    JZB Scale Uncertainty\n";
116	>	if(uPDF>0) datacard << "PDF      lnN    " << 1+uPDF << "          -          -    uncertainty from PDFs\n";
117	>	datacard.close();
118	>	}
119	>
120	>	void prepare_datadriven_datacard(string RunDirectory, TFile *f, float uJSU, float uPDF) {
121	>	TH1F *dataob = (TH1F*)f->Get("data_obs");
122	>	TH1F *signal = (TH1F*)f->Get("signal");
123	>	TH1F *Tbackground = (TH1F*)f->Get("TTbarBackground");
124	>	TH1F *Zbackground  = (TH1F*)f->Get("ZJetsBackground");
125	>
126	>	assert(dataob);
127	>	assert(signal);
128	>	assert(Tbackground);
129	>	assert(Zbackground);
130
131	<	ofstream datacard;
132	<	write_warning(__FUNCTION__,"Need to rethink systematics we want to use !");
133	<	datacard.open ((RunDirectory+"/susydatacard.txt").c_str());
134	<	datacard << "Writing this to a file.\n";
135	<	datacard << "imax 1\n"; // number of channels
136	<	datacard << "jmax 2\n"; // number of backgrounds (Z+Jets prediction and TTbar prediction)
137	<	datacard << "kmax *\n"; // number of nuisance parameters (sources of systematic uncertainties)
138	<	datacard << "---------------\n";
139	<	datacard << "shapes * * limitfile.root $PROCESS $PROCESS_$SYSTEMATIC\n";
140	<	datacard << "---------------\n";
141	<	datacard << "bin 1\n";
142	<	datacard << "observation "<<dataob->Integral()<<"\n";
143	<	datacard << "------------------------------\n";
144	<	datacard << "bin             1          1          1\n";
145	<	datacard << "process         signal     TTbarBackground     ZJetsBackground\n";
146	<	datacard << "process         0          1          2\n";
147	<	datacard << "rate            "<<signal->Integral()<<"         "<<Tbackground->Integral()<<"         "<<Zbackground->Integral()<<"\n";
148	<	datacard << "--------------------------------\n";
149	<	datacard << "lumi     lnN    " << 1+ PlottingSetup::lumiuncert << "       -       -    luminosity uncertainty\n"; // only affects MC -> signal!
150	<	datacard << "Stat   shape    -          1          -    statistical uncertainty (ttbar)\n";
151	<	datacard << "Stat   shape    -          -          1    statistical uncertainty (zjets)\n";
152	<	datacard << "Sys    shape    -          1          -    systematic uncertainty on ttbar\n";
153	<	datacard << "Sys    shape    -          -          1    systematic uncertainty on zjets\n";
154	<	datacard << "Stat   shape    1          -          -    statistical uncertainty (signal)\n";
155	<	datacard << "JES    shape    1          -          -    uncertainty on Jet Energy Scale\n";
156	<	datacard << "JSU      lnN    " << 1+uJSU << "          -          -    JZB Scale Uncertainty\n";
157	<	if(uPDF>0) datacard << "PDF      lnN    " << 1+uPDF << "          -          -    uncertainty from PDFs\n";
131	>
132	>	ofstream datacard;
133	>	write_warning(__FUNCTION__,"Need to rethink systematics we want to use !");
134	>	datacard.open ((RunDirectory+"/susydatacard.txt").c_str());
135	>	datacard << "Writing this to a file.\n";
136	>	datacard << "imax " << dataob->GetNbinsX() << "\n"; // number of channels
137	>	datacard << "jmax 2\n"; // number of backgrounds (Z+Jets prediction and TTbar prediction)
138	>	datacard << "kmax *\n"; // number of nuisance parameters (sources of systematic uncertainties)
139	>	datacard << "---------------\n";
140	>	datacard << "bin\t";
141	>	for(int i=1;i<=dataob->GetNbinsX();i++) datacard << " " << dataob->GetBinLowEdge(i) << "to" << dataob->GetBinLowEdge(i)+dataob->GetBinWidth(i) << " \t";
142	>	datacard << "\n";
143	>	datacard << "observation\t";
144	>	for(int i=1;i<=dataob->GetNbinsX();i++) datacard << " " << dataob->GetBinContent(i) << " \t";
145	>	datacard<<"\n";
146	>	datacard << "------------------------------\n";
147	>	datacard << "bin\t";
148	>	for(int i=1;i<=dataob->GetNbinsX();i++) {
149	>	datacard << " " << dataob->GetBinLowEdge(i) << "to" << dataob->GetBinLowEdge(i)+dataob->GetBinWidth(i) << "\t" <<
150	>	" " << dataob->GetBinLowEdge(i) << "to" << dataob->GetBinLowEdge(i)+dataob->GetBinWidth(i) << "\t" <<
151	>	" " << dataob->GetBinLowEdge(i) << "to" << dataob->GetBinLowEdge(i)+dataob->GetBinWidth(i) << "\t";
152	>	}
153	>	datacard << "\n";
154	>	datacard << "process\t";
155	>	for(int i=1;i<=dataob->GetNbinsX();i++) datacard << "\t signal \t TTbarBackground \t ZJetsBackground";
156	>	datacard << "\n";
157	>	datacard << "process\t";
158	>	for(int i=1;i<=dataob->GetNbinsX();i++) datacard << "\t 0 \t 1 \t 2";
159	>	datacard << "\n";
160	>
161	>
162	>	datacard << "rate\t";
163	>	for(int i=1;i<=dataob->GetNbinsX();i++) datacard << "\t " << signal->GetBinContent(i) << " \t " << Tbackground->GetBinContent(i) << " \t " << Zbackground->GetBinContent(i) << "\t";
164	>	datacard<<"\n";
165	>
166	>	datacard << "--------------------------------\n";
167	>
168	>
169	>	datacard << "lumi     lnN    \t";
170	>	for(int i=1;i<=dataob->GetNbinsX();i++) datacard << " " << 1+ PlottingSetup::lumiuncert << "\t - \t -";
171	>	datacard << "       luminosity uncertainty\n"; // only affects MC -> signal!
172	>
173	>	// Statistical uncertainty
174	>	datacard << "Stat     lnN    \t";
175	>	for(int i=1;i<=dataob->GetNbinsX();i++) {
176	>	//Signal
177	>	float central = signal->GetBinContent(i);
178	>	float up      = ((TH1F*)f->Get("signal_StatDown"))->GetBinContent(i);
179	>	float down    = ((TH1F*)f->Get("signal_StatUp"))->GetBinContent(i);
180	>	float suncert=0;
181	>	if(central>0) {
182	>	if(abs(up-central)>abs(down-central)) suncert=abs(up-central)/central;
183	>	else suncert=abs(central-down)/central;
184	>	}
185	>
186	>	//TTbar
187	>	central = Tbackground->GetBinContent(i);
188	>	up      = ((TH1F*)f->Get("TTbarBackground_StatUp"))->GetBinContent(i);
189	>	down    = ((TH1F*)f->Get("TTbarBackground_StatDown"))->GetBinContent(i);
190	>	float tuncert=0;
191	>	if(central>0) {
192	>	if(abs(up-central)>abs(down-central)) tuncert=abs(up-central)/central;
193	>	else tuncert=abs(central-down)/central;
194	>	}
195	>	//ZJets
196	>	central = Zbackground->GetBinContent(i);
197	>	up      = ((TH1F*)f->Get("ZJetsBackground_StatUp"))->GetBinContent(i);
198	>	down    = ((TH1F*)f->Get("ZJetsBackground_StatDown"))->GetBinContent(i);
199	>	float zuncert=0;
200	>	if(central>0) {
201	>	if(abs(up-central)>abs(down-central)) zuncert=abs(up-central)/central;
202	>	else zuncert=abs(central-down)/central;
203	>	}
204	>	datacard << " " << 1+suncert << " \t " << 1+tuncert << "\t " << 1+zuncert << " \t ";
205	>	}
206	>
207	>	datacard << "       statistical uncertainty\n";
208	>
209	>	// Statistical uncertainty
210	>	datacard << "Sys     lnN    \t";
211	>	for(int i=1;i<=dataob->GetNbinsX();i++) {
212	>	float central = Tbackground->GetBinContent(i);
213	>	float up      = ((TH1F*)f->Get("TTbarBackground_SysUp"))->GetBinContent(i);
214	>	float down    = ((TH1F*)f->Get("TTbarBackground_SysDown"))->GetBinContent(i);
215	>	float tuncert=0;
216	>	if(central>0) {
217	>	if(abs(up-central)>abs(down-central)) tuncert=abs(up-central)/central;
218	>	else tuncert=abs(central-down)/central;
219	>	}
220	>	central = Zbackground->GetBinContent(i);
221	>	up      = ((TH1F*)f->Get("ZJetsBackground_SysUp"))->GetBinContent(i);
222	>	down    = ((TH1F*)f->Get("ZJetsBackground_SysDown"))->GetBinContent(i);
223	>	float zuncert=0;
224	>	if(central>0) {
225	>	if(abs(up-central)>abs(down-central)) zuncert=abs(up-central)/central;
226	>	else zuncert=abs(central-down)/central;
227	>	}
228	>	datacard << " - \t " << 1+tuncert << "\t " << 1+zuncert << " \t ";
229	>	}
230	>
231	>	datacard << "       systematic uncertainty\n";
232	>
233	>
234	>	float JESup = ((TH1F*)f->Get("signal_JESUp"))->Integral();
235	>	float JESdn = ((TH1F*)f->Get("signal_JESDown"))->Integral();
236	>	float central = signal->Integral();
237	>	float uJES=0;
238	>	if(abs(JESup-central)>abs(JESdn-central)) uJES=abs(JESup-central)/central;
239	>	else uJES=abs(central-JESdn)/central;
240	>	datacard << "JES    lnN";
241	>	for(int i=1;i<=dataob->GetNbinsX();i++) datacard << " " << 1+uJES << "\t - \t  - \t";
242	>	datacard << "uncertainty on Jet Energy Scale\n";
243	>
244	>	datacard << "JSU      lnN    ";
245	>	for(int i=1;i<=dataob->GetNbinsX();i++) datacard << " " << 1+uJSU << "\t - \t  - \t";
246	>	datacard << "JZB Scale Uncertainty\n";
247	>
248	>	if(uPDF>0) {
249	>	datacard << "PDF      lnN    ";
250	>	for(int i=1;i<=dataob->GetNbinsX();i++) datacard << " " << 1+uPDF << "\t - \t - \t";
251	>	datacard << "uncertainty from PDFs\n";
252	>	}
253	>
254	>	datacard.close();
255	>	}
256	>
257	>	void prepare_limit_datacard(string RunDirectory, TFile *f, float uJSU, float uPDF) {
258
259	<	// datacard << "peak   shape    1          1    uncertainty on signal resolution.
260	<	datacard.close();
259	>	if(FullMCAnalysis) {
260	>	//dealing with MC analysis - need to use shapes.
261	>	prepare_MC_datacard(RunDirectory,f,uJSU,uPDF);
262	>	} else {
263	>	//doing mutibin analysis
264	>	prepare_datadriven_datacard(RunDirectory,f,uJSU,uPDF);
265	>	}
266	>
267		}
268
269		float QuickDrawNevents=0;
#	Line 117 | Line 296 | TH1F* QuickDraw(TTree *events, string hn
296		return histo;
297		}
298
120	–	/*void do_stat_up(TH1F *h, float sign=1.0) {
121	–	for(int i=1;i<=h->GetNbinsX();i++) {
122	–	h->SetBinContent(i,h->GetBinContent(i)+sign*h->GetBinError(i));
123	–	}
124	–	h->Write();
125	–	}
126	–
127	–	void do_stat_dn(TH1F *h) {
128	–	do_stat_up(h,-1.0);
129	–	}*/
299
300		void SQRT(TH1F *h) {
301		for (int i=1;i<=h->GetNbinsX();i++) {
#	Line 142 | Line 311 | TH1F* Multiply(TH1F *h1, TH1F *h2) {
311		return h;
312		}
313
314	+
315	+	string ReduceNumericHistoName(string origname) {
316	+	int pos = origname.find("__h_");
317	+	if(pos == -1) return origname;
318	+	return origname.substr(0,pos);
319	+	}
320	+
321	+
322	+
323	+	void generate_mc_shapes(bool dotree, TTree *signalevents, string mcjzb,string datajzb,vector<float> binning, TCut weight, TCut JetCut, string addcut, string identifier, map <  pair<float, float>, map<string, float>  > &xsec) {
324	+	//weight can be varied to take care of efficiency errors (weightEffUp, weightEffDown)
325	+	vector<float> mbinning;
326	+	mbinning.push_back(-8000);
327	+	for(int i=binning.size()-1;i>=0;i--) mbinning.push_back(-binning[i]);
328	+	mbinning.push_back(0);
329	+	for(int i=0;i<binning.size();i++) mbinning.push_back(binning[i]);
330	+	mbinning.push_back(8000);
331	+
332	+
333	+	TCanvas *shapecan = new TCanvas("shapecan","shapecan");
334	+	TH1F* Signal;
335	+
336	+	stringstream sInverseCutWeight;
337	+	sInverseCutWeight << "(1.0/" << (const char*)cutWeight<<")";
338	+	TCut InverseCutWeight(sInverseCutWeight.str().c_str());
339	+	if(weight==cutWeight) InverseCutWeight = TCut("1.0");
340	+	if(!dotree) {
341	+	//dealing with ALL MC samples (when we save the standard file)
342	+	THStack McPredicted = allsamples.DrawStack("McPred",mcjzb,mbinning, "JZB [GeV]", "events", ((cutmass&&cutOSSF&&JetCut)*(weight*InverseCutWeight)),mc,luminosity);
343	+	int nhists = McPredicted.GetHists()->GetEntries();
344	+	for (int i = 0; i < nhists; i++) {
345	+	TH1* hist = static_cast<TH1*>(McPredicted.GetHists()->At(i));
346	+	//          cout << hist->GetName() << " : " << hist->Integral() << endl;
347	+	//          cout << "     " << ReduceNumericHistoName(hist->GetName()) << endl;
348	+	if(identifier=="StatUp") {
349	+	float appliedweight = hist->Integral() / hist->GetEntries();
350	+	for(int ibin=1;ibin<=hist->GetNbinsX();ibin++) hist->SetBinContent(ibin,hist->GetBinContent(ibin)+appliedweight*TMath::Sqrt(hist->GetBinContent(ibin)/appliedweight));
351	+	}
352	+	if(identifier=="StatDown") {
353	+	float appliedweight = hist->Integral() / hist->GetEntries();
354	+	for(int ibin=1;ibin<=hist->GetNbinsX();ibin++) hist->SetBinContent(ibin,hist->GetBinContent(ibin)-appliedweight*TMath::Sqrt(hist->GetBinContent(ibin)/appliedweight));
355	+	}
356	+	hist->SetName(("mc_"+ReduceNumericHistoName(hist->GetName())+"_"+identifier).c_str());
357	+	hist->Write();
358	+	}
359	+	} else {
360	+	string histoname="mc_signal";
361	+	if(identifier!="") histoname=histoname+"_"+identifier;
362	+	Signal = QuickDraw(signalevents,histoname,mcjzb,binning,"JZB", "events",((cutmass&&cutOSSF&&JetCut&&basiccut))*(weight*InverseCutWeight),addcut,mc,luminosity,xsec);
363	+	if(identifier=="StatUp") {
364	+	float appliedweight = Signal->Integral() / Signal->GetEntries();
365	+	for(int ibin=1;ibin<=Signal->GetNbinsX();ibin++) Signal->SetBinContent(ibin,Signal->GetBinContent(ibin)+appliedweight*TMath::Sqrt(Signal->GetBinContent(ibin)/appliedweight));
366	+	}
367	+	if(identifier=="StatDown") {
368	+	float appliedweight = Signal->Integral() / Signal->GetEntries();
369	+	for(int ibin=1;ibin<=Signal->GetNbinsX();ibin++) Signal->SetBinContent(ibin,Signal->GetBinContent(ibin)-appliedweight*TMath::Sqrt(Signal->GetBinContent(ibin)/appliedweight));
370	+	}
371	+	Signal->Write();
372	+	}
373	+
374	+	/*
375	+	*
376	+	* Jet energy scale (easy, use different JetCut)
377	+	* JZB Scale Uncertainty (will be a number - signal only)
378	+	* Efficiency (will be weightEffUp, weightEffDown)
379	+	* PDFs (signal only) -> Will be a number yet again, computed in the traditional way
380	+	*/
381	+	delete shapecan;
382	+
383	+	}
384	+
385	+
386		void generate_shapes_for_systematic(bool signalonly, bool dataonly,TFile *limfile, TTree *signalevents, string identifier, string mcjzb, string datajzb, int JES,vector<float> binning, TCanvas *limcan, string addcut, map <  pair<float, float>, map<string, float>  > &xsec) {
387
388		binning.push_back(8000);
#	Line 173 | Line 414 | void generate_shapes_for_systematic(bool
414		if(signalonly) {
415		dout << "Processing a signal with mcjzb: " << mcjzb << " (identifier: '" << identifier << "')" << endl;
416		TH1F *ZOSSFP = QuickDraw(signalevents,"ZOSSFP",mcjzb,binning,     "JZB", "events",cutmass&&cutOSSF&&limitnJetcut&&basiccut,addcut,mc,luminosity,xsec);
176	–	TH1F *ZOSOFP = QuickDraw(signalevents,"ZOSOFP",mcjzb,binning,     "JZB", "events",cutmass&&cutOSOF&&limitnJetcut&&basiccut,addcut,mc,luminosity,xsec);
417		TH1F *ZOSSFN = QuickDraw(signalevents,"ZOSSFN","-"+mcjzb,binning, "JZB", "events",cutmass&&cutOSSF&&limitnJetcut&&basiccut,addcut,mc,luminosity,xsec);
418	<	TH1F *ZOSOFN = QuickDraw(signalevents,"ZOSOFN","-"+mcjzb,binning, "JZB", "events",cutmass&&cutOSOF&&limitnJetcut&&basiccut,addcut,mc,luminosity,xsec);
418	>	TH1F *ZOSOFP;
419	>	TH1F *ZOSOFN;
420	>
421	>	if(!PlottingSetup::FullMCAnalysis) {
422	>	ZOSOFP = QuickDraw(signalevents,"ZOSOFP",mcjzb,binning,     "JZB", "events",cutmass&&cutOSOF&&limitnJetcut&&basiccut,addcut,mc,luminosity,xsec);
423	>	ZOSOFN = QuickDraw(signalevents,"ZOSOFN","-"+mcjzb,binning, "JZB", "events",cutmass&&cutOSOF&&limitnJetcut&&basiccut,addcut,mc,luminosity,xsec);
424	>	}
425
426		TH1F *SBOSSFP;
427		TH1F *SBOSOFP;
428		TH1F *SBOSSFN;
429		TH1F *SBOSOFN;
430
431	<	if(PlottingSetup::RestrictToMassPeak) {
431	>	if(PlottingSetup::RestrictToMassPeak&&!PlottingSetup::FullMCAnalysis) {
432		SBOSSFP = QuickDraw(signalevents,"SBOSSFP",mcjzb,binning, "JZB", "events",cutOSSF&&limitnJetcut&&basiccut&&sidebandcut,addcut,mc,luminosity,xsec);
433		SBOSOFP = QuickDraw(signalevents,"SBOSOFP",mcjzb,binning, "JZB", "events",cutOSOF&&limitnJetcut&&basiccut&&sidebandcut,addcut,mc,luminosity,xsec);
434		SBOSSFN = QuickDraw(signalevents,"SBOSSFN","-"+mcjzb,binning, "JZB", "events",cutOSSF&&limitnJetcut&&basiccut&&sidebandcut,addcut,mc,luminosity,xsec);
#	Line 201 | Line 447 | void generate_shapes_for_systematic(bool
447		TH1F *flippedLpred = new TH1F("flippedLpred","flippedLpred",binning.size()-1,&binning[0]);
448
449		Lobs->Add(ZOSSFP);
450	<	Lpred->Add(ZOSSFN);
450	>	if(!PlottingSetup::FullMCAnalysis) Lpred->Add(ZOSSFN);
451
452		dout << "SITUATION FOR SIGNAL: " << endl;
453	<	dout << " OSSF     JZB> 0 : " << ZOSSFP->Integral() << "     JZB < 0 :" << ZOSSFN->Integral() << endl;
454	<	dout << " OSOF     JZB> 0 : " << ZOSOFP->Integral() << "     JZB < 0 :" << ZOSOFN->Integral() << endl;
455	<	if(PlottingSetup::RestrictToMassPeak) {
456	<	dout << " OSSF SB  JZB> 0 : " << SBOSSFP->Integral() << "     JZB < 0 :" << SBOSSFN->Integral() << endl;
457	<	dout << " OSOF SB  JZB> 0 : " << SBOSOFP->Integral() << "     JZB < 0 :" << SBOSOFN->Integral() << endl;
453	>
454	>
455	>	if(PlottingSetup::FullMCAnalysis) {
456	>	dout << " OSSF     JZB> 0 : " << ZOSSFP->Integral() << endl;
457	>	} else {
458	>	dout << " OSSF     JZB> 0 : " << ZOSSFP->Integral() << "     JZB < 0 :" << ZOSSFN->Integral() << endl;
459	>	dout << " OSOF     JZB> 0 : " << ZOSOFP->Integral() << "     JZB < 0 :" << ZOSOFN->Integral() << endl;
460	>	if(PlottingSetup::RestrictToMassPeak&&!PlottingSetup::FullMCAnalysis) {
461	>	dout << " OSSF SB  JZB> 0 : " << SBOSSFP->Integral() << "     JZB < 0 :" << SBOSSFN->Integral() << endl;
462	>	dout << " OSOF SB  JZB> 0 : " << SBOSOFP->Integral() << "     JZB < 0 :" << SBOSOFN->Integral() << endl;
463	>	}
464		}
465
466
467		flippedLobs->Add(ZOSSFN);
468	<	flippedLpred->Add(ZOSSFP);
468	>	if(!PlottingSetup::FullMCAnalysis) flippedLpred->Add(ZOSSFP);
469
470	<	if(PlottingSetup::RestrictToMassPeak) {
471	<	Lpred->Add(ZOSOFP,1.0/3);
472	<	Lpred->Add(ZOSOFN,-1.0/3);
473	<	Lpred->Add(SBOSSFP,1.0/3);
474	<	Lpred->Add(SBOSSFN,-1.0/3);
475	<	Lpred->Add(SBOSOFP,1.0/3);
476	<	Lpred->Add(SBOSOFN,-1.0/3);
477	<
478	<	//flipped prediction
479	<	flippedLpred->Add(ZOSOFP,-1.0/3);
480	<	flippedLpred->Add(ZOSOFN,1.0/3);
481	<	flippedLpred->Add(SBOSSFP,-1.0/3);
482	<	flippedLpred->Add(SBOSSFN,1.0/3);
483	<	flippedLpred->Add(SBOSOFP,-1.0/3);
484	<	flippedLpred->Add(SBOSOFN,1.0/3);
485	<
486	<	} else {
487	<	Lpred->Add(ZOSOFP,1.0);
488	<	Lpred->Add(ZOSOFN,-1.0);
489	<
490	<	//flipped prediction
491	<	flippedLpred->Add(ZOSOFP,-1.0);
492	<	flippedLpred->Add(ZOSOFN,1.0);
470	>	if(!PlottingSetup::FullMCAnalysis) {
471	>	if(PlottingSetup::RestrictToMassPeak) {
472	>	Lpred->Add(ZOSOFP,1.0/3);
473	>	Lpred->Add(ZOSOFN,-1.0/3);
474	>	Lpred->Add(SBOSSFP,1.0/3);
475	>	Lpred->Add(SBOSSFN,-1.0/3);
476	>	Lpred->Add(SBOSOFP,1.0/3);
477	>	Lpred->Add(SBOSOFN,-1.0/3);
478	>
479	>	//flipped prediction
480	>	flippedLpred->Add(ZOSOFP,-1.0/3);
481	>	flippedLpred->Add(ZOSOFN,1.0/3);
482	>	flippedLpred->Add(SBOSSFP,-1.0/3);
483	>	flippedLpred->Add(SBOSSFN,1.0/3);
484	>	flippedLpred->Add(SBOSOFP,-1.0/3);
485	>	flippedLpred->Add(SBOSOFN,1.0/3);
486	>	} else {
487	>	Lpred->Add(ZOSOFP,1.0);
488	>	Lpred->Add(ZOSOFN,-1.0);
489	>
490	>	//flipped prediction
491	>	flippedLpred->Add(ZOSOFP,-1.0);
492	>	flippedLpred->Add(ZOSOFN,1.0);
493	>	}
494		}
495
496		TH1F *signal = (TH1F*)Lobs->Clone("signal");
497	<	signal->Add(Lpred,-1);
497	>	if(!PlottingSetup::FullMCAnalysis) signal->Add(Lpred,-1);
498		signal->SetName(signalname.c_str());
499		signal->SetTitle(signalname.c_str());
500		signal->Write();
501
502		TH1F *flippedsignal = (TH1F*)flippedLobs->Clone();
503	<	flippedsignal->Add(flippedLpred,-1);
503	>	if(!PlottingSetup::FullMCAnalysis) flippedsignal->Add(flippedLpred,-1);
504		flippedsignal->SetName(("flipped_"+signalname).c_str());
505		flippedsignal->Write();
506
#	Line 261 | Line 514 | void generate_shapes_for_systematic(bool
514		signalStatDn->SetTitle(((string)signal->GetTitle()+"_StatDown").c_str());
515
516		for(int i=1;i<=signalStatDn->GetNbinsX();i++) {
517	<	float staterr = TMath::Sqrt(Lpred->GetBinContent(i) + Lobs->GetBinContent(i));
518	<	dout << "Stat err in bin " << i << " : " << staterr << endl;
519	<	dout << "    prediction: " << Lpred->GetBinContent(i) << " , observation: " << Lobs->GetBinContent(i) << " --> signal: " << signal->GetBinContent(i) << endl;
520	<	dout << "    we obtain : " << signal->GetBinContent(i)-staterr << " , " << signal->GetBinContent(i)+staterr << endl;
517	>	float staterr;
518	>	if(!PlottingSetup::FullMCAnalysis) staterr = TMath::Sqrt(Lpred->GetBinContent(i) + Lobs->GetBinContent(i));
519	>	else staterr = TMath::Sqrt(Lobs->GetBinContent(i));
520	>
521	>	if(!PlottingSetup::FullMCAnalysis) {
522	>	dout << "Stat err in bin " << i << " : " << staterr << endl;
523	>	dout << "    prediction: " << Lpred->GetBinContent(i) << " , observation: " << Lobs->GetBinContent(i) << " --> signal: " << signal->GetBinContent(i) << endl;
524	>	dout << "    we obtain : " << signal->GetBinContent(i)-staterr << " , " << signal->GetBinContent(i)+staterr << endl;
525	>	}
526		if(signal->GetBinContent(i)-staterr>0) signalStatDn->SetBinContent(i,signal->GetBinContent(i)-staterr);
527		else signalStatDn->SetBinContent(i,0);
528		signalStatUp->SetBinContent(i,signal->GetBinContent(i)+staterr);
#	Line 286 | Line 544 | void generate_shapes_for_systematic(bool
544		flippedsignalStatDn->SetTitle(((string)flippedsignal->GetTitle()+"_StatDown").c_str());
545
546		for(int i=1;i<=flippedsignalStatDn->GetNbinsX();i++) {
547	<	float staterr = TMath::Sqrt(flippedLpred->GetBinContent(i) + flippedLobs->GetBinContent(i));
547	>	float staterr;
548	>	if(!PlottingSetup::FullMCAnalysis) staterr = TMath::Sqrt(flippedLpred->GetBinContent(i) + flippedLobs->GetBinContent(i));
549	>	else staterr = TMath::Sqrt(flippedLobs->GetBinContent(i));
550		if(flippedsignal->GetBinContent(i)-staterr>0) flippedsignalStatDn->SetBinContent(i,flippedsignal->GetBinContent(i)-staterr);
551		else flippedsignalStatDn->SetBinContent(i,0);
552		flippedsignalStatUp->SetBinContent(i,flippedsignal->GetBinContent(i)+staterr);
#	Line 686 | Line 946 | void generate_shapes_for_systematic(bool
946
947		}
948
949	<	ShapeDroplet LimitsFromShapes(bool asymptotic, float firstGuess, TTree *events, string addcut, string name, string mcjzb, string datajzb, vector<float> jzbbins, float jzbpeakerrormc, float jzbpeakerrordata) {
949	>	void DoFullCls(string RunDirectory,float firstGuess) {
950	>	vector<float> epoints = get_expected_points(firstGuess,25,RunDirectory,true);//get 20 points;
951	>	ofstream RealScript;
952	>	dout << "Going to write the full CLs script to " << RunDirectory << "/LimitScript.sh" << endl;
953	>	RealScript.open ((RunDirectory+"/LimitScript.sh").c_str());
954	>	RealScript << "#!/bin/bash \n";
955	>	RealScript << "\n";
956	>	RealScript << "RunDirectory=" << RunDirectory << "\n";
957	>	RealScript << "CMSSWDirectory=" << PlottingSetup::CMSSW_dir << "\n";
958	>	RealScript << "ORIGTMP=$TMP\n";
959	>	RealScript << "ORIGTMPDIR=$TMPDIR\n";
960	>	RealScript << "export TMP=" << RunDirectory << "\n";
961	>	RealScript << "export TMPDIR=" << RunDirectory << "\n";
962	>	RealScript << "\n";
963	>	RealScript << "echo \"Going to compute limit in $RunDirectory\"\n";
964	>	RealScript << "ORIGDIR=`pwd`\n";
965	>	RealScript << "\n";
966	>	RealScript << "pushd $CMSSWDirectory\n";
967	>	RealScript << "cd ~/final_production_2011/CMSSW_4_2_8/src/HiggsAnalysis/\n";
968	>	RealScript << "origscramarch=$SCRAM_ARCH\n";
969	>	RealScript << "origbase=$CMSSW_BASE\n";
970	>	RealScript << "export SCRAM_ARCH=slc5_amd64_gcc434\n";
971	>	RealScript << "eval `scram ru -sh`\n";
972	>	RealScript << "\n";
973	>	RealScript << "mkdir -p $RunDirectory\n";
974	>	RealScript << "cd $RunDirectory\n";
975	>	RealScript << "echo `pwd`\n";
976	>	RealScript << "mkdir -p FullCLs\n";
977	>	RealScript << "cd FullCLs\n";
978	>	RealScript << "\n";
979	>	RealScript << "combineEXE=\"${CMSSW_BASE}/bin/${SCRAM_ARCH}/combine\"\n";
980	>	RealScript << "\n";
981	>	RealScript << "dobreak=0\n";
982	>	RealScript << "npoints=0\n";
983	>	RealScript << "time for i in {";
984	>	RealScript << epoints[0]/4 << ","; // adding a VERY VERY low point just to be totally safe
985	>	RealScript << epoints[0]/2 << ","; // adding a VERY low point just to be safe
986	>	for(int i=0;i<epoints.size();i++) RealScript << epoints[i] << ",";
987	>	RealScript << 2*epoints[epoints.size()-1] << ","; // adding a VERY high point just to be safe
988	>	RealScript << 4*epoints[epoints.size()-1]; // adding a VERY VERY high point just to be totally safe
989	>	RealScript << "}; do \n";
990	>	RealScript << "let \"npoints = npoints +1\"\n";
991	>	RealScript << "if [[ $dobreak -gt 0 ]]; then \n";
992	>	RealScript << "    continue \n";
993	>	RealScript << "fi \n";
994	>	RealScript << "echo -e \" Going to compute CLs for $i \\n\"\n";
995	>	RealScript << "echo $(date +%d%b%Y,%H:%M:%S)\n";
996	>	RealScript << "time " << PlottingSetup::cbafbasedir << "/DistributedModelCalculations/Utilities/TimeOut \"$combineEXE ../susydatacard.txt -M HybridNew --freq -T 500 --clsAcc 0 -v 0 -n TestPoint --saveHybridResult --saveToys -s -1 -i 8 --singlePoint $i >/dev/null 2>/dev/null\"\n";
997	>	//  RealScript << "$combineEXE ../susydatacard.txt -M HybridNew --freq -T 500 --clsAcc 0 -v 0 -n TestPoint --saveHybridResult --saveToys -s -1 -i 8 --singlePoint $i >/dev/null 2>/dev/null \n";
998	>	RealScript << "errorsize=$?\n";
999	>	RealScript << "rm " << RunDirectory << "/rstat* 2>/dev/null \n";
1000	>	RealScript << "if [[ $errorsize -gt 0 ]]; then \n";
1001	>	RealScript << "    echo \"Apparently something went wrong (had to be aborted)\"\n";
1002	>	RealScript << "    if [[ $npoints -gt 10 ]]; then \n";
1003	>	RealScript << "        dobreak=1 \n";
1004	>	RealScript << "    fi \n";
1005	>	RealScript << "fi \n";
1006	>	RealScript << "done\n";
1007	>
1008	>	RealScript << "\n";
1009	>	RealScript << "hadd -f FullCLs.root higgsCombine*.root\n";
1010	>	RealScript << "mkdir originalFiles\n";
1011	>	RealScript << "mv higgsCombine* originalFiles\n";
1012	>	RealScript << "$combineEXE ../susydatacard.txt -M HybridNew --freq --grid=FullCLs.root\n";
1013	>	RealScript << "$combineEXE ../susydatacard.txt -M HybridNew --freq --grid=FullCLs.root --expectedFromGrid 0.5\n";
1014	>	RealScript << "$combineEXE ../susydatacard.txt -M HybridNew --freq --grid=FullCLs.root --expectedFromGrid 0.84\n";
1015	>	RealScript << "$combineEXE ../susydatacard.txt -M HybridNew --freq --grid=FullCLs.root --expectedFromGrid 0.16\n";
1016	>	RealScript << "$combineEXE ../susydatacard.txt -M HybridNew --freq --grid=FullCLs.root --expectedFromGrid 0.975\n";
1017	>	RealScript << "$combineEXE ../susydatacard.txt -M HybridNew --freq --grid=FullCLs.root --expectedFromGrid 0.025\n";
1018	>	RealScript << "\n";
1019	>	RealScript << "hadd FinalResult.root higgsCombineTest.HybridNew*\n";
1020	>	RealScript << "\n";
1021	>	RealScript << "g++ " << PlottingSetup::cbafbasedir << "/DistributedModelCalculations/ShapeLimits/ReadAndSave.C -o ReadAndSave.exec `root-config --glibs --cflags` \n";
1022	>	RealScript << "./ReadAndSave.exec FinalResult.root " << RunDirectory << "/FullCLsShapeDropletResult.txt \n";
1023	>	RealScript << "\n";
1024	>	RealScript << "rm " << RunDirectory << "/rstat* \n";
1025	>	RealScript << "\n";
1026	>	RealScript << "#####\n";
1027	>	RealScript << "echo \"Finalizing ...\"\n";
1028	>	RealScript << "cd $ORIGDIR\n";
1029	>	RealScript << "echo $ORIGDIR\n";
1030	>	RealScript << "ls -ltrh ../SandBox/ | grep combine\n";
1031	>	RealScript << "export SCRAM_ARCH=$origscramarch\n";
1032	>	RealScript << "export TMP=$ORIGTMP\n";
1033	>	RealScript << "export TMPDIR=$ORIGTMPDIR\n";
1034	>	RealScript << "exit 0\n";
1035	>	RealScript << "\n";
1036	>	RealScript.close();
1037	>	gSystem->Exec(("bash "+RunDirectory+"/LimitScript.sh").c_str());
1038	>	}
1039	>
1040	>
1041	>	ShapeDroplet LimitsFromShapes(float low_fullCLs, float high_CLs, TTree *events, string addcut, string name, string mcjzb, string datajzb, vector<float> jzbbins, float jzbpeakerrormc, float jzbpeakerrordata) {
1042		map <  pair<float, float>, map<string, float>  >  xsec;
1043		if(SUSYScanSpace::SUSYscantype==mSUGRA) xsec=getXsec(PlottingSetup::cbafbasedir+"/"+mSUGRAxsFile);
1044
#	Line 723 | Line 1075 | ShapeDroplet LimitsFromShapes(bool asymp
1075		bool dataonly=false;
1076
1077		generate_shapes_for_systematic(signalonly,dataonly,limfile,events,"",mcjzb,datajzb,noJES,jzbbins,limcan,addcut,xsec);
726	–	//  generate_shapes_for_systematic(signalonly,dataonly,limfile,events,"peakUp",newjzbexpression(mcjzb,jzbpeakerrormc),newjzbexpression(datajzb,jzbpeakerrordata),noJES,jzbbins,limcan,addcut,xsec);
727	–	//  generate_shapes_for_systematic(signalonly,dataonly,limfile,events,"peakDown",newjzbexpression(mcjzb,-jzbpeakerrormc),newjzbexpression(datajzb,-jzbpeakerrordata),noJES,jzbbins,limcan,addcut,xsec);
1078		generate_shapes_for_systematic(signalonly,dataonly,limfile,events,"JESUp",mcjzb,datajzb,JESup,jzbbins,limcan,addcut,xsec);
1079		generate_shapes_for_systematic(signalonly,dataonly,limfile,events,"JESDown",mcjzb,datajzb,JESdown,jzbbins,limcan,addcut,xsec);
1080
#	Line 793 | Line 1143 | ShapeDroplet LimitsFromShapes(bool asymp
1143		limfile->Close();
1144		dout << "Info: Shape root file and datacard have been generated in " << RunDirectory << endl;
1145		stringstream command;
1146	<	if(asymptotic) {
1147	<	if(firstGuess>0) command << "bash CreateModel.sh " << RunDirectory << " susydatacard.txt" << " 1 0 0 " << firstGuess; // ASYMPTOTIC LIMITS WITH FIRST GUESS
1148	<	else command << "bash CreateModel.sh " << RunDirectory << " susydatacard.txt" << " 0 0 0 0"; // ASYMPTOTIC LIMITS
1149	<	}
1150	<	else command << "bash CreateModel.sh " << RunDirectory << " susydatacard.txt" << " 2 0 " << int(0.5 * firstGuess) << " " << int(2*firstGuess); // ASYMPTOTIC LIMITS
1146	>	string DropletLocation;
1147	>	int CreatedModelFileExitCode;
1148	>	//----------------------------------------------------------------
1149	>	//do an asymptotic limit first
1150	>	command << "bash " << PlottingSetup::cbafbasedir<< "/DistributedModelCalculations/ShapeLimits/CreateModel.sh " << RunDirectory << " susydatacard.txt" << " 0 0 0 0"; // ASYMPTOTIC LIMITS
1151		dout <<"Going to run : " << command.str() << endl;
1152	<	int CreatedModelFileExitCode = gSystem->Exec(command.str().c_str());
1152	>	CreatedModelFileExitCode = gSystem->Exec(command.str().c_str());
1153		dout << "exit code of limit algorithm (CreateModel.sh) : " << CreatedModelFileExitCode << endl;
1154		if(!(CreatedModelFileExitCode==0)) {
1155	<	write_warning(__FUNCTION__,"Something bad happened. It looks like a shape analysis is not the way to go. ");
1156	<	ShapeDroplet alpha;
1157	<	alpha.observed=-12345; // this is the flag to say watch out something went wrong with the signal ...
1158	<	alpha.SignalIntegral=1;
1159	<	return alpha;
1155	>	write_warning(__FUNCTION__,"Something bad happened. It looks like a shape analysis is not the way to go. ");
1156	>	ShapeDroplet beta;
1157	>	beta.observed=-12345; // this is the flag to say watch out something went wrong with the signal ...
1158	>	beta.SignalIntegral=1;
1159	>	return beta;
1160	>	}
1161	>	DropletLocation=RunDirectory+"/ShapeDropletResult.txt";
1162	>	ShapeDroplet beta;
1163	>	beta.readDroplet(DropletLocation);
1164	>	float obsLimit = beta.observed/SignalIntegral;
1165	>	dout << "Checking if the obtained limit (" << beta.observed << " / " << SignalIntegral << " = " << beta.observed/SignalIntegral << " is in [" << low_fullCLs << " , " << high_CLs << "]" << endl;
1166	>	if(obsLimit<high_CLs && obsLimit>low_fullCLs) {
1167	>	//if(PlottingSetup::scantype!=mSUGRA||(obsLimit<high_CLs && obsLimit>low_fullCLs)) {
1168	>	dout << " It is! Full CLs ahead!" << endl;
1169	>	DoFullCls(RunDirectory,beta.observed);
1170	>	DropletLocation=RunDirectory+"/FullCLsShapeDropletResult.txt";
1171	>	} else {
1172	>	dout << " It is not ... happy with asymptotic limits." << endl;
1173		}
1174	+	//----------------------------------------------------------------
1175		ShapeDroplet alpha;
1176	<	alpha.readDroplet(RunDirectory+"/ShapeDropletResult.txt");
1176	>	alpha.readDroplet(DropletLocation);
1177		alpha.PDF=PDFuncert;
1178		alpha.JSU=JZBscale;
1179		alpha.SignalIntegral=SignalIntegral;
1180		dout << "Done reading limit information from droplet" << endl;
1181		dout << alpha << endl;
1182	<
1182	>
1183		dout << "Everything is saved in " << RunDirectory << endl;
1184		dout << "Cleaning up ... " << std::flush;
1185	<	/*  dout << "   1) Make sure models directory exists ... " << std::flush;
1186	<	gSystem->Exec("mkdir -p models/");
823	<	dout << " ok!" << endl;
824	<	dout << "   2) Deleting any previous model files with the same name ... " << std::flush;
825	<	stringstream delcommand;
826	<	delcommand << "rm models/model_" << name << ".root";
827	<	gSystem->Exec(delcommand.str().c_str());
828	<	stringstream delcommand2;
829	<	delcommand2 << "rm models/model_" << name << ".txt";
830	<	gSystem->Exec(delcommand2.str().c_str());
831	<	dout << " ok!" << endl;
832	<	dout << "   3) Transfering the new model files ... " << std::flush;
833	<	stringstream copycommand;
834	<	copycommand << "cp " << RunDirectory << "/model.root models/model_" << name << ".root";
835	<	gSystem->Exec(copycommand.str().c_str());
836	<	stringstream copycommand2;
837	<	copycommand2 << "cp " << RunDirectory << "/susydatacard.txt models/model_" << name << ".txt";
838	<	gSystem->Exec(copycommand2.str().c_str());
839	<	stringstream copycommand3;
840	<	copycommand3 << "cp " << RunDirectory << "/limitfile.root models/model_" << name << "_histo.root";
841	<	gSystem->Exec(copycommand3.str().c_str());
842	<	dout << " ok!" << endl;
843	<	dout << "   4) Removing original working directory (" << RunDirectory << ") ... " << std::flush;*/
844	<	gSystem->Exec(("rm -r "+RunDirectory).c_str());
1185	>	write_warning(__FUNCTION__,"NOT CLEANING UP");
1186	>	gSystem->Exec(("rm -rf "+RunDirectory).c_str());
1187		dout << " ok!" << endl;
1188		delete limcan;
1189		return alpha;
1190		}
1191
1192	<	void PrepareDataShapes(string datajzb, vector<float> jzbbins, float jzbpeakerrordata) {
1193	<	//  dout << "Generating data shape templates ... " << endl;
1192	>	void PrepareDataShapes(string mcjzb, string datajzb, vector<float> jzbbins, float jzbpeakerrordata) {
1193	>
1194		map <  pair<float, float>, map<string, float>  >  xsec;
1195		TFile *datafile = new TFile("../DistributedModelCalculations/StoredShapes.root","RECREATE");
1196		TCanvas *limcan = new TCanvas("limcan","Canvas for calculating limits");
1197		TTree *faketree;
1198		bool dataonly=true;
1199		bool signalonly=false;
1200	<	string mcjzb="JZBforMCinPrepareDataShapes";//this string is not used.
859	<	float jzbpeakerrormc=0;
1200	>
1201		generate_shapes_for_systematic(signalonly,dataonly,datafile,faketree,"",mcjzb,datajzb,noJES,jzbbins,limcan,"",xsec);
1202	<	// don't need these effects for obs & pred, only for signal!
1203	<	//  generate_shapes_for_systematic(signalonly,dataonly,datafile,faketree,"peakUp",newjzbexpression(mcjzb,jzbpeakerrormc),newjzbexpression(datajzb,jzbpeakerrordata),noJES,jzbbins,limcan,"",xsec);
1204	<	//  generate_shapes_for_systematic(signalonly,dataonly,datafile,faketree,"peakDown",newjzbexpression(mcjzb,-jzbpeakerrormc),newjzbexpression(datajzb,-jzbpeakerrordata),noJES,jzbbins,limcan,"",xsec);
1205	<	//  generate_shapes_for_systematic(signalonly,dataonly,datafile,faketree,"JESUp",mcjzb,datajzb,JESup,jzbbins,limcan,"",xsec);
1206	<	//  generate_shapes_for_systematic(signalonly,dataonly,datafile,faketree,"JESDown",mcjzb,datajzb,JESdown,jzbbins,limcan,"",xsec);
1202	>
1203	>	generate_mc_shapes(false,faketree,mcjzb,datajzb,jzbbins,cutWeight,cutnJets,"","",xsec);
1204	>	generate_mc_shapes(false,faketree,mcjzb,datajzb,jzbbins,cutWeight,cutnJetsJESup,"","JESUp",xsec);
1205	>	generate_mc_shapes(false,faketree,mcjzb,datajzb,jzbbins,cutWeight,cutnJetsJESdown,"","JESDown",xsec);
1206	>	generate_mc_shapes(false,faketree,mcjzb,datajzb,jzbbins,TCut("weightEffUp"),cutnJets,"","EffUp",xsec);
1207	>	generate_mc_shapes(false,faketree,mcjzb,datajzb,jzbbins,TCut("weightEffDown"),cutnJets,"","EffDown",xsec);
1208	>	generate_mc_shapes(false,faketree,mcjzb,datajzb,jzbbins,cutWeight,cutnJets,"","StatUp",xsec);
1209	>	generate_mc_shapes(false,faketree,mcjzb,datajzb,jzbbins,cutWeight,cutnJets,"","StatDown",xsec);
1210	>	/*
1211		datafile->Close();
1212	+	*/
1213		}
1214

Diff Legend

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