C++ TGraphErrors::SetLineStyle方法代码示例

void make1DLimit(TString combine_dir,TString type= "WH",bool blind=true){
  //TString combine_dir = "test_runSusyHgg/signalInj_sms_ChiWH_0_175/";
  //WH
  wh_limits.push_back(wh_125);
  wh_limits.push_back(wh_150);
  wh_limits.push_back(wh_175);
  wh_limits.push_back(wh_200);
  //HH
  hh_limits.push_back(hh_125);
  hh_limits.push_back(hh_150);
  hh_limits.push_back(hh_175);
  hh_limits.push_back(hh_200);

  TGraph obser( (200-125)/25 );
  TGraph graph( (200-125)/25 );
  TGraphAsymmErrors error( (200-125)/25 );
  TGraphAsymmErrors error2S( (200-125)/25 );

  TGraph obser_r( (200-125)/25 );
  TGraph graph_r( (200-125)/25 );
  TGraphAsymmErrors error_r( (200-125)/25 );
  TGraphAsymmErrors error_r2S( (200-125)/25 );
  
  TGraphErrors* theo = 0;
  if(type=="WH") theo = getTheoXSec("xsecs/CharginoNeutralino.txt");
  else theo = getTheoXSec("xsecs/Higgsino_ElectroHiggs.txt");
  //else theo = getTheoXSec("/home/amott/HggApp/SusyHgg/xsecs/Higgsino.txt");

  for(int m=125;m<=200;m+=25) {
    int i=(m-125)/25;
    TFile limit_file(Form("%s/higgsCombineChi%s_0_%d.Asymptotic.mH120.root",combine_dir.Data(),type.Data(),m) );
    TTree *limit_tree = (TTree*)limit_file.Get("limit");
    TTreeFormula limit_form("get_limit","limit",limit_tree);
    
    float down_2s = -1;
    float down = -1;
    float exp = -1;
    float up = -1;
    float up_2s = -1;
    float obs = -1;

    if( type == "WH" )
      {
	down_2s = wh_limits.at(i)[0];
	down = wh_limits.at(i)[1];
	exp = wh_limits.at(i)[2];
	up = wh_limits.at(i)[3];
	up_2s = wh_limits.at(i)[4];
	obs = wh_limits.at(i)[5];
      }
    else if ( type == "HH")
      {
	down_2s = hh_limits.at(i)[0];
	down = hh_limits.at(i)[1];
	exp = hh_limits.at(i)[2];
	up = hh_limits.at(i)[3];
	up_2s = hh_limits.at(i)[4];
	obs = hh_limits.at(i)[5];
      }
    else
      {
	std::cerr << "UNRECOGNIZED OPTION!!! QUITTING" << std::endl;
      }
    
    if(i==0) m+=5; //first point is actually at m=130
    graph.SetPoint(i,float(m), exp);
    error.SetPoint(i,float(m), exp);
    error2S.SetPoint(i, float(m), exp);
    error.SetPointError(i, 0, 0, exp-down, up-exp);
    error2S.SetPointError(i, 0 , 0 , exp-down_2s, up_2s-exp);
  
    graph_r.SetPoint(i,float(m),exp/theo->Eval(m));
    error_r.SetPoint(i,float(m),exp/theo->Eval(m));
    error_r2S.SetPoint(i,float(m),exp/theo->Eval(m));
    error_r.SetPointError(i,0,0,(exp-down)/theo->Eval(m),(up-exp)/theo->Eval(m));
    error_r2S.SetPointError(i, 0, 0, (exp-down_2s)/theo->Eval(m), (up_2s-exp)/theo->Eval(m) );
    
  
    obser.SetPoint(i,float(m),obs);
    obser_r.SetPoint(i,float(m),obs/theo->Eval(m));
    if(i==0) m-=5;
  }


    TCanvas cv;
    cv.SetLogy();
    cv.SetGrid(1,1);
    theo->SetMaximum(1e2);
    theo->SetMinimum(1e-2);
    theo->GetYaxis()->SetLabelSize(0.05);
    theo->GetYaxis()->SetTitleSize(0.06);
    theo->GetYaxis()->SetTitleOffset(0.8);
    theo->GetYaxis()->SetTitle("95% CL #sigma upper limit (pb)");
    theo->GetXaxis()->SetTitle("m_{chargino} (GeV)");
    if(type=="HH") theo->GetXaxis()->SetTitle("m_{neutralino} (GeV)");
    theo->SetFillColor(kBlue);
    theo->SetLineStyle(kDotted);
    theo->SetLineWidth(2);
    
    error.SetMaximum(1e2);
//.........这里部分代码省略.........

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

  //cout << "at bottom contrib plot" << endl;
  TCanvas* c1 = new TCanvas("c1","Bottom Contribution",150,0,1150,1000);
  TGraphErrors *gr0     = new TGraphErrors(plotCount0-1,ptOFF1,Rb0,dx,eb0);
  TGraphErrors *gr2     = new TGraphErrors(plotCount2-1,ptOFF2,Rb2,dx,eb2);
  TGraphErrors *grC     = new TGraphErrors(plotCount-1,pT,RbC,dx,ebC);
  TGraphErrors *grF     = new TGraphErrors(l-1,xF,yF);
  TGraphErrors *grFmax  = new TGraphErrors(l-1,xF,maxF);
  TGraphErrors *grFmin  = new TGraphErrors(l-1,xF,minF);
  TGraphErrors *grP     = new TGraphErrors(p-1,xP,yP,0,dyP);
  TGraphErrors *grPr    = new TGraphErrors(2,pTP,RbP,0,EbP);
  TGraphErrors *grPPr   = new TGraphErrors(2,pTPP,RbPP,0,EbPP);


  c1->cd(1);

  grP->SetTitle("");
  grP->GetXaxis()->SetTitle("NPE p_{T} (GeV/c)");
  grP->GetYaxis()->SetTitle("r_{B}");
  gr0->SetMarkerStyle(20);
  gr0->SetMarkerSize(1);
  gr0->SetLineColor(kBlue);
  gr0->SetMarkerColor(kBlue);
  gr2->SetMarkerStyle(22);
  gr2->SetMarkerSize(1);
  gr2->SetLineColor(kRed);
  gr2->SetMarkerColor(kRed);
  grC->SetMarkerStyle(21);
  grC->SetMarkerSize(1);
  grC->SetLineColor(kRed);
  grC->SetMarkerColor(kRed);
  grP->GetXaxis()->SetLimits(0,10);
  grP->GetYaxis()->SetRangeUser(0,1);
  grF->SetLineStyle(1);
  grFmax->SetLineStyle(2);
  grFmin->SetLineStyle(2);
  grP->SetMarkerStyle(33);
  grP->SetMarkerColor(kBlack);
  grPr->SetMarkerStyle(29);
  grPr->SetMarkerColor(9);
  grPr->SetLineColor(9);
  grPPr->SetMarkerStyle(29);
  grPPr->SetMarkerColor(49);
  grPPr->SetLineColor(49);


  grP->Draw("AP");
  //grC->Draw("same P");
  gr2->Draw("same P");
  grF->Draw("same");
  grFmax->Draw("same");
  grFmin->Draw("same");
  gr0->Draw("same P");
 // grPr->Draw("same P");
  //grPPr->Draw("same P");

  TLegend* leg2 = new TLegend(0.15,0.68,0.48,0.85);
  leg2->AddEntry(gr0,"STAR Run 12 - Low p_{T} Analysis","pe");
  leg2->AddEntry(gr2,"STAR Run 12 - High p_{T} Analysis","pe");
  //leg2->AddEntry(grC,"Combined Trigs","pe");
  leg2->AddEntry(grP,"STAR Run 6 Analysis (Stat. Uncertainty)","pe");
//  leg2->AddEntry(grPr,"Run 12 Data, Run 5/6 Templates)","pe");
  //leg2->AddEntry(grPPr,"Run 5/6 Refit (prev Template)","pe");
  leg2->AddEntry(grF,"FONLL Calculation","l");
  leg2->Draw("same");

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 make1DLimitHH(TString combine_dir,bool blind=true){
  //TString combine_dir = "test_runSusyHgg/signalInj_sms_ChiHH_0_175/";

  TGraph obser( (275-125)/25 );
  TGraph graph( (275-125)/25 );
  TGraphAsymmErrors error( (275-125)/25 );

    for(int m=125;m<501;m+=25) {
      int i=(m-125)/25;
      TFile limit_file(Form("%s/higgsCombineChiHH_0_%d.Asymptotic.mH120.root",combine_dir.Data(),m) );
      TTree *limit_tree = (TTree*)limit_file.Get("limit");
      TTreeFormula limit_form("get_limit","limit",limit_tree);

      limit_tree->GetEntry(1);
      float down = limit_form.EvalInstance();
      limit_tree->GetEntry(2);
      float exp = limit_form.EvalInstance();
      limit_tree->GetEntry(3);
      float up = limit_form.EvalInstance();
      limit_tree->GetEntry(5);
      float obs = limit_form.EvalInstance();

      graph.SetPoint(i,float(m),exp);
      error.SetPoint(i,float(m),exp);
      error.SetPointError(i,0,0,exp-down,up-exp);

      obser.SetPoint(i,float(m),obs);
    }

    TGraphErrors* theo = getTheoXSec("/home/amott/HggApp/SusyHgg/xsecs/Higgsino.txt");

    TCanvas cv;
    cv.SetLogy();

    theo->SetMaximum(1e2);
    theo->SetMinimum(1e-2);
    theo->GetYaxis()->SetTitle("95% CL #sigma upper limit (pb)");
    theo->GetXaxis()->SetTitle("m_{chargino}");
    theo->SetFillColor(kBlue);
    theo->SetLineStyle(kDotted);
    theo->SetLineWidth(2.0);
    
    error.SetMaximum(1e2);
    error.SetMinimum(1e-2);
    error.GetYaxis()->SetTitle("95% CL #sigma upper limit (pb)");
    error.GetXaxis()->SetTitle("m_{chargino}");
    error.SetFillColor(kGreen);
    error.Draw("A3");

    theo->Draw("3C");

    graph.SetLineStyle(kDashed);
    graph.SetLineWidth(2);
    graph.Draw("C");

    obser.SetLineStyle(1);
    obser.SetLineWidth(2);
    if(!blind) obser.Draw("C");

    TLegend leg(0.7,0.7,0.85,0.85);
    leg.SetFillColor(0);
    leg.SetBorderSize(0);
    leg.AddEntry(&graph,"expected","l");
    leg.AddEntry(&error,"expected #pm1#sigma","F");
    leg.AddEntry(theo,"theoretical","f");
    if(!blind)     leg.AddEntry(&obser,"observed","l");

    leg.Draw("SAME");

    TLatex prelim(0.65,0.96,"CMS Preliminary");
    prelim.SetNDC();
    prelim.SetTextSize(0.045);
    prelim.Draw();

    TLatex lbl(0.5,0.86,"#sqrt{s} = 8 TeV  #int L dt = 19.78 fb^{-1}");
    lbl.SetNDC();
    lbl.SetTextSize(0.045);
    lbl.Draw();

    cv.SaveAs(combine_dir+"expected_exclusion_HH_1D.png");


}

//.........这里部分代码省略.........
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0};
   Double_t Graph3_fey1002[18] = {
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0,
   0};
   gre = new TGraphErrors(18,Graph3_fx1002,Graph3_fy1002,Graph3_fex1002,Graph3_fey1002);
   gre->SetName("Graph3");
   gre->SetTitle("Graph");
   gre->SetLineStyle(7);
   gre->SetLineWidth(3);
   gre->SetMarkerSize(0);
   
   TH1F *Graph_Graph1002 = new TH1F("Graph_Graph1002","Graph",100,9,5000);
   Graph_Graph1002->SetMinimum(1e-38);
   Graph_Graph1002->SetMaximum(5e-36);
   Graph_Graph1002->SetDirectory(0);
   Graph_Graph1002->SetStats(0);

   ci = TColor::GetColor("#000099");
   Graph_Graph1002->SetLineColor(ci);
   Graph_Graph1002->GetXaxis()->SetLabelFont(42);
   Graph_Graph1002->GetXaxis()->SetLabelSize(0.035);
   Graph_Graph1002->GetXaxis()->SetTitleSize(0.035);
   Graph_Graph1002->GetXaxis()->SetTitleFont(42);
   Graph_Graph1002->GetYaxis()->SetLabelFont(42);
   Graph_Graph1002->GetYaxis()->SetLabelSize(0.035);
   Graph_Graph1002->GetYaxis()->SetTitleSize(0.035);
   Graph_Graph1002->GetYaxis()->SetTitleFont(42);
   Graph_Graph1002->GetZaxis()->SetLabelFont(42);
   Graph_Graph1002->GetZaxis()->SetLabelSize(0.035);
   Graph_Graph1002->GetZaxis()->SetTitleSize(0.035);
   Graph_Graph1002->GetZaxis()->SetTitleFont(42);
   gre->SetHistogram(Graph_Graph1002);
   
   gre->Draw("pc");
   
   Double_t Graph4_fx1003[18] = {
   20,
   30,
   40,
   50,

//.........这里部分代码省略.........
					}
					std::string fileName = filepath + material[mat] + interaction[inter] + std::string("_QE_") + std::string(pdeChar) + std::string("_SPTR_80_ASIC_30_DT300_histos.root"); std::string(Form("%d", 1));
					std::cout << fileName << std::endl;
					TFile *fIn = new TFile(fileName.c_str(), "read");

					std::string baseName = "AVG.";
					std::string smearName = "";
					if(smear==1){
						smearName = "Smear";
					}

					//NPE
					for(int npe=0;npe<10;npe++){
						std::string histName = baseName + smearName + std::string("DTOF2") + std::string("N") + std::string(Form("%d", npes[npe]));
						std::cout << histName << std::endl;

						TH1F *h1 = (TH1F*) fIn->Get(histName.c_str());
						TF1* gauF1 = new TF1("gauF1","gaus",-200,200);
						h1->Fit("gauF1","","e",-200,200);
						h1->Draw();

						std::string path = "/home/jmbenlloch/next/petalo/work/histo/phys/5ps/jitter_avg/npe/";
						std::string histo = path + material[mat] + interaction[inter] + std::string("_QE_") + std::string(pdeChar) + std::string("_") + histName + std::string(".png");
						std::cout << histo << std::endl;

						c1->Print(histo.c_str());

						// FWHM
						sigmasNPE[mat][inter][smear][pde][npe] = 2.35 * h1->GetFunction("gauF1")->GetParameter(2);
						errorsNPE[mat][inter][smear][pde][npe] = 2.35 * h1->GetFunction("gauF1")->GetParError(2);
					}

					fIn->Close();
				}
			}
		}
	}
	c1->Close();

	double npesD[10] = {1,2,3,4,5,6,7,8,9,10};
	//EColor colors[5] = {kRed,kBlue,kGreen+3,kCyan};
	EColor colors[5] = {kRed,kBlack,kGreen+3,kBlue};

	for(int mat=0; mat<2;mat++){
		if (mat==0){
			int pdes[5] = {10,15,20,100};
		}else{
			int pdes[5] = {30,50,70,100};
		}

		for(int inter=1; inter<2;inter++){
			for(int smear=0; smear<1;smear++){
				std::string smearName = "";
				if(smear==1){
					smearName = "_smear";
				}
				//NPE
				//TCanvas *c2 = new TCanvas("c2","multipads",900,700);
				TCanvas *c2 = tcanvset("c");

				//TLegend *leg = new TLegend(0.7, 0.7, 0.9, 0.9);
				TLegend *leg = tlegset();
				leg->SetFillColor(0);

				for(int pde=0; pde<4;pde++){
					TGraphErrors *gLxe = new TGraphErrors(10, npesD, sigmasNPE[mat][inter][smear][pde], 0, errorsNPE[mat][inter][smear][pde]);
					gLxe->SetLineColor(colors[pde]);
					gLxe->SetLineWidth(2);
					gLxe->SetLineStyle(4-pde);
					grerrplot(*gLxe);

					leg->AddEntry(gLxe,Form("PDE %d%%", pdes[pde]), "lp");
					std::cout << "pde: " << pde << std::endl;

					if(pde==0){
						gLxe->SetTitle("");
						gLxe->GetXaxis()->SetTitle("Number of photoelectrons");
						gLxe->GetXaxis()->SetLimits(0,10);
						gLxe->GetYaxis()->SetTitle("CRT (ps)");
						gLxe->SetMinimum(0.);
						gLxe->SetMaximum(250.);
						gLxe->Draw();

					}else{
						gLxe->Draw("same");
					}
					if(pde==3){
						leg->Draw("same");
					}
				}
				std::string path = "/home/jmbenlloch/next/petalo/PESOA/PETALOS/root/paper/plots/";
				std::string histo = path + material[mat] + interaction[inter] + std::string("_avg_npe") + smearName + std::string(".png");

				c2->Print(histo.c_str());
				c2->Close();

			}
		}
	}
}

//.........这里部分代码省略.........
    p++;
  }
  fp1.close();

  //cout << "at bottom contrib plot" << endl;
  TCanvas* c1 = new TCanvas("c1","Bottom Contribution",150,0,1150,1000);
  TGraphErrors *gr0     = new TGraphErrors(plotCount0-1,ptOFF1,Rb0,dx,eb0);
  TGraphErrors *gr2     = new TGraphErrors(plotCount2-1,ptOFF2,Rb2,dx,eb2);
  TGraphErrors *grC     = new TGraphErrors(plotCount-1,pT,RbC,dx,ebC);
  TGraphErrors *grF     = new TGraphErrors(l-1,xF,yF);
  TGraphErrors *grFmax  = new TGraphErrors(l-1,xF,maxF);
  TGraphErrors *grFmin  = new TGraphErrors(l-1,xF,minF);
  TGraphErrors *grP     = new TGraphErrors(p-1,xP,yP,0,dyP);

  c1->cd(1);

  gr0->SetTitle("Bottom Contribution");
  gr0->GetXaxis()->SetTitle("p_{T,e}");
  gr0->GetYaxis()->SetTitle("#frac{r_{B}}{(r_{B}+r_{C})}");
  gr0->SetMarkerStyle(20);
  gr0->SetMarkerSize(1);
  gr0->SetLineColor(kBlue);
  gr0->SetMarkerColor(kBlue);
  gr2->SetMarkerStyle(22);
  gr2->SetMarkerSize(1);
  gr2->SetLineColor(kRed);
  gr2->SetMarkerColor(kRed);
  grC->SetMarkerStyle(21);
  grC->SetMarkerSize(1);
  grC->SetLineColor(kRed);
  grC->SetMarkerColor(kRed);
  gr0->GetXaxis()->SetLimits(1,14);
  gr0->GetYaxis()->SetRangeUser(0,1);
  grF->SetLineStyle(1);
  grFmax->SetLineStyle(2);
  grFmin->SetLineStyle(2);
  grP->SetMarkerStyle(33);
  grP->SetMarkerColor(kBlack);
  
  
  gr0->Draw("AP");
  // grC->Draw("same P");
  gr2->Draw("same P");
  grF->Draw("same");
  grFmax->Draw("same");
  grFmin->Draw("same");
  grP->Draw("same P");

  TLegend* leg2 = new TLegend(0.15,0.68,0.4,0.85);
  leg2->AddEntry(gr0,"High Tower 0 Trigs","pe");
  leg2->AddEntry(gr2,"High Tower 2 Trigs","pe");
  // leg2->AddEntry(grC,"Combined Trigs","pe");
  leg2->AddEntry(grF,"FONLL (Uncertainty: Scale Only)","l");
  leg2->AddEntry(grP,"Run 5/6 Analysis (Stat Uncertainty)","pe");
  leg2->Draw("same");
  
   // Make PDF with output canvases
  if(makePDF)
    {
      //Set front page
      TCanvas* fp = new TCanvas("fp","Front Page",100,0,1000,900);
      fp->cd();
      TBox *bLabel = new TBox(0.01, 0.88, 0.99, 0.99);
      bLabel->SetFillColor(38);
      bLabel->Draw();
      TLatex tl;

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

	(*tfyv).operator()(pp) = (*fyv)[k];
	(*tfyev).operator()(pp) = (*fyev)[k];
	(*tyv).operator()(pp) = (*yv)[k];
	(*tyev).operator()(pp) = (*yev)[k];

	// Fill the full ratio arrays
	(*ratio).ResizeTo(pp+1);
	(*eratio).ResizeTo(pp+1);

	(*ratio).operator()(pp)=((*fyv)[k]/al_fy_from_pass0((*yv)[k]));
	//std::cout<<ratio(k)<<std::endl;
	(*eratio).operator()(pp)=((*fyev)[k]/al_fy_from_pass0((*yv)[k]));
	std::cout<<" @ y = "<<(*yv)[k]<<" (ratio, eratio) = "<<"("<<(*ratio)[pp]<<","<<(*eratio)[pp]<<")"<<std::endl;


	pp++;

      }
    }
  }  
  


  TLegend * leg= new TLegend(0.2,0.6,0.4,0.8);

  myfunc();

  TGraphErrors *gp = new TGraphErrors(*tyv,*tfyv,*tyev,*tfyev);
  c->cd(1);
  gPad->SetGridy();
  gPad->SetGridx();
  gp->Draw("AEP");
  gp->SetMarkerStyle(4);
  gp->SetMarkerColor(kBlue);
  gp->SetLineStyle(1);
  gp->SetLineColor(kBlue);
  gp->Fit("myfunc","EBMR");
  TF1 *tfit = gp->GetFunction("myfunc");
  tfit->SetLineColor(kRed);


  gp->GetYaxis()->SetRangeUser(0,4);
  gp->GetXaxis()->SetRangeUser(-1.0,0.2);
  gp->GetXaxis()->SetTitle("y(GeV/c)");
  gp->GetYaxis()->SetTitle("F(y)(GeV/c)^{-1}");
  gp->SetTitle("");
  leg->AddEntry(gp,"data","p");
  leg->AddEntry(tfit,"F(y) from pass1","l");

  
  // From model:
  myfunc();
  TF1 *ft = new TF1("ft"," al_fy_from_pass0(x)",-1.0,0.1); 
  ft->SetLineColor(kBlack);
  ft->SetLineStyle(1);
 
  // draw the model
  ft->Draw("same");
  leg->AddEntry(ft,"F(y) model","l");
  leg->Draw();

  // Draw in log scale
  c->cd(2);
  gPad->SetGridy();
  gPad->SetGridx();
  gp->Draw("AEP");
  ft->Draw("same");
  gPad->Modified();
  gPad->Update();
  gPad->SetLogy();
  gp->GetYaxis()->SetRangeUser(1e-6,10);
  gp->GetXaxis()->SetRangeUser(-1.0,0.2);
  gp->SetTitle(Form("%s Al, Y scaling with DIS subtraction;y(GeV/c);F(y,q)(GeV/c)^{-1}",tgt_pos.Data()));
  leg->Draw();

 
  TCanvas * c1 = new TCanvas("c1", "",0,0,700,500);  
  c1->Draw();
  gPad->SetGridy();
  gPad->SetGridx();
  // Fill the ratio graphs
  TGraphErrors * grr = new TGraphErrors(*tyv,*ratio,*tyev,*eratio); 
  grr->SetMarkerStyle(2);
  grr->SetMarkerColor(2);
  grr->SetLineStyle(1);
  grr->SetLineColor(2);
  grr->Draw("AEP");
  //grr->GetYaxis()->SetRangeUser(-100);
  grr->GetXaxis()->SetRangeUser(-1.0,0.2);
  grr->GetXaxis()->SetTitle("y(GeV/c)");
  grr->GetYaxis()->SetTitle("data_{F(y)}/model_{F(y)}");


 //close the text file
  fclose(fp);

  // // // c->Print(Form("xsec_rad_exp_mod_ratio_x_geq_%2.1f_and_x_leq_%2.1f_%s_dummy_%s_%s.png",xlow,xup,tgt_pos.Data(),angle.Data(),fit.Data()));
 
}

//.........这里部分代码省略.........
  for(int i=0;i<Nbins; i++)
  {
  double temp = vy_significance[i] + vy_significance_erry[i]; 
  if(temp > highest_sig){highest_sig = temp; hightest_bin = i;}

  }
  cout<<endl;
  cout<< "Zprime_mass[k]: " << Zprime_mass[k] << endl;
  cout<<"hightest_bin : "<< hightest_bin << endl;  
  cout<<"highest_sig : "<<highest_sig << endl;
  cout<<"vx[hightest_bin] : "<< vx[hightest_bin] <<endl;
  cout<<"sig: "<< vy_significance[hightest_bin] <<" +/- "<< vy_significance_erry[hightest_bin] << endl;
  cout<<endl;

  // draw the TGraphError

  c1[k]->cd();

  TGraphErrors *gr 		  = new TGraphErrors(Nbins , vx , vy_significance, vy_significance_errx, vy_significance_erry );
  TGraphErrors *gr_eff 		  = new TGraphErrors(Nbins , vx , vy_sig_eff, vy_sig_eff_errx, vy_sig_eff_erry );
  TGraphErrors *gr_background_eff = new TGraphErrors(Nbins , vx , vy_background_eff, vy_background_eff_errx, vy_background_eff_erry );


  // draw setting
  gr->GetYaxis()->SetRangeUser(0, 1.1);
  gr->SetTitle(c1_title);
//  gr->GetXaxis()->SetTitle("Minimum p_{T} threshold of leading electron [ GeV]");
  gr->GetXaxis()->SetTitle("| delta eta |");

  gr->SetLineColor(2);
  gr_eff->SetLineColor(kBlue);
  gr_background_eff->SetLineColor(kGreen+2);

  gr->SetLineStyle(1);
  gr_eff->SetLineStyle(1);
  gr_background_eff->SetLineStyle(1);

  gr->SetLineWidth(3);
  gr_eff->SetLineWidth(3);
  gr_background_eff->SetLineWidth(3);

  gr->SetMarkerColor(2);
  gr_eff->SetMarkerColor(kBlue);
  gr_background_eff->SetMarkerColor(kGreen+2);

  gr->SetMarkerStyle(4);
  gr_eff->SetMarkerStyle(5);
  gr_background_eff->SetMarkerStyle(8);

  gr->SetMarkerSize(1);
  gr_eff->SetMarkerSize(1);
  gr_background_eff->SetMarkerSize(1);

  gr->Draw("ALP");
  gr_eff->Draw("same lp");
  gr_background_eff->Draw("same lp");

  leg->AddEntry( gr , leg_name1 ,"pl");
  leg->AddEntry( gr_eff , leg_name2 ,"pl");
  leg->AddEntry( gr_background_eff , leg_name3 ,"pl");
  leg->Draw();

  TString path_name = "/Users/index0192/13_TeV_Task/Zprime_Electron_channel/plot_SubLedaing_ele_pt_optimize/" ;
  c1_name = path_name + c1_name + ".png";
  c1[k]->SaveAs(c1_name);

//.........这里部分代码省略.........
   0.25,
   0.2};
   Double_t Graph0_fehy3021[10] = {
   5.988963,
   6.218694,
   6.307709,
   6.368316,
   5.847823,
   5.655323,
   5.083818,
   4.561256,
   3.722445,
   2.542957};
   TGraphAsymmErrors *grae = new TGraphAsymmErrors(10,Graph0_fx3021,Graph0_fy3021,Graph0_felx3021,Graph0_fehx3021,Graph0_fely3021,Graph0_fehy3021);
   grae->SetName("Graph0");
   grae->SetTitle("Graph");

   Int_t ci;      // for color index setting
   TColor *color; // for color definition with alpha
   ci = TColor::GetColor("#ffff00");
   grae->SetFillColor(ci);

   ci = TColor::GetColor("#ff0000");
   grae->SetLineColor(ci);
   grae->SetLineWidth(4);
   grae->SetMarkerStyle(20);
   grae->SetMarkerSize(0);
   
   TH1F *Graph_Graph3021 = new TH1F("Graph_Graph3021","Graph",100,-2.5,2.5);
   Graph_Graph3021->SetMinimum(0);
   Graph_Graph3021->SetMaximum(140);
   Graph_Graph3021->SetDirectory(0);
   Graph_Graph3021->SetStats(0);
   Graph_Graph3021->SetLineStyle(0);
   Graph_Graph3021->SetMarkerStyle(20);
   Graph_Graph3021->GetXaxis()->SetNdivisions(505);
   Graph_Graph3021->GetXaxis()->SetLabelFont(42);
   Graph_Graph3021->GetXaxis()->SetLabelOffset(0.007);
   Graph_Graph3021->GetXaxis()->SetLabelSize(0);
   Graph_Graph3021->GetXaxis()->SetTitleSize(0.07142857);
   Graph_Graph3021->GetXaxis()->SetTitleOffset(1.1);
   Graph_Graph3021->GetXaxis()->SetTitleFont(42);
   Graph_Graph3021->GetYaxis()->SetTitle("d#sigma (W^{+}#rightarrow#font[12]{l}^{+}#nu) / d#eta_{lab} [nb]");
   Graph_Graph3021->GetYaxis()->SetLabelFont(42);
   Graph_Graph3021->GetYaxis()->SetLabelOffset(0.007);
   Graph_Graph3021->GetYaxis()->SetLabelSize(0.07142857);
   Graph_Graph3021->GetYaxis()->SetTitleSize(0.07142857);
   Graph_Graph3021->GetYaxis()->SetTitleOffset(1.05);
   Graph_Graph3021->GetYaxis()->SetTitleFont(42);
   Graph_Graph3021->GetZaxis()->SetLabelFont(42);
   Graph_Graph3021->GetZaxis()->SetLabelOffset(0.007);
   Graph_Graph3021->GetZaxis()->SetLabelSize(0.05);
   Graph_Graph3021->GetZaxis()->SetTitleSize(0.06);
   Graph_Graph3021->GetZaxis()->SetTitleFont(42);
   grae->SetHistogram(Graph_Graph3021);
   
   grae->Draw("a2");
   
   Double_t Graph1_fx3022[10] = {
   2.2,
   1.75,
   1.25,
   0.75,
   0.25,
   -0.25,
   -0.75,

//.........这里部分代码省略.........
	  for(unsigned int i = 0; i < count[a]; ++i)
	  {
		  tree->GetEntry(i);
		  TAxis *xaxis = h->GetXaxis();
		  int bin1 = xaxis->FindBin(p1.true_E)-1;
		  int bin2 = xaxis->FindBin(p2.true_E)-1;

		  if(keep_photon(p1) && bin1 >= 0 && bin1<n)
		  {
		      elems[bin1]++;
		  }
		  
		  if(keep_photon(p2) && bin2 >= 0 && bin2<n)
		  {
		      elems[bin2]++;
		  }
	  }
   }
  
  for(int k = 0; k < n; ++k)
  {
      err[k] = new double[elems[k]];
      index[k] = 0;
  }
  
  
  for (unsigned a = 0; a < N; ++a)
  {
	  TTree *tree = trees[a];
	  for(unsigned int i = 0; i < count[a]; ++i)
	  {
		  tree->GetEntry(i);
		  TAxis *xaxis = h->GetXaxis();
		  int bin1 = xaxis->FindBin(p1.true_E)-1;
		  int bin2 = xaxis->FindBin(p2.true_E)-1;
		  
		  

		  if(keep_photon(p1) && bin1 >= 0 &&  bin1 < n)
		  {
		      err[bin1][index[bin1]++] = p1.E-p1.true_E;
		  }
		  
		  if(keep_photon(p2) && bin2 >= 0 &&  bin2 < n)
		  {
		      err[bin2][index[bin2]++] = p2.E-p2.true_E;
		  }
	  }
   }
  
  for(int k = 0; k < n; ++k)
  {
      double mean = 0;
      for(int j = 0; j < elems[k]; ++j)
      {
         mean += err[k][j];
      }
      mean /= double(elems[k]);
      
      double tot = 0;
      //mean = 0;
      for(int j = 0; j < elems[k]; ++j)
      {
         //if(k == 2) printf("err: %.3f \n", err[k][j]);
         tot += (err[k][j]-mean) * (err[k][j]-mean);
      }
      tot /= double(elems[k]);
      double E = float(E_min + float(k)/float(n) * (E_max-E_min));
      /*printf("%.2f %.6f %.6f %.6f %% %.7f %% %.7f %d\n", E, float(mean), float(sqrt(tot)), 100 * float(sqrt(tot))/float(E), 100 * float(sqrt(tot))/float(E)/sqrt(index[k]-1), float(1)/float(n) * (E_max-E_min)/2.0, index[k]);*/
      printf("%.2f %.6f %% %.7f %% %.7f %d\n", E, 100 * float(sqrt(tot))/float(E), 100 * float(sqrt(tot))/float(E)/sqrt(index[k]-1), float(1)/float(n) * (E_max-E_min)/2.0, index[k]);
      en[k] = E;
      den[k] = 100 * float(sqrt(tot))/float(E);
      yerr[k] = 100 * float(sqrt(tot))/float(E)/sqrt(index[k]-1);
      xerr[k] = float(1)/float(n) * (E_max-E_min)/2.0;
  }
  
   TF1 *energy = new TF1("energy", " [0]/sqrt(x) + [1]", 100, 1000);
   TGraphErrors *gr = new TGraphErrors(n,en,den,xerr,yerr);
   gr->Fit("energy");
   gr->SetTitle("\\mbox{Resolution en energie};E\\mbox{ (GeV) };\\sigma_{E} / E \\mbox{ \\% }");
   gr->SetLineStyle(0);
   gr->SetLineWidth(0);
   /* h->Sumw2();
  h2->Sumw2();
  
  h->Fit("f1");

  TH1F *r = (TH1F *)h->Clone();
  r->Divide(h, h2, 1., 1., "B");
  //r->Draw();
  h->Draw();*/
  gr->Draw("AP");
  

  printf("tight ratio: %.6f\n", float(keptEntries)/float(totalEntries));
    /*h->Multiply(fa);
    h->Draw("E");*/
 // h->Draw(); // plot the histogram
  //h2->Draw();
}

void PlotPhCal(const char *filename = "phCalibration_C0.dat", int col, int row)
{

  Init();

  int value, a, b, n;
  double x[10], y[10];
  char string[200];
  double aovergain, gain, par1;


  TF1 *pol2Fit = new TF1("pol2Fit", "pol2");
  pol2Fit->SetRange(50,500);
  TF1 *linFit = new TF1("linFit", "pol1");
  linFit->SetRange(150,600);
  TF1 *ursFit = new TF1("ursFit", Fitfcn, 40., 1500., 4);
  ursFit->SetParameter(0,0.00382);
  ursFit->SetParameter(1,0.886);
  ursFit->SetParameter(2,112.7);
  ursFit->SetParameter(3,113.0);


  TGraphErrors *graph = new TGraphErrors();
  TGraph *gSlope = new TGraph();

  FILE *inputFile = fopen(filename, "r");

  for (int i = 0; i < 4; i++) fgets(string, 200, inputFile);
  
  for (int icol = 0; icol < 52; ++icol)
    {
      for (int irow = 0; irow < 80; ++irow)
	{
	  n = 0;
	  for (int point = 0; point < 10; point++)
	    {
	      fscanf(inputFile, "%s", string);
	      
	      if (strcmp(string, "N/A") == 0)  value = 7777;
	      else 
		{
		  value  = atoi(string);
		}		      
	      if ((icol == col) && (irow == row))
		{

		  if (point == 0) x[n] = 50;
		  else if (point == 1) x[n] = 100;
		  else if (point == 2) x[n] = 150;
		  else if (point == 3) x[n] = 200;
		  else if (point == 4) x[n] = 250;
		  else if (point == 5) x[n] = 210;
		  else if (point == 6) x[n] = 350;
		  else if (point == 7) x[n] = 490;
		  else if (point == 8) x[n] = 560;
		  else if (point == 9) x[n] = 1400;
		  
		  if (value != 7777)
		    {
		      //  y[n] = value;  // in testboard adc range
		      y[n] = (value + 400) / 7;  // in 8 bit experiment adc range
		      graph->SetPoint(n, x[n], y[n]);	
		      //  graph->SetPointError(n, 2, 10);   // in testboard adc range
		      graph->SetPointError(n, 2, 10/6);   // in 8 bit experiment adc range

		      n++;
		    }
		}
	    }
	  fscanf(inputFile, "%s %2i %2i", string, &a, &b);  //comment
	}
    }
  graph->GetXaxis()->SetTitle("Vcal [low range DAC units]");
  graph->GetYaxis()->SetTitle("PH [ADC units]");
  graph->GetYaxis()->SetTitleOffset(1.2);


  graph->GetXaxis()->SetTitleSize(0.055);
  graph->GetYaxis()->SetTitleSize(0.055);
  graph->GetXaxis()->SetLabelSize(0.04);
  graph->GetYaxis()->SetLabelSize(0.05);
  graph->GetXaxis()->SetTitleOffset(1.15);
  graph->GetYaxis()->SetTitleOffset(1.2);
  graph->SetTitle("");
  graph->SetLineStyle(1);
  graph->SetLineWidth(2);

  //   graph->Fit("linFit","R");
//   gain = linFit->GetParameter(1);


  graph->Fit("ursFit","R");
  graph->Draw("A*");
  par1 = ursFit->GetParameter(1);

  
  NewLatex()->DrawLatex(0.55, 0.5, Form("p_{1} = %.1f ", par1));
  //  NewLatex()->DrawLatex(0.45, 0.4, "NON-LINEAR pixel");       // PlotPhCal("/home/l_tester/ptr/moduleDB/M0567-070115.08:23/T-10a/phCalibration_C12.dat",7,23)
  NewLatex()->DrawLatex(0.5, 0.4, "LINEAR pixel");                // PlotPhCal("/home/l_tester/ptr/moduleDB/M0493-061206.10:45/T-10a/phCalibration_C1.dat",7,23)
}

//.........这里部分代码省略.........
  TH2D* hist = new TH2D("h","h",100,m0min,m0max,100,120,m12max);
  hist->Draw();  
  hist->GetXaxis()->SetTitle("m_{0} (GeV/c^{2})");
  hist->GetYaxis()->SetTitle("m_{1/2} (GeV/c^{2})");
  hist->GetYaxis()->SetTitleOffset(1.);
  hist->GetXaxis()->SetNdivisions(506);
  //  if (tanBeta_ == 50)  hist->GetXaxis()->SetNdivisions(504);
  hist->GetYaxis()->SetNdivisions(506);

  //int col[]={2,3,4};

  
  //TFile *f = TFile::Open("exclusion_Spring11_CLs.root");        
  //TFile *f = TFile::Open("exclusion_Fall10_tcmet_JPT.root"); 
  //TFile *f = TFile::Open("exclusion_Fall10_pfmet_pfjets.root"); 
  //TFile *f = new TFile("exclusion_Fall10_pfmet_pfjets_CLs.root");
  
  //TH2F* h = (TH2F*) f->Get("hexcl_NLO_obs");
  //TH2F* h = (TH2F*) f->Get("hexcl_NLO_exp");
  //TH2F* h = (TH2F*) f->Get("hexcl_NLO_expp1");
  //TH2F* h = (TH2F*) f->Get("hexcl_NLO_expm1");

  //h->SetMaximum(3);
  //h->Draw("samecolz");
  
  TSpline3 *sFirst = new TSpline3("sFirst",First);
  sFirst->SetLineColor(kRed);
  sFirst->SetLineWidth(3);
  First->SetLineColor(kRed);
  First->SetLineWidth(3);

  TSpline3 *sSecond = new TSpline3("sSecond",Second);
  sSecond->SetLineColor(kBlue);
  sSecond->SetLineStyle(2);
  sSecond->SetLineWidth(3);
  Second->SetLineColor(kBlue);
  Second->SetLineStyle(2);
  Second->SetLineWidth(3);

  TSpline3 *sSecond_up = new TSpline3("sSecond_up",Second_up);
  sSecond_up->SetLineColor(kCyan);
  sSecond_up->SetLineStyle(1);
  sSecond_up->SetLineWidth(3);
  Second_up->SetLineColor(kBlue);
  //Second_up->SetLineColor(1);
  Second_up->SetLineWidth(2);

  TSpline3 *sSecond_low = new TSpline3("sSecond_low",Second_low);
  sSecond_low->SetLineColor(kCyan);
  sSecond_low->SetLineStyle(1);
  sSecond_low->SetLineWidth(3);
  Second_low->SetLineColor(kBlue);
  //Second_low->SetLineColor(1);
  Second_low->SetLineWidth(2);

  Third->SetLineColor(kCyan);
  Third->SetLineWidth(30);

  
  // TSpline3 *sThird = new TSpline3("sThird",Third);
  // sThird->SetLineColor(kGreen+2);
  // sThird->SetLineStyle(4);
  // sThird->SetLineWidth(3);
  // Third->SetLineColor(kGreen+2);
  // Third->SetLineStyle(4);
  // Third->SetLineWidth(3);