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

TGraph *MarcosExclusionLine(TH2F *exclusionshape, int scantype) {
//  write_warning(__FUNCTION__,"USING MARIAS ALGORITHM...");
//  return get_exclusion_line(exclusionshape);
  
  TH2F *fakehisto = flipth2f(exclusionshape);
  TGraph *fakegraph = get_mSUGRA_exclusion_line(fakehisto, scantype);
  TGraph *realgraph = new TGraph(fakegraph->GetN());
  double x,y;
  float last_x=0;
  float last_y=0;
  int counter=0;
  for(int i=0;i<fakegraph->GetN();i++) {
	fakegraph->GetPoint(i,x,y);
	if(scantype==PlottingSetup::SMS) {
	  if(y-x<75) {
	    realgraph->SetPoint(counter,last_x,last_y);
	    counter++;
	    continue;
	  }
	}
	realgraph->SetPoint(counter,y,x);
	last_x=y;
	last_y=x;
	counter++;
  }
  realgraph->SetLineColor(TColor::GetColor("#151515")); //nice black
  realgraph->SetLineWidth(2);
  //realgraph->SetLineWidth(4);//paper style
  
  delete fakegraph;

  return realgraph;
}

void
plotSignificance(const char* filename, const char* channel, double min_=0., double max_=7., bool log_=false, std::string dataset_="CMS Preliminary,  H#rightarrow#tau#tau,  4.9 fb^{-1} at 7 TeV, 19.4 fb^{-1} at 8 TeV", std::string xaxis_="m_{H} [GeV]", std::string yaxis_="Significance", bool mssm_=false)
{
  TFile* file = TFile::Open(filename);
  // retrieve TGraphs from file
  TGraph* expected = (TGraph*)file->Get(std::string(channel).append("/expected").c_str());
  TGraph* observed = (TGraph*)file->Get(std::string(channel).append("/observed").c_str());
  TGraphAsymmErrors* innerBand = (TGraphAsymmErrors*)file->Get(std::string(channel).append("/innerBand").c_str());
  TGraphAsymmErrors* outerBand = (TGraphAsymmErrors*)file->Get(std::string(channel).append("/outerBand").c_str());
  // set up styles
  SetStyle();
  // create the unit lines
  TGraph* unit3 = new TGraph();
  TGraph* unit5 = new TGraph();
  for(int imass=0, ipoint=0; imass<expected->GetN(); ++imass){
    unit3->SetPoint(ipoint, expected->GetX()[imass], 3.);
    unit5->SetPoint(ipoint, expected->GetX()[imass], 5.);
    ++ipoint;
  }
  // set proper maximum
  float max = maximum(expected, max_);
  // do the plotting 
  TCanvas canv = TCanvas("canv", "Limits", 600, 600);
  // do the plotting 
  plottingSignificance(canv, innerBand, outerBand, expected, observed, unit3, unit5, xaxis_, yaxis_, min_, max, log_, mssm_);
  /// setup the CMS Preliminary
  CMSPrelim(dataset_.c_str(), "", 0.145, 0.835);
  // write results to files
  canv.Print(std::string(channel).append("_significance").append(".png").c_str());
  canv.Print(std::string(channel).append("_significance").append(".pdf").c_str()); 
  canv.Print(std::string(channel).append("_significance").append(".eps").c_str()); 
  return;
}

void ButkevichDrawer::adjust_graph(TGraph& graph) {
  for (int iPoint=0;iPoint<graph.GetN();iPoint++) {
    double x,y;
    graph.GetPoint(iPoint,x,y);
    graph.SetPoint(iPoint,x,y*1.0e-38);
  }
}

void getStartFinFrom(const TGraph& gr,
                     Double_t& start, Double_t& fin,
                     const Bool_t filterZero=kTRUE) {
   // get a nice time axis for the domain of the graph
   start = fin = 0;
   const Int_t n = gr.GetN();
   Bool_t first = kTRUE;
   const Double_t* y = gr.GetX();
   for (Int_t i=0; i<n; ++i, ++y) {
      if (filterZero && ((*y)==0)) {
         continue;
      } else {
         if (first) {
            start = fin = *y;
            first = kFALSE;
         } else {
            if ( (*y) < start ) {
               start = *y;
            }
            if ( (*y) > fin ) {
               fin = *y;
            }
         }
      }
   }
}

void
plotMaxLikelihood(const char* filename, const char* channel, double min_=0., double max_=-1., bool log_=false, std::string dataset_="CMS Preliminary,  H#rightarrow#tau#tau,  4.9 fb^{-1} at 7 TeV, 19.4 fb^{-1} at 8 TeV", std::string xaxis_="m_{H} [GeV]", std::string yaxis_="best fit for #sigma/#sigma_{SM}", bool mssm_=false)
{
  TFile* file = TFile::Open(filename);
  // retrieve TGraphs from file
  TGraph* expected = (TGraph*)file->Get(std::string(channel).append("/expected").c_str());
  TGraphAsymmErrors* innerBand = (TGraphAsymmErrors*)file->Get(std::string(channel).append("/innerBand").c_str());
  // set up styles
  SetStyle();
  // create the unit line
  TGraph* unit = 0;
  if(!mssm_){
    unit = new TGraph();
    for(int imass=0, ipoint=0; imass<expected->GetN(); ++imass){
      unit->SetPoint(ipoint, expected->GetX()[imass], 1.); ++ipoint;
    }
  }
  // set proper maximum
  float max = maximum(expected, max_);
  // do the plotting 
  TCanvas canv = TCanvas("canv", "Limits", 600, 600);
  plottingLimit(canv, innerBand, 0, expected, 0, unit, xaxis_, yaxis_, min_, max, log_, "BESTFIT", "", mssm_);
  // setup CMS Preliminary
  CMSPrelim(dataset_.c_str(), "", 0.145, 0.835);
  // write results to files
  canv.Print(std::string(channel).append("_").append("bestfit").append(".png").c_str());
  canv.Print(std::string(channel).append("_").append("bestfit").append(".pdf").c_str()); 
  canv.Print(std::string(channel).append("_").append("bestfit").append(".eps").c_str()); 
  return;
}

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 Calculation(){
  gStyle->SetOptStat("neRMI");
  Int_t CalibrationNumber  = 0; 

  TFile* tf  = new TFile(Form("CalibrationData/Calibration_4231_%d.root",CalibrationNumber));
  const int nCSI = 2716;  
  TH1D*  hisCalFactor[nCSI];
  int    nCal[nCSI];
  double CalFactor[nCSI];
  double CalRMS[nCSI];

  TH1D* hisCalDistrib = new TH1D("hisCalDistrib","",200,0,2);
  TH1D* hisRMSDistrib = new TH1D("hisRMSDistrib","",100,0,0.1);
  TH1D* hisHitDistrib = new TH1D("hisHitDistrib","",100,0,100);
  TGraph* grChannel   = new TGraph();
  std::ofstream ofs(Form("Data/CalibrationFactor%d.txt",CalibrationNumber+1));

  for( int i = 0; i< nCSI; i++){
    hisCalFactor[i] = (TH1D*)tf->Get(Form("his_Calibration_%04d",i));
    nCal[i]         = hisCalFactor[i]->Integral();
    if( hisCalFactor[i]->Integral() < 20 ){
      continue;
    }

    hisCalFactor[i]->Fit("gaus","Q","");
    TF1* calFunction= hisCalFactor[i]->GetFunction("gaus");
    CalFactor[i]    = hisCalFactor[i]->GetMean();
    CalRMS[i]       = hisCalFactor[i]->GetRMS();
    //CalFactor[i] = calFunction->GetParameter(1);
    //CalRMS[i] = calFunction->GetParameter(2);
    hisCalDistrib->Fill(CalFactor[i]);
    hisHitDistrib->Fill(nCal[i]);
    grChannel->SetPoint(grChannel->GetN(),i,nCal[i]);
    if( nCal[i] != 0){
      hisRMSDistrib->Fill(CalRMS[i]/CalFactor[i]);      
    }
    ofs << i << "\t" << CalFactor[i] << std::endl;    
  }
  ofs.close();
  
  TCanvas* can = new TCanvas("can","",1200,1200);
  can->Divide(2,2);
  can->cd(1);
  hisCalDistrib->Draw();
  can->cd(2);
  hisRMSDistrib->Draw();
  can->cd(3);
  hisHitDistrib->Draw();
  can->cd(4);
  grChannel->SetMarkerStyle(4);
  grChannel->Draw("AP");

}

void KVIDentifier::CopyGraph(const TGraph& graph)
{
   // Copy coordinates of points from the TGraph

   Double_t x, y;
   //change number of points
   Set(graph.GetN());
   for (int i = 0; i < GetN(); i++) {
      graph.GetPoint(i, x, y);
      SetPoint(i, x, y);
   }
}

TGraph* ContourGraph( TH2F* hist,double xmin=16, double xmax=90) {

    //temporary canvas
    TCanvas* MOO = new TCanvas( TString::Format("dummy_canvas_%s", hist->GetName()), "A scan of m_{0} versus m_{12}", 0, 0, 650,640);
    MOO->cd();

    TGraph* gr0 = new TGraph();
    TH2F* h = (TH2F*)hist->Clone();
    TGraph* gr = (TGraph*)gr0->Clone(TString::Format("gr_%s", h->GetName()));

    cout << "==> Will dumb histogram: " << h->GetName() << " into a graph" <<endl;

    h->SetContour( 1 );
    //h->GetXaxis()->SetRangeUser(250,1200);
    h->GetXaxis()->SetRangeUser(xmin, xmax);
    //h->GetYaxis()->SetRangeUser(2,50);

    double pval = CombinationGlob::cl_percent[1];
    std::cout << pval << std::endl; 
    double signif = TMath::NormQuantile(1-pval);
    h->SetContourLevel( 0, signif );
    h->Draw("CONT LIST");
    h->SetDirectory(0);
    gPad->Update();

    TObjArray *contours = (TObjArray*) gROOT->GetListOfSpecials()->FindObject("contours");
    Int_t ncontours     = contours->GetSize();
    cout << "Found " << ncontours << " contours " << endl;

    TList *list = (TList*)contours->At(0);
    contours->Print("v");
    if(!list) return NULL;

    gr = (TGraph*)list->First();
    if(!gr) return NULL;

    gr->SetName(TString::Format("gr_%s", hist->GetName()));
    //gr->SetName(hist->GetName());
    int N = gr->GetN();
    double x0, y0;
    for(int j=0; j<N; j++) {
        gr->GetPoint(j,x0,y0);
        cout << j << " : " << x0 << " : "<<y0 << endl;
    }
    //  //  gr->SetMarkerSize(2.0);    
    //gr->Draw("ALP");

    delete MOO;

    cout << "Generated graph " << gr << " with name " << gr->GetName() << endl;
    return gr;
}

TGraph *PeakFit::GetComponentVsAngle(unsigned int state, double scale, const Float_t barAngles[][2]) {
	std::vector< int > barMap = GetBarMap(barAngles);
	TGraph *gr = new TGraph();

	for (unsigned int barIndex = 0; barIndex < GetNumBars(); barIndex++) {
		int bar = barMap[barIndex];
		gr->SetPoint(gr->GetN(), barAngles[bar][0], components_[bar][state] * scale);
	}
	gr->SetLineColor(kRed);
	gr->SetLineWidth(2);

	return gr;
}

void
plotTanb(const char* filename, const char* channel, bool draw_injected_=false, double min_=0., double max_=60., bool log_=false, std::string dataset_="CMS Preliminary,  H#rightarrow#tau#tau,  4.9 fb^{-1} at 7 TeV, 19.8 fb^{-1} at 8 TeV", std::string xaxis_="m_{A} [GeV]", std::string yaxis_="#bf{tan#beta}")
{
  TFile* file = TFile::Open(filename);
  // retrieve TGraphs from file
  TGraph* expected = (TGraph*)file->Get(std::string(channel).append("/expected").c_str());
  TGraph* expected_low = (TGraph*)file->Get(std::string(channel).append("/expected_low").c_str());
  TGraph* observed = (TGraph*)file->Get(std::string(channel).append("/observed").c_str());
  TGraph* observed_low = (TGraph*)file->Get(std::string(channel).append("/observed_low").c_str());
  TGraph* upperLEP = (TGraph*)file->Get(std::string(channel).append("/upperLEP").c_str());
  TGraph* lowerLEP = (TGraph*)file->Get(std::string(channel).append("/lowerLEP").c_str());
  TGraphAsymmErrors* innerBand = (TGraphAsymmErrors*)file->Get(std::string(channel).append("/innerBand").c_str());
  TGraphAsymmErrors* innerBand_low = (TGraphAsymmErrors*)file->Get(std::string(channel).append("/innerBand_low").c_str());
  TGraphAsymmErrors* outerBand = (TGraphAsymmErrors*)file->Get(std::string(channel).append("/outerBand").c_str());
  TGraphAsymmErrors* outerBand_low = (TGraphAsymmErrors*)file->Get(std::string(channel).append("/outerBand_low").c_str()); 
  TGraphAsymmErrors* plain = (TGraphAsymmErrors*)file->Get(std::string(channel).append("/plain").c_str());
  TGraphAsymmErrors* plain_low = (TGraphAsymmErrors*)file->Get(std::string(channel).append("/plain_low").c_str());

  // this is new for injected plot 
  TGraph* injected = 0;;
  if(draw_injected_) {injected = (TGraph*)file->Get("injected/observed");}

  if(draw_injected_){
    for( int i=0; i<expected->GetN(); i++){
      double shift = injected->GetY()[i]-expected->GetY()[i];
      innerBand->SetPoint(i, innerBand->GetX()[i], innerBand->GetY()[i]+shift);
      outerBand->SetPoint(i, outerBand->GetX()[i], outerBand->GetY()[i]+shift);
    }
  }

  // this functionality is not yet supported
  std::map<double, TGraphAsymmErrors*> higgsBands;
  // this functionality is not yet supported
  std::map<std::string, TGraph*> comparisons;

  // set up styles
  SetStyle();
  // do the plotting 
  TCanvas canv = TCanvas("canv", "Limits", 600, 600);
  // do the plotting 
  plottingTanb(canv, plain, plain_low, innerBand, innerBand_low, outerBand, outerBand_low, expected, expected_low, observed, observed_low, lowerLEP, upperLEP, higgsBands, comparisons, xaxis_, yaxis_, injected, min_, max_, log_);
  /// setup the CMS Preliminary
  CMSPrelim(dataset_.c_str(), "", 0.145, 0.835);
  // write results to files
  canv.Print(std::string(channel).append("_tanb").append(".png").c_str());
  canv.Print(std::string(channel).append("_tanb").append(".pdf").c_str()); 
  canv.Print(std::string(channel).append("_tanb").append(".eps").c_str()); 
  return;
}

void BoundaryFinder::StoreDataFromGraph(const TGraph &gr)
{
    Double_t *x = gr.GetX();
    Double_t *y = gr.GetY();

    ComputeCenter(gr);

    Double_t xr=0, yr=0;

    for (UInt_t i=0; i<gr.GetN(); i++)
    {
        xr = x[i]-fCenter.first;
        yr = y[i]-fCenter.second;
        point_t p = {::sqrt(xr*xr + yr*yr), i};
        fSortedVals.insert(std::make_pair(::atan2(yr,xr),p));
    }
}

void BoundaryFinder::ComputeCenter(const TGraph& gr)
{
    Double_t *x = gr.GetX();
    Double_t *y = gr.GetY();
    UInt_t n = gr.GetN();

    Double_t sum_x = 0;
    Double_t sum_y = 0;
    for (UInt_t i=0; i<n; i++)
    {
        sum_x += x[i];
        sum_y += y[i];
    }

    fCenter.first = sum_x/n;
    fCenter.second = sum_y/n;
}

void drawtext()
{
   Int_t i,n;
   Double_t x,y;
   TLatex *l;

   TGraph *g = (TGraph*)gPad->GetListOfPrimitives()->FindObject("Graph");
   n = g->GetN();
   for (i=1; i<n; i++) {
      g->GetPoint(i,x,y);
      l = new TLatex(x,y+0.2,Form("%4.2f",y));
      l->SetTextSize(0.025);
      l->SetTextFont(42);
      l->SetTextAlign(21);
      l->Paint();
   }
}

    TGraph*
ContourGraph( TH2F* hist)
{
    TGraph* gr0 = new TGraph();
    TH2F* h = (TH2F*)hist->Clone();
    gr = (TGraph*)gr0->Clone(h->GetName());
    //  cout << "==> Will dumb histogram: " << h->GetName() << " into a graph" <<endl;
    h->SetContour( 1 );
    double pval = CombinationGlob::cl_percent[1];
    double signif = TMath::NormQuantile(1-pval);
    h->SetContourLevel( 0, signif );
    h->Draw("CONT LIST");
    h->SetDirectory(0);
    gPad->Update();
    TObjArray *contours = gROOT->GetListOfSpecials()->FindObject("contours");
    Int_t ncontours     = contours->GetSize();
    TList *list = (TList*)contours->At(0);
    Int_t number_of_lists = list->GetSize();
    gr = (TGraph*)list->At(0);
    TGraph* grTmp = new TGraph();
    for (int k = 0 ; k<number_of_lists ; k++){
        grTmp = (TGraph*)list->At(k);
        Int_t N = gr->GetN();
        Int_t N_tmp = grTmp->GetN();
        if(N < N_tmp) gr = grTmp;
        //    mg->Add((TGraph*)list->At(k));
    }

    gr->SetName(hist->GetName());
    int N = gr->GetN();
    double x0, y0;

    //  for(int j=0; j<N; j++) {
    //    gr->GetPoint(j,x0,y0);
    //    cout << j << " : " << x0 << " : "<<y0 << endl;
    //  }
    //  //  gr->SetMarkerSize(2.0);
    //  gr->SetMarkerSize(2.0);
    //  gr->SetMarkerStyle(21);
    //  gr->Draw("LP");
    //  cout << "Generated graph " << gr << " with name " << gr->GetName() << endl;
    return gr;
}