C++ TAxis::GetBinWidth方法代码示例

//______________________________________________________________________________
TEveCaloDataVec* MakeVecData(Int_t ncells=0)
{
   // Example how to fill data when bins can be iregular.
   // If ncells = 0 (default) whole histogram is taken,
   // otherwise just ncells cells around the maximum.

   TFile::SetCacheFileDir(".");
   TFile* hf = TFile::Open(histFile, "CACHEREAD");
   TH2F* h1 = (TH2F*)hf->Get("ecalLego");
   TH2F* h2 = (TH2F*)hf->Get("hcalLego");

   TEveCaloDataVec* data = new TEveCaloDataVec(2);
   data->RefSliceInfo(0).Setup("ECAL", 0.3, kRed);
   data->RefSliceInfo(1).Setup("HCAL", 0.1, kBlue);

   TAxis *ax =  h1->GetXaxis();
   TAxis *ay =  h1->GetYaxis();

   Int_t xm = 1, xM = ax->GetNbins();
   Int_t ym = 1, yM = ay->GetNbins();
   if (ncells != 0)
   {
      Int_t cx, cy, cz;
      h1->GetMaximumBin(cx, cy, cz);
      xm = TMath::Max(xm, cx-ncells);
      xM = TMath::Min(xM, cx+ncells);
      ym = TMath::Max(ym, cy-ncells);
      yM = TMath::Min(yM, cy+ncells);
   }

   // Take every second cell and set a random size.
   for(Int_t i=xm; i<=xM; i+=2)
   {
      for(Int_t j=ym; j<=yM; j+=2)
      {
         if ( (i+j) % 3)
         {
            data->AddTower(ax->GetBinLowEdge(i), ax->GetBinUpEdge(i),
                           ay->GetBinLowEdge(j), ay->GetBinUpEdge(j));
            data->FillSlice(0, h1->GetBinContent(i, j));
            data->FillSlice(1, h2->GetBinContent(i, j));
         }
         else
         {
            data->AddTower(ax->GetBinLowEdge(i),
                           2 * ax->GetBinWidth(i) + ax->GetBinLowEdge(i),
                           ay->GetBinLowEdge(j),
                           2 * ay->GetBinWidth(j) + ay->GetBinLowEdge(j));
            data->FillSlice(0, h2->GetBinContent(i, j));
            data->FillSlice(1, h2->GetBinContent(i, j));
         }
      }
   }

   data->SetEtaBins(ax);
   data->SetPhiBins(ay);
   data->DataChanged();
   return data;
}

void GetInt(TH1F *hzj,float xmin,float xmax,double& integralfull, double& integral)
{
  TAxis *axis = hzj->GetXaxis();
  int bmin = axis->FindBin(xmin); 
  int bmax = axis->FindBin(xmax); 
  integralfull = hzj->Integral();
  integral = hzj->Integral(bmin,bmax);
  integral -= hzj->GetBinContent(bmin)*(xmin-axis->GetBinLowEdge(bmin))/axis->GetBinWidth(bmin);
  integral -= hzj->GetBinContent(bmax)*(axis->GetBinUpEdge(bmax)-xmax)/axis->GetBinWidth(bmax);
  fout<<"Full Integral: "<<integralfull<<endl;
}

double integrate(TH1F* h, double xmin, double xmax){
  TAxis *axis = h->GetXaxis();
  int bmin = axis->FindBin(xmin); //in your case xmin=-1.5
  int bmax = axis->FindBin(xmax); //in your case xmax=0.8
  double integral = h->Integral(bmin,bmax);
  integral -= h->GetBinContent(bmin)*(xmin-axis->GetBinLowEdge(bmin))/axis->GetBinWidth(bmin);
  integral -= h->GetBinContent(bmax)*(axis->GetBinUpEdge(bmax)-xmax)/axis->GetBinWidth(bmax);

  return integral;

}

double integralInRange( TH1F * var_h, float lower, float upper ){
  TAxis *axis = var_h->GetXaxis();
  int bmin = axis->FindBin(lower);
  int bmax = axis->FindBin(upper);
  double integr = var_h->Integral(bmin,bmax);
  integr -= ( (var_h->GetBinContent(bmin)) * (lower - (axis->GetBinLowEdge(bmin))) )/
    axis->GetBinWidth(bmin);
  integr -= ( (var_h->GetBinContent(bmax)) * ( (axis->GetBinUpEdge(bmax)) - upper) )/
    axis->GetBinWidth(bmax);
  return integr;
}

void GetIntg(TH1F *hzj,float xmin,float xmax,double& count1, double& count2)
{
 
  TAxis *axis = hzj->GetXaxis();
  int bmin = axis->FindBin(xmin); 
  int bmax = axis->FindBin(xmax); 
  double integralfull = hzj->Integral();
  double integral = hzj->Integral(bmin,bmax);
  integral -= hzj->GetBinContent(bmin)*(xmin-axis->GetBinLowEdge(bmin))/axis->GetBinWidth(bmin);
  integral -= hzj->GetBinContent(bmax)*(axis->GetBinUpEdge(bmax)-xmax)/axis->GetBinWidth(bmax);
  count1=count1+integralfull;
  count2=count2+integral;
  // cout<<count1<<'\t'<<count2<<endl;

}

TH1 *PullHisto(const TList *list, const char *name, Int_t min, Int_t max,
               Double_t &mean)
{
  THnSparse *hs = list->FindObject(name);
  if (!hs) return 0;
  TAxis *atmp = hs->GetAxis(1);
  atmp->SetRange(min, max);
  // !!!!!!!!!!!!!!!!!!!!
  hs->GetAxis(2)->SetRangeUser(-0.5, 0.5);
  TH1 *hfin = hs->Projection(0);
  hfin->SetTitle(Form("p_{t} #in (%4.2f, %4.2f) GeV/c",
                      atmp->GetBinLowEdge(min),
                      atmp->GetBinLowEdge(max) + atmp->GetBinWidth(max)));
  mean = atmp->GetBinLowEdge(min) +
    (atmp->GetBinLowEdge(max) + atmp->GetBinWidth(max) -
     atmp->GetBinLowEdge(min))/2.0;
  // !!!!!!!!!!!!!!!!!!!!
  return hfin;//->Rebin();
}

void divideByBinWidth(TH1* histogram)
{
  if ( !histogram ) return;
  TAxis* xAxis = histogram->GetXaxis();
  int numBins = xAxis->GetNbins();
  for ( int iBin = 1; iBin <= numBins; ++iBin ) {
    double binContent = histogram->GetBinContent(iBin);
    double binError = histogram->GetBinError(iBin);
    double binWidth = xAxis->GetBinWidth(iBin);
    histogram->SetBinContent(iBin, binContent/binWidth);
    histogram->SetBinError(iBin, binError/binWidth);
  }
}

TH1* divideHistogramByBinWidth(TH1* histogram)
{
  std::string histogramDensityName = Form("%s_density", histogram->GetName());
  TH1* histogramDensity = (TH1*)histogram->Clone(histogramDensityName.data());
  TAxis* xAxis = histogram->GetXaxis();
  int numBins = xAxis->GetNbins();
  for ( int iBin = 1; iBin <= numBins; ++iBin ) {
    double binContent = histogram->GetBinContent(iBin);
    double binError = histogram->GetBinError(iBin);
    double binWidth = xAxis->GetBinWidth(iBin);
    histogramDensity->SetBinContent(iBin, binContent/binWidth);
    histogramDensity->SetBinError(iBin, binError/binWidth);
  }
  return histogramDensity;
}

//________________________________________________________________________
void CorrectFromTheWidth(TH1D *h1) {
  //
  // Correct from the width of the bins --> dN/dp_{T} (GeV/c)^{-1}
  //

  TAxis *axis = h1->GetXaxis();
  Int_t nbinX = h1->GetNbinsX();

  for(Int_t i = 1; i <= nbinX; i++) {

    Double_t width = axis->GetBinWidth(i);
    Double_t content = h1->GetBinContent(i);
    Double_t error = h1->GetBinError(i); 
    h1->SetBinContent(i,content/width); 
    h1->SetBinError(i,error/width);
  }

}

void rebin() {
   //create a fix bin histogram
   TH1F *h = new TH1F("h","test rebin",100,-3,3);
   Int_t nentries = 1000;
   h->FillRandom("gaus",nentries);
   Double_t xbins[1001];
   Int_t k=0;
   TAxis *axis = h->GetXaxis();
   for (Int_t i=1;i<=100;i++) {
      Int_t y = (Int_t)h->GetBinContent(i);
      if (y <=0) continue;
      Double_t dx = axis->GetBinWidth(i)/y;
      Double_t xmin = axis->GetBinLowEdge(i);
      for (Int_t j=0;j<y;j++) {
         xbins[k] = xmin +j*dx;
         k++;
      }
   }
   xbins[k] = axis->GetXmax();
   //create a variable bin-width histogram out of fix bin histogram
   //new rebinned histogram should have about 10 entries per bin
   TH1F *hnew = new TH1F("hnew","rebinned",k,xbins);
   hnew->FillRandom("gaus",10*nentries);

   //rebin hnew keeping only 50% of the bins
   Double_t xbins2[501];
   Int_t kk=0;
   for (Int_t j=0;j<k;j+=2) {
      xbins2[kk] = xbins[j];
      kk++;
   }
   xbins2[kk] = xbins[k];
   TH1F *hnew2 = (TH1F*)hnew->Rebin(kk,"hnew2",xbins2);

   //draw the 3 histograms
   TCanvas *c1 = new TCanvas("c1","c1",800,1000);
   c1->Divide(1,3);
   c1->cd(1);
   h->Draw();
   c1->cd(2);
   hnew->Draw();
   c1->cd(3);
   hnew2->Draw();
}

void th22mama(TH2* m, const char* filename, const char* comment="none")
{
    char tmp[128];
    time_t now = time(0);
    strftime(tmp, sizeof(tmp), "%Y-%m-%d %T", localtime(&now)); // not the 23-Mar-07 16:02:34 format
    
    TAxis *xax = m->GetXaxis(), *yax = m->GetYaxis();
    
    const int nx = std::min(xax->GetNbins(), 2048);
    const int ny = yax->GetNbins();
    std::cout << "matrix is " << nx << 'x' << ny << "; comment='" << comment << "'" << std::endl;
    
    std::ofstream mama(filename);
    mama << "!FILE=Disk \n"
    "!KIND=Matrix \n"
    "!LABORATORY=Oslo Cyclotron Laboratory (OCL) \n"
    "!EXPERIMENT=siri2root \n"
    "!COMMENT=" << comment << "\n"
    "!TIME=" << tmp << "\n"
    "!CALIBRATION EkeV=6";
    const float cal[6] = { xax->GetBinLowEdge(1), xax->GetBinWidth(1), 0,
        yax->GetBinLowEdge(1), yax->GetBinWidth(1), 0 };
    for(int i=0; i<6; ++i) {
        snprintf(tmp, sizeof(tmp), ",%13.6E", cal[i]);
        mama << tmp;
    }
    mama << "\n!PRECISION=16\n";
    snprintf(tmp, sizeof(tmp), "!DIMENSION=2,0:%4d,0:%4d\n", nx-1, ny-1);
    mama << tmp;
    snprintf(tmp, sizeof(tmp), "!CHANNEL=(0:%4d,0:%4d)\n", nx-1, ny-1);
    mama << tmp;
    
    for(int iy=1; iy<=ny; ++iy) {
        for(int ix=1; ix<=nx; ++ix)
            mama << m->GetBinContent(ix, iy) << ' ';
        mama << "\n";
    }
    mama << "!IDEND=" << endl << endl;
    
    mama.close();
}

void DrawMeasurement( TH1 *h, Int_t nbin, Double_t val, Double_t err )
{
   TAxis *axis   = h->GetYaxis();
   
   Double_t dy;
   Double_t x1,y1,x2,y2, ymin;
   
   Int_t i = nbin;
   
   ymin = axis->GetBinCenter(i)+1;
   dy   = axis->GetBinWidth(i);
   x1   = val - err;
   y1   = ymin - 0.05*dy;
   x2   = val + err;
   y2   = ymin + 0.05*dy;
  
   TLine* line = new TLine;
   TLine* lbar = new TLine;
   TLine* rbar = new TLine;
   TMarker* m2 = new TMarker( val, ymin, 20 );
   Int_t color = 1;
   line->SetLineColor( color );
   line->SetLineWidth( 2 );
   lbar->SetLineColor( color );
   lbar->SetLineWidth( 2 );
   rbar->SetLineColor( color );
   rbar->SetLineWidth( 2 );
   m2->SetMarkerColor( color );
   m2->SetMarkerSize( 1.2 );
   m2->SetMarkerStyle( 21 );
   
   line->DrawLine( x1, ymin, x2, ymin );
   lbar->DrawLine( x1, y1, x1, y2 );
   rbar->DrawLine( x2, y1, x2, y2 );
   m2->Draw();
}

Double_t effSigma(TH1 *hist )
{

  TAxis *xaxis = hist->GetXaxis();
  Int_t nb = xaxis->GetNbins();
  if(nb < 10) {
    cout << "effsigma: Not a valid histo. nbins = " << nb << endl;
    return 0.;
  }
  
  Double_t bwid = xaxis->GetBinWidth(1);
  if(bwid == 0) {
    cout << "effsigma: Not a valid histo. bwid = " << bwid << endl;
    return 0.;
  }
  Double_t xmin = xaxis->GetXmin();
  Double_t ave = hist->GetMean();
  Double_t rms = hist->GetRMS();

  Double_t total=0.;
  for(Int_t i=0; i<nb+2; i++) {
    total+=hist->GetBinContent(i);
  }
  if(total < 100.) {
    cout << "effsigma: Too few entries " << total << endl;
    return 0.;
  }
  Int_t ierr=0;
  Int_t ismin=999;
  
  Double_t rlim=0.683*total;
  Int_t nrms=rms/(bwid);    // Set scan size to +/- rms
  if(nrms > nb/10) nrms=nb/10; // Could be tuned...

  Double_t widmin=9999999.;
  for(Int_t iscan=-nrms;iscan<nrms+1;iscan++) { // Scan window centre
    Int_t ibm=(ave-xmin)/bwid+1+iscan;
    Double_t x=(ibm-0.5)*bwid+xmin;
    Double_t xj=x;
    Double_t xk=x;
    Int_t jbm=ibm;
    Int_t kbm=ibm;
    Double_t bin=hist->GetBinContent(ibm);
    total=bin;
    for(Int_t j=1;j<nb;j++){
      if(jbm < nb) {
        jbm++;
        xj+=bwid;
        bin=hist->GetBinContent(jbm);
        total+=bin;
        if(total > rlim) break;
      }
      else ierr=1;
      if(kbm > 0) {
        kbm--;
        xk-=bwid;
        bin=hist->GetBinContent(kbm);
        total+=bin;
        if(total > rlim) break;
      }
      else ierr=1;
    }
    Double_t dxf=(total-rlim)*bwid/bin;
    Double_t wid=(xj-xk+bwid-dxf)*0.5;
    if(wid < widmin) {
      widmin=wid;
      ismin=iscan;
    }   
  }
  if(ismin == nrms || ismin == -nrms) ierr=3;
  if(ierr != 0) cout << "effsigma: Error of type " << ierr << endl;
  
  return widmin;
  
}

void make( string charge = "p", bool lines = false, int iCen = 0, string cs = "Pi", string hFile = "histograms.root" ){

	Bichsel bgen;
	gStyle->SetPadLeftMargin(0.16);
	gStyle->SetPadBottomMargin(0.12);
	gStyle->SetPadRightMargin(0.12);

	set_plot_style();

	RooPlotLib rpl;
	string rpName = "rp_" + charge + "_dedx_beta_lines.pdf";
	if ( false == lines )
		rpName = "rp_" + charge + "_dedx_beta.pdf";

	string chargeName = "Positive Tracks : ";
	if ( "n" == charge )
		chargeName = "Negative Tracks : ";

	Reporter rp( rpName, 600, 800 );

	TFile * f = new TFile( hFile.c_str(), "READ" );

	TH1D * nlBeta = (TH1D*)f->Get( "nlBeta" );


	for ( int iPt = 0; iPt < 100; iPt ++ ){
		if ( iPt != 16 ) continue;


		TAxis * x = nlBeta->GetXaxis();
		double lpT = x->GetBinLowEdge( iPt + 1 );
		double hpT = x->GetBinLowEdge( iPt + 1 ) + x->GetBinWidth( iPt + 1 );
		double avgP = (lpT + hpT) / 2.0;

		string name = "dedx_tof/dedx_tof_" + cs + "_" + charge + "_" + ts(iCen) + "_" + ts( iPt);
		TH2D * h2 = (TH2D*)f->Get( name.c_str() );

		if ( !h2 )
			break;
		if ( iPt == 5 ) // what went wrong?
			continue;


		double x1 = h2->GetXaxis()->GetBinLowEdge( 1 );
		double x2 = h2->GetXaxis()->GetBinLowEdge( h2->GetNbinsX() ) + h2->GetXaxis()->GetBinWidth( h2->GetNbinsX() );
		double y1 = h2->GetYaxis()->GetBinLowEdge( 1 );
		double y2 = h2->GetYaxis()->GetBinLowEdge( h2->GetNbinsY() ) + h2->GetYaxis()->GetBinWidth( h2->GetNbinsY() );

		  

		rpl.style( h2 ).set( "draw", "col" )
			.set( "optstat", 0 ).set( "logz", 1 )
			.set( "title", " ; ln(dE/dx) - ln(dE/dx)_{#pi}  ; #beta^{-1} - #beta^{-1}_{#pi} " )
			.set( "ts", 0.16 )
			.set( "xts", 0.06 )
			.set( "xto", 0.75 )
			.set( "yts", 0.07 )
			.set( "yto", 0.85 )
			.set( "yoffset", 0.01 ).draw(); 



		if ( 16 == iPt ){
			gStyle->SetTitleFontSize( 0.45 );
			rpl.style(h2)
				.set( "xr", -0.3, 0.95 )
				.set( "yr", -0.15, 0.6 ).draw();
		}
		TLatex *text = new TLatex( -0.5, 0.64, (chargeName + dts(lpT) + " < p_{T} [GeV/c] < " + dts(hpT)).c_str() );
		text->SetTextSize(0.055);
		text->Draw("");


		if ( lines  ){
			double bKaon = one_beta( mK, avgP ) - one_beta( mPi, avgP );
			double bProton = one_beta( mP, avgP ) - one_beta( mPi, avgP );

			double dKaon = bgen.meanLog( avgP, mK, -1, 1000 ) - bgen.meanLog( avgP, mPi, -1, 1000 );
			double dProton = bgen.meanLog( avgP, mP, -1, 1000 ) - bgen.meanLog( avgP, mPi, -1, 1000 );

			drawTofLines( 0, kRed+1, x1, x2 );		
			drawTpcLines( 0, kRed + 1, y1, y2);
			// around Kaons
			drawTofLines( bKaon, kOrange+1, x1, x2 );	
			drawTpcLines( dKaon, kOrange + 1, y1, y2);
			// around Protons
			drawTofLines( bProton, kBlack, x1, x2 );	
			drawTpcLines( dProton, kBlack, y1, y2);
		}

		// double dElec = bgen.meanLog( avgP, mE, -1, 1000 ) - bgen.meanLog( avgP, mPi, -1, 1000 );
		// double bElec = one_beta( mE, avgP ) - one_beta( mPi, avgP );
		// //drawTpcLines( dMerged, kBlack, y1, y2);
		// TEllipse * ell = new TEllipse( dElec,bElec, sigTpc * 3, sigTof * 3 );
		// ell->SetFillColorAlpha( kBlack, 0 );
		// ell->Draw();
		// 
		
		gPad->SetRightMargin( 0.01 );
		gPad->SetBottomMargin( 0.10 );
//.........这里部分代码省略.........