C++ TCanvas::Clear方法代码示例

//.........这里部分代码省略.........
	dataV.plotOn(w1frame, MarkerColor(kBlue),MarkerStyle(27),MarkerSize(1.9),XErrorSize(0), Rescale(rescale*.95823), DataError(RooAbsData::None));
	myPDFV->plotOn(w1frame, LineColor(kBlue),LineWidth(2), Normalization(rescale*.95823));
      }
      if ( rescale != 1. )
      w1frame->GetYaxis()->SetRangeUser(0, ymax_w1 / 1000. * 1.5);

      // 
      //  wminus
      // 
      RooPlot* w2frame =  wminusmass->frame(50);

      w2frame->GetXaxis()->CenterTitle();
      w2frame->GetYaxis()->CenterTitle();
      w2frame->GetYaxis()->SetTitle(" ");
      
      double ymax_w2;
      TH1F *w2a = new TH1F("w2a", "w2a", 50, 1e-09, 120);
      tinA->Project("w2a", "wminusmass");
      ymax_w2 = w2a->GetMaximum();
      
      TH1F *w2_minus = new TH1F("w2_minus", "w2_minus", 50, 1e-09, 120);
      tinV->Project("w2_minus", "wminusmass")
      ymax_w2 = w2_minus->GetMaximum() > ymax_w2 ? w2_minus->GetMaximum() : ymax_w2;
      
      if ( drawa ) {
	dataA.plotOn(w2frame, MarkerColor(kRed),MarkerStyle(4),MarkerSize(1.5),LineWidth(0),XErrorSize(0), Rescale(rescale), DataError(RooAbsData::None));
	myPDFA->plotOn(w2frame, LineColor(kRed),LineWidth(2), Normalization(rescale));
      }
      if ( drawv ) {
	//dataV.plotOn(w2frame, MarkerColor(kBlue),MarkerStyle(27),MarkerSize(1.9),XErrorSize(0), Rescale(rescale), DataError(RooAbsData::None));
	//myPDFV->plotOn(w2frame, LineColor(kBlue),LineWidth(2), Normalization(rescale));
	dataV.plotOn(w2frame, MarkerColor(kBlue),MarkerStyle(27),MarkerSize(1.9),XErrorSize(0), Rescale(rescale*.95823), DataError(RooAbsData::None));
	myPDFV->plotOn(w2frame, LineColor(kBlue),LineWidth(2), Normalization(rescale*.95823));
      }
      if ( rescale != 1. ) 
	w2frame->GetYaxis()->SetRangeUser(0, ymax_w2 / 1000. * 1.5);
    }
    if ( drawpaper ) {
      TCanvas* can =new TCanvas("can","can",600,600);

      if ( offshell ) {
	w1frame->GetXaxis()->SetTitle("m_{l#nu} [GeV]");
	w1frame->Draw();
	can->Print(Form("paperplots/wplusmass_%.0fGeV_spin1_2in1_ww.eps", mH));
	can->SaveAs(Form("paperplots/wplusmass_%.0fGeV_spin1_2in1_ww.C", mH));
      }
      
      can->Clear();
      hsframe->Draw();
      can->Print(Form("paperplots/costhetastar_%.0fGeV_spin1_2in1_ww.eps", mH));      
      can->SaveAs(Form("paperplots/costhetastar_%.0fGeV_spin1_2in1_ww.C", mH));      
      
      can->Clear();
      Phi1frame->Draw();
      can->Print(Form("paperplots/phistar1_%.0fGeV_spin1_2in1_ww.eps", mH));      
      can->SaveAs(Form("paperplots/phistar1_%.0fGeV_spin1_2in1_ww.C", mH));      

      can->Clear();
      h1frame->GetXaxis()->SetTitle("cos#theta_{1} or cos#theta_{2}");
      h1frame->Draw();
      can->Print(Form("paperplots/costheta1_%.0fGeV_spin1_2in1_ww.eps", mH));
      can->SaveAs(Form("paperplots/costheta1_%.0fGeV_spin1_2in1_ww.C", mH));

      can->Clear();
      Phiframe->Draw();
      can->Print(Form("paperplots/phi_%.0fGeV_spin1_2in1_ww.eps", mH));      
      can->SaveAs(Form("paperplots/phi_%.0fGeV_spin1_2in1_ww.C", mH));      


    }

    else {
      
      TCanvas* cww = new TCanvas( "cww", "cww", 1000, 600 );
      cww->Divide(4,2);
      if ( offshell ) {
	cww->cd(1);
	w1frame->Draw();
	cww->cd(2);
	w2frame->Draw();
      }
      cww->cd(3);
      hsframe->Draw();
      cww->cd(4);
      box3->Draw();
      cww->cd(5);
      Phi1frame->Draw();
      cww->cd(6);
      h1frame->Draw();
      cww->cd(7);
      h2frame->Draw();
      cww->cd(8);
      Phiframe->Draw();
      
      cww->Print(Form("epsfiles/angles_VWW%.0f_JHU_7D.eps", mH));
      cww->Print(Form("pngfiles/angles_VWW%.0f_JHU_7D.png", mH));
      delete cww;
    }

}

//.........这里部分代码省略.........
  pltWpEN->AddHist1D(h1_WpT_WpToEN_MCNom_Log,"H",3,1,1001,2,0);
  pltWmEN->AddHist1D(h1_WpT_WmToEN_MCNom_Log,"H",3,1,1001,2,0);
  pltWpEN->AddHist1D(h1_WpT_WpToEN_RDNom,"E",kBlack,1,0,21,1);
  pltWmEN->AddHist1D(h1_WpT_WmToEN_RDNom,"E",kBlack,1,0,21,1);

  pltWpMN->SetLegend(0.53,0.7,0.7,0.8);
  pltWmMN->SetLegend(0.53,0.7,0.7,0.8);
  pltWpEN->SetLegend(0.53,0.7,0.7,0.8);
  pltWmEN->SetLegend(0.53,0.7,0.7,0.8);
  //cout << "check4" << endl;

  pltWpMN->GetLegend()->SetTextSize(0.03);
  pltWmMN->GetLegend()->SetTextSize(0.03);
  pltWpEN->GetLegend()->SetTextSize(0.03);
  pltWmEN->GetLegend()->SetTextSize(0.03);

  pltWpMN->GetLegend()->AddEntry(h1_WpT_WpToMN_MCNom_Log,"Nominal Truth Post","l");
  pltWmMN->GetLegend()->AddEntry(h1_WpT_WmToMN_MCNom_Log,"Nominal Truth Post","l");
  pltWpMN->GetLegend()->AddEntry(h1_WpT_WpToMN_RDNom,"DetUnf Fit Data ","p");
  pltWmMN->GetLegend()->AddEntry(h1_WpT_WmToMN_RDNom,"DetUnf Fit Data ","p");

  pltWpEN->GetLegend()->AddEntry(h1_WpT_WpToEN_MCNom_Log,"Nominal Truth Post","l");
  pltWmEN->GetLegend()->AddEntry(h1_WpT_WmToEN_MCNom_Log,"Nominal Truth Post","l");
  pltWpEN->GetLegend()->AddEntry(h1_WpT_WpToEN_RDNom,"DetUnf Fit Data","p");
  pltWmEN->GetLegend()->AddEntry(h1_WpT_WmToEN_RDNom,"DetUnf Fit Data","p");

  pltWpMN->SetLogx();
  pltWmMN->SetLogx();
  pltWpEN->SetLogx();
  pltWmEN->SetLogx();
  
  
  pltWpMN->Draw(myCan,kTRUE,format);
  pltWmMN->Draw(myCan,kTRUE,format);
  pltWpEN->Draw(myCan,kTRUE,format);
  pltWmEN->Draw(myCan,kTRUE,format);

  sprintf(tmpName,"WpT_WpMN_ratio");
  pltWpMN_ratio= new CPlot(tmpName,"",tmpName,"ratio/bin");
  sprintf(tmpName,"WpT_WmMN_ratio");
  pltWmMN_ratio= new CPlot(tmpName,"",tmpName,"ratio/bin");

  sprintf(tmpName,"WpT_WpEN_ratio");
  pltWpEN_ratio= new CPlot(tmpName,"",tmpName,"ratio/bin");
  sprintf(tmpName,"WpT_WmEN_ratio");
  pltWmEN_ratio= new CPlot(tmpName,"",tmpName,"ratio/bin");

  pltWpMN_ratio->setOutDir(resultDir);
  pltWmMN_ratio->setOutDir(resultDir);
  pltWpEN_ratio->setOutDir(resultDir);
  pltWmEN_ratio->setOutDir(resultDir);

  pltWpMN_ratio->AddHist1D(h1_WpT_WpToMN_ratio,"E",1,1,1001,21,1);
  pltWmMN_ratio->AddHist1D(h1_WpT_WmToMN_ratio,"E",1,1,1001,21,1);
  pltWpEN_ratio->AddHist1D(h1_WpT_WpToEN_ratio,"E",1,1,1001,21,1);
  pltWmEN_ratio->AddHist1D(h1_WpT_WmToEN_ratio,"E",1,1,1001,21,1);

  pltWpMN_ratio->SetLegend(0.7,0.8,0.9,0.9);
  pltWmMN_ratio->SetLegend(0.7,0.8,0.9,0.9);
  pltWpEN_ratio->SetLegend(0.7,0.8,0.9,0.9);
  pltWmEN_ratio->SetLegend(0.7,0.8,0.9,0.9);

  pltWpMN_ratio->GetLegend()->SetTextSize(0.03);
  pltWmMN_ratio->GetLegend()->SetTextSize(0.03);
  pltWpEN_ratio->GetLegend()->SetTextSize(0.03);
  pltWmEN_ratio->GetLegend()->SetTextSize(0.03);

  //pltWpMN_ratio->GetLegend()->AddEntry(h1_WpT_WpToMN_ratio,"WpToMN PostFSR Data/MC TruthPost");
  //pltWmMN_ratio->GetLegend()->AddEntry(h1_WpT_WmToMN_ratio,"WmToMN PostFSR Data/MC TruthPost");
  //pltWpEN_ratio->GetLegend()->AddEntry(h1_WpT_WpToEN_ratio,"WpToEN PostFSR Data/MC TruthPost");
  //pltWmEN_ratio->GetLegend()->AddEntry(h1_WpT_WmToEN_ratio,"WmToEN PostFSR Data/MC TruthPost");

  pltWpMN_ratio->GetLegend()->AddEntry(h1_WpT_WpToMN_ratio,"WpToMN DetUnfFit Data/MC TruthPost");
  pltWmMN_ratio->GetLegend()->AddEntry(h1_WpT_WmToMN_ratio,"WmToMN DetUnfFit Data/MC TruthPost");
  pltWpEN_ratio->GetLegend()->AddEntry(h1_WpT_WpToEN_ratio,"WpToEN DetUnfFit Data/MC TruthPost");
  pltWmEN_ratio->GetLegend()->AddEntry(h1_WpT_WmToEN_ratio,"WmToEN DetUnfFit Data/MC TruthPost");
  
  pltWpMN_ratio->SetLogx();
  pltWmMN_ratio->SetLogx();
  pltWpEN_ratio->SetLogx();
  pltWmEN_ratio->SetLogx();
  
  pltWpMN_ratio->Draw(myCan,kTRUE,format);
  pltWmMN_ratio->Draw(myCan,kTRUE,format);
  pltWpEN_ratio->Draw(myCan,kTRUE,format);
  pltWmEN_ratio->Draw(myCan,kTRUE,format);

  //TFile *RootOut= new TFile("./WpT_PostFSRData_MC_ratio.root","recreate");
  TFile *RootOut= new TFile("../RstSelection/WpT_FitDetUnfData_MC_ratio.root","recreate");
  h1_WpT_WpToMN_ratio->Write();
  h1_WpT_WmToMN_ratio->Write();
  h1_WpT_WpToEN_ratio->Write();
  h1_WpT_WmToEN_ratio->Write();
  RootOut->Close();
  
  Fout.close();
  myCan->Clear();

  return 0;
}

//.........这里部分代码省略.........
	IntCharMap cuts;
	cuts.insert(ICPair(protE, "protE"));
	cuts.insert(ICPair(copl, "copl"));
	cuts.insert(ICPair(balance, "balance"));
	cuts.insert(ICPair(dAlpha, "dAlpha"));
	cuts.insert(ICPair(missM, "missM"));
	cuts.insert(ICPair(invM, "invM"));
	cuts.insert(ICPair(copl_balance, "copl_balance"));
	cuts.insert(ICPair(balance_missM, "balance_missM"));
	cuts.insert(ICPair(protE_copl, "protE_copl"));
	cuts.insert(ICPair(copl_missM, "copl_missM"));
	cuts.insert(ICPair(missM_invM, "missM_invM"));
	cuts.insert(ICPair(balance_dAlpha, "balance_dAlpha"));
	cuts.insert(ICPair(copl_dAlpha, "copl_dAlpha"));
	cuts.insert(ICPair(protE_copl_balance, "protE_copl_balance"));
	cuts.insert(ICPair(copl_balance_dAlpha, "copl_balance_dAlpha"));
	cuts.insert(ICPair(all, "allCuts"));

	std::vector<const char*> props;
	props.push_back("invM");
	props.push_back("missM");

	const char* windows[2] = {"prompt", "random"};

	// connect to the file to read from
	TFile f(file, "READ");
	if (!f.IsOpen()) {
		printf("Error opening file %s!\n", file);
		return 1;
	}
	printf("Read ");
	f.Print();
	// load the desired folder from file (not needed in case the EndFileMacro was used)
	//TFolder *t = (TFolder*)f.Get("ROOT Memory");

	// Get the Prompt-Random-Ratio from the analysed root file
	TH1F *h_ratio = (TH1F*)f.FindObjectAny("PHYS_promptRandomRatio");
	const double ratio = h_ratio->GetMean();
	std::cout << "[INFO] Prompt-Random-Ratio: " << ratio << std::endl;

	// declare some temporary used variables
	TH1F* h[N_WINDOWS];
	char buffer[50];
	// create an empty canvas for histogram drawing
	TCanvas *c = new TCanvas("c", "plot", 20, 10, 700, 500);
	c->SetFillColor(0);
	c->SetBorderMode(0);
	c->SetBorderSize(2);
	//c->SetLogy();
	c->SetFrameBorderMode(0);
	c->SetLeftMargin(.14);
	c->SetRightMargin(.1);
	c->SetBottomMargin(.125);
	c->SetTopMargin(.05);

	for (std::vector<const char*>::iterator it_prop = props.begin(); it_prop != props.end(); ++it_prop) {
		for (ICIter it_cut = cuts.begin(); it_cut != cuts.end(); ++it_cut) {
			// skip the current combination if the cut contains the current variable, e. g. invM with some cut involving invM
			if (strstr(it_cut->second, *it_prop))
				continue;
			// get prompt and random histograms for every property cut combination
			for (unsigned int i = 0; i < N_WINDOWS; i++) {
				sprintf(buffer, "%s_%s_%s", *it_prop, it_cut->second, windows[i]);
				h[i] = (TH1F*)f.FindObjectAny(buffer);
				if (!h[i]) {
					std::cout << "[ERROR] Couldn't find histogram " << buffer << std::endl;
					exit(1);
				}
			}
			c->Clear();
			// subtract to chosen window size normalized random-background from events inside the prompt-window
			h[PROMPT]->Add(h[RANDOM], -1*ratio);
			prepare_hist(h[PROMPT]);
			/* Rebin histogram */
			int n = 50;
			h[PROMPT]->Rebin(n);  // combine n bins => n MeV binning
			/*char* label;
			if (label)
				free(label);
			char* label = malloc(strlen("#Events / ") + 5 + strlen(" MeV^{2}") + 1);
			//sprintf(label, "%s%.2f%s", "#Events / ", n/1000., " GeV^{2}");
			sprintf(label, "%s%d%s", "#Events / ", n, " MeV^{2}");
			h->GetYaxis()->SetTitle(label);*/
			/* set axes according to plotted values and prepare histogram for saving */
			sprintf(buffer, "%s%d%s", "#Events / ", n, " MeV");
			h[PROMPT]->GetYaxis()->SetTitle(buffer);
			h[PROMPT]->SetLineColor(color);
			sprintf(buffer, "%s [MeV]", *it_prop);
			h[PROMPT]->GetXaxis()->SetTitle(buffer);
			h[PROMPT]->Draw();
			sprintf(buffer, "%s_%s.%s", *it_prop, it_cut->second, suffix);
			c->Print(buffer);
		}
		std::cout << "[INFO] Created all random-subtracted plots for " << *it_prop << std::endl;
	}

	f.Close();

	return EXIT_SUCCESS;
}

//.........这里部分代码省略.........
    p->SetLogy();
  }
  p->SetLogz();
  h->Draw("colz");

  TObjArray* fits = new TObjArray;
  h->FitSlicesY(0, 1, -1, 0, "QN", fits);

  TF1* mean = new TF1("mean", "pol1");
  TF1* var  = new TF1("var", "pol1");
  // mean->FixParameter(0, 0);
  // var->FixParameter(0, 0);
  for (Int_t i = 0; i < 3; i++) { 
    p = c->cd(2+i);
    if (spd) { 
      p->SetLogx();
      p->SetLogy();
    }
    TH1* hh = static_cast<TH1*>(fits->At(i));
    hh->Draw();

    if (i == 0) continue;
    
    hh->Fit((i == 1? mean : var), "+Q");
    
  }

  TGraphErrors* g1 = new TGraphErrors(h->GetNbinsX());
  g1->SetFillColor(kBlue-10);
  g1->SetFillStyle(3001);
  g1->SetLineStyle(1);
  TGraph* u1 = new TGraph(h->GetNbinsX());
  TGraph* l1 = new TGraph(h->GetNbinsX());
  u1->SetLineColor(kBlue+1);
  l1->SetLineColor(kBlue+1);
  u1->SetName("u1");
  l1->SetName("l1");
  TGraphErrors* g2 = new TGraphErrors(h->GetNbinsX());
  g2->SetFillColor(kRed-10);
  g2->SetFillStyle(3001);
  g2->SetLineStyle(2);
  TGraph* u2 = new TGraph(h->GetNbinsX());
  TGraph* l2 = new TGraph(h->GetNbinsX());
  u2->SetLineColor(kRed+1);
  l2->SetLineColor(kRed+1);
  u2->SetName("u2");
  l2->SetName("l2");
  for (Int_t i = 1; i <= h->GetNbinsX(); i++) {
    Double_t x  = hh->GetXaxis()->GetBinCenter(i);
    Double_t y  = mean->Eval(x);
    Double_t e  = var->Eval(y);
    Double_t e1 = nVar * e;
    if (spd) e1 *= TMath::Log10(e);
    // Printf("%10e -> %10e +/- %10e", x, y, ee);
    g1->SetPoint(i-1, x, y);
    g1->SetPointError(i-1, 0, e1);
    u1->SetPoint(i-1, x, y+e1);
    l1->SetPoint(i-1, x, y-e1);
    // Printf("%3d: %f -> %f +/- %f", i, x, y, ee);

    Double_t e2 = nVar*0.05*x;
    g2->SetPoint(i-1, x, x);
    g2->SetPointError(i-1, 0, e2);
    u2->SetPoint(i-1, x, x+e2);
    l2->SetPoint(i-1, x, x-e2);
  }

  p = c->cd(1);
  c->Clear();
  c->cd();
  c->SetLogz();
  h->Draw("colz");
  g1->Draw("3 same");
  u1->Draw("l same");
  l1->Draw("l same");
  g2->Draw("3 same");
  u2->Draw("l same");
  l2->Draw("l same");

  Double_t ly = 0.9;
  Double_t dy = 0.06;
  TLatex* ltx = new TLatex(0.15, ly, Form("#LTy#GT = %f + %f x",
					   mean->GetParameter(0),
					   mean->GetParameter(1)));
  ltx->SetNDC();
  ltx->SetTextSize(dy);
  ltx->SetTextAlign(13);
  ltx->Draw();

  ly -= dy + 0.01;
  ltx->DrawLatex(0.15, ly, Form("#sigma_{y} = %f + %f x", 
				var->GetParameter(0),
				var->GetParameter(1)));
  
  ly -= dy + 0.01;
  ltx->DrawLatex(0.15, ly, Form("#delta = %f #sigma %s", 
				nVar, (spd ? "log_{10}(#sigma" : "")));
	    
					   
}

void startCorrelation(int run,int eventNumber){

  //int entry=65000;
  char runRootDirectory[FILENAME_MAX];
  sprintf(runRootDirectory,"/unix/anita3/flight0809/root/");

  char headName[FILENAME_MAX];
  char eventName[FILENAME_MAX];
  char gpsName[FILENAME_MAX];
  RawAnitaEvent *evPtr = 0;
  RawAnitaHeader *hdPtr = 0;
  Adu5Pat *patPtr =0;

  sprintf(headName,"%s/run%d/headFile%d.root",runRootDirectory,run,run);
  sprintf(eventName,"%s/run%d/eventFile%d.root",runRootDirectory,run,run);
  sprintf(gpsName,"%s/run%d/gpsFile%d.root",runRootDirectory,run,run);

  TFile *eventFile = new TFile(eventName);
  TFile *headFile = new TFile(headName);
  TFile *fpGps = new TFile(gpsName);

  if(!fpGps){
    std::cout << "no GPS file\n";
    return;
  }

  TTree *adu5PatTree = (TTree*) fpGps->Get("adu5PatTree");
  TTree *eventTree = (TTree*)eventFile->Get("eventTree");
  TTree *headTree = (TTree*)headFile->Get("headTree");

  eventTree->SetBranchAddress("event",&evPtr);
  headTree->SetBranchAddress("header",&hdPtr);
  headTree->BuildIndex("eventNumber");
  adu5PatTree->SetBranchAddress("pat",&patPtr);
  adu5PatTree->BuildIndex("realTime");

  int entry=headTree->GetEntryNumberWithBestIndex(eventNumber);

  eventTree->GetEntry(entry);
  headTree->GetEntry(entry);



  /*
  //for(int phi=0;phi<16;phi++){
    while(hdPtr->triggerTimeNs<349.99e6 || hdPtr->triggerTimeNs>350.005e6 || hdPtr->l3TrigPattern==0){
      if(lastEvent==hdPtr->eventNumber){
	std::cout << "no more entries in run" << std::endl;
	return;
      }
      std::cout << "opened entry " << entry << " (event " << hdPtr->eventNumber << ") with triggerTimeNs " << hdPtr->triggerTimeNs << std::endl;
      lastEvent=hdPtr->eventNumber;
      entry++;
      eventTree->GetEntry(entry);
      headTree->GetEntry(entry);
      adu5PatTree->GetEntry(entry);
    }
    //}
    */

  std::cout << "opened entry " << entry << " (event " << hdPtr->eventNumber << ") with triggerTimeNs " << hdPtr->triggerTimeNs << std::endl;

  //get the pat ptr that corresponds to the timing of the event
  Int_t patEntry;
  patEntry = adu5PatTree->GetEntryNumberWithBestIndex(hdPtr->realTime);
  adu5PatTree->GetEntry(patEntry);

  std::cout << patPtr->realTime << " " << hdPtr->realTime << std::endl;

  if(patPtr->realTime != hdPtr->realTime){
    std:: cout << "pat time doesn't match head time, pat realTime: " << patPtr->realTime << " head realTime: " << hdPtr->realTime << std::endl;
    //return;
  }


  for(int chan=0;chan<90;chan++){
    if(evPtr->xMax[chan]>130 || evPtr->xMin[chan]<-130){
      saturatedChannel[chan]=1;
      std::cout << "saturated " << chan << std::endl;
    }
    else saturatedChannel[chan]=0;
  }


  //entry++;
  TH2D *crossCorrelation = crossCorrelate(evPtr,hdPtr,patPtr);
 
  sprintf(headName,"crossCorrCan");
  TCanvas *crossCorrCan = (TCanvas*)gROOT->FindObject(headName);
  if(!crossCorrCan)
    crossCorrCan = new TCanvas(headName,headName,800,400);
  crossCorrCan->Clear();
  //sumCrossCorrs->Draw("aitoff");
  crossCorrelation->Draw("colz");
  TAxis *phiSectors = new TAxis(16,1,16);
  phiSectors->Draw();

  //std::cout << "event " << hdPtr->eventNumber << " time " << 

  double theta,phi;
//.........这里部分代码省略.........

//.........这里部分代码省略.........
	hSum[10]->Add(h[38]);
	hSum[10]->Add(h[41]);
	hSum[11]->Add(h[34]);
	hSum[11]->Add(h[37]);
	hSum[11]->Add(h[40]);
	hRatio[6]->Divide(hSum[11], hSum[10]);
	cv->Divide(2, 1);

	pd = cv->cd(1);
	pd->Divide(1, 2);
	pd->cd(1);
	hSum[0]->SetTitle("Up - Middle - Down;MeV;mHz");
	hSum[0]->SetLineColor(kRed);
	hSum[2]->SetLineColor(kGreen);
	hSum[1]->SetLineColor(kBlue);
	for (i=0; i<4; i++) hSum[i]->SetLineWidth(2);
	hSum[0]->Draw();
	hSum[2]->Draw("same");
	hSum[1]->Draw("same");
	val = hSum[0]->IntegralAndError(1, 35, err);
	sprintf(str, "Up = %5.2f+-%4.2f", val, err);
	txt->SetTextColor(kRed);
	txt->DrawTextNDC(0.4, 0.85, str);
	val = hSum[2]->IntegralAndError(1, 35, err);
	sprintf(str, "Mid = %5.2f+-%4.2f", val, err);
	txt->SetTextColor(kGreen);
	txt->DrawTextNDC(0.4, 0.8, str);
	val = hSum[1]->IntegralAndError(1, 35, err);
	sprintf(str, "Down = %5.2f+-%4.2f", val, err);
	txt->SetTextColor(kBlue);
	txt->DrawTextNDC(0.4, 0.75, str);
	
	pd->cd(2);
	hSum[3]->SetTitle("Average of all positions;MeV;mHz");
	hSum[3]->SetLineColor(kBlack);
	hSum[3]->Draw();
	val = hSum[3]->IntegralAndError(1, 35, err);
	sprintf(str, "Avr = %5.2f+-%4.2f", val, err);
	txt->SetTextColor(kBlack);
	txt->DrawTextNDC(0.4, 0.8, str);

	pd = cv->cd(2);
	pd->Divide(1, 4);
	pd->cd(1);
	hRatio[0]->SetTitle("Down/Up all by all;MeV;");
	hRatio[0]->SetLineColor(kBlack);
	hRatio[0]->SetLineWidth(2);
	hRatio[0]->SetMinimum(0.5);
	hRatio[0]->SetMaximum(1.0);
	hRatio[0]->Fit("pol0", "", "", 1, 7);
	
	pd->cd(2);
	hRatio[2]->SetTitle("Down/Up pair by pair;MeV;");
	hRatio[2]->SetLineColor(kBlack);
	hRatio[2]->SetLineWidth(2);
	hRatio[2]->SetMinimum(0.5);
	hRatio[2]->SetMaximum(1.0);
	hRatio[2]->Fit("pol0", "", "", 1, 7);

	pd->cd(3);
	hRatio[3]->SetTitle("Mid/Up all by all;MeV;");
	hRatio[3]->SetLineColor(kBlack);
	hRatio[3]->SetLineWidth(2);
	hRatio[3]->SetMinimum(0.65);
	hRatio[3]->SetMaximum(1.15);
	hRatio[3]->Fit("pol0", "", "", 1, 7);

	pd->cd(4);
	hRatio[1]->SetTitle("First half/Last half;MeV;");
	hRatio[1]->SetLineColor(kBlack);
	hRatio[1]->SetLineWidth(2);
	hRatio[1]->SetMinimum(0.75);
	hRatio[1]->SetMaximum(1.25);
	hRatio[1]->Fit("pol0", "", "", 1, 7);

	cv->Update();
	cv->Print(pname);

	cv->Clear();
	cv->Divide(2, 2);
	for (i=0; i<3; i++) {
		cv->cd(i+1);
		sprintf(str, "Down/Up all by all for %s;MeV;", periods[i]);
		hRatio[4+i]->SetTitle(str);
		hRatio[4+i]->SetLineColor(kBlack);
		hRatio[4+i]->SetLineWidth(2);
		hRatio[4+i]->SetMinimum(0.5);
		hRatio[4+i]->SetMaximum(1.0);
		hRatio[4+i]->Fit("pol0", "", "", 1, 7);
	}
	cv->Update();
	cv->Print(pname);

	sprintf(str, "%s]", pname);
	cv->Print(str);
	for (i=0; i<12; i++) hSum[i]->Write();
	for (i=0; i<7; i++) hRatio[i]->Write();
	f->Close();
	delete cv;
}

TObjArray getDistributionFromPlotter(TString plot,TString baseURL="~/scratch0/top-nosyst/plotter.root")
{
  using namespace std;

  setStyle();
  
  //the relevant processes
  TString procs[]={"Di-bosons","Single top", "W+jets", "Z-#gamma^{*}+jets#rightarrow ll", "other t#bar{t}", "t#bar{t} dileptons", "data"};
  const size_t nprocs=sizeof(procs)/sizeof(TString);

  //containers for the histos
  TList *data     = new TList;
  TH1F *totalData=0;
  TList *mc       = new TList;
  TH1F *totalMC=0;
  TH1F *ttbarMC=0;
  TList *spimpose = new TList;
  
  //open file with histograms
  TFile *f=TFile::Open(baseURL);
  for(size_t iproc=0; iproc<nprocs; iproc++)
    {
      TH1F *histo = (TH1F *) f->Get(procs[iproc]+"/"+plot);
      if(histo==0) { cout << "[Warning] " << plot << " not found for " << procs[iproc] << " ...skipping" << endl; continue; }
      histo->SetDirectory(0);
      if(procs[iproc].Contains("data")) 
	{
	  data->Add(histo);
	  if(totalData==0) { totalData = (TH1F *) histo->Clone(plot+"totaldata"); totalData->SetDirectory(0); }
	  else             { totalData->Add(histo); }
	}
      else
	{
	  mc->Add(histo);
	  if(totalMC==0) { totalMC = (TH1F *) histo->Clone(plot+"totalmc"); totalMC->SetDirectory(0); }
	  else           { totalMC->Add(histo); }
	  if(procs[iproc]=="t#bar{t} dileptons") { ttbarMC = (TH1F *)histo->Clone(plot+"ttbar"); ttbarMC->SetDirectory(0); }
	}
    }
  f->Close();

  cout << "Drawing now" << endl;
  
  //draw
  TCanvas *cnv = getNewCanvas(plot+"c",plot+"c",false);
  cnv->Clear();
  cnv->SetWindowSize(600,600);
  cnv->cd(); 
  TLegend *leg=showPlotsAndMCtoDataComparison(cnv,*mc,*spimpose,*data,false);
  formatForCmsPublic(cnv,leg,"CMS preliminary", 4);
  cnv->SaveAs(plot+".C");

  cout << "Stat comparison for " << plot << endl
       << "Sample \t Average \t\t RMS " << endl
       << "----------------------------------------" << endl;
  if(totalData) cout << "Data \t " << totalData->GetMean() << " +/- " << totalData->GetMeanError() << " \t " << totalData->GetRMS() << " +/- " << totalData->GetRMSError() << endl;
  if(totalMC)   cout << "MC \t "   << totalMC->GetMean() << " +/- " << totalMC->GetMeanError() << " \t " << totalMC->GetRMS() << " +/- " << totalMC->GetRMSError() << endl;
  if(ttbarMC)   cout << "Signal \t "   << ttbarMC->GetMean() << " +/- " << ttbarMC->GetMeanError() << " \t " << ttbarMC->GetRMS() << " +/- " << ttbarMC->GetRMSError() << endl;

  TObjArray res;
  res.Add(totalData);
  res.Add(totalMC);
  res.Add(ttbarMC);
  res.Add(mc);
  return res;
}

//.........这里部分代码省略.........
      for(std::vector<TChain *>::const_iterator mc_itr = mc_vec.begin();
	  mc_itr != mc_vec.end();
	  mc_itr++){
	TChain *mc = *mc_itr;
 	mc->SetBranchStatus("*",0);
 	for(std::set<TString>::const_iterator itr = branchList.begin();
 	    itr != branchList.end();
 	    itr++){
 	  //std::cout << "[STATUS] Enabling branch: " << *itr << std::endl;
 	  mc->SetBranchStatus(*itr, 1);
 	}
	mc->SetBranchStatus(branchname, 1);
	mc->SetBranchStatus("smearEle", 1);
	mc->SetBranchStatus("puWeight",1);
	mc->SetBranchStatus("mcGenWeight", 1);
	mc->SetBranchStatus("r9Weight", 1);
 	

	TString mcHistName;  mcHistName+="hist_"; mcHistName+=mc_itr-mc_vec.begin();//better for .C generation
	//TString mcHistName; mcHistName+=mc_itr-mc_vec.begin(); mcHistName+="_hist";//better for .C generation
	
	//decide this for MC
	//TString weights="mcGenWeight";
	//if(pdfIndex!="") weights+="*(pdfWeights_cteq66["+pdfIndex+"]/pdfWeights_cteq66[0])";
	//if(usePU) weights+="*puWeight";
	//if(useR9Weight) weights+="*r9Weight";
	//std::cout<<"Complete selection for MC is "<<selection_MC *weights.Data()<<std::endl;
	//mc->Draw(branchNameMC+">>"+mcHistName+binning, selection_MC *weights.Data());
	mc->Draw(branchNameMC+">>"+mcHistName+binning, selection_MC);
      }
    }


  c->Clear();
  TLegend *leg = new TLegend(0.5,0.7,0.7,0.85);
  leg->SetBorderSize(1);
  leg->SetFillColor(0);
  leg->SetTextSize(0.04);
//   if(dataLabel !="" && mcLabel !="") leg->Draw();
//   //c->GetListOfPrimitives()->Add(leg,"");


  TH1F *d = (TH1F *) gROOT->FindObject("data_hist");
  if(dataLabel !="") leg->AddEntry(d,dataLabel,"p");
  d->SetStats(0);
  d->SetTitle("");

  d->SetMarkerStyle(20);
  d->SetMarkerSize(1);

  if(d->GetEntries()==0 || d->Integral()==0){
    //d=(TH1F *) gROOT->FindObject("0_hist");
    d=(TH1F *) gROOT->FindObject("hist_0");
    d->SetMarkerSize(0);
  }
  //d->SaveAs("tmp/d_hist.root");
  //s->SaveAs("tmp/s_hist.root");

  yLabel.Form("Events /(%.2f %s)", d->GetBinWidth(2), yLabelUnit.Data());
  
  std::cout << "nEvents data: " << d->Integral() << "\t" << d->GetEntries() << std::endl;
  std::cout << "Now normalizing data to 1: "<<std::endl;
  std::cout << "80 is in bin "<<d->FindBin(80)<<std::endl;
  std::cout << "100 is in bin "<<d->FindBin(100)<<std::endl;
  std::cout << "d integral is "<<d->Integral()<<std::endl;
  std::cout << "d integral(1,100) is "<<d->Integral(1,100)<<std::endl;

//.........这里部分代码省略.........
     h_Nsig->SetLineColor(c_sig); h_Nsig_GJ->SetLineColor(c_sig); h_Nsig_QCD->SetLineColor(c_sig);
     h_Nsig->SetLineStyle(ls_cmb); h_Nsig_GJ->SetLineStyle(ls_gj); h_Nsig_QCD->SetLineStyle(ls_qcd);
     h_Nsig->SetLineWidth(2); h_Nsig_GJ->SetLineWidth(2); h_Nsig_QCD->SetLineWidth(2);
     //Background purity
     TH1F* h_Nbkg_GJ  = (TH1F*)infile_GJ->Get("h_Nbkg_"+isovar+"_"+ptrange+"_"+cut);
     TH1F* h_Nbkg_QCD = (TH1F*)infile_QCD->Get("h_Nbkg_"+isovar+"_"+ptrange+"_"+cut);
     h_Nbkg_GJ->Rebin(3); h_Nbkg_QCD->Rebin(3);
             //TH1F* h_Nbkg  = (TH1F*)h_Nbkg_QCD->Clone("h_Nbkg");
     TH1F* h_Nbkg  = (TH1F*)h_Nbkg_GJ->Clone("h_Nbkg");
     h_Nbkg->Add(h_Nbkg_QCD);
     h_Nbkg->SetLineColor(c_bkg); h_Nbkg_GJ->SetLineColor(c_bkg); h_Nbkg_QCD->SetLineColor(c_bkg);
     h_Nbkg->SetLineStyle(ls_cmb); h_Nbkg_GJ->SetLineStyle(ls_gj); h_Nbkg_QCD->SetLineStyle(ls_qcd);
     h_Nbkg->SetLineWidth(2); h_Nbkg_GJ->SetLineWidth(2); h_Nbkg_QCD->SetLineWidth(2);

     //Denominator
     TH1F* h_Deno_GJ  = (TH1F*)h_Nsig_GJ->Clone("h_Deno_GJ");
     h_Deno_GJ->Add(h_Nbkg_GJ);

     TH1F* h_Deno_QCD  = (TH1F*)h_Nsig_QCD->Clone("h_Deno_QCD");
     h_Deno_QCD->Add(h_Nbkg_QCD);

//     TH1F* h_Deno_GJ  = (TH1F*)infile_GJ->Get("h_Deno_"+isovar+"_"+ptrange+"_"+cut);
//     TH1F* h_Deno_QCD = (TH1F*)infile_QCD->Get("h_Deno_"+isovar+"_"+ptrange+"_"+cut);
//     h_Deno_GJ->Rebin(3); h_Deno_QCD->Rebin(3);
//             //TH1F* h_Deno  = (TH1F*)h_Deno_QCD->Clone("h_Deno");
     TH1F* h_Deno  = (TH1F*)h_Deno_GJ->Clone("h_Deno");
     h_Deno->Add(h_Deno_QCD);

     h_Deno->SetLineColor(c_den); h_Deno_GJ->SetLineColor(c_den); h_Deno_QCD->SetLineColor(c_den);
     h_Deno->SetLineStyle(ls_cmb); h_Deno_GJ->SetLineStyle(ls_gj); h_Deno_QCD->SetLineStyle(ls_qcd);
     h_Deno->SetLineWidth(1); h_Deno_GJ->SetLineWidth(1); h_Deno_QCD->SetLineWidth(1);

     // should probably add these as last bin..
     h_Nsig->ClearUnderflowAndOverflow();
     h_Nbkg->ClearUnderflowAndOverflow();
     h_Deno->ClearUnderflowAndOverflow();
     
     //Define purity
     TGraphAsymmErrors *purity_signal     = new TGraphAsymmErrors(h_Nsig,h_Deno,"n");
     purity_signal->SetLineWidth(2);
     purity_signal->SetLineColor(4);
     purity_signal->SetMarkerColor(4);
     purity_signal->SetMarkerStyle(20);
     
     TGraphAsymmErrors *purity_background = new TGraphAsymmErrors(h_Nbkg,h_Deno,"n");
     purity_background->SetLineWidth(2);
     purity_background->SetLineColor(2);
     purity_background->SetMarkerColor(2);
     purity_background->SetMarkerStyle(24);

     //-----------------------------------------------------
     // Calculations
     log<<boost::format("%8s %15s %7s %7s ") % ptrange % cut % isovar % rebin  ;

     // Start lower bound calculation
     Int_t nbins;
     nbins = purity_background->GetN();
     log<<boost::format("(%2i bins) \n") % nbins ; 

     std::vector<Double_t> thresholds, xsatthresh, ysatthresh, xsatqt, avgpur;
     std::vector<Int_t> batthresh, batqt;

     thresholds.clear(); // purity thresholds
     xsatthresh.clear(); // x st. f(x) > thresh
     ysatthresh.clear(); // y = f(x)x at thresh
     batthresh.clear();  // bin corresponding to x

//.........这里部分代码省略.........
    mc->SetBranchAddress("HLTNames", &HLTNames);
    mc->SetBranchAddress("HLTResults", &HLTResults);
    
    mc->GetEntries();
    for(Long64_t jentry=0; jentry < mc->GetEntriesFast(); jentry++){
      mc->GetEntry(jentry);
      std::vector<Bool_t>::const_iterator hltRes_itr=HLTResults->begin();
      for(std::vector<std::string>::const_iterator hlt_itr=HLTNames->begin();
	  hlt_itr!=HLTNames->end();
	  hlt_itr++,hltRes_itr++){
	std::cout << *hlt_itr << std::endl;
	hlt_mc->Fill(hlt_itr-HLTNames->begin(),*hltRes_itr);
      }
    }

    for(std::vector<std::string>::const_iterator hlt_itr=HLTNames->begin();
	hlt_itr!=HLTNames->end();
	hlt_itr++){
      int index=hlt_itr-HLTNames->begin()+1;
      if(hlt_mc->GetBinContent(index)>0)
      hlt_mc->GetXaxis()->SetBinLabel(index,(*hlt_itr).c_str());
    }
  } else {

  std::cout<<"qui"<<std::endl;
    data->Draw(branchNameData+">>data_hist"+binning, selection_data);
    if(mc!=NULL){
      if(usePU)  mc->Draw(branchNameMC+">>mc_hist"+binning, selection_MC *"puWeight");
    else  mc->Draw(branchNameMC+">>mc_hist"+binning, selection_MC);
  std::cout<<"qui"<<std::endl;
    }
  }

  c->Clear();
  TH1F *d = (TH1F *) gROOT->FindObject("data_hist");
  TH1F *s = (TH1F *) gROOT->FindObject("mc_hist");
  if(s==NULL) s=d;
  //d->SaveAs("tmp/d_hist.root");
  s->SaveAs("tmp/s_hist.root");
  std::cout<<"qui"<<std::endl;

  yLabel.Form("Events /(%.2f %s)", s->GetBinWidth(2), yLabelUnit.Data());
  float max = 1.1 * std::max(
			     d->GetMaximum(),///d->Integral(),
			     s->GetMaximum() ///s->Integral()
			     );
  max=1.1*d->GetMaximum();
  std::cout << "max = " << max << std::endl;
  std::cout << "nEvents data: " << d->Integral() << "\t" << d->GetEntries() << std::endl;
  std::cout << "nEvents signal: " << s->Integral() << "\t" << s->GetEntries() << std::endl;
  if(logy){
    max*=10;
    d->GetYaxis()->SetRangeUser(0.1,max);
    s->GetYaxis()->SetRangeUser(0.1,max);
    c->SetLogy();
  } else {
    d->GetYaxis()->SetRangeUser(0,max);
    s->GetYaxis()->SetRangeUser(0,max);
  }
  s->GetYaxis()->SetTitle(yLabel);
  s->GetXaxis()->SetTitle(xLabel);
  d->GetYaxis()->SetTitle(yLabel);
  d->GetXaxis()->SetTitle(xLabel);


  d->SetMarkerStyle(20);

void MakePlots(TString filename, TString energy="8TeV", TString lumi=""){
  TString outDir=filename; outDir.ReplaceAll("fitres","img");
  outDir="tmp/";
  //std::map<TString, TH2F *> deltaNLL_map;
  
  /*------------------------------ Plotto */
  TCanvas *c = new TCanvas("c","c");
  
  TFile f_in(filename, "read");
  if(f_in.IsZombie()){
    std::cerr << "File opening error: " << filename << std::endl;
    return;
  }
  
  TList *KeyList = f_in.GetListOfKeys();
  std::cout << KeyList->GetEntries() << std::endl;
  for(int i =0; i <  KeyList->GetEntries(); i++){
    c->Clear();
    TKey *key = (TKey *)KeyList->At(i);
    if(TString(key->GetClassName())!="RooDataSet") continue;
    RooDataSet *dataset = (RooDataSet *) key->ReadObj();
    
    TString constTermName = dataset->GetName();
    TString alphaName=constTermName; alphaName.ReplaceAll("constTerm","alpha");

    TTree *tree = dataset2tree(dataset);
    TGraphErrors bestFit_ = bestFit(tree, alphaName, constTermName);
    TH2F *hist = prof2d(tree, alphaName, constTermName, "nll", "(12,-0.0005,0.0115,29,-0.0025,0.1425)",true);

//     //    deltaNLL_map.insert(std::pair <TString, TH2F *>(keyName,hist));
    hist->SaveAs(outDir+"/deltaNLL-"+constTermName+".root");
    hist->Draw("colz");
    bestFit_.Draw("P same");
    bestFit_.SetMarkerSize(2);

    Int_t iBinX, iBinY;
    Double_t x,y;
    hist->GetBinWithContent2(0,iBinX,iBinY);
    x= hist->GetXaxis()->GetBinCenter(iBinX);
    y= hist->GetYaxis()->GetBinCenter(iBinY);
    TGraph nllBestFit(1,&x,&y);

    nllBestFit.SetMarkerStyle(3);
    nllBestFit.SetMarkerColor(kRed);
    TList* contour68 = contourFromTH2(hist, 0.68);

    hist->Draw("colz");
    hist->GetZaxis()->SetRangeUser(0,50);
    bestFit_.Draw("P same");
    nllBestFit.Draw("P same");
    //contour68->Draw("same");
    c->SaveAs(outDir+"/deltaNLL-"+constTermName+".png");
    hist->SaveAs("tmp/hist-"+constTermName+".root");
    nllBestFit.SaveAs("tmp/nllBestFit.root");
    contour68->SaveAs("tmp/contour68.root");
    delete hist;
    hist = prof2d(tree, alphaName, constTermName, "nll", "(12,-0.0005,0.0115,29,-0.0025,0.1425)");
    RooHistPdf *histPdf = nllToL(hist);
    delete hist;
    RooDataSet *gen_dataset=histPdf->generate(*histPdf->getVariables(),1000000,kTRUE,kFALSE);
    TTree *genTree = dataset2tree(gen_dataset);
    genTree->SaveAs("tmp/genTree-"+constTermName+".root");
    delete gen_dataset;
    delete histPdf;
    
    TGraphErrors toyGraph = g(genTree, constTermName);
    TGraphErrors bestFitGraph = g(tree,alphaName, constTermName);
    TGraphErrors bestFitScanGraph = g(y, x);
    delete genTree;
    delete tree;
    toyGraph.SetFillColor(kGreen);
    toyGraph.SetLineColor(kBlue);
    toyGraph.SetLineStyle(2);
    bestFitGraph.SetLineColor(kBlack);
    bestFitScanGraph.SetLineColor(kRed);
    bestFitScanGraph.SetLineWidth(2);


    
    TMultiGraph g_multi("multigraph","");
    g_multi.Add(&toyGraph,"L3");
    g_multi.Add(&toyGraph,"L");
    g_multi.Add(&bestFitGraph, "L");
    g_multi.Add(&bestFitScanGraph, "L");
   
    g_multi.Draw("A");

    c->Clear();
    g_multi.Draw("A");
    c->SaveAs(outDir+"/smearing_vs_energy-"+constTermName+".png");

    //    TPaveText *pv = new TPaveText(0.7,0.7,1, 0.8);    
//     TLegend *legend = new TLegend(0.7,0.8,0.95,0.92);
//     legend->SetFillStyle(3001);
//     legend->SetFillColor(1);
//     legend->SetTextFont(22); // 132
//     legend->SetTextSize(0.04); // l'ho preso mettendo i punti con l'editor e poi ho ricavato il valore con il metodo GetTextSize()
//   //  legend->SetFillColor(0); // colore di riempimento bianco
//     legend->SetMargin(0.4);  // percentuale della larghezza del simbolo
    //    SetLegendStyle(legend);
//.........这里部分代码省略.........

   void miplots2sig_18fb( bool doLogy=false ) {

       int fillStyle1sig = 1 ;
       int fillStyle2sig = 1 ;
       int fillColor1sig = kBlue-9 ;
       int fillColor2sig = kBlue-10 ;

       gStyle->SetOptStat(0) ;
       gStyle->SetLabelSize(0.12,"x") ;
       gStyle->SetLabelSize(0.08,"y") ;
       gStyle->SetTitleW(0.8) ;
       gStyle->SetTitleH(0.08) ;
       gStyle->SetTitleAlign(13) ;
       gStyle->SetPadLeftMargin(0.14) ;
       gStyle->SetPadRightMargin(0.02) ;
       gStyle->SetLabelOffset(0.02,"y") ;
       gStyle->SetPadGridY(1) ;

       if ( doLogy ) {
          gStyle->SetOptLogy(1) ;
       } else {
          gStyle->SetOptLogy(0) ;
       }

       gDirectory->Delete("h*") ;

       double xvals[4] = { 1.05, 2.05, 3.05, 4.05 } ;
       double xel[4] = {0.5, 0.5, 0.5, 0.5} ;
       double xeh[4] = {0.5, 0.5, 0.5, 0.5} ;


       TCanvas* cmip = (TCanvas*) gDirectory->FindObject("cmip") ;
       if ( cmip == 0x0 ) {
          cmip = new TCanvas("cmip", "MIP", 800, 600) ;
       }
       cmip->Clear() ;

       cmip->Divide(4,1) ;



    //--------

       cmip->cd(1) ;

       TH1F* h_data_nb2_m4 = new TH1F( "h_data_nb2_m4", "Data, nB=2, MET4", 4, 0.5, 4.5 ) ;

       h_data_nb2_m4 -> GetXaxis() -> SetBinLabel( 1, "HT1" ) ;
       h_data_nb2_m4 -> GetXaxis() -> SetBinLabel( 2, "HT2" ) ;
       h_data_nb2_m4 -> GetXaxis() -> SetBinLabel( 3, "HT3" ) ;
       h_data_nb2_m4 -> GetXaxis() -> SetBinLabel( 4, "HT4" ) ;

       h_data_nb2_m4 -> SetBinContent( 1, 0 ) ;
       h_data_nb2_m4 -> SetBinContent( 2, 57 ) ;
       h_data_nb2_m4 -> SetBinContent( 3, 17 ) ;
       h_data_nb2_m4 -> SetBinContent( 4, 19 ) ;

       h_data_nb2_m4 -> SetMarkerStyle(20) ;
       h_data_nb2_m4 -> SetMarkerSize(1.2) ;
       h_data_nb2_m4 -> SetLineWidth(2) ;
       if ( doLogy) h_data_nb2_m4 -> SetMinimum( 10 ) ;
       h_data_nb2_m4 -> SetMaximum( 100 ) ;


       h_data_nb2_m4 -> Draw("e") ;


       double nb2_m4_yvals[4] = { -99., 72.2, 20.0, 18.4 } ;
       double nb2_m4_yerrh[4] = { 0.0, 10.0, 4.4, 4.7 } ;
       double nb2_m4_yerrl[4] = { 0.0, 9.0, 3.8, 3.8 } ;
       double nb2_m4_yerrh_2s[4] = { 0.0, 21.3, 9.6, 10.5 } ;
       double nb2_m4_yerrl_2s[4] = { 0.0, 17.0, 7.0, 6.9 } ;

       TGraphAsymmErrors* gr_nb2_m4 = new TGraphAsymmErrors( 4, xvals, nb2_m4_yvals, xel, xeh,
            nb2_m4_yerrl, nb2_m4_yerrh ) ;

       TGraphAsymmErrors* gr_nb2_m4_2s = new TGraphAsymmErrors( 4, xvals, nb2_m4_yvals, xel, xeh,
            nb2_m4_yerrl_2s, nb2_m4_yerrh_2s ) ;

       gr_nb2_m4 -> SetFillStyle(fillStyle1sig) ;
       gr_nb2_m4 -> SetFillColor(fillColor1sig) ;
       gr_nb2_m4 -> SetLineColor(4) ;
       gr_nb2_m4 -> SetLineWidth(1) ;
       gr_nb2_m4 -> SetMarkerColor(4) ;
       //gr_nb2_m4 -> SetMarkerStyle(24) ;

       gr_nb2_m4_2s -> SetFillStyle(fillStyle2sig) ;
       gr_nb2_m4_2s -> SetFillColor(fillColor2sig) ;
       gr_nb2_m4_2s -> SetLineColor(4) ;
       gr_nb2_m4_2s -> SetMarkerColor(4) ;
       //gr_nb2_m4_2s -> SetMarkerStyle(24) ;

       gr_nb2_m4_2s -> Draw("p2") ;
       gr_nb2_m4_2s -> Draw("p") ;
       gr_nb2_m4 -> Draw("p2") ;
       gr_nb2_m4 -> Draw("p") ;
       h_data_nb2_m4 -> Draw("same e") ;



//.........这里部分代码省略.........

//.........这里部分代码省略.........
	    
      ttbar_Nom[h].hist[ch]->SetMaximum(MaxHisto*maxh);
      
      //ttbar_Nom[h].hist[ch]->SetMinimum(0.1);
      
      ttbar_Nom[h].hist[ch]->Draw("hist");
      ttbar_Up[h].hist[ch]->Draw("histSAME");
      ttbar_Down[h].hist[ch]->Draw("histSAME");
      

      //ttbar_1[h].hist[ch]->Draw("histoSAME");

      /***********************
             Legends
      ***********************/
      TLegend *leg;
      float legx1=0.70;
      float legy1=0.75;
      float legx2=0.93;
      float legy2=0.93;
      leg = new TLegend(legx1,legy1,legx2,legy2);
      leg->SetFillColor(0);
      leg->SetLineColor(1);
      leg->SetTextFont(62);
      leg->SetTextSize(0.03);

      leg->AddEntry((TObject*)0, cat + ": " + nSyst,"");
      leg->AddEntry((TObject*)0,"","");
      leg->AddEntry(ttbar_Nom[h].hist[ch],"Nom","l");
      leg->AddEntry(ttbar_Up[h].hist[ch],"Up","l");
      leg->AddEntry(ttbar_Down[h].hist[ch],"Down","l");
      //leg->AddEntry((TObject*)0,"","");
      leg->Draw("SAME");

      //-------------------------------------------------------
      // CMS Legend
      //-------------------------------------------------------
      TString htitleCMSChannel[3];
      htitleCMSChannel[0] = "#mu^{#pm}+jets channel";
      htitleCMSChannel[1] = "e^{#pm}+jets channel";
      htitleCMSChannel[2] = "l^{#pm}+jets channel";
      
      TLatex *titlePr;
      titlePr  = new TLatex(-20.,50.,"Preliminary");
      titlePr->SetNDC();
      titlePr->SetTextAlign(12);
      titlePr->SetX(0.25);
      titlePr->SetY(0.97);
      titlePr->SetTextColor(2);
      titlePr->SetTextFont(42);
      titlePr->SetTextSize(0.05);
      titlePr->SetTextSizePixels(24);
      titlePr->Draw("SAME");

      TLatex *title;
      title  = new TLatex(-20.,50.,"CMS #sqrt{s} = 13TeV, L = 36.5 fb^{-1}");
      title->SetNDC();
      title->SetTextAlign(12);
      title->SetX(0.20);
      title->SetY(0.80);
      title->SetTextFont(42);
      title->SetTextSize(0.03);
      title->SetTextSizePixels(24);
      title->Draw("SAME");
      
      TLatex *chtitle;
      chtitle  = new TLatex(-20.,50.,htitleCMSChannel[ch]);
      chtitle->SetNDC();
      chtitle->SetTextAlign(12);
      chtitle->SetX(0.20);
      chtitle->SetY(0.86);
      chtitle->SetTextFont(42);
      chtitle->SetTextSize(0.03);
      chtitle->SetTextSizePixels(24);
      chtitle->Draw("SAME");

            
      /***********************
            Save Histos
      ***********************/    
      TString dirfigname_log;
      if(LogScale) dirfigname_log = "_log";
      else dirfigname_log = "";
      TString dirfigname_pdf;
      TString dirfigname_png;
      dirfigname_pdf = "../" + dirnameIn + "figuresNomUpDown_" + fl + "/ttbbNomUpDown_" + cat + "/pdf" + dirfigname_log + "/";
      //dirfigname_png = dirnameIn + "figuresSystComp_" + fl + "/ttbb/png" + dirfigname_log + "/";
      // make a dir if it does not exist!!
      gSystem->mkdir(dirfigname_pdf,       kTRUE);
      histocanvas->SaveAs(dirfigname_pdf + nSyst + "_" + ttbar_Nom[h].hist[ch]->GetName() + ".pdf");
      //gSystem->mkdir(dirfigname_png,       kTRUE);
      //histocanvas->SaveAs(dirfigname_png + WJets[h].hist[ch]->GetName() + ".png");
      
      // clear Canvas
      histocanvas->Clear();    

    }
  }
  
} //end Plots.C

//.........这里部分代码省略.........
      printf("\n\n Simple average 0lep/1lep = %5.3f +/- %5.3f\n\n", simpleAveR_0over1, simpleAveR_0over1_err ) ;


      TH1F* hscalefactor_ttwj_0over1ratio_1b = (TH1F*) hmctruth_ttwj_0over1ratio_1b->Clone("hscalefactor_ttwj_0lep_1b") ;
      TH1F* hscalefactor_ttwj_0over1ratio_2b = (TH1F*) hmctruth_ttwj_0over1ratio_2b->Clone("hscalefactor_ttwj_0lep_2b") ;
      TH1F* hscalefactor_ttwj_0over1ratio_3b = (TH1F*) hmctruth_ttwj_0over1ratio_3b->Clone("hscalefactor_ttwj_0lep_3b") ;

      hscalefactor_ttwj_0over1ratio_1b->Scale(1./simpleAveR_0over1) ;
      hscalefactor_ttwj_0over1ratio_2b->Scale(1./simpleAveR_0over1) ;
      hscalefactor_ttwj_0over1ratio_3b->Scale(1./simpleAveR_0over1) ;


      hmctruth_ttwj_0over1ratio_1b->SetMinimum(0.) ;
      hmctruth_ttwj_0over1ratio_2b->SetMinimum(0.) ;
      hmctruth_ttwj_0over1ratio_3b->SetMinimum(0.) ;
      hmctruth_ttwj_0over1ratio_1b->SetMaximum(3.) ;
      hmctruth_ttwj_0over1ratio_2b->SetMaximum(3.) ;
      hmctruth_ttwj_0over1ratio_3b->SetMaximum(3.) ;

      hscalefactor_ttwj_0over1ratio_1b->SetMinimum(0.) ;
      hscalefactor_ttwj_0over1ratio_2b->SetMinimum(0.) ;
      hscalefactor_ttwj_0over1ratio_3b->SetMinimum(0.) ;
      hscalefactor_ttwj_0over1ratio_1b->SetMaximum(2.5) ;
      hscalefactor_ttwj_0over1ratio_2b->SetMaximum(2.5) ;
      hscalefactor_ttwj_0over1ratio_3b->SetMaximum(2.5) ;



      gStyle->SetPadRightMargin(0.08) ;
      TCanvas* cttwj = (TCanvas*) gDirectory->FindObject("cttwj") ;
      if ( cttwj == 0x0 ) {
         cttwj = new TCanvas("cttwj","ttwj closure", 700, 950) ;
      }
      cttwj->Clear() ;
      cttwj->Divide(1,2) ;

      TLegend* legend_ttwj = new TLegend( 0.91, 0.77,  0.99, 0.93 ) ;
      legend_ttwj->SetFillColor(kWhite) ;
      legend_ttwj->AddEntry( hmctruth_ttwj_0over1ratio_1b, "=1b") ;
      legend_ttwj->AddEntry( hmctruth_ttwj_0over1ratio_2b, "=2b") ;
      legend_ttwj->AddEntry( hmctruth_ttwj_0over1ratio_3b, ">=3b") ;

      cttwj->cd(1) ;
      hmctruth_ttwj_0over1ratio_1b->SetTitle("ttwj: 0 Lepton / 1 Lepton, Ratio") ;
      hmctruth_ttwj_0over1ratio_1b->Draw() ;
      hmctruth_ttwj_0over1ratio_2b->Draw("same") ;
      hmctruth_ttwj_0over1ratio_3b->Draw("same") ;
      legend_ttwj->Draw() ;

      cttwj->cd(2) ;
      hscalefactor_ttwj_0over1ratio_1b->SetTitle("ttwj: 0 Lepton / 1 Lepton, Scale Factor") ;
      hscalefactor_ttwj_0over1ratio_1b->Draw() ;
      hscalefactor_ttwj_0over1ratio_2b->Draw("same") ;
      hscalefactor_ttwj_0over1ratio_3b->Draw("same") ;
      line->DrawLine(0.5,1.,hscalefactor_ttwj_0over1ratio_1b->GetNbinsX()+0.5,1.) ;
      legend_ttwj->Draw() ;

      TString outttwj( infileStr ) ;
      outttwj.ReplaceAll("rootfiles","outputfiles") ;
      outttwj.ReplaceAll(".root","-mcclosure-ttwj1.pdf") ;
      cttwj->SaveAs( outttwj ) ;




   //----

void background_draw(const char *rootname = "background_plots2.root")
{
	char strs[128];
	char strl[1024];
	const char titles[NHISTS][16] = {"gtDiff", "R1", "R2", "RZ", "PX", "PY", "PZ", 
		"NX", "NY", "NZ", "NE", "NH", "PH", "AH", "AE", "AM", "AMO", "P2AZ", "AH1", "AE1", "P2AZ1"};
	const char suffix[4][10] = {"A-rand", "A-diff", "B-diff", "C-diff"};
	const Color_t color[4] = {kGreen+2, kBlue, kRed, kOrange};
	const int marker[4] = {kOpenCircle, kFullCircle, kOpenSquare, kOpenStar};
	TH1D *h[NHISTS][4];
	int i, j;
	double hMax;
	int iMax;
	char pdfname[1024];
	char *ptr;

	strcpy(pdfname, rootname);
	ptr = strstr(pdfname, ".root");
	if (ptr) {
		strcpy(ptr, ".pdf");
	} else {
		strcat(pdfname, ".pdf");
	}
	
	gROOT->SetStyle("Plain");
	gStyle->SetOptStat(0);
	gStyle->SetOptFit(1);
	
	gStyle->SetTitleXSize(0.05);
	gStyle->SetTitleYSize(0.05);
	gStyle->SetLabelSize(0.05);
	gStyle->SetPadLeftMargin(0.16);
	gStyle->SetPadBottomMargin(0.16);
	
	TFile *fRoot = new TFile(rootname);
	if (!fRoot->IsOpen()) {
		printf("root-file %s not found  - run background_calc() first!\n", rootname);
		return;
	}
	for (i=0; i<NHISTS; i++) {
		for (j=0; j<4; j++) {
			sprintf(strs, "h%s%s", titles[i], suffix[j]);
			h[i][j] = (TH1D*) fRoot->Get(strs);
			if (!h[i][j]) {
				printf("%s not found  - rerun background_calc() to create all hists!\n", strs);
				fRoot->Close();
				return;
			}
			h[i][j]->SetLineWidth(2);
			h[i][j]->SetLineColor(color[j]);
			h[i][j]->SetMarkerColor(color[j]);
			h[i][j]->SetMarkerStyle(marker[j]);
			h[i][j]->GetYaxis()->SetLabelSize(0.05);
			h[i][j]->SetMinimum(0);
			h[i][j]->GetYaxis()->SetTitle("");
		}
	}
	
	TCanvas *cv = new TCanvas("CV", "Background plots", 1200, 900);
	TLegend *lg = new TLegend(0.7, 0.8, 0.95, 0.95);
	lg->AddEntry(h[0][0], "Random", "LP");
	lg->AddEntry(h[0][1], "Neutrino", "LP");
	lg->AddEntry(h[0][2], "Cosmic-A", "LP");
	lg->AddEntry(h[0][3], "Cosmic-B", "LP");
	
	sprintf(strl, "%s[", pdfname);
	cv->SaveAs(strl);
	
	for (i=0; i<NHISTS; i++) {
		cv->Clear();

		hMax = 0;
		iMax = 0;
		for (j=0; j<3; j++) if (h[i][j]->GetMaximum() > hMax) {
			hMax = h[i][j]->GetMaximum();
			iMax = j;
		}
		h[i][iMax]->Draw();
		for (j=0; j<4; j++) if (j != iMax) h[i][j]->Draw("same");
		lg->Draw();
		cv->SaveAs(pdfname);
	}
	
	sprintf(strl, "%s]", pdfname);
	cv->SaveAs(strl);
	fRoot->Close();
}