StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
makeBkgdFilesForRatio.C
1 // ********************************************************
2 // This code makes the necessary root histos for the
3 // background subtraction and systematic uncertainty calculations
4 // Set charge == 1 for postive charges
5 // Set charge == -1 for negative charges
6 // Set charge == 0 for the sum of both
7 //
8 // Set two_or_four == 2 for 2 GeV wide bin histograms
9 // Set two_or_four == 4 for 4 GeV wide bin hisrograms
10 //
11 // Set |eta| boundries from 0 to 1 by default (bins are 0.01 wide in eta)
12 // Set etaLow == 0 for |eta| = 0
13 // Set etaHigh == 100 for |eta| = 1
14 //
15 // This code uses the W MC sample to set the signal in the normalization
16 // region. It is less sensitive to low statistics fluctuations.
17 //
18 // ********************************************************
19 
20 void makeBkgdFilesForRatio(int charge, int two_or_four, int etaLow=0, int etaHigh=100) {
21 
22  // ******************************************************
23  // Read in the data set and all the histograms needed
24  // to do the actual calculation
25  // ******************************************************
26 
27  string dataDir="/star/data01/pwg/stevens4/wAnalysis/xSecPaper/sl11b/data/";
28  string dataName="run9setABCD.wana.hist.root";
29 
30  gStyle->SetOptDate(0);
31  TFile *f1 = new TFile(Form("%s%s",dataDir,dataName));
32 
33  // get the signal w/ and w/o the EEMC in the algo
34  if (charge == 1) {
35  TH1F *signal = (TH1F*)((TH2F*)f1->Get("pos_muclustpTbal_wE_etaBin2"))->ProjectionY("signal_py",etaLow,etaHigh);
36  TH1F *signal_wo_eemc = (TH1F*)((TH2F*)f1->Get("pos_muclustpTbal_noE_etaBin2"))->ProjectionY("signal_no_eemc_py",etaLow,etaHigh);
37  } else if (charge == -1) {
38  TH1F *signal = (TH1F*)((TH2F*)f1->Get("neg_muclustpTbal_wE_etaBin2"))->ProjectionY("signal_py",etaLow,etaHigh);
39  TH1F *signal_wo_eemc = (TH1F*)((TH2F*)f1->Get("neg_muclustpTbal_noE_etaBin2"))->ProjectionY("signal_no_eemc_py",etaLow,etaHigh);
40  } else if (charge == 0) { //not eta dependent yet
41  TH1F *signal = (TH1F*)f1->Get("muclustPtBal");
42  TH1F *signal_wo_eemc = (TH1F*)f1->Get("muclustPtBalnoE");
43  }
44 
45  // ******************************************************
46  // Read in all the MC data sets for background
47  // subtractions and studies
48  // ******************************************************
49 
50  TFile *MC_fs[6];
51  string mcDir="/star/u/stevens4/wAnalysis/efficXsec/outEmb/gainUp2/";
52  string tpcHV="";
53  //string tpcHV="HighTpcHV";
54 
55  MC_fs[0] = new TFile(Form("%sWplus%s.wana.hist.root",mcDir,tpcHV)); // W+ -> e++nu
56  MC_fs[1] = new TFile(Form("%sWminus%s.wana.hist.root",mcDir,tpcHV)); // W- -> e-+nu
57  MC_fs[2] = new TFile(Form("%sWtau%s.wana.hist.root",mcDir,tpcHV)); // W -> tau+nu
58  MC_fs[3] = new TFile(Form("%sZany%s.wana.hist.root",mcDir,tpcHV)); // Z -> any
59  MC_fs[4] = new TFile(Form("%sZe+e-Interf%s.wana.hist.root",mcDir,tpcHV)); // Z/gamma* -> e+ e-
60 
61  //get eta dependent signal and background from MC
62  TH1F *MC_dists_raw[5][3];
63  for (int i=0; i<5; i++) {
64  if (charge == 1) {
65  MC_dists_raw[i][0] = (TH1F*)((TH2F*)MC_fs[i]->Get("pos_muclustpTbal_wE_etaBin2"))->ProjectionY(Form("%d_0_py",i),etaLow,etaHigh);
66  MC_dists_raw[i][1] = (TH1F*)((TH2F*)MC_fs[i]->Get("pos_muclustpTbal_noE_etaBin2"))->ProjectionY(Form("%d_1_py",i),etaLow,etaHigh);
67  MC_dists_raw[i][2] = (TH1F*)((TH2F*)MC_fs[i]->Get("pos_muclustpTbal_back_etaBin2"))->ProjectionY(Form("%d_2_py",i),etaLow,etaHigh);
68  //cout<<i<<" "<<MC_dists_raw[i][0]->Integral()<<endl;
69  } else if (charge == -1) {
70  MC_dists_raw[i][0] = (TH1F*)((TH2F*)MC_fs[i]->Get("neg_muclustpTbal_wE_etaBin2"))->ProjectionY(Form("%d_0_py",i),etaLow,etaHigh);
71  MC_dists_raw[i][1] = (TH1F*)((TH2F*)MC_fs[i]->Get("neg_muclustpTbal_noE_etaBin2"))->ProjectionY(Form("%d_1_py",i),etaLow,etaHigh);
72  MC_dists_raw[i][2] = (TH1F*)((TH2F*)MC_fs[i]->Get("neg_muclustpTbal_back_etaBin2"))->ProjectionY(Form("%d_2_py",i),etaLow,etaHigh);
73  } else if (charge == 0) { //not eta dependent yet
74  MC_dists_raw[i][0] = (TH1F*)MC_fs[i]->Get("muclustPtBal");
75  MC_dists_raw[i][1] = (TH1F*)MC_fs[i]->Get("muclustPtBalnoE");
76  MC_dists_raw[i][2] = (TH1F*)MC_fs[i]->Get("muclustPtBal_bckgrd");
77  }
78  }
79 
80  //pure lumi from pythia to scale MC
81  float lumi[5] = {128.6,385.0,96.2,53.1,531.9};
82 
83  //scale MC to match lumi of the dataset
84  float lumi_fact[6];
85  for (int i=0; i<5; i++) {lumi_fact[i] = 13.18/lumi[i];}
86  for (int i=0; i<5; i++) {
87  for (int j=0; j<3; j++) {
88  MC_dists_raw[i][j]->Scale(lumi_fact[i]);
89  //cout<<MC_dists_raw[i][j]->Integral()<<endl;
90  }
91  }
92 
93  // Repack the MC histograms to mesh with the odd staggered binning that
94  // is being used
95  char str[200];
96  TH1F *MC_dists_repack[5][3];
97  for (int i=0; i<5; i++) {
98  sprintf(str,"mcclustPtBal_%d",i);
99  MC_dists_repack[i][0] = new TH1F(str,str,49,1.,99.);
100  sprintf(str,"mcclustPtBalnoE_%d",i);
101  MC_dists_repack[i][1] = new TH1F(str,str,49,1.,99.);
102  sprintf(str,"mcclustPtBal_bckgrd_%d",i);
103  MC_dists_repack[i][2] = new TH1F(str,str,49,1.,99.);
104  }
105 
106  for (int i=0; i<5; i++) {
107  for (int j=0; j<3; j++) {
108  for (int k=1; k<=49; k++) {
109  MC_dists_repack[i][j]->SetBinContent(k,MC_dists_raw[i][j]->GetBinContent(2*k)+MC_dists_raw[i][j]->GetBinContent(2*k+1));
110  }
111  }
112  }
113 
114  // Make the EEMC background the real one of interest now
115  for (int i=0; i<5; i++) {
116  MC_dists_repack[i][1]->Add(MC_dists_repack[i][0],-1);
117  MC_dists_raw[i][1]->Add(MC_dists_raw[i][0],-1.);
118  }
119 
120  // **********************************************
121  // Do some other stuff
122  // **********************************************
123 
124  TH1F *eemc_bkgd = signal_wo_eemc->Clone();
125  eemc_bkgd->Add(signal,-1.); //this is the '2nd endcap' bkgd
126  eemc_bkgd->Add(MC_dists_raw[4][1],-1.);//remove known Z from '2nd endcap' bkgd
127 
128  TH1F *signal_final = signal->Clone();
129  signal_final->Add(eemc_bkgd,-1);
130 
131  //all histos rebinned by to match binning starting at 25 GeV
132  TH1F *eemc_bkgd2 = new TH1F("eemc_bkgd2","eemc_bkgd2",49,1,99);
133  TH1F *zsig_bkgd2 = new TH1F("zsig_bkgd2","zsig_bkgd2",49,1,99);
134  TH1F *zeemc_bkgd2 = new TH1F("zeemc_bgkd2","zeemc_bkgd2",49,1,99);
135  TH1F *zback_bkgd2 = new TH1F("zback_bkgd2","zback_bkgd2",49,1,99);
136 
137  TH1F *zanysig_bkgd2 = new TH1F("zanysig_bkgd2","zanysig_bkgd2",49,1,99);
138  TH1F *zanyeemc_bkgd2 = new TH1F("zanyeemc_bgkd2","zanyeemc_bkgd2",49,1,99);
139  TH1F *zanyback_bkgd2 = new TH1F("zanyback_bkgd2","zanyback_bkgd2",49,1,99);
140 
141  TH1F *zsig = MC_dists_raw[4][0]->Clone();
142  TH1F *zeemc = MC_dists_raw[4][1]->Clone();
143  TH1F *zback = MC_dists_raw[4][2]->Clone();
144 
145  //subtract off MC Z contribution to signal
146  signal_final->Add(zsig,-1.);
147 
148  // First get the "nominal" QCD background shape (charge summed)
149  TH1F *bkgd_shape1 = (TH1F*)((TH2F*)f1->Get("pos_muclustpTbal_back_etaBin2"))->ProjectionY("back1_py",etaLow,etaHigh);
150  bkgd_shape1->Add((TH1F*)((TH2F*)f1->Get("neg_muclustpTbal_back_etaBin2"))->ProjectionY("back2_py",etaLow,etaHigh));
151 
152  //rebin background shape
153  TH1F *bkgd_shape_nom = new TH1F("bkgd_shape","bkgd_shape",49,1,99);
154  for (int i=1; i<=49; i++) {
155  bkgd_shape_nom->SetBinContent(i,bkgd_shape1->GetBinContent(2*i)+
156  bkgd_shape1->GetBinContent(2*i+1));
157  }
158  TH1F *bkgd_shape_nom2 = bkgd_shape_nom->Clone();
159 
160  TH1F *signal_final2 = new TH1F("signal_final2","signal_final2",49,1,99);
161  signal_final2->SetLineColor(2);
162  signal_final2->SetLineWidth(2.*signal_final2->GetLineWidth());
163  TH1F *signal2 = new TH1F("signal2","signal2",49,1,99);
164  for (int i=1; i<=49; i++) { //repack all distributions
165  signal_final2->SetBinContent(i,signal_final->GetBinContent(2*i)+
166  signal_final->GetBinContent(2*i+1));
167  signal2->SetBinContent(i,signal->GetBinContent(2*i)+
168  signal->GetBinContent(2*i+1));
169  eemc_bkgd2->SetBinContent(i,eemc_bkgd->GetBinContent(2*i)+
170  eemc_bkgd->GetBinContent(2*i+1));
171  zsig_bkgd2->SetBinContent(i,zsig->GetBinContent(2*i)+
172  zsig->GetBinContent(2*i+1));
173  zeemc_bkgd2->SetBinContent(i,zeemc->GetBinContent(2*i)+
174  zeemc->GetBinContent(2*i+1));
175  zback_bkgd2->SetBinContent(i,zback->GetBinContent(2*i)+
176  zback->GetBinContent(2*i+1));
177  zanysig_bkgd2->SetBinContent(i,MC_dists_raw[3][0]->GetBinContent(2*i)+
178  MC_dists_raw[3][0]->GetBinContent(2*i+1));
179  zanyeemc_bkgd2->SetBinContent(i,MC_dists_raw[3][1]->GetBinContent(2*i)+
180  MC_dists_raw[3][1]->GetBinContent(2*i+1));
181  zanyback_bkgd2->SetBinContent(i,MC_dists_raw[3][2]->GetBinContent(2*i)+
182  MC_dists_raw[3][2]->GetBinContent(2*i+1));
183  }
184 
185  TCanvas *can2 = new TCanvas("can2","can2",0,0,600,400);
186  signal2->Draw();
187  signal_final2->Draw("same");
188  //can2->Print("signal_w_eemc.eps");
189  //can2->Print("signal_w_eemc.png");
190 
191 
192  // **********************************************
193  // plot all signals and backgrounds from data and simu
194  // **********************************************
195  eemc_bkgd->SetLineColor(8);
196  eemc_bkgd->SetLineWidth(2.*eemc_bkgd->GetLineWidth());
197  bkgd_shape_nom2->SetLineColor(4);
198  bkgd_shape_nom2->SetLineWidth(2.*bkgd_shape_nom2->GetLineWidth());
199 
200  TCanvas *can8 = new TCanvas("can8","can8",0,0,3000,1800);
201  can8->Divide(5,3);
202  for (int i=0; i<5; i++) {
203  for (int j=0; j<3; j++) {
204  can8->cd(5*(j)+i+1);
205  //cout << 5*(j)+i+1 << endl;
206  gPad->SetGridx(0);
207  gPad->SetGridy(0);
208  if (j == 0) {
209  signal_final2->Draw();
210  } else if (j == 1) {
211  gPad->SetLogy();
212  eemc_bkgd->Draw();
213  } else if (j == 2) {
214  gPad->SetLogy();
215  bkgd_shape_nom2->Draw();
216  }
217  MC_dists_repack[i][j]->Draw("same");
218  }
219  }
220  //can8->Print("phys_bkgds.eps");
221  //can8->Print("phys_bkgds.png");
222 
223  for (int i=0; i<5; i++) {
224  float signal_sum = 0.;
225  float eemc_sum = 0.;
226  float back_sum = 0.;
227  for (int j=13; j<=49; j++) {
228  signal_sum += MC_dists_repack[i][0]->GetBinContent(j);
229  eemc_sum += MC_dists_repack[i][1]->GetBinContent(j);
230  back_sum += MC_dists_repack[i][2]->GetBinContent(j);
231  }
232  }
233 
234 
235  float signal_in_normMC[3];
236  if (charge == 1) {
237  signal_in_normMC[0]=MC_dists_repack[0][0]->Integral(8,8);
238  signal_in_normMC[1]=MC_dists_repack[0][0]->Integral(9,9);
239  signal_in_normMC[2]=MC_dists_repack[0][0]->Integral(10,10);
240  }
241  if (charge == -1) {
242  signal_in_normMC[0]=MC_dists_repack[1][0]->Integral(8,8);
243  signal_in_normMC[1]=MC_dists_repack[1][0]->Integral(9,9);
244  signal_in_normMC[2]=MC_dists_repack[1][0]->Integral(10,10);
245  }
246 
247 
248  //Tau needs to be scale due to polarized tau decay (note from Carl)
249  float taufrac=1.5;
250  // subtract off the tau background
251  signal_final2->Add(MC_dists_repack[2][0],-1.*taufrac);
252 
253  // ******************************************************
254  // Do the iterative normalization of the W signal to
255  // to the background for the nominal shape
256  // ******************************************************
257 
258  // Now lets try a new background function
259  // Fit the peak with a line from 23<ET<39
260  TF1 *func1 = new TF1("func1","[0]+[1]*x",23,39);
261  func1->SetParameter(0,0.);
262  func1->SetParameter(1,0.);
263 
264  TCanvas *can4 = new TCanvas("can4","can4",0,0,600,400);
265  signal_final2->Draw();
266 
267  // Calculate the background for the nominal signal
268  // by doing the 20 iterations over the signal peak
269  float signal_in_norm[50];
270  TH1F *bkgd_shape_unnorm[20];
271  TH1F *signal_for_new[20];
272  for (int i=0; i<20; i++) {
273  bkgd_shape_unnorm[i] = (TH1F*)bkgd_shape_nom->Clone();
274  signal_for_new[i] = (TH1F*)signal_final2->Clone();
275 
276  // subtract off the remaining Z signal from the background
277  bkgd_shape_unnorm[i]->Add(zback_bkgd2,-1.);
278 
279  // calculate the W signal in the normalization bins
280  signal_in_norm[8] = func1->Integral(15,17);
281  signal_in_norm[9] = func1->Integral(17,19);
282  signal_in_norm[10] = func1->Integral(19,21);
283  //cout << "Integrals = " << signal_in_norm[8] << " " << signal_in_norm[9] << " " << signal_in_norm[10] << endl;
284  for (int j=8; j<=10; j++) {
285  if (signal_in_norm[j] < 0) {signal_in_norm[j] = 0.;}
286  }
287 
288  //use MC for signal in normalization window
289  signal_in_norm[8]=signal_in_normMC[0];
290  signal_in_norm[9]=signal_in_normMC[1];
291  signal_in_norm[10]=signal_in_normMC[2];
292 
293  // calculate the normalization factor using the
294  // signal subtraction normalization bins
295  float normt = 0, normb = 0.;
296  for (int k=8; k<10; k++) {
297  if (bkgd_shape_unnorm[i]->GetBinContent(k) > 0) {
298  normt += signal_final2->GetBinContent(k)-signal_in_norm[k];
299  normb += bkgd_shape_unnorm[i]->GetBinContent(k);
300  }
301  }
302  if (normb > 0 && normt > 0) {
303  float norm = normt/normb;
304  bkgd_shape_unnorm[i]->Scale(norm);
305  bkgd_shape_unnorm[i]->Draw("same");
306  //cout << "norm " << i << " = " << norm << endl;
307  }
308  else if(normb > 0 && normt < 0){
309  float norm = 0;
310  bkgd_shape_unnorm[i]->Scale(norm);
311  bkgd_shape_unnorm[i]->Draw("same");
312  }
313 
314  // With the new signal estimate, calculate the normalization
315  // factor
316  for (int j=1; j<=49; j++) {
317  if (bkgd_shape_unnorm[i]->GetBinContent(j) < 0) {bkgd_shape_unnorm[i]->SetBinContent(j,0.);}
318  }
319  signal_for_new[i]->Add(bkgd_shape_unnorm[i],-1.);
320  signal_for_new[i]->Fit(func1,"RQ");
321  }
322  //can4->Print("new_norm.eps");
323  //can4->Print("new_norm.png");
324 
325  TH1F *signal_in_norm_region = new TH1F("signal_in_norm_region","signal_in_norm_region",49,1.,99.);
326 
327  signal_in_norm_region->SetBinContent(8,signal_in_norm[8]);
328  signal_in_norm_region->SetBinContent(9,signal_in_norm[9]);
329 
330  TH1F *new_bkgd = new TH1F("new_bkgd","new_bkgd",49,1.,99.);
331  new_bkgd = (TH1F*)bkgd_shape_unnorm[19]->Clone();
332  new_bkgd->SetName("new_bkgd");
333 
334  TCanvas *can5 = new TCanvas("can5","can5",0,0,600,400);
335  signal_final2->Draw();
336  new_bkgd->Draw("same");
337  signal_in_norm_region->Draw("same");
338  //can5->Print("signal_new.eps");
339  //can5->Print("signal_new.png");
340 
341  // ******************************************************
342  // Calculate all the 60 shapes for the background
343  // systematic shape study
344  // ******************************************************
345 
346  // get the various background shapes made by
347  // using different parameters for sPtBal cut
348  TH1F *bkgd_hists_from_file[21];
349  TH1F *bkgd_hists_from_file2[21];
350  char str[200];
351  for (int i=0; i<=20; i++) {
352 
353  //use projections from eta dependent histos
354  sprintf(str,"pos_failsPtBal_sPtBal_bin_%d_py",i);
355  bkgd_hists_from_file[i] = (TH1F*)((TH2F*)f1->Get(Form("pos_failsPtBal_sPtBal_bin_%d_etaBin2",i)))->ProjectionY(str,etaLow,etaHigh);
356  sprintf(str,"neg_failsPtBal_sPtBal_bin_%d_py",i);
357  bkgd_hists_from_file2[i] = (TH1F*)((TH2F*)f1->Get(Form("neg_failsPtBal_sPtBal_bin_%d_etaBin2",i)))->ProjectionY(str,etaLow,etaHigh);
358  bkgd_hists_from_file[i]->Add(bkgd_hists_from_file2[i]);
359 
360  }
361 
362  // Now do the rebinning
363  TH1F *bkgd_hists1[21];
364  TH1F *bkgd_hists2[21];
365  TH1F *bkgd_hists3[21];
366  for (int i=0; i<=20; i++) {
367  sprintf(str,"bkgd_hist1_%d",i);
368  bkgd_hists1[i] = new TH1F(str,str,49,1,99);
369  sprintf(str,"bkgd_hist2_%d",i);
370  bkgd_hists2[i] = new TH1F(str,str,49,1,99);
371  sprintf(str,"bkgd_hist3_%d",i);
372  bkgd_hists3[i] = new TH1F(str,str,49,1,99);
373  for (int k=1; k<=49; k++) {
374  bkgd_hists1[i]->SetBinContent(k,bkgd_hists_from_file[i]->GetBinContent(2*k)+bkgd_hists_from_file[i]->GetBinContent(2*k+1));
375  bkgd_hists2[i]->SetBinContent(k,bkgd_hists_from_file[i]->GetBinContent(2*k)+bkgd_hists_from_file[i]->GetBinContent(2*k+1));
376  bkgd_hists3[i]->SetBinContent(k,bkgd_hists_from_file[i]->GetBinContent(2*k)+bkgd_hists_from_file[i]->GetBinContent(2*k+1));
377  }
378  }
379 
380  // initiaize the iterative fit functions
381  TF1 *func2 = new TF1("func2","[0]+[1]*x",23,39);
382  func2->SetParameter(0,0.);
383  func2->SetParameter(1,0.);
384  TF1 *func3 = new TF1("func3","[0]+[1]*x",23,39);
385  func3->SetParameter(0,0.);
386  func3->SetParameter(1,0.);
387 
388  // Now loop over the the 60 possibilities (20 loops and 3 normalization regions)
389  float final_sig_in_norm[21][3];
390  float final_chisquare[21];
391  float signal_in_norm1[50];
392  float signal_in_norm2[50];
393  float signal_in_norm3[50];
394  TH1F *new_bkgd_hists1[21];
395  TH1F *new_bkgd_hists2[21];
396  TH1F *new_bkgd_hists3[21];
397  TH1F *bkgd_shape_unnorm1[20];
398  TH1F *bkgd_shape_unnorm2[20];
399  TH1F *bkgd_shape_unnorm3[20];
400  TH1F *signal_for_new1[20];
401  TH1F *signal_for_new2[20];
402  TH1F *signal_for_new3[20];
403  for (int i=0; i<=20; i++) { //loop over possible sPtBalance cuts
404 
405  func1->SetParameter(0,0.);
406  func1->SetParameter(1,0.);
407  func2->SetParameter(0,0.);
408  func2->SetParameter(1,0.);
409  func3->SetParameter(0,0.);
410  func3->SetParameter(1,0.);
411 
412  for (int l=0; l<20; l++) { //loop over 20 iterations
413  bkgd_shape_unnorm1[l] = (TH1F*)bkgd_hists1[i]->Clone();
414  bkgd_shape_unnorm2[l] = (TH1F*)bkgd_hists2[i]->Clone();
415  bkgd_shape_unnorm3[l] = (TH1F*)bkgd_hists3[i]->Clone();
416  signal_for_new1[l] = (TH1F*)signal_final2->Clone();
417  signal_for_new2[l] = (TH1F*)signal_final2->Clone();
418  signal_for_new3[l] = (TH1F*)signal_final2->Clone();
419 
420  // subtract off the Z contamination
421  bkgd_shape_unnorm1[l]->Add(zback_bkgd2,-1.);
422  bkgd_shape_unnorm2[l]->Add(zback_bkgd2,-1.);
423  bkgd_shape_unnorm3[l]->Add(zback_bkgd2,-1.);
424 
425  // calculate the signal in the normalization bins
426  signal_in_norm1[8] = func1->Integral(15,17);
427  signal_in_norm1[9] = func1->Integral(17,19);
428  signal_in_norm1[10] = func1->Integral(19,21);
429  signal_in_norm2[8] = func2->Integral(15,17);
430  signal_in_norm2[9] = func2->Integral(17,19);
431  signal_in_norm2[10] = func2->Integral(19,21);
432  signal_in_norm3[8] = func3->Integral(15,17);
433  signal_in_norm3[9] = func3->Integral(17,19);
434  signal_in_norm3[10] = func3->Integral(19,21);
435 
436  for (int m=8; m<=10; m++) {
437  if (signal_in_norm1[m] < 0) {signal_in_norm1[m] = 0.;}
438  if (signal_in_norm2[m] < 0) {signal_in_norm2[m] = 0.;}
439  if (signal_in_norm3[m] < 0) {signal_in_norm3[m] = 0.;}
440  }
441 
442  //use MC for signal in normalization window
443  signal_in_norm1[8]=signal_in_normMC[0];
444  signal_in_norm1[9]=signal_in_normMC[1];
445  signal_in_norm1[10]=signal_in_normMC[2];
446  signal_in_norm2[8]=signal_in_normMC[0];
447  signal_in_norm2[9]=signal_in_normMC[1];
448  signal_in_norm2[10]=signal_in_normMC[2];
449  signal_in_norm3[8]=signal_in_normMC[0];
450  signal_in_norm3[9]=signal_in_normMC[1];
451  signal_in_norm3[10]=signal_in_normMC[2];
452 
453  // calculate the normalization factor for 1 bin
454  float normt = 0, normb = 0.;
455  for (int k=8; k<=8; k++) {
456  if (bkgd_shape_unnorm1[l]->GetBinContent(k) > 0) {
457  normt += signal_final2->GetBinContent(k)-signal_in_norm1[k];
458  normb += bkgd_shape_unnorm1[l]->GetBinContent(k);
459  }
460  }
461  if (normb > 0 && normt > 0) {
462  float norm = normt/normb;
463  bkgd_shape_unnorm1[l]->Scale(norm);
464  }
465  else if(normb > 0 && normt < 0){
466  float norm = 0;
467  bkgd_shape_unnorm1[l]->Scale(norm);
468  }
469 
470  for (int m=1; m<=49; m++) {
471  if (bkgd_shape_unnorm1[l]->GetBinContent(m) < 0) {bkgd_shape_unnorm1[l]->SetBinContent(m,0.);}
472  }
473  signal_for_new1[l]->Add(bkgd_shape_unnorm1[l],-1.);
474  signal_for_new1[l]->Fit(func1,"RQ");
475 
476  // calculate the normalization factor for 2 bins
477  normt = 0.; normb = 0.;
478  for (int k=8; k<=9; k++) {
479  if (bkgd_shape_unnorm2[l]->GetBinContent(k) > 0) {
480  normt += signal_final2->GetBinContent(k)-signal_in_norm2[k];
481  normb += bkgd_shape_unnorm2[l]->GetBinContent(k);
482  }
483  }
484  if (normb > 0 && normt > 0) {
485  float norm = normt/normb;
486  bkgd_shape_unnorm2[l]->Scale(norm);
487  //cout << "norm " << i << " " << l << " = " << norm << endl;
488  }
489  else if(normb > 0 && normt < 0){
490  float norm = 0;
491  bkgd_shape_unnorm2[l]->Scale(norm);
492  }
493 
494  for (int m=1; m<=49; m++) {
495  if (bkgd_shape_unnorm2[l]->GetBinContent(m) < 0) {bkgd_shape_unnorm2[l]->SetBinContent(m,0.);}
496  }
497  signal_for_new2[l]->Add(bkgd_shape_unnorm2[l],-1.);
498  signal_for_new2[l]->Fit(func2,"RQ");
499 
500  // calculate the normalization factor for 3 bins
501  normt = 0.; normb = 0.;
502  for (int k=8; k<=10; k++) {
503  if (bkgd_shape_unnorm3[l]->GetBinContent(k) > 0) {
504  normt += signal_final2->GetBinContent(k)-signal_in_norm3[k];
505  normb += bkgd_shape_unnorm3[l]->GetBinContent(k);
506  }
507  }
508  if (normb > 0 && normt > 0) {
509  float norm = normt/normb;
510  bkgd_shape_unnorm3[l]->Scale(norm);
511  }
512  else if(normb > 0 && normt < 0){
513  float norm = 0;
514  bkgd_shape_unnorm3[l]->Scale(norm);
515  }
516 
517  for (int m=1; m<=49; m++) {
518  if (bkgd_shape_unnorm3[l]->GetBinContent(m) < 0) {bkgd_shape_unnorm3[l]->SetBinContent(m,0.);}
519  }
520  signal_for_new3[l]->Add(bkgd_shape_unnorm3[l],-1.);
521  signal_for_new3[l]->Fit(func3,"RQ");
522  } // end of for loop over l
523 
524  // Save the last iteration as the background histogram
525  new_bkgd_hists1[i] = (TH1F*)bkgd_shape_unnorm1[19]->Clone();
526  new_bkgd_hists2[i] = (TH1F*)bkgd_shape_unnorm2[19]->Clone();
527  new_bkgd_hists3[i] = (TH1F*)bkgd_shape_unnorm3[19]->Clone();
528 
529  }
530 
531  // plot all the 60 histograms
532  gStyle->SetTitleBorderSize(0);
533  TCanvas *can6 = new TCanvas("can6","can6",0,0,600,400);
534  signal_final2->SetStats(kFALSE);
535  if (charge == 1) {
536  signal_final2->SetTitle("W+ Background Shapes");
537  } else if (charge == -1) {
538  signal_final2->SetTitle("W- Background Shapes");
539  }
540  signal_final2->GetXaxis()->SetRangeUser(0.,70.);
541  signal_final2->GetXaxis()->SetTitle("2x2 Cluster E_{T} (GeV)");
542  signal_final2->Draw();
543  for (int i=0; i<=20; i++) {
544  new_bkgd_hists1[i]->Draw("same");
545  new_bkgd_hists2[i]->Draw("same");
546  new_bkgd_hists3[i]->Draw("same");
547  }
548  new_bkgd->SetLineColor(4);
549  //new_bkgd->SetLineWidth(4.*new_bkgd->GetLineWidth());
550  new_bkgd->Draw("same");
551  if (charge == 1) {
552  //can6->Print("Wplus_bkgd_shapes.eps");
553  //can6->Print("Wplus_bkgd_shapes.png");
554  } else if (charge == -1) {
555  //can6->Print("Wminus_bkgd_shapes.eps");
556  //can6->Print("Wminus_bkgd_shapes.png");
557  }
558 
559  TH1F *chi2s = new TH1F("chi2s","chi2s",50,0.,10.);
560  for (int i=0; i<=20; i++) {
561  chi2s->Fill(final_chisquare[i]);
562  }
563 
564  TCanvas *can7 = new TCanvas("can7","can7",0,0,600,400);
565  chi2s->Draw();
566  //can7->Print("chi2s.eps");
567  //can7->Print("chi2s.png");
568 
569  // ************************************************
570  // Now calculate all background numbers and their
571  // systematic uncertainties (and spit them out to histograms)
572  // ************************************************
573 
574  // First get the simple numbers (backgrounds and their statistical uncertainties)
575  TH1F *tauhist = MC_dists_repack[2][0]->Clone();
576  tauhist->Scale(taufrac);
577  float tau_norm = lumi_fact[2];
578  float Z_norm = lumi_fact[4];
579  float bkgd_sum = 0.;
580  float signal_sum = 0.;
581  float raw_sum = 0.;
582  float QCD_sum = 0., tau_sum = 0., eemc_sum = 0.;
583  float QCD_raw_sum = 0.;
584  float Zany_bkgd_sum = 0.;
585  float Zany_eemc_sum = 0.;
586  float zsig_sum = 0., zeemc_sum = 0.,zback_sum = 0.;
587  float zanysig_sum = 0., zanyeemc_sum = 0.,zanyback_sum = 0.;
588  for (int i=13; i<=49; i++) { //get counts in ET>25 for diff contrib.
589  bkgd_sum += new_bkgd->GetBinContent(i);
590  bkgd_sum += tauhist->GetBinContent(i);
591  bkgd_sum += eemc_bkgd2->GetBinContent(i);
592  QCD_raw_sum += bkgd_shape_nom->GetBinContent(i);
593  QCD_sum += new_bkgd->GetBinContent(i);
594  tau_sum += tauhist->GetBinContent(i);
595  eemc_sum += eemc_bkgd2->GetBinContent(i);
596  signal_sum += signal_final2->GetBinContent(i);
597  raw_sum += signal2->GetBinContent(i);
598  zsig_sum += zsig_bkgd2->GetBinContent(i);
599  zeemc_sum += zeemc_bkgd2->GetBinContent(i);
600  zback_sum += zback_bkgd2->GetBinContent(i);
601  zanysig_sum += zanysig_bkgd2->GetBinContent(i);
602  zanyeemc_sum += zanyeemc_bkgd2->GetBinContent(i);
603  zanyback_sum += zanyback_bkgd2->GetBinContent(i);
604  }
605  cout << "The total background for ET>25 is " << bkgd_sum+zsig_sum << endl;
606  cout << "QCD = " << QCD_sum << ", tau = " << tau_sum << ", eemc = " << eemc_sum << ", and Z = " << zsig_sum << endl;
607  cout << "Raw = " << raw_sum << endl;
608  cout << "W Signal (w/o tau) = " << signal_sum-QCD_sum << endl;
609  cout << "Z in sig = " << zsig_sum << endl;
610  cout << "Z in eemc = " << zeemc_sum << endl;
611  cout << "Z in back = " << zback_sum << endl;
612  cout << "Zany in sig = " << zanysig_sum << endl;
613  cout << "Zany in eemc = " << zanyeemc_sum << endl;
614  cout << "Zany in back = " << zanyback_sum << endl;
615  cout << "QCD raw in back = " << QCD_raw_sum << endl;
616  cout << "The QCD stat unc. is " << sqrt(norm*QCD_sum) << endl;
617  float tau_stat = tau_norm*taufrac*sqrt(tau_sum/(tau_norm*taufrac));
618  float tau_syst = 0.13*tau_sum;
619  cout << "The tau stat unc. is " << tau_stat << " syst is " << tau_syst << " and total is " << sqrt(tau_stat*tau_stat + tau_syst*tau_syst) << endl;
620  float eemc_stat = sqrt(eemc_sum);
621  float eemc_syst = 0.13*zeemc_sum;
622  cout << "The eemc stat unc. is " << eemc_stat << " syst is " << eemc_syst << " and total is " << sqrt(eemc_stat*eemc_stat + eemc_syst*eemc_syst) << endl;
623  float Z_stat = Z_norm*sqrt(zsig_sum/Z_norm);
624  float Z_syst = 0.13*zsig_sum;
625  cout << "The Z stat unc. is " << Z_stat << " syst is " << Z_syst << " and total is " << sqrt(Z_stat*Z_stat + Z_syst*Z_syst) << endl;
626 
627  //for use with asymmetry analysis
628  //cout << "f_tau = " << tau_sum/raw_sum << endl;
629  //cout << "f_QCD = " << QCD_sum/raw_sum << endl;
630  //cout << "f_EEMC = " << eemc_sum/raw_sum << endl;
631  //cout << "f_Z = " << zsig_sum/raw_sum << endl;
632 
633  // Set up some histograms to hold all the errors that are calculated
634  TH1F *raw_stat_err2 = new TH1F("raw_stat_err2","raw_stat_err2",49,1.,99.);
635  TH1F *QCD_stat_err2 = new TH1F("QCD_stat_err2","QCD_stat_err2",49,1.,99.);
636  TH1F *eemc_stat_err2 = new TH1F("eemc_stat_err2","eemc_stat_err2",49,1.,99.);
637  TH1F *tau_stat_err2 = new TH1F("tau_stat_err2","tau_stat_err2",49,1.,99.);
638  TH1F *QCD_syst_high_err = new TH1F("QCD_syst_high_err","QCD_syst_high_err",49,1.,99.);
639  TH1F *QCD_syst_low_err = new TH1F("QCD_syst_low_err","QCD_syst_low_err",49,1.,99.);
640  TH1F *zsig_stat_err2 = new TH1F("zsig_stat_err2","zsig_stat_err2",49,1.,99.);
641  TH1F *zback_stat_err2 = new TH1F("zback_stat_err2","zback_stat_err2",49,1.,99.);
642  TH1F *zeemc_stat_err2 = new TH1F("zeemc_stat_err2","zeemc_stat_err2",49,1.,99.);
643 
644 
645  for (int i=1; i<=49; i++) {
646  raw_stat_err2->SetBinContent(i,signal2->GetBinContent(i));
647  QCD_stat_err2->SetBinContent(i,fabs(norm*new_bkgd->GetBinContent(i)));
648  eemc_stat_err2->SetBinContent(i,max(eemc_bkgd2->GetBinContent(i),0)); //remove negative error
649  tau_stat_err2->SetBinContent(i,tau_norm*taufrac*tauhist->GetBinContent(i));
650  zsig_stat_err2->SetBinContent(i,Z_norm*zsig_bkgd2->GetBinContent(i));
651  zback_stat_err2->SetBinContent(i,norm*norm*Z_norm*zback_bkgd2->GetBinContent(i));
652  zeemc_stat_err2->SetBinContent(i,Z_norm*zeemc_bkgd2->GetBinContent(i));
653 
654  //cout << "Error " << i << " " << raw_stat_err2->GetBinCenter(i) << " " << raw_stat_err2->GetBinContent(i) << " " << QCD_stat_err2->GetBinContent(i) << " " << eemc_stat_err2->GetBinContent(i) << " " << tau_stat_err2->GetBinContent(i) << endl;
655  }
656 
657  // Now go through all the 60 background shapes and find the
658  // high and low in each ET bin to give the maximum extent uncertainty
659  // for each ET bin (and also sum the low and high to give the
660  // overall number used in the xsec analysis)
661 
662  TH1F *low_bkgd = new_bkgd->Clone();
663  low_bkgd->SetName("low_bkgd");
664  TH1F *high_bkgd = new_bkgd->Clone();
665  high_bkgd->SetName("high_bkgd");
666 
667  float low_sum = 0.;
668  float high_sum = 0.;
669  for (int i=7; i<=49; i++) {
670  float high = 0.;
671  float low = 10000.;
672  for (int j=0; j<=20; j++) {
673  if (new_bkgd_hists1[j]->GetBinContent(i) < low) {
674  if (new_bkgd_hists1[j]->GetBinContent(i) >= 0) {
675  low = new_bkgd_hists1[j]->GetBinContent(i);
676  }
677  }
678  if (new_bkgd_hists1[j]->GetBinContent(i) > high) {
679  high = new_bkgd_hists1[j]->GetBinContent(i);
680  }
681 
682  if (new_bkgd_hists2[j]->GetBinContent(i) < low) {
683  if (new_bkgd_hists2[j]->GetBinContent(i) >= 0) {
684  low = new_bkgd_hists2[j]->GetBinContent(i);
685  }
686  }
687  if (new_bkgd_hists2[j]->GetBinContent(i) > high) {
688  high = new_bkgd_hists2[j]->GetBinContent(i);
689  }
690 
691  if (new_bkgd_hists3[j]->GetBinContent(i) < low) {
692  if (new_bkgd_hists3[j]->GetBinContent(i) >= 0) {
693  low = new_bkgd_hists3[j]->GetBinContent(i);
694  }
695  }
696  if (new_bkgd_hists3[j]->GetBinContent(i) > high) {
697  high = new_bkgd_hists3[j]->GetBinContent(i);
698  }
699  //cout << i << " low = " << low << " high = " << high << endl;
700  //} // end of k-loop
701  } // end of j-loop
702 
703  // calculate the sum
704  low_bkgd->SetBinContent(i,0.);
705  if ((low != 10000) && (new_bkgd->GetBinContent(i)-low > 0)) {
706  if ((i >= 13)&&(i<=49)) {low_sum += low;}
707  low_bkgd->SetBinContent(i,low);
708  }
709  if ((i >= 13)&&(i<=49)) {high_sum += high;}
710  high_bkgd->SetBinContent(i,high);
711  //cout << i << " " << low_sum << " " << high_sum << endl;
712  //cout << QCD_syst_low_err->GetBinCenter(i) << " low = " << low << " high = " << high << " nom = " << new_bkgd->GetBinContent(i) << endl;
713  // set the bin-by-bin error too
714  if ((low != 10000) && (new_bkgd->GetBinContent(i)-low > 0)) {
715  QCD_syst_low_err->SetBinContent(i,new_bkgd->GetBinContent(i)-low);
716  } else {
717  QCD_syst_low_err->SetBinContent(i,0.);
718  }
719  QCD_syst_high_err->SetBinContent(i,high-new_bkgd->GetBinContent(i));
720 
721  } // end of i=loop
722 
723  cout << "QCD shape sys. unc. calc************" << endl;
724  cout << "The QCD low sum = " << low_sum << endl;
725  cout << "The QCD high sum = " << high_sum << endl;
726 
727  cout << "The QCD low error = " << QCD_sum-low_sum << endl;
728  cout << "The QCD high error = " << high_sum-QCD_sum << endl;
729 
730  //totals for ET > 25
731  float tot_stat = sqrt(tau_stat*tau_stat+eemc_stat*eemc_stat+Z_stat*Z_stat+norm*QCD_sum);
732  cout << "total stat unc. is " << tot_stat << endl;
733  float tot_syst_low = sqrt(tau_syst*tau_syst+eemc_syst*eemc_syst+Z_syst*Z_syst+(QCD_sum-low_sum)*(QCD_sum-low_sum));
734  float tot_syst_high = sqrt(tau_syst*tau_syst+eemc_syst*eemc_syst+Z_syst*Z_syst+(QCD_sum-high_sum)*(QCD_sum-high_sum));
735  cout << "total syst unc. is low: " << tot_syst_low << " and high: " << tot_syst_high << endl;
736  cout << "total unc is low: " << sqrt(tot_syst_low*tot_syst_low+tot_stat*tot_stat) << " high: " << sqrt(tot_syst_high*tot_syst_high+tot_stat*tot_stat) <<endl;
737 
738  //final signal histo
739  TH1F *signal_final3 = new TH1F("signal_final3","signal_final3",49,1.,99.);
740  for (int i=1; i<=49; i++) {
741  signal_final3->SetBinContent(i,signal_final2->GetBinContent(i));
742  }
743  signal_final3->Add(new_bkgd,-1.);
744 
745 
746  //need to get 4 GeV wide bins starting at ET=25 to correct yields
747 #if 1
748  cout<<"ET min, ET max, N_obs, N_obs-Nbkgd, QCD stat, EEMC stat, Tau stat, Zee stat, QCD syst+, QCD syst-, total stat, MC bkgd tot"<<endl;
749  for (int i=7; i<=49; i=i+2) {
750  cout.setf(ios::fixed);
751  cout.precision(2);
752  float Qcd_stat_err2=QCD_stat_err2->GetBinContent(i)+QCD_stat_err2->GetBinContent(i+1);
753  float Eemc_stat_err2=eemc_stat_err2->GetBinContent(i)+eemc_stat_err2->GetBinContent(i+1);
754  float Tau_stat_err2=tau_stat_err2->GetBinContent(i)+tau_stat_err2->GetBinContent(i+1);
755  float Z_stat_err2=zsig_stat_err2->GetBinContent(i)+zsig_stat_err2->GetBinContent(i+1);
756  cout << signal2->GetBinCenter(i)-1 << "," << signal2->GetBinCenter(i)+3 << "," << signal2->GetBinContent(i)+signal2->GetBinContent(i+1) << "," << signal_final3->GetBinContent(i)+signal_final3->GetBinContent(i+1) << "," << sqrt(Qcd_stat_err2) << "," << sqrt(Eemc_stat_err2) << "," << sqrt(Tau_stat_err2) << "," << sqrt(Z_stat_err2) << "," << QCD_syst_high_err->GetBinContent(i)+QCD_syst_high_err->GetBinContent(i+1) << "," << QCD_syst_low_err->GetBinContent(i)+QCD_syst_low_err->GetBinContent(i+1) ;
757  //total statistical uncertainty
758  cout<< "," <<sqrt(Tau_stat_err2+Z_stat_err2+Qcd_stat_err2+Eemc_stat_err2);
759  //total MC syst uncertainty
760  cout<< "," << zsig_bkgd2->GetBinContent(i)+zsig_bkgd2->GetBinContent(i+1) + tauhist->GetBinContent(i)+tauhist->GetBinContent(i+1)<<endl;
761 
762  }
763 #endif
764 
765  // ******************************************************
766  // Write out the histograms of interest to a
767  // file
768  // ******************************************************
769 
770  if (charge == 1) {
771  TFile *f2 = new TFile("bkgd_histos_pos_final.root","recreate");
772  } else if (charge == -1) {
773  TFile *f2 = new TFile("bkgd_histos_neg_final.root","recreate");
774  } else if (charge == 0) {
775  TFile *f2 = new TFile("bkgd_histos_sum.root","recreate");
776  }
777 
778  f2->cd();
779 
780  if (two_or_four == 2) {
781 
782  TH1F *signal_final3 = new TH1F("signal_final3","signal_final3",49,1.,99.);
783  for (int i=1; i<=49; i++) {
784  signal_final3->SetBinContent(i,signal_final2->GetBinContent(i));
785  }
786  signal_final3->Add(new_bkgd,-1.);
787 
788  tauhist->Write();
789  zsig_bkgd2->Write();
790  signal2->Write();
791  signal_final2->Write();
792  signal_final3->Write();
793  eemc_bkgd2->Write();
794  new_bkgd->Write();
795  low_bkgd->Write();
796  high_bkgd->Write();
797  signal_in_norm_region->Write();
798 
799  raw_stat_err2->Write();
800  QCD_stat_err2->Write();
801  eemc_stat_err2->Write();
802  tau_stat_err2->Write();
803  QCD_syst_high_err->Write();
804  QCD_syst_low_err->Write();
805  zsig_stat_err2->Write();
806  zback_stat_err2->Write();
807  zeemc_stat_err2->Write();
808 
809 #if 0
810  for (int i=1; i<=49; i++) {
811  cout.setf(ios::fixed);
812  cout.precision(2);
813  cout << " " << signal2->GetBinCenter(i) << " & " << signal2->GetBinContent(i) << " & " << signal_final3->GetBinContent(i) << " & " << QCD_stat_err2->GetBinContent(i) << " & " << eemc_stat_err2->GetBinContent(i) << " & " << tau_stat_err2->GetBinContent(i) << " & " << zsig_stat_err2->GetBinContent(i) << " & " << zback_stat_err2->GetBinContent(i) << " & " << QCD_syst_high_err->GetBinContent(i) << " & " << QCD_syst_low_err->GetBinContent(i) << " \\\\" << endl;
814 
815  }
816 #endif
817 
818  f2->Close();
819 
820  } else if (two_or_four == 4) {
821 
822  TH1F *tauhist_repack = new TH1F("tauhist_r","tauhist_r",24,3.,99.);
823  TH1F *zsig_bkgd2_repack = new TH1F("zsig_bkgd2_r","zsig_bkgd2_r",24,3,99.);
824  TH1F *signal2_repack = new TH1F("signal2_r","signal2_r",24,3.,99.);
825  TH1F *signal_final2_repack = new TH1F("signal_final2_r","signal_final2_r",24,3.,99.);
826  TH1F *signal_final3_repack = new TH1F("signal_final3_r","signal_final3_r",24,3.,99.);
827  TH1F *eemc_bkgd2_repack = new TH1F("eemc_bkgd2_r","eemc_bkgd2_r",24,3.,99.);
828  TH1F *new_bkgd_repack = new TH1F("new_bkgd_r","new_bkgd_r",24,3.,99.);
829  TH1F *low_bkgd_repack = new TH1F("low_bkgd_r","low_bkgd_r",24,3.,99.);
830  TH1F *high_bkgd_repack = new TH1F("high_bkgd_r","high_bkgd_r",24,3.,99.);
831  TH1F *signal_in_norm_region_repack = new TH1F("signal_in_norm_region_r","signal_in_norm_region_r",24,3.,99.);
832 
833  TH1F *raw_stat_err2_repack = new TH1F("raw_stat_err2_r","raw_stat_err2_r",24,3.,99.);
834  TH1F *QCD_stat_err2_repack = new TH1F("QCD_stat_err2_r","QCD_stat_err2_r",24,3.,99.);
835  TH1F *eemc_stat_err2_repack = new TH1F("eemc_stat_err2_r","eemc_stat_err2_r",24,3.,99.);
836  TH1F *tau_stat_err2_repack = new TH1F("tau_stat_err2_r","tau_stat_err2_r",24,3.,99.);
837  TH1F *QCD_syst_high_err_repack = new TH1F("QCD_syst_high_err_r","QCD_syst_high_err_r",24,3.,99.);
838  TH1F *QCD_syst_low_err_repack = new TH1F("QCD_syst_low_err_r","QCD_syst_low_err_r",24,3.,99.);
839  TH1F *zsig_stat_err2_repack = new TH1F("zsig_stat_err2_r","zsig_stat_err2_r",24,3.,99.);
840  TH1F *zback_stat_err2_repack = new TH1F("zback_stat_err2_r","zback_stat_err2_r",24,3.,99.);
841  TH1F *zeemc_stat_err2_repack = new TH1F("zeemc_stat_err2_r","zeemc_stat_err2_r",24,3.,99.);
842 
843  for (int i=1; i<=24; i++) {
844  tauhist_repack->SetBinContent(i,tauhist->GetBinContent(2*i)+
845  tauhist->GetBinContent(2*i+1));
846  zsig_bkgd2_repack->SetBinContent(i,zsig_bkgd2->GetBinContent(2*i)+
847  zsig_bkgd2->GetBinContent(2*i+1));;
848  signal2_repack->SetBinContent(i,signal2->GetBinContent(2*i)+
849  signal2->GetBinContent(2*i+1));
850  signal_final2_repack->SetBinContent(i,signal_final2->GetBinContent(2*i)+
851  signal_final2->GetBinContent(2*i+1));
852  signal_final3_repack->SetBinContent(i,signal_final3->GetBinContent(2*i)+
853  signal_final3->GetBinContent(2*i+1));
854  eemc_bkgd2_repack->SetBinContent(i,eemc_bkgd2->GetBinContent(2*i)+
855  eemc_bkgd2->GetBinContent(2*i+1));
856  new_bkgd_repack->SetBinContent(i,new_bkgd->GetBinContent(2*i)+
857  new_bkgd->GetBinContent(2*i+1));
858  low_bkgd_repack->SetBinContent(i,low_bkgd->GetBinContent(2*i)+
859  low_bkgd->GetBinContent(2*i+1));
860  high_bkgd_repack->SetBinContent(i,high_bkgd->GetBinContent(2*i)+
861  high_bkgd->GetBinContent(2*i+1));
862  signal_in_norm_region_repack->SetBinContent(i,signal_in_norm_region->GetBinContent(2*i)+signal_in_norm_region->GetBinContent(2*i+1));
863 
864  raw_stat_err2_repack->SetBinContent(i,raw_stat_err2->GetBinContent(2*i)+
865  raw_stat_err2->GetBinContent(2*i+1));
866  QCD_stat_err2_repack->SetBinContent(i,QCD_stat_err2->GetBinContent(2*i)+
867  QCD_stat_err2->GetBinContent(2*i+1));
868  eemc_stat_err2_repack->SetBinContent(i,eemc_stat_err2->GetBinContent(2*i)+
869  eemc_stat_err2->GetBinContent(2*i+1));
870  tau_stat_err2_repack->SetBinContent(i,tau_stat_err2->GetBinContent(2*i)+
871  tau_stat_err2->GetBinContent(2*i+1));
872  QCD_syst_high_err_repack->SetBinContent(i,QCD_syst_high_err->GetBinContent(2*i)+QCD_syst_high_err->GetBinContent(2*i+1));
873  QCD_syst_low_err_repack->SetBinContent(i,QCD_syst_low_err->GetBinContent(2*i)+QCD_syst_low_err->GetBinContent(2*i+1));
874  zsig_stat_err2_repack->SetBinContent(i,zsig_stat_err2->GetBinContent(2*i)+
875  zback_stat_err2->GetBinContent(2*i+1));
876  zback_stat_err2_repack->SetBinContent(i,zback_stat_err2->GetBinContent(2*i)+
877  zback_stat_err2->GetBinContent(2*i+1));
878  zeemc_stat_err2_repack->SetBinContent(i,zeemc_stat_err2->GetBinContent(2*i)+
879  zeemc_stat_err2->GetBinContent(2*i+1));
880  }
881 
882  tauhist_repack->Write();
883  zsig_bkgd2_repack->Write();
884  signal2_repack->Write();
885  signal_final2_repack->Write();
886  signal_final3_repack->Write();
887  eemc_bkgd2_repack->Write();
888  new_bkgd_repack->Write();
889  low_bkgd_repack->Write();
890  high_bkgd_repack->Write();
891  signal_in_norm_region_repack->Write();
892 
893  raw_stat_err2_repack->Write();
894  QCD_stat_err2_repack->Write();
895  eemc_stat_err2_repack->Write();
896  tau_stat_err2_repack->Write();
897  QCD_syst_high_err_repack->Write();
898  QCD_syst_low_err_repack->Write();
899  zsig_stat_err2_repack->Write();
900  zback_stat_err2_repack->Write();
901  zeemc_stat_err2_repack->Write();
902 
903 #if 0
904  for (int i=1; i<=24; i++) {
905  cout.setf(ios::fixed);
906  cout.precision(2);
907  cout << " " << signal2_repack->GetBinCenter(i) << " & " << signal2_repack->GetBinContent(i) << " & " << signal_final3_repack->GetBinContent(i) << " & " << QCD_stat_err2_repack->GetBinContent(i) << " & " << eemc_stat_err2_repack->GetBinContent(i) << " & " << tau_stat_err2_repack->GetBinContent(i) << " & " << zsig_stat_err2_repack->GetBinContent(i) << " & " << zback_stat_err2_repack->GetBinContent(i) << " & " << QCD_syst_high_err_repack->GetBinContent(i) << " & " << QCD_syst_low_err_repack->GetBinContent(i) << " \\\\" << endl;
908  }
909 #endif
910 
911  f2->Close();
912 
913  }
914 
915 
916 } // end of macro
Generated by   doxygen 1.8.5