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

//.........这里部分代码省略.........
   29.36996,
   29.50732,
   29.64004,
   29.76811,
   29.89155,
   30.01034,
   30.1245,
   30.23401,
   30.33889,
   30.43912,
   30.53472,
   30.62567,
   30.71198,
   30.79366,
   30.87069,
   30.94308,
   31.01084,
   31.07395,
   31.13242,
   31.18625,
   31.23544,
   31.28,
   31.31991,
   31.35518,
   31.38581,
   31.4118,
   31.43315,
   31.43315};
   graph = new TGraph(102,bkg_fx2,bkg_fy2);
   graph->SetName("bkg");
   graph->SetTitle("Projection of totalPdf");
   graph->SetFillColor(1);
   graph->SetLineColor(2);
   graph->SetLineStyle(2);
   graph->SetLineWidth(3);
   graph->SetMarkerStyle(20);
   
   TH1F *Graph_bkg2 = new TH1F("Graph_bkg2","Projection of totalPdf",102,1080,1152);
   Graph_bkg2->SetMinimum(3.818479);
   Graph_bkg2->SetMaximum(33.94357);
   Graph_bkg2->SetDirectory(0);
   Graph_bkg2->SetStats(0);
   Graph_bkg2->SetLineWidth(2);
   Graph_bkg2->SetMarkerStyle(20);
   Graph_bkg2->GetXaxis()->SetNdivisions(505);
   Graph_bkg2->GetXaxis()->SetLabelFont(132);
   Graph_bkg2->GetXaxis()->SetLabelOffset(0.01);
   Graph_bkg2->GetXaxis()->SetLabelSize(0.06);
   Graph_bkg2->GetXaxis()->SetTitleSize(0.072);
   Graph_bkg2->GetXaxis()->SetTitleOffset(0.95);
   Graph_bkg2->GetXaxis()->SetTitleFont(132);
   Graph_bkg2->GetYaxis()->SetLabelFont(132);
   Graph_bkg2->GetYaxis()->SetLabelOffset(0.01);
   Graph_bkg2->GetYaxis()->SetLabelSize(0.06);
   Graph_bkg2->GetYaxis()->SetTitleSize(0.072);
   Graph_bkg2->GetYaxis()->SetTitleOffset(0.95);
   Graph_bkg2->GetYaxis()->SetTitleFont(132);
   Graph_bkg2->GetZaxis()->SetLabelFont(132);
   Graph_bkg2->GetZaxis()->SetLabelSize(0.06);
   Graph_bkg2->GetZaxis()->SetTitleSize(0.072);
   Graph_bkg2->GetZaxis()->SetTitleOffset(1.2);
   Graph_bkg2->GetZaxis()->SetTitleFont(132);
   graph->SetHistogram(Graph_bkg2);
   
   graph->Draw("l");
   

void PlotLimit ( string LimitFile , string filePrefix, string LimTitle , bool DoObsLim , bool DoExpLim ) {

    setTDRStyle();

    vector<float> vMass           ;
    vector<float> vObsLimit       ; 
    vector<float> vMeanExpLimit   ; 
    vector<float> vMedianExpLimit ; 
    vector<float> vExpLim68Down   ; 
    vector<float> vExpLim68Up     ; 
    vector<float> vExpLim95Down   ; 
    vector<float> vExpLim95Up     ;

    ifstream indata;
    indata.open(LimitFile.c_str());
    if(!indata) { // file couldn't be opened
        cerr << "Error: file could not be opened" << endl;
        return;
    }
    float Mass           ;
    float ObsLimit       ; 
    float MeanExpLimit   ; 
    float MedianExpLimit ; 
    float ExpLim68Down   ; 
    float ExpLim68Up     ; 
    float ExpLim95Down   ; 
    float ExpLim95Up     ;

    while ( indata >> Mass >> ObsLimit >> MeanExpLimit >> MedianExpLimit >> ExpLim95Down >> ExpLim68Down >> ExpLim68Up >> ExpLim95Up )  {
        cout << Mass << " " << MeanExpLimit  << " " << MedianExpLimit <<" "<< ExpLim68Down <<" "<< ExpLim68Up <<" "<< ExpLim95Down <<" "<< ExpLim95Up << endl;

        vMass           .push_back(Mass           );
        vObsLimit       .push_back(ObsLimit       ); 
        vMeanExpLimit   .push_back(MeanExpLimit   ); 
        vMedianExpLimit .push_back(MedianExpLimit ); 
        vExpLim68Down   .push_back(ExpLim68Down   ); 
        vExpLim68Up     .push_back(ExpLim68Up     ); 
        vExpLim95Down   .push_back(ExpLim95Down   ); 
        vExpLim95Up     .push_back(ExpLim95Up     );
    }

    TCanvas* cLimit = new TCanvas("c1","c1",900,600);
    cLimit->cd();

    float x1 = vMass.at(0) - 5. ;
    float x2 = vMass.at(vMass.size()-1) + 5. ; 

    // Expected Limit
    TGraph* ExpLim = NULL ;
    TGraphAsymmErrors* ExpBand68 = NULL ;
    TGraphAsymmErrors* ExpBand95 = NULL ;
    float min = 999999., max = 0;
    float minZoom = 999999., maxZoom = 0;
    if ( DoExpLim ) {
        float x[100];
        float ex[100];
        float y[100];
        float yu68[100];
        float yd68[100];
        float yu95[100];
        float yd95[100]; 
        for ( int i = 0 ; i < (signed) vMass.size() ; ++i ) {
            x[i] = vMass.at(i) ; ex[i] = 0 ; 
            y[i] = vMedianExpLimit.at(i) ;   
	    if(y[i]    > max) max = y[i]   ; if(y[i]    < min) min = y[i]   ;
	    if(x[i]<=250){ if(y[i]    > maxZoom) maxZoom = y[i]   ; if(y[i]    < minZoom) minZoom = y[i]   ;}
            yu68[i] = vExpLim68Up.at(i)   -y[i];   if(yu68[i] > max) max = yu68[i]; if(yu68[i] < min) min = yu68[i];
            yd68[i] = y[i] - vExpLim68Down.at(i);  if(yd68[i] > max) max = yd68[i]; if(yd68[i] < min) min = yd68[i];
            yu95[i] = vExpLim95Up.at(i)   -y[i];   if(yu95[i] > max) max = yu95[i]; if(yu95[i] < min) min = yu95[i];
            yd95[i] = y[i] - vExpLim95Down.at(i);  if(yd95[i] > max) max = yd95[i]; if(yd95[i] < min) min = yd95[i];
        }
        ExpBand95 = new TGraphAsymmErrors((signed) vMass.size(),x,y,ex,ex,yd95,yu95);
        ExpBand95->SetFillColor(kYellow); 
        ExpBand95->GetYaxis()->SetRangeUser(0.,50);
        ExpBand95->GetXaxis()->SetRangeUser(x1,x2);
        ExpBand95->GetXaxis()->SetTitle("Higgs mass [GeV/c^{2}]");
        ExpBand95->GetYaxis()->SetTitle("95% Limit on #sigma/#sigma_{SM} ");
        ExpBand95->Draw("A3");
        ExpBand95->GetYaxis()->SetRangeUser(0.,50);
        ExpBand95->Draw("A3");
        ExpBand68 = new TGraphAsymmErrors((signed) vMass.size(),x,y,ex,ex,yd68,yu68);
        ExpBand68->SetFillColor(kGreen); 
        ExpBand68->Draw("3");

        ExpLim = new TGraph((signed) vMass.size(),x,y);    
        ExpLim->SetLineWidth(2);
        ExpLim->SetLineStyle(1);
        ExpLim->Draw("l");
    }

    // Observed Limit
    TGraph* ObsLim = NULL ;
    if ( DoObsLim ) {
        float x[100];
        float y[100];    
        for ( int i = 0 ; i < (signed) vMass.size() ; ++i ) { 
            x[i] = vMass.at(i) ; 
            y[i] = vObsLimit.at(i) ; 
	    if(y[i] > max) max = y[i]; if(y[i] < min) min = y[i];
	    if(x[i]<=250){ if(y[i]    > maxZoom) maxZoom = y[i]   ; if(y[i]    < minZoom) minZoom = y[i]   ;}
//.........这里部分代码省略.........

//Draw an error envelop from a TGraphAsymmErrors
void drawGraphEnv(TGraphAsymmErrors* g, float scale = 1, bool doAxis = false,
                  int linestyle = 1, int linewidth = 2, int linecolor = 4,
                  bool fill = false,
                  TH1** ppGl = 0, TH1** ppGh = 0){
  //                  const char* opth = "", const char* optl = ""){
  char buffer[256];
  int n = g->GetN();
  double* x = new double[n+2];
  double* yh = new double[n+2];
  double* yl = new double[n+2];

  for(int i = 0; i < n; ++i){
    x[i+1]  = g->GetX()[i];
    yh[i+1] = g->GetY()[i] + g->GetEYhigh()[i] * scale;
    yl[i+1] = g->GetY()[i] - g->GetEYlow()[i] * scale;
  }

  x[0] = g->GetX()[0] - g->GetEXlow()[0];
  yl[0] = yl[1] - (yl[2]-yl[1])/(x[2]-x[1])*g->GetEXlow()[0];
  yh[0] = yh[1] - (yh[2]-yh[1])/(x[2]-x[1])*g->GetEXlow()[0];

  x[n+1] = g->GetX()[n-1] + g->GetEXhigh()[n-1];
  yl[n+1] = yl[n] + (yl[n]-yl[n-1])/(x[n]-x[n-1])*g->GetEXhigh()[n-1];
  yh[n+1] = yh[n] + (yh[n]-yh[n-1])/(x[n]-x[n-1])*g->GetEXhigh()[n-1];

  TGraph* gl = new TGraph(n+2, x, yl);
  gl->SetTitle(g->GetTitle());
  TGraph* gh = new TGraph(n+2, x, yh);
  gh->SetTitle(g->GetTitle());
  const char* opt = "l,same";
  
  sprintf(buffer, "%s_eh", g->GetName());
  gh->SetName(buffer);
  
  sprintf(buffer, "%s_el", g->GetName());
  gl->SetName(buffer);

  if(doAxis){
    g->Draw("AP");
    TH1* h = g->GetHistogram();
    gPad->Clear();
    h->Draw();
  }

  if(fill){
    gl->SetFillStyle(linestyle);
    gl->SetFillColor(linecolor);
    gl->SetLineColor(linecolor);
  
    gh->SetFillStyle(linestyle);
    gh->SetFillColor(linecolor);
    gh->SetLineColor(linecolor);
  } else{
    gl->SetLineStyle(linestyle);
    gl->SetLineWidth(linewidth);
    gl->SetLineColor(linecolor);
    
    gh->SetLineStyle(linestyle);
    gh->SetLineWidth(linewidth);
    gh->SetLineColor(linecolor);
  }
  
  gh->Draw(opt);
  gl->Draw(opt);

  if(ppGh) *ppGh = gh;
  if(ppGl) *ppGl = gl;
}

//.........这里部分代码省略.........
   92.26923,
   57.33169,
   37.75782,
   26.04724,
   13.78295,
   8.096796,
   3.441144,
   1.753948,
   1.006317,
   0.6274052,
   0.415992,
   0.2891995,
   0.2087556,
   0.1553637,
   0.118593,
   0.0924765,
   0.07343641,
   0.05924083,
   0.04844944,
   0.04010506,
   0.03355513,
   0.02834456,
   0.02414957,
   0.02073562,
   0.01793004,
   0.01560389,
   0.01365953,
   0.01365953};
   graph = new TGraph(128,sigL_fx2,sigL_fy2);
   graph->SetName("sigL");
   graph->SetTitle("Projection of totalPdf");
   graph->SetFillColor(1);
   graph->SetLineColor(2);
   graph->SetLineStyle(2);
   graph->SetLineWidth(3);
   graph->SetMarkerStyle(20);
   
   TH1F *Graph_sigL2 = new TH1F("Graph_sigL2","Projection of totalPdf",128,4126,6214);
   Graph_sigL2->SetMinimum(1.186233e-20);
   Graph_sigL2->SetMaximum(3895.701);
   Graph_sigL2->SetDirectory(0);
   Graph_sigL2->SetStats(0);
   Graph_sigL2->SetLineWidth(2);
   Graph_sigL2->SetMarkerStyle(20);
   Graph_sigL2->GetXaxis()->SetNdivisions(505);
   Graph_sigL2->GetXaxis()->SetLabelFont(132);
   Graph_sigL2->GetXaxis()->SetLabelOffset(0.01);
   Graph_sigL2->GetXaxis()->SetLabelSize(0.06);
   Graph_sigL2->GetXaxis()->SetTitleSize(0.072);
   Graph_sigL2->GetXaxis()->SetTitleOffset(0.95);
   Graph_sigL2->GetXaxis()->SetTitleFont(132);
   Graph_sigL2->GetYaxis()->SetLabelFont(132);
   Graph_sigL2->GetYaxis()->SetLabelOffset(0.01);
   Graph_sigL2->GetYaxis()->SetLabelSize(0.06);
   Graph_sigL2->GetYaxis()->SetTitleSize(0.072);
   Graph_sigL2->GetYaxis()->SetTitleOffset(0.95);
   Graph_sigL2->GetYaxis()->SetTitleFont(132);
   Graph_sigL2->GetZaxis()->SetLabelFont(132);
   Graph_sigL2->GetZaxis()->SetLabelSize(0.06);
   Graph_sigL2->GetZaxis()->SetTitleSize(0.072);
   Graph_sigL2->GetZaxis()->SetTitleOffset(1.2);
   Graph_sigL2->GetZaxis()->SetTitleFont(132);
   graph->SetHistogram(Graph_sigL2);
   
   graph->Draw("l");
   

//.........这里部分代码省略.........
                }
                new(graphsV2D[nGraphs]) TGraph(*gr) ;
                nGraphs++;
            }
        }
    }

    // Ion probability
    hIonProb2D->GetZaxis()->SetRangeUser(0.00501,80);

    TH2F *hIonProb2Dc = (TH2F*) hIonProb2D->Clone("hIonProb2Dc");
    const Int_t NcontI = 4;
    Double_t contI[NcontI] = {0.01,0.1,1.0,10.0};
    hIonProb2Dc->SetContour(NcontI, contI);
    hIonProb2Dc->Draw("cont list");

    c->Update();
    TObjArray *contsI2D = (TObjArray*) gROOT->GetListOfSpecials()->FindObject("contours");
    TClonesArray graphsI2D("TGraph",NcontI);
    {
        Int_t ncontours = contsI2D->GetSize();
        TList* clist = NULL;
        Int_t nGraphs = 0;
        TGraph *gr = NULL;
        for(Int_t i = 0; i < ncontours; i++) {
            clist = (TList*) contsI2D->At(i);

            for(Int_t j = 0 ; j < clist->GetSize(); j++) {
                gr = (TGraph*) clist->At(j);
                if(!gr) continue;

                if( !(i%2) ) {
                    gr->SetLineWidth(1);
                    gr->SetLineStyle(2);
                    gr->SetLineColor(kOrange-3);
                } else {
                    gr->SetLineWidth(1);
                    gr->SetLineStyle(1);
                    gr->SetLineColor(kOrange-3);
                }

                new(graphsI2D[nGraphs]) TGraph(*gr) ;
                nGraphs++;
            }
        }
    }




    // "Axis range" in Osiris units:
    Double_t ylow  = hDen2D[0]->GetYaxis()->GetBinLowEdge(FirstyBin);
    Double_t yup = hDen2D[0]->GetYaxis()->GetBinUpEdge(LastyBin);
    Double_t xmin = hDen2D[0]->GetXaxis()->GetXmin();
    Double_t xmax = hDen2D[0]->GetXaxis()->GetXmax();

    TLine *lineYzero = new TLine(xmin,0.0,xmax,0.0);
    lineYzero->SetLineColor(kGray+2);
    lineYzero->SetLineStyle(2);

    TLine *lineYup = new TLine(xmin,yup,xmax,yup);
    lineYup->SetLineColor(kGray+1);
    lineYup->SetLineStyle(2);

    TLine *lineYdown = new TLine(xmin,ylow,xmax,ylow);
    lineYdown->SetLineColor(kGray+1);

//.........这里部分代码省略.........
ofstream out_efficiency;
ofstream out_discrimination;
out_efficiency.open(trigger[set_type]+dir_name+"efficiency.txt"); 


Float_t MC_data[2] = {1.,1.};////////////for efficiency 
//Float_t MC_data[2] = {1./0.708,1./.9};////////////TO BE CHANGED
//Float_t MC_data[2] = {1./0.704,1./.899};////////////for trigger logic
//Float_t MC_data[2] = {1./0.703,1./.964};////////////for blike single btag logic
Float_t lumi_qcd=lumi/MC_data[set_type];
do{
	
   TFile *file_initial = new TFile(file_names[files]);
//	file_initial->ls();   
	TH1F *histos[100];
	for (int hist=0;hist<nhistos;++hist){
		if (region_type==2) ((TH1F*)file_initial->Get(hist_names[hist]))->Rebin(2);
		histos[hist] = (TH1F*)file_initial->Get(hist_names[hist])->Clone("h");
		if (files==0) data_histos[hist] = (TH1F*)file_initial->Get(hist_names[hist])->Clone("data");
		if ((files>0)&&(files<7))histos[hist]->Scale(lumi);
		if (files>=7)histos[hist]->Scale(lumi_qcd);
		
		if (files==1) 	sum_histos[hist] = (TH1F*)histos[hist]->Clone(sum_histos_names[hist]);
		if (files>1)	sum_histos[hist]->Add(histos[hist]); 

		if (files>0) histos[hist]->Sumw2(kFALSE);
		if (hist==4) cout<<files<<"   "<<histos[4]->Integral() <<endl;

		if (files==1) {
			signal_histos[hist] = (TH1F*)file_initial->Get(hist_names[hist])->Clone(hist_names_sum[hist]+"newhist");
			signal_histos[hist]->Scale(lumi);
			signal_histos[hist]->Sumw2(kFALSE);
			signal_histos[hist]->SetLineColor(LINECOLOR[files]);
			signal_histos[hist]->SetLineStyle(LINESTYLE[files]);
			signal_histos2[hist]=(TH1F*)signal_histos[hist]->Clone("signalHist2");
		}
		if (files==2) {
			tthbb_histos[hist] = (TH1F*)file_initial->Get(hist_names[hist])->Clone(hist_names_sum[hist]+"tthb");
			tthbb_histos[hist]->Scale(lumi);
			tthbb_histos[hist]->Sumw2(kFALSE);
			tthbb_histos[hist]->SetLineColor(LINECOLOR[files]);
			tthbb_histos[hist]->SetLineStyle(LINESTYLE[files]);
		}
		if (files==3) {
			tthnbb_histos[hist] = (TH1F*)file_initial->Get(hist_names[hist])->Clone(hist_names_sum[hist]+"tthnb");
			tthnbb_histos[hist]->Scale(lumi);
			tthnbb_histos[hist]->Sumw2(kFALSE);
			tthnbb_histos[hist]->SetLineColor(LINECOLOR[files]);
			tthnbb_histos[hist]->SetLineStyle(LINESTYLE[files]);
		}
		histos[hist]->SetLineColor(LINECOLOR[files]);
		histos[hist]->SetLineStyle(LINESTYLE[files]);
		histos[hist]->SetLineWidth(LINEWIDTH[files]);
		histos[hist]->SetFillStyle(FILLSTYLE[files]);
		if ((files!=0)) histos[hist]->SetFillColor(FILLCOLOR[files]);
		
		if (files==0) {
			histos[hist]->SetMarkerStyle(20);
			data_histos[hist]->SetLineColor(1);
			data_histos[hist]->SetMarkerStyle(20);
			data_histos[hist]->SetMarkerSize(.8);
		}
	 	if (files>=4) stacks[hist]->Add(histos[hist]);
		if (hist==0) histos_for_legened[files] = (TH1F*)histos[0]->Clone("newd");
		if (files==4)	discr_histos[hist] = (TH1F*)file_initial->Get(hist_names[hist])->Clone("discr");
		if (files>4)	discr_histos[hist]->Add(histos[hist]); 

void showGraphs(double canvasSizeX, double canvasSizeY,
		const std::vector<std::string>& entriesToPlot,
		std::map<std::string, TGraph*>& graphs,
		std::map<std::string, std::string>& legendEntries,
		int colors[], int lineWidths[], int lineStyles[],
		double legendTextSize, double legendPosX, double legendPosY, double legendSizeX, double legendSizeY, 
		double xMin, double xMax, const std::string& xAxisTitle, double xAxisOffset,
		double yMin, double yMax, const std::string& yAxisTitle, double yAxisOffset,
		const std::string& outputFileName)
{
  TCanvas* canvas = new TCanvas("canvas", "canvas", canvasSizeX, canvasSizeY);
  canvas->SetFillColor(10);
  canvas->SetBorderSize(2);
  
  canvas->SetLeftMargin(0.14);
  canvas->SetBottomMargin(0.12);

  TH1* dummyHistogram = new TH1D("dummyHistogram", "dummyHistogram", 100, xMin, xMax);
  dummyHistogram->SetTitle("");
  dummyHistogram->SetStats(false);
  dummyHistogram->SetMinimum(yMin);
  dummyHistogram->SetMaximum(yMax);

  TAxis* xAxis = dummyHistogram->GetXaxis();
  xAxis->SetTitle(xAxisTitle.data());
  xAxis->SetTitleSize(0.045);
  xAxis->SetTitleOffset(xAxisOffset);

  TAxis* yAxis = dummyHistogram->GetYaxis();
  yAxis->SetTitle(yAxisTitle.data());
  yAxis->SetTitleSize(0.045);
  yAxis->SetTitleOffset(yAxisOffset);

  dummyHistogram->Draw("axis");

  TLegend* legend = new TLegend(legendPosX, legendPosY, legendPosX + legendSizeX, legendPosY + legendSizeY, "", "brNDC"); 
  legend->SetBorderSize(0);
  legend->SetFillColor(0);
  legend->SetTextSize(legendTextSize);

  size_t numEntriesToPlot = entriesToPlot.size();
  for ( size_t iEntryToPlot = 0; iEntryToPlot < numEntriesToPlot; ++iEntryToPlot ) {
    const std::string& entryToPlot = entriesToPlot[iEntryToPlot];
    TGraph* graph = graphs[entryToPlot];
    assert(graph);
    graph->SetLineColor(colors[iEntryToPlot]);
    graph->SetLineWidth(lineWidths[iEntryToPlot]);
    graph->SetLineStyle(lineStyles[iEntryToPlot]);
    graph->SetMarkerColor(colors[iEntryToPlot]);
    graph->Draw("L");

    legend->AddEntry(graph, legendEntries[entryToPlot].data(), "l");
  }

  legend->Draw();

  canvas->Update();
  std::string outputFileName_plot = "plots/";
  size_t idx = outputFileName.find_last_of('.');
  outputFileName_plot.append(std::string(outputFileName, 0, idx));
  if ( idx != std::string::npos ) canvas->Print(std::string(outputFileName_plot).append(std::string(outputFileName, idx)).data());
  canvas->Print(std::string(outputFileName_plot).append(".png").data());
  canvas->Print(std::string(outputFileName_plot).append(".pdf").data());
  
  delete dummyHistogram;
  delete legend;
  delete canvas;  
}

//.........这里部分代码省略.........
       gMATUp->SetPoint(i, e/GeV, (e-eVarMATUp)/eVarMATUp); 
       gMATDown->SetPoint(i, e/GeV, (e-eVarMATDown)/eVarMATDown); 
       
       gPSUp->SetPoint(i, e/GeV, (e-eVarPSUp)/eVarPSUp); 
       gPSDown->SetPoint(i, e/GeV, (e-eVarPSDown)/eVarPSDown); 
       
       gLowUp->SetPoint(i, e/GeV, (e-eVarLowUp)/eVarLowUp); 
       gLowDown->SetPoint(i, e/GeV, (e-eVarLowDown)/eVarLowDown); 
       
     } else {
       
       gNom->SetPoint(i, e/GeV, (e-eNominal)/eNominal - alpha); 
     
       gZeeAllUp->SetPoint(i, e/GeV, (e-eVarZeeAllUp)/eVarZeeAllUp - alpha); 
       gZeeAllDown->SetPoint(i, e/GeV, (e-eVarZeeAllDown)/eVarZeeAllDown - alpha); 
       
       gZeeUp->SetPoint(i, e/GeV, (e-eVarZeeUp)/eVarZeeUp - alpha); 
       gZeeDown->SetPoint(i, e/GeV, (e-eVarZeeDown)/eVarZeeDown - alpha); 
       
       gZeeGenUp->SetPoint(i, e/GeV, (e-eVarGenUp)/eVarGenUp - alpha); 
       gZeeGenDown->SetPoint(i, e/GeV, (e-eVarGenDown)/eVarGenDown - alpha); 
       
       gMATUp->SetPoint(i, e/GeV, (e-eVarMATUp)/eVarMATUp - alpha); 
       gMATDown->SetPoint(i, e/GeV, (e-eVarMATDown)/eVarMATDown - alpha); 
       
       gPSUp->SetPoint(i, e/GeV, (e-eVarPSUp)/eVarPSUp - alpha); 
       gPSDown->SetPoint(i, e/GeV, (e-eVarPSDown)/eVarPSDown - alpha); 
       
       gLowUp->SetPoint(i, e/GeV, (e-eVarLowUp)/eVarLowUp - alpha); 
       gLowDown->SetPoint(i, e/GeV, (e-eVarLowDown)/eVarLowDown - alpha); 
       
     }

     i++;

     if( logx )
       e *= exp(eStep);
     else
       e += eStep;

   }
 
   TCanvas* canvas = new TCanvas();
   if( logx )
     canvas->SetLogx();

   if( !subtractAverage ) {
     gMATUp->SetMaximum(alpha+.005);
     gMATUp->SetMinimum(alpha+-.007);
   } else {
     gMATUp->SetMaximum(.005);
     gMATUp->SetMinimum(-.007);
   }

   char grafname[99];
   if( ptype==0 )
     sprintf(grafname,"ES  Linearity (%s), #eta = %4.2f, Electrons", year.c_str(), eta);
   else if( ptype==1 )
     sprintf(grafname,"ES  Linearity (%s), #eta = %4.2f, Unconverted photons", year.c_str(), eta);
   else if( ptype==2 )
     sprintf(grafname,"ES  Linearity (%s), #eta = %4.2f, Converted photons", year.c_str(), eta);
   gMATUp->SetTitle(grafname);
   
   gMATUp->GetXaxis()->SetTitleOffset(1.2*gNom->GetXaxis()->GetTitleOffset());
   gMATUp->GetYaxis()->SetTitleOffset(1.2*gNom->GetYaxis()->GetTitleOffset());
   gMATUp->GetXaxis()->SetTitle("E [GeV]");
   gMATUp->GetYaxis()->SetTitle("#alpha(E)");
   gMATUp->SetLineWidth(2);
   gMATUp->SetLineStyle(2);
   gMATUp->SetLineColor(4);
   gMATUp->Draw("AL");

   gMATDown->SetLineWidth(2);
   gMATDown->SetLineStyle(3);
   gMATDown->SetLineColor(4);
   gMATDown->Draw("L");

   gNom->SetLineWidth(2);
   gNom->SetLineColor(4);
   gNom->Draw("L");

   gZee->SetMarkerStyle(20);
   gZee->Draw("P");

   TLegend* leg1 = new TLegend(.7,.15,.89,.3);
   leg1->SetBorderSize(0);
   leg1->SetFillColor(0);
   leg1->AddEntry(gZee,"#alpha, Z peak", "p");
   leg1->AddEntry(gNom,"PS, central", "l");
   leg1->AddEntry(gMATUp,"PS, up", "l");
   leg1->AddEntry(gMATDown,"PS, down", "l");

   leg1->Draw("same");

   // char plotname[99];
   // sprintf(plotname,"linearity_%d_%2.1f.png",ptype,eta);
   //   canvas->Print("overlay.ps");
  
   return canvas;
}

void Polarization(){

    TGraph *gObs;
    TGraph *gExp;
    TGraph *gObsL;
    TGraph *gExpL;
    TGraph *gObsR;
    TGraph *gExpR;
    
    TFile *fc = TFile::Open("Limits2DHistograms_T2tt_postfit.root");
    TFile *fl = TFile::Open("Limits2DHistograms_T2tt_lefthanded_postfit.root");
    TFile *fr = TFile::Open("Limits2DHistograms_T2tt_righthanded_postfit.root");
    
    TGraph *g;
    fc->cd();
    g = (TGraph*)fc->Get("gObs");
    gObs = (TGraph*)g->Clone("Obs");
    g = (TGraph*)fc->Get("gExp");
    gExp = (TGraph*)g->Clone("Exp");
    fl->cd();
    g = (TGraph*)fl->Get("gObs_2");
    gObsL = (TGraph*)g->Clone("ObsL");
    g = (TGraph*)fl->Get("gExp_2");
    gExpL = (TGraph*)g->Clone("ExpL");
    fr->cd();
    g = (TGraph*)fr->Get("gObs");
    gObsR = (TGraph*)g->Clone("ObsR");
    g = (TGraph*)fr->Get("gExp");
    gExpR = (TGraph*)g->Clone("ExpR");
    
    gObs->SetLineWidth(4);
    gExp->SetLineWidth(4);
    gObsL->SetLineWidth(2);
    gExpL->SetLineWidth(2);
    gObsR->SetLineWidth(2);
    gExpR->SetLineWidth(2);
    gObs->SetLineStyle(1);
    gExp->SetLineStyle(7);
    gObsL->SetLineStyle(1);
    gExpL->SetLineStyle(7);
    gObsR->SetLineStyle(1);
    gExpR->SetLineStyle(7);
    gObs->SetLineColor(kBlack);
    gExp->SetLineColor(kBlack);
    gObsL->SetLineColor(kBlue);
    gExpL->SetLineColor(kBlue);
    gObsR->SetLineColor(kRed);
    gExpR->SetLineColor(kRed);
    if(killlowdiag){
    for( int i = gObs->GetN()-1; i>=0;--i){
        double x,y;
        gObs->GetPoint(i,x,y);
        if(x-y<172.5) gObs->RemovePoint(i);
    }
    for( int i = gExp->GetN()-1; i>=0;--i){
        double x,y;
        gExp->GetPoint(i,x,y);
        if(x-y<172.5) gExp->RemovePoint(i);
    }
    }
    for( int i = gObsL->GetN()-1; i>=0;--i){
        double x,y;
        gObsL->GetPoint(i,x,y);
        if(x-y<172.5) gObsL->RemovePoint(i);
    }
    for( int i = gObsR->GetN()-1; i>=0;--i){
        double x,y;
        gObsR->GetPoint(i,x,y);
        if(x-y<172.5) gObsR->RemovePoint(i);
    }
    for( int i = gExpL->GetN()-1; i>=0;--i){
        double x,y;
        gExpL->GetPoint(i,x,y);
        if(x-y<172.5) gExpL->RemovePoint(i);
    }
    for( int i = gExpR->GetN()-1; i>=0;--i){
        double x,y;
        gExpR->GetPoint(i,x,y);
        if(x-y<172.5) gExpR->RemovePoint(i);
    }

    
   TCanvas *c1 = new TCanvas("c1", "c1",50,50,600,600);
   gStyle->SetOptFit(1);
   gStyle->SetOptStat(0);
   gStyle->SetOptTitle(0);
   gStyle->SetErrorX(0.5); 
   //c1->Range(-6.311689,-1.891383,28.75325,4.56342);
   c1->SetFillColor(0);
   c1->SetBorderMode(0);
   c1->SetBorderSize(2);
   //c1->SetLogy();
   c1->SetTickx(1);
   c1->SetTicky(1);
   c1->SetLeftMargin(0.15);
   c1->SetRightMargin(0.05);
   c1->SetTopMargin(0.07);
   c1->SetBottomMargin(0.15);
   c1->SetFrameFillStyle(0);
   c1->SetFrameBorderMode(0);
//.........这里部分代码省略.........

void SignfificanceT2tt(bool pval, int exp=0){


    TFile *f = TFile::Open("Significances2DHistograms_T2tt.root");
    TH2F *h;
    if(pval){
        if(exp==0) h = (TH2F*)f->Get("hpObs");
        else if(exp==1) h = (TH2F*)f->Get("hpExpPosteriori");
        else if(exp==2) h = (TH2F*)f->Get("hpExpPriori");
    }
    else {
        if(exp==0) h = (TH2F*)f->Get("hObs");
        else if(exp==1) h = (TH2F*)f->Get("hExpPosteriori");
        else if(exp==2) h = (TH2F*)f->Get("hExpPriori");
    }
    h->GetXaxis()->SetTitle("m_{#tilde{t}} [GeV]");
    h->GetXaxis()->SetLabelFont(42);
    h->GetXaxis()->SetLabelSize(0.035);
    h->GetXaxis()->SetTitleSize(0.05);
    h->GetXaxis()->SetTitleOffset(1.2);
    h->GetXaxis()->SetTitleFont(42);
    h->GetYaxis()->SetTitle("m_{#tilde{#chi}_{1}^{0}} [GeV]");
    h->GetYaxis()->SetLabelFont(42);
    h->GetYaxis()->SetLabelSize(0.035);
    h->GetYaxis()->SetTitleSize(0.05);
    h->GetYaxis()->SetTitleOffset(1.35);
    h->GetYaxis()->SetTitleFont(42);
    double maximum = h->GetMaximum();
    double minimum = h->GetMinimum();
    double sigmin = 99; int sigminx=-1; int sigminy=-1; if(pval) sigmin = -99;
    h->GetZaxis()->SetRangeUser(minimum,maximum);
    for(int x = 1; x<=h->GetNbinsX();++x){
        for(int y = 1; y<=h->GetNbinsX();++y){
            if(!pval&&h->GetBinContent(x,y)<sigmin){ sigmin =h->GetBinContent(x,y); sigminx = x; sigminy = y; }
            if( pval&&h->GetBinContent(x,y)>sigmin){ sigmin =h->GetBinContent(x,y); sigminx = x; sigminy = y; }
            if(!pval&&h->GetXaxis()->GetBinLowEdge(x)<h->GetYaxis()->GetBinLowEdge(y)+75) h->SetBinContent(x,y,-999);
            if(h->GetXaxis()->GetBinLowEdge(x)<374) continue;
            if(h->GetXaxis()->GetBinLowEdge(x)>424) continue;
            if(h->GetYaxis()->GetBinLowEdge(y)<199) continue;
            if(h->GetYaxis()->GetBinLowEdge(y)>249) continue;
            if(!pval&&h->GetBinContent(x,y)>0.3) h->SetBinContent(x,y,0.05);
        }
    }
    h->GetZaxis()->SetRangeUser(minimum,maximum);
    if(!pval) cout << "minimal significance " << sigmin << " at " << h->GetXaxis()->GetBinLowEdge(sigminx) << "-" << h->GetYaxis()->GetBinLowEdge(sigminy) << endl;
    else cout << "maximal p- value " << sigmin << " at " << h->GetXaxis()->GetBinLowEdge(sigminx) << "-" << h->GetYaxis()->GetBinLowEdge(sigminy) << endl;
    
   TCanvas *c1 = new TCanvas("c1", "c1",50,50,600,600);
   gStyle->SetOptFit(1);
   gStyle->SetOptStat(0);
   gStyle->SetOptTitle(0);
   gStyle->SetErrorX(0.5); 
   //c1->Range(-6.311689,-1.891383,28.75325,4.56342);
   c1->SetFillColor(0);
   c1->SetBorderMode(0);
   c1->SetBorderSize(2);
   //c1->SetLogy();
   c1->SetTickx(1);
   c1->SetTicky(1);
   c1->SetLeftMargin(0.15);
//   c1->SetRightMargin(0.05);
    c1->SetRightMargin(0.15);
   c1->SetTopMargin(0.07);
   c1->SetBottomMargin(0.15);
   c1->SetFrameFillStyle(0);
   c1->SetFrameBorderMode(0);
   c1->SetFrameFillStyle(0);
   c1->SetFrameBorderMode(0);
   gStyle->SetHatchesLineWidth(0);

   TH1F *hSum = new TH1F("hSum","",10,100,900);
   hSum->SetMinimum(0.);
   hSum->SetMaximum(550);
   hSum->SetDirectory(0);
   hSum->SetStats(0);
    hSum->Draw();
    hSum->GetYaxis()->SetRangeUser(0,500);
    hSum->GetXaxis()->SetRangeUser(100,900);

   Int_t ci;   // for color index setting
   ci = TColor::GetColor("#000099");
   hSum->SetLineColor(ci);
   hSum->SetLineStyle(0);
   hSum->SetMarkerStyle(20);
   hSum->GetXaxis()->SetTitle("m_{#tilde{t}} [GeV]");
   //hSum->GetXaxis()->SetBit(TAxis::kLabelsVert);
   hSum->GetXaxis()->SetLabelFont(42);
   //hSum->GetXaxis()->SetLabelOffset(0.005);
   hSum->GetXaxis()->SetLabelSize(0.035);
   hSum->GetXaxis()->SetTitleSize(0.06);
   hSum->GetXaxis()->SetTitleOffset(1.2);
   hSum->GetXaxis()->SetTitleFont(42);
   hSum->GetYaxis()->SetTitle("m_{#tilde{#chi}}_{1}^{0} [GeV]");
   hSum->GetYaxis()->SetLabelFont(42);
   //hSum->GetYaxis()->SetLabelOffset(0.007);
   hSum->GetYaxis()->SetLabelSize(0.035);
   hSum->GetYaxis()->SetTitleSize(0.05);
   hSum->GetYaxis()->SetTitleOffset(1.3);
   hSum->GetYaxis()->SetTitleFont(42);
    
//.........这里部分代码省略.........

void makePlots()
{
 // **********************************************
 // *            Input parameters                *
 // **********************************************
 // switch to include/exclude sytematics uncertainties 
  bool plot7TeV = false;
  bool plot10TeV = true;

 // array of signal efficiencies
 Double_t S_eff[5] = {0.36, 0.43, 0.54, 0.59, 0.66};
 
 // array of N_background for 100 pb^-1
 Double_t N_bkg_100[5] = {10.35, 6.37, 2.83, 1.27, 1.27};
 Double_t N_bkg_100_7TeV[5] = {4.068, 1.937, 0.816, 0.3841, 0.3841};
 
 // arrays of leptoquark masses
 string mass[5] = {"m250", "m300", "m400", "m500", "m600"};
 Double_t m[5] = {250, 300, 400, 500, 600};
 
 
 // array of theoretical cross-sections for different leptoquark mass
 Double_t xsTh[5] = {9.64045, 3.6081, 0.691125, 0.176217, 0.0538333};
 Double_t xs7Th[10] = {3.48, 1.21, 0.206, 0.0463, 0.0124, 0.00676, 0.00377, 0.00124, 0.000437, 0.00016};
 
 // arrays of beta^2 ranges that contain 95% CL exclusion of different leptoquark mass hypotheses
 Double_t rangeMin[5] = {0.018, 0.03, 0.09, 0.3,  0.9};
 Double_t rangeMax[5] = { 0.03, 0.08,  0.2, 0.55, 1.5};
 Double_t rangeMin_7TeV[5] = {0.03,  0.08,  0.35,  1.0, 3.5};
 Double_t rangeMax_7TeV[5] = {0.08, 0.25,  0.55, 1.5, 5.0};

 // number of points used for interpolation
 Int_t NPts = 5;
  
 // filename for the final plot (NB: changing the name extension changes the file format)
 string fileName = "beta_vs_m_excl.eps";
 if (plot7TeV && !plot10TeV) fileName = "7TeV_beta_vs_m_excl.eps";
 if (!plot7TeV && plot10TeV) fileName = "10TeV_beta_vs_m_excl.eps";
  
 // axes labels for the final plot 
 string title = ";m (GeV/c^{2});#beta";

 TH2F *bg = new TH2F("bg",title.c_str(), 100, 200., 650., 100, 0., 1.);
 bg->SetStats(kFALSE);
 bg->SetTitleOffset(1.,"X");
 bg->SetTitleOffset(1.5,"Y");

 //TF1 *fit1 = new TF1("fit1","[0]*exp([1]*x)",200,600);
 TF1 *fit1 = new TF1("fit1","[0]+[1]*x+[2]*pow(x,2)",200,600);

 //TF1 *fit2 = new TF1("fit2","[0]*exp([1]*x)",200,600);
 TF1 *fit2 = new TF1("fit2","[0]+[1]*x+[2]*pow(x,2)",200,600);
 
 // region excluded by Tevatron limits
 Double_t x_shaded[14] = {200,205,210,215,220,225,230,235,240,245,250,255,256,200};  // 250 pb-1
 Double_t y_shaded[14] = {sqrt(0.0267547),sqrt(0.0380903),sqrt(0.0531225),sqrt(0.0736446),sqrt(0.10221),sqrt(0.142373),sqrt(0.194058),sqrt(0.265091),sqrt(0.385415),sqrt(0.520114),sqrt(0.721552),sqrt(0.945015),1,1};
 Double_t NEW_x_shaded[14] = {216,220,225,231,235,254,268,276,284,289,293,296,297,299};
 Double_t NEW_y_shaded[14] = {0.02,0.04,0.06,0.08,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1};  // 1 fb-1

 // **********************************************
 // *  Don't change anything below this point!   *
 // **********************************************
 
  // turn on/off batch mode
 gROOT->SetBatch(kTRUE);
 // set ROOT style
 style();
 
 TCanvas *c = new TCanvas("c","",800,800);
 c->cd();

 bg->Draw();
 
 TGraph *grexcl = new TGraph(14,x_shaded,y_shaded);
 grexcl->SetFillColor(kGray+1);
 grexcl->Draw("f");
 TGraph *NEWexcl = new TGraph(14,NEW_x_shaded,NEW_y_shaded);
 NEWexcl->SetLineStyle(6);
 NEWexcl->SetLineWidth(3);
 NEWexcl->Draw("C");

 Int_t size = sizeof(S_eff)/sizeof(*S_eff);

 // array of minimum beta^2 values for 95% CL exclusion
//  Double_t beta[size];
 Double_t beta[5] = {0.156865, 0.214541, 0.382863, 0.655572, 1.11471};
 Double_t beta_ln[5] = {log(0.156865), log(0.214541), log(0.382863), log(0.655572), log(1.11471)};
 Double_t beta_7TeV[5] = {0.221188, 0.307225, 0.59993, 1.14052, 2.08235};
 Double_t beta_7TeV_ln[5] = {log(0.221188), log(0.307225), log(0.59993), log(1.14052), log(2.08235)};

//  for(Int_t i = 0; i < size; i++){
  
//   beta_7TeV[i] = beta_for_exclusion(NPts, xs7Th[i], rangeMin_7TeV[i], rangeMax_7TeV[i], 0, S_eff[i], 0, N_bkg_100_7TeV[i], 0, mass[i]);
//   beta_7TeV_ln[i] = log(beta_7TeV[i]);
//  }

 TGraph *beta_vs_m = new TGraph(size, m, beta);
 TGraph *beta_vs_m_fit = new TGraph(size, m, beta_ln);
 beta_vs_m->SetLineWidth(2);
 beta_vs_m->SetLineStyle(7);
//.........这里部分代码省略.........

//.........这里部分代码省略.........
  TGraph *gEzSlope = new TGraph(Nsim,Current,EzSlope);
  TGraph *gEzSlopeMaxFocus = new TGraph(Nsim,Current,EzSlopeMaxFocus);
  TGraph *gEzSlopeLO = new TGraph(Nsim,Current,EzSlopeLO);

  TGraph *gTransRatio = new TGraph(Nsim,Current,TransRatio);
  TGraph *gTransRatioFocus = new TGraph(Nsim,Current,TransRatioFocus);

  // Int_t color1 = kMagenta-8;
  // Int_t color2 = kMagenta-2;
  // Int_t color3 = kMagenta-4;
  // Int_t color4 = kAzure+1;
  // Int_t color5 = kGray+2;
  // Int_t markersty = 20;
  // Float_t markersiz = 1.0;
  // Int_t linewit = 2;

  Int_t color1 = kGray+3;
  Int_t color2 = kGray+2;
  Int_t color3 = kAzure+1;
  Int_t color4 = kRed;
  Int_t color5 = kGray+1;
  Int_t markersty = 20;
  Float_t markersiz = 1.0;
  Int_t linewit = 1;

  gRadiusX->SetLineColor(color1);
  gRadiusX->SetLineWidth(linewit);
  gRadiusX->SetMarkerColor(color1);
  gRadiusX->SetMarkerStyle(markersty);
  gRadiusX->SetMarkerSize(markersiz);

  gRadiusLU->SetLineColor(color1);
  gRadiusLU->SetLineWidth(1);
  gRadiusLU->SetLineStyle(4);
  gRadiusLU->SetMarkerColor(color1);
  gRadiusLU->SetMarkerStyle(markersty);
  gRadiusLU->SetMarkerSize(0.0);

  gRadiusZ->SetLineColor(color2);
  gRadiusZ->SetLineWidth(linewit);
  gRadiusZ->SetMarkerColor(color2);
  gRadiusZ->SetMarkerStyle(markersty);
  gRadiusZ->SetMarkerSize(markersiz);

  gRadiusLO->SetLineColor(color2);
  gRadiusLO->SetLineWidth(1);
  gRadiusLO->SetLineStyle(3);
  gRadiusLO->SetMarkerColor(color2);
  gRadiusLO->SetMarkerStyle(markersty);
  gRadiusLO->SetMarkerSize(0.0);

  gRadiusPsi->SetLineColor(color3);
  gRadiusPsi->SetLineWidth(linewit);
  gRadiusPsi->SetMarkerColor(color3);
  gRadiusPsi->SetMarkerStyle(markersty);
  gRadiusPsi->SetMarkerSize(markersiz);

  gEzMax->SetLineColor(color2);
  gEzMax->SetLineWidth(linewit);
  gEzMax->SetMarkerColor(color2);
  gEzMax->SetMarkerStyle(markersty);
  gEzMax->SetMarkerSize(markersiz);

  gEzMaxUse->SetLineColor(color2);
  gEzMaxUse->SetLineWidth(linewit);
  gEzMaxUse->SetMarkerColor(color2);

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

    hFrame[i]->GetXaxis()->SetLabelSize(0.07);
    hFrame[i]->GetXaxis()->SetLabelOffset(0.01/vfactor[i]);

    hFrame[i]->GetXaxis()->SetTitleSize(0.07);
    hFrame[i]->GetXaxis()->SetTitleOffset(1.1);///vfactor[i]);
    
    hFrame[i]->Draw("axis");    

    if(i==2) {
      Float_t yMin = 0.0001; // minYmean[0] - (maxYmean[0]-minYmean[0])*0.1;
      Float_t yMax = maxYmean[0] + (maxYmean[0]-minYmean[0])*0.1;
      hFrame[i]->GetYaxis()->SetRangeUser(yMin,yMax);
      hFrame[i]->GetXaxis()->SetTitle("propagation length [mm]");
      hFrame[i]->GetYaxis()->SetTitle("p_{z} [GeV/c]");
      //    PlasmaGlob::SetH1LabelSize(hFrame[i]);
      
      gEneRms->SetFillColor(kGray);
      gEneRms->Draw("3");
      gYmean[0]->Draw("L");  

      // Charge on right axis:
      Float_t yrMin = 0.00001; // minCharge - (maxCharge-minCharge)*0.1;
      Float_t yrMax = maxCharge + (maxCharge-minCharge)*0.1;
      Float_t slope = (yMax-yMin)/(yrMax-yrMin);
      
      Double_t *x = gCharge->GetX();
      Double_t *y = gCharge->GetY();
      for(Int_t j=0;j<gCharge->GetN();j++) {
      	gCharge->SetPoint(j,x[j],slope*(y[j]-yrMin)+yMin);
      }

      //hFrame[i]->GetYaxis()->SetRangeUser(yrMin,yrMax);
      gCharge->SetLineStyle(2);
      gCharge->SetLineWidth(3);
      gCharge->SetLineColor(kAzure-8);      
      gCharge->Draw("L");

      pad[i]->Update();

    } else if(i==1) {
      Float_t yMin,yMax;
      if(minXrms[1]<minXrms[0])
	yMin = minXrms[1];
      else
	yMin = minXrms[0];
      if(maxXrms[1]>maxXrms[0])
	yMax = maxXrms[1];
      else
	yMax = maxXrms[0];

      Float_t yDist = yMax - yMin;
      yMin -= 0.1*yDist; 
      yMax += 0.1*yDist;
      hFrame[i]->GetYaxis()->SetRangeUser(yMin,yMax);
      hFrame[i]->GetXaxis()->SetTitle("propagation length [mm]");
      hFrame[i]->GetYaxis()->SetTitle("rms size [#mum]");
      //    PlasmaGlob::SetH1LabelSize(hFrame[i]);

      gXrms[0]->SetLineStyle(1);
      gXrms[0]->SetLineWidth(3);
      gXrms[0]->SetLineColor(kOrange+10);      
      gXrms[0]->Draw("L");
      gXrms[1]->SetLineStyle(2);
      gXrms[1]->SetLineWidth(3);
      gXrms[1]->SetLineColor(kAzure-8);

void plotConvergence(int step = 1){

  int count = 0;
  int tempcount=0;

  ifstream myInput(Form("converge_steps_%i.txt", step));

//    ifstream myInput("lindaMacros/ANITAsymmetric.txt");

  double meanH, rmsH, gradH, meanV, rmsV, gradV;

  double sumMeanH[1000];
  double sumRmsH[1000];
  double sumGradH[1000];
  double sumMeanV[1000];
  double sumRmsV[1000];
  double sumGradV[1000];
  double x[1000];

  if (myInput.is_open()){
    while(myInput.good()){
      
      myInput >> meanH >> rmsH >> gradH >> meanV >> rmsV >> gradV;
      
      if (tempcount%100==0){
	sumMeanH[count] = meanH;
	sumRmsH[count] = rmsH;
	sumGradH[count] = gradH;
	sumMeanV[count] = meanV;
	sumRmsV[count] = rmsV;
	sumGradV[count] = gradV;
	x[count] = tempcount*1.;
	count++;
	cout << meanH << " " << rmsH << " " << gradH << endl;
      }
      
      tempcount++;
    }
  }
  
  
  count--;

  TGraph *gMeanH = new TGraph(count, x, sumMeanH);
  TGraph *gRmsH  = new TGraph(count, x, sumRmsH);
  TGraph *gGradH = new TGraph(count, x, sumGradH);
  TGraph *gMeanV = new TGraph(count, x, sumMeanV);
  TGraph *gRmsV  = new TGraph(count, x, sumRmsV);
  TGraph *gGradV = new TGraph(count, x, sumGradV);

  gMeanH->SetLineColor(kCyan+1);
  gMeanH->SetLineWidth(2);

  gRmsH->SetLineColor(kMagenta+3);
  gRmsH->SetLineStyle(7);
  gRmsH->SetLineWidth(2);

  gGradH->SetLineColor(kGreen+2);
  gGradH->SetLineWidth(2);
  gGradH->SetLineStyle(8);


  gMeanV->SetLineColor(kCyan+1);
  gMeanV->SetLineWidth(2);

  gRmsV->SetLineColor(kMagenta+3);
  gRmsV->SetLineStyle(7);
  gRmsV->SetLineWidth(2);

  gGradV->SetLineColor(kGreen+2);
  gGradV->SetLineWidth(2);
  gGradV->SetLineStyle(8);


  TCanvas *c1 = new TCanvas("c1");
 

  gMeanH->Draw("Al");
  gMeanH->SetTitle("Convergence:horizontal antenna pairs;Iteration;Value");

  gMeanH->GetYaxis()->SetRangeUser(0.01,10000);

  c1->SetLogy();
  gMeanH->Draw("al");

  gRmsH->Draw("l");
  gGradH->Draw("l");

  TLegend *leg = new TLegend(0.89, 0.89, 0.55, 0.7);
  leg->SetFillColor(kWhite);
  leg->AddEntry(gMeanH, "#Sigma MEAN", "l");
  leg->AddEntry(gRmsH,  "#Sigma RMS" , "l");
  leg->AddEntry(gGradH, "#Sigma GRAD", "l");

  leg->Draw();

  c1->Print(Form("Convergence_HorizontalAntennaPairs_steps_%i.png", step));
  c1->Print(Form("Convergence_HorizontalAntennaPairs_steps_%i.pdf", step));


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

//.........这里部分代码省略.........
  //this only apply to NarrowResonnance case
  for(unsigned int i=0;i<N;i++){
    if(strengthLimit){
    ObsLimit[i]  /= XSecScaleFactor;
    ExpLimitm2[i]/= XSecScaleFactor;
    ExpLimitm1[i]/= XSecScaleFactor;
    ExpLimit  [i]/= XSecScaleFactor;
    ExpLimitp1[i]/= XSecScaleFactor;
    ExpLimitp2[i]/= XSecScaleFactor;
    }
    ThXSec[i]    *= XSecScaleFactor;
  }

    
  //limits in terms of signal strength
  TCanvas* c = new TCanvas("c", "c",600,600);
  TH1F* framework = new TH1F("Graph","Graph",1,90,160);
  framework->SetStats(false);
  framework->SetTitle("");
  framework->GetXaxis()->SetTitle("Higgs boson mass [GeV]");
  framework->GetYaxis()->SetTitleOffset(1.70);
  if(strengthLimit){
  framework->GetYaxis()->SetTitle("#mu = #sigma_{95%} / #sigma_{th}");
  framework->GetYaxis()->SetRangeUser(9E-1,1E3);
  c->SetLogy(true);
  }else{
  framework->GetYaxis()->SetTitle("#sigma_{95%} (fb)");
  framework->GetYaxis()->SetRangeUser(1E-1,1E3);
  c->SetLogy(true);
  }
  framework->Draw();

  TGraph* TGObsLimit   = new TGraph(N,MassAxis,ObsLimit);  TGObsLimit->SetLineWidth(2);
  TGraph* TGExpLimit   = new TGraph(N,MassAxis,ExpLimit);  TGExpLimit->SetLineWidth(2); TGExpLimit->SetLineStyle(2);
  TCutG* TGExpLimit1S  = GetErrorBand("1S", N, MassAxis, ExpLimitm1, ExpLimitp1);  
  TCutG* TGExpLimit2S  = GetErrorBand("2S", N, MassAxis, ExpLimitm2, ExpLimitp2);  TGExpLimit2S->SetFillColor(5);
  TGraph* THXSec        = new TGraph(N,MassAxis,ThXSec); THXSec->SetLineWidth(2); THXSec->SetLineStyle(1); THXSec->SetLineColor(4);


  TGExpLimit->SetLineColor(2);  TGExpLimit->SetLineStyle(1);
  TGObsLimit->SetLineWidth(2);  TGObsLimit->SetMarkerStyle(20);
  TGExpLimit2S->Draw("fc same");
  TGExpLimit1S->Draw("fc same");
  if(!blind) TGObsLimit->Draw("same CP");
  TGExpLimit->Draw("same C");

  if(strengthLimit){
     TLine* SMLine = new TLine(framework->GetXaxis()->GetXmin(),1.0,framework->GetXaxis()->GetXmax(),1.0);
     SMLine->SetLineWidth(2); SMLine->SetLineStyle(1); SMLine->SetLineColor(4);      
     SMLine->Draw("same C");
  }else{
     THXSec->Draw("same C");
  }


  TPaveText *pave = new TPaveText(0.1,0.96,0.99,0.99,"NDC");
  char LumiLabel[1024];
  if(energy<9){  sprintf(LumiLabel,"CMS preliminary,  #sqrt{s}=%.0f TeV, #scale[0.5]{#int} L=%6.1ffb^{-1} - %20s",energy, luminosity,legendName.Data());
  }else{         sprintf(LumiLabel,"CMS preliminary,  #sqrt{s}=%.0f TeV #scale[0.5]{#int} L=%6.1ffb^{-1}, #sqrt{s}=%.0f TeV #scale[0.5]{#int} L=%6.1ffb^{-1}",7.0,5.0,8.0,19.7);
  }
  pave->SetBorderSize(0);
  pave->SetFillStyle(0);
  pave->SetTextFont(42);
  TObjArray* tokens = (TString(LumiLabel)).Tokenize("\\\\");
  int nt = tokens->GetEntries();
  for(int it=0; it<nt; ++it){