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

TGraph* autogain152(TH1 *hist) {

   hist->GetXaxis()->SetRangeUser(200.,16000.);
   TSpectrum *s = new TSpectrum();
   Int_t nfound = s->Search(hist,6,"",0.08); //This will be dependent on the source used.
   printf("Found %d candidate peaks to fit\n",nfound);
   if(nfound > 6)
      nfound = 6;

   std::vector<float> vec;
   for(int x=0;x<nfound;x++)
      vec.push_back(s->GetPositionX()[x]);

   std::sort(vec.begin(),vec.end());

   Float_t energies[] = {121.7830, 244.6920, 344.276, 778.903, 964.131, 1408.011};
   TGraph* slopefit = new TGraph(nfound, &(vec[0]), energies);

   printf("Now fitting: Be patient\n");
   slopefit->Fit("pol1");
   if(slopefit->GetFunction("pol1")->GetChisquare() > 10.) {
      slopefit->RemovePoint(slopefit->GetN()-1);
      slopefit->Fit("pol1");
   }
   TChannel *chan = 0;
   slopefit->Draw("AC*");

   return slopefit;
}

void test5()
{
   TCanvas *c1 = new TCanvas("c1", "c1", 0,   0, 600, 300);
   TCanvas *c2 = new TCanvas("c2", "c2", 605, 0, 600, 300);

   const Int_t n = 20;
   Double_t x[n], y[n];
   for (Int_t i = 0; i < n; i++) {
      x[i] = i;
      y[i] = i;
   }

   TGraph *gr = new TGraph(n, x, y);
   gr->SetMarkerStyle(20);
   c1->cd();
   c1->SetGrid();
   gr->Draw("APC");
   gr->SetHighlight();
   TGraph::SetHighlightPad(c2);

   TH1F *h[n];
   for (Int_t i = 0; i < n; i++) {
      h[i] = new TH1F(TString::Format("h_%02d", i), "", 100, -5.0, 5.0);
      h[i]->SetTitle(h[i]->GetName());
      h[i]->FillRandom("gaus");
      gr->AddHighlight(i, h[i]);
   }

   gr->GetListOfHighlights()->SetOwner(kTRUE);
   gr->RemovePoint(5);
   gr->RemovePoint(9); // point with x = 10 (after remove previous point)
   delete h[14];
   delete h[17];
}

TGraph* getContour(TH2D* h, TString name) {

  TGraph* graph = new TGraph(100);
  graph->SetName(name);
  int ip = 0;
  int nx = h->GetXaxis()->GetNbins();
  int ny = h->GetYaxis()->GetNbins();

  // for y>x
  int ix = -1;
  for(int j=ny;true; j--) {
    int k = -1;
    for(int i=2; i<nx-1; i++) {
      if(h->GetBinContent(i,j) < 0) {
	std::cout << "i,j,z : " << i << ", " << j << ", " << h->GetBinContent(i,j) << std::endl;
        k = i;
        break;
      }
    }// for i
    if(k<0) continue;
    double y = h->GetYaxis()->GetBinCenter(j);
    double x1 = h->GetXaxis()->GetBinCenter(k-1);
    double x2 = h->GetXaxis()->GetBinCenter(k);
    double z1 = h->GetBinContent(k-1,j);
    double z2 = h->GetBinContent(k,j);
    double x = x1 + (x2-x1)*fabs(z1)/fabs(z2-z1);
    std::cout << "y, x1, x2, z1, z2, x : " << y << ", " << x1 << ", " << x2 << ", " << x << ", " << z1 << ", " << z2 << std::endl;
    graph->SetPoint(ip++,x,y);

    if(h->GetYaxis()->GetBinCenter(j) < h->GetXaxis()->GetBinCenter(k)) {
      ix = k;
      break;
    }
  }// for j

  if(ix < 0) std::cout << "Something wrong...." << std::endl;

  // for y<x
  for(int i=ix; i<=nx; i++) {
    int k = -1;
    for(int j=2; j<ny-1; j++) {
      if(h->GetBinContent(i,j) < 0) {
        k = j;
        break;
      }
    }// for j
    if(k<0) continue;
    double x = h->GetXaxis()->GetBinCenter(i);
    double y1 = h->GetYaxis()->GetBinCenter(k-1);
    double y2 = h->GetYaxis()->GetBinCenter(k);
    double z1 = h->GetBinContent(i,k-1);
    double z2 = h->GetBinContent(i,k);
    double y = y1 + (y2-y1)*fabs(z1)/fabs(z2-z1);
    std::cout << "x, y1, y2, z1, z2, y : " << x << ", " << y1 << ", " << y2 << ", " << y << ", " << z1 << ", " << z2 << std::endl;
    graph->SetPoint(ip++,x,y);
  }// for i

  ip = graph->GetN()-1;
  while(1) {
    double x, y;
    graph->GetPoint(ip,x,y);
    if(x>1) break;
    else graph->RemovePoint(ip);
    ip--;
  }

  return graph;
}

//-------------------------------------------------------------------------------
//-------------------------------------------------------------------------------
TGraph* compMVAcut(const TH2* histogramMVAoutput_vs_Pt, const TH1* histogramPt, double Efficiency_or_FakeRate)
{
  const TAxis* xAxis = histogramMVAoutput_vs_Pt->GetXaxis();
  int numBinsX = xAxis->GetNbins();

  const TAxis* yAxis = histogramMVAoutput_vs_Pt->GetYaxis();
  int numBinsY = yAxis->GetNbins();

  TGraph* graph = new TGraphAsymmErrors(numBinsX);
  std::string graphName = Form("%s_graph", histogramMVAoutput_vs_Pt->GetName());
  graph->SetName(graphName.data());

  int numPoints = 0;

  for ( int iBinX = 1; iBinX <= numBinsX; ++iBinX ) {
    double ptMin = xAxis->GetBinLowEdge(iBinX);
    double ptMax = xAxis->GetBinUpEdge(iBinX);

    int binLowIndex = const_cast<TH1*>(histogramPt)->FindBin(ptMin);
    int binUpIndex  = const_cast<TH1*>(histogramPt)->FindBin(ptMax);
    //std::cout << "ptMin = " << ptMin << ", ptMax = " << ptMax << ": binLowIndex = " << binLowIndex << ", binUpIndex = " << binUpIndex << std::endl;
    histogramPt->GetXaxis()->SetRange(binLowIndex, binUpIndex);

    // CV: skip bins of low statistics
    if ( histogramPt->GetEntries() < 100 ) {
      std::cout << "Warning: bin @ x = " << xAxis->GetBinCenter(iBinX) << " has low statistics (#entries = " << histogramPt->GetEntries() << ") --> skipping !!" << std::endl;
      continue;
    }

    double x = histogramPt->GetMean();

    //std::cout << "iBinX = " << iBinX << ": xMean = " << x << std::endl;

    double normalization = 0.;
    for ( int iBinY = numBinsY; iBinY >= 1; --iBinY ) {
      normalization += histogramMVAoutput_vs_Pt->GetBinContent(iBinX, iBinY);
    }
    // CV: skip bins of low statistics
    //if ( normalization < 100 ) {
    //  std::cout << "Warning: bin @ x = " << xAxis->GetBinCenter(iBinX) << " has low statistics (normalization = " << normalization << ") --> skipping !!" << std::endl;
    //  continue;
    //}
    
    double y = -1.;

    double runningSum = 0.;
    for ( int iBinY = numBinsY; iBinY >= 1; --iBinY ) {
      double binContent_normalized = histogramMVAoutput_vs_Pt->GetBinContent(iBinX, iBinY)/normalization;
      if ( (runningSum + binContent_normalized) > Efficiency_or_FakeRate ) {
	y = yAxis->GetBinUpEdge(iBinY) - ((Efficiency_or_FakeRate - runningSum)/binContent_normalized)*yAxis->GetBinWidth(iBinY);
	//std::cout << "iBinY = " << iBinY << " (yCenter = " << yAxis->GetBinCenter(iBinY) << "): binContent = " << binContent_normalized << std::endl;
	//std::cout << "--> setting y = " << y << std::endl;
	assert(y >= yAxis->GetBinLowEdge(iBinY));
	break;
      } else {
	runningSum += binContent_normalized;
	//std::cout << "iBinY = " << iBinY << " (yCenter = " << yAxis->GetBinCenter(iBinY) << "): runningSum = " << runningSum << std::endl;
      }
    }
    
    graph->SetPoint(numPoints, x, y);
    ++numPoints;
  }

  for ( int iPoint = numBinsX; iPoint > numPoints; --iPoint ) {
    //std::cout << "removing point #" << iPoint << std::endl;
    graph->RemovePoint(graph->GetN() - 1);
  }

  // reset x-axis range selection 
  histogramPt->GetXaxis()->SetRange(1., 0.);

  return graph;
}

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);
//.........这里部分代码省略.........