C++ TF1::SetParName方法代码示例

TF1* Mathematics::landau_func()
{
	//  Want plotting consistancy so have defined my own landau function
	TF1 *land = new TF1( "mylandau", "[0]*TMath::Landau( x, [1], [2] )" );
	land->SetParName( 0, "Constant" );
	land->SetParName( 1, "Mean" );
	land->SetParName( 2, "Sigma" );
	return land;
}

	  TF1 * getDoubleFit(TH1 * hist, bool autorange, double inner, double outer, int which=1){
		  /**
		   * @brief performes a double gaussian Fit
		   */

		  Double_t parameters[6] = {0,0,0,0,0,0};
		  Double_t parameters_new[6] = {0,0,0,0,0,0};


		  double center = hist->GetMean();

		  if (autorange){
			  double inner = fabs(center - hist->GetXaxis()->GetXmax()) /10;
			  double outer = fabs(center - hist->GetXaxis()->GetXmax()) /10*8;
		  }

		  TF1 * fit1 = new TF1("fit1", "gaus", center-inner, center+inner );

		  hist->Fit(fit1, "Q0R");
		  fit1->GetParameters(&parameters[0]);

		  TF1 * fit2 = new TF1("fit2", "gaus", center-outer, center+outer );

		  hist->Fit(fit2, "Q0R");
		  fit2->GetParameters(&parameters[3]);

		  TF1 * fitproper = new TF1("fitproper", "gaus(0)+gaus(3)", center-outer, center+outer);
		  fitproper->SetParameters(parameters);
		  fitproper->SetParName(0, "Const.(inner)");
		  fitproper->SetParName(1, "Mean (inner)");
		  fitproper->SetParName(2, "Sigma (inner)");
		  fitproper->SetParName(3, "Const.(outer)");
		  fitproper->SetParName(4, "Mean (outer)");
		  fitproper->SetParName(5, "Sigma (outer)");

		  hist->Fit(fitproper, "Q0R");
		  fitproper->GetParameters(&parameters_new[0]);


		  fitproper->SetLineColor(hist->GetLineColor());
		  fitproper->SetLineWidth(hist->GetLineWidth());
		  fitproper->SetLineStyle(2);

		  if (which==1){
		  	  TF1 * Gausinner = new TF1("Gausinner", "gaus(0)", center-inner, center+inner);
			  Gausinner->SetParameters(parameters_new);
		  	  return Gausinner;
		  }
		  else{
			  TF1 * Gausouter = new TF1("Gausouter", "gaus(3)", center-outer, center+outer);
			  Gausouter->SetParameters(parameters_new);
			  return Gausouter;
		  }

	  }

TF1* Mathematics::raw_gamma_func( int OutputLevel )
{
	(void) OutputLevel;
	TF1* output = new TF1( "gammaf", "[0]*TMath::GammaDist( x, [1], [2], [3] )" );

	output->SetParName( 0, "Const" );
	output->SetParName( 1, "#gamma" );
	output->SetParName( 2, "#mu" );
	output->SetParName( 3, "#beta" );

	return output;
}

TF1* GetFit(const char *address, int lyr, double xinit) {
  const char *bgr = Form("[6]*TMath::Exp([7]*(x-%f))",xinit);
  const char *la1 = "(1-[3]-[4]-[5])*TMath::Landau(x,[1],[2],1)";
  const char *la2 = "[3]*TMath::Landau(x,2*[1]+1.4*[2],2.0*[2],1)";
  const char *la3 = "[4]*TMath::Landau(x,3*[1]+3.3*[2],3.0*[2],1)";
  const char *la4 = "[5]*TMath::Landau(x,4*[1]+5.6*[2],4.0*[2],1)";
  TF1 *ret = new TF1("fit_H1",
                     Form("[0]*(%s+%s+%s+%s)+%s",la1,la2,la3,la4,bgr) );
  ret->SetParName(0,"A");
  ret->SetParName(1,"lambda");
  ret->SetParName(2,"sigma");
  ret->SetParName(3,"f2");
  ret->SetParameter(0,1e4);  ret->SetParLimits(0,1e2,1e7);
  ret->SetParameter(1,21);   ret->SetParLimits(1,13,25.0);
  ret->SetParameter(2,3.0);  ret->SetParLimits(2,1.5,5.0);
  if(lyr>1) {
    ret->SetParameter(1,15);   ret->SetParLimits(1,2,30);
    ret->SetParameter(2,5.0);  ret->SetParLimits(2,0.8,10.0);
  }
  ret->SetParameter(3,0.20); ret->SetParLimits(3,0,0.40);
  ret->SetParameter(4,0.0);
  ret->SetParameter(5,0.0);
  ret->SetParameter(6,500); ret->SetParLimits(6,1e2,1e7);
  ret->SetParameter(7,-1);  ret->SetParLimits(7,-10,-0.1);

  ret->SetLineColor(kRed-3);

  ifstream infit;
  infit.open( address );
  double tmp;
  bool found=false;
  for(int n=0; n!=ret->GetNpar(); ++n) {
    infit >> tmp;
    if(!infit.good()) break;
    ret->SetParameter(n,tmp);
    infit >> tmp;
    found = true;
  }
  infit.close();
  if(found) {
    cout << " Previous fit results found" << endl;
    ret->SetParLimits(4,0,0.40);
    ret->SetParLimits(5,0,0.40);
  } else {
    ret->SetParLimits(4,+1,-1);
    ret->SetParLimits(5,+1,-1);
  }
  return ret;
}

int main(int argc, char const *argv[]) {
	TFile *myFile = new TFile("output.root","RECREATE");

	Float_t E1,Px1,Py1,Pz1,E2,Px2,Py2,Pz2;
	Double_t mass = 0;

	TH1D *m12 = new TH1D("m12","m12",60,0.99,1.08);
	TF1 *myBW = new TF1("myBW","TMath::BreitWigner(x,[0],[1])", 1.01, 1.03);

	myBW->SetParName(0,"Mass");
	myBW->SetParName(1,"#Gamma");

	TFile *inputFile = new TFile("data.root");
	TTree *data_tree = (TTree*)inputFile->Get("data");

	data_tree->SetBranchAddress("E1", &E1);
	data_tree->SetBranchAddress("Px1", &Px1);
	data_tree->SetBranchAddress("Py1", &Py1);
	data_tree->SetBranchAddress("Pz1", &Pz1);
	data_tree->SetBranchAddress("E2", &E2);
	data_tree->SetBranchAddress("Px2", &Px2);
	data_tree->SetBranchAddress("Py2", &Py2);
	data_tree->SetBranchAddress("Pz2", &Pz2);

	TLorentzVector vec1(0.0,0.0,0.0,0.0);
	TLorentzVector vec2(0.0,0.0,0.0,0.0);

	int count = data_tree->GetEntries();

	for (int i = 0; i < count; i++) {	
		data_tree->GetEntry(i);
		vec1.SetPxPyPzE(Px1,Py1,Pz1,E1);
		vec2.SetPxPyPzE(Px2,Py2,Pz2,E2);
		mass = (vec1 + vec2).M();

		m12->Fill(mass);
	}
	myBW->SetParameters(0,1.02,0.05);
	myBW->SetParLimits(0,1.01,1.03);
	m12->Fit("myBW","","",0,2);
	m12->GetXaxis()->SetTitle("Mass (GeV)");
	gStyle->SetOptFit(1111);

	myFile->cd();
	myFile->Write();
	myFile->Close();

	return 0;
}

//------------------------------------------------------------------------------
std::pair<double,double> fitskwlin( char* hs, double xl=-0.1, double xr=0.1 ) {
  TH1 *h = (TH1*)gDirectory->Get(hs);

  if( h == NULL ) {
    cout << hs << " does not exist\n";
    return std::make_pair(0.,0.);
  }

  int nb = h->GetNbinsX();
  double x1 = h->GetBinCenter(1); // first
  double x9 = h->GetBinCenter(nb); // last

  if( xl > x1 && xl < x9 ) x1 = xl; // left
  if( xr > x1 && xr < x9 ) x9 = xr; // right

  // create a TF1 with the range from x1 to x9 and 3 parameters
  TF1 *tanhFcn = new TF1( "tanhFcn", fitSkw, x1, x9, 2 );
  tanhFcn->SetParName( 0, "dyoff" );
  tanhFcn->SetParName( 1, "dyslp" );

  // set start values:
  tanhFcn->SetParameter( 0, 0 ); // dy off [um]
  tanhFcn->SetParameter( 1, 99 ); // dy slope [um/skw]

  h->Fit( "tanhFcn", "R Q", "p" );// R = range from tanhFcn

  return std::make_pair(tanhFcn->GetParameter(0),tanhFcn->GetParameter(1));
}

	  TF1 * tripleGaussFit(TH1 * hist, bool autorange, double inner=0.01, double outer=0.1){
		  /**
		   * @brief performes a double gaussian Fit
		   */
		  Double_t parameters[9] = {0,0,0,0,0,0,0,0,0};


		  double center = hist->GetMean();


		  if (autorange){
			  double inner = fabs(center - hist->GetXaxis()->GetXmax()) /10;
			  double outer = fabs(center - hist->GetXaxis()->GetXmax()) /10*8;
		  }
		  double middle = (outer - inner)/2;

		  TF1 * fit1 = new TF1("fit1", "gaus", center-inner, center+inner );

		  hist->Fit(fit1, "0R");
		  fit1->GetParameters(&parameters[0]);

		  TF1 * fitmid = new TF1("fitmid", "gaus", center-middle, center+middle );
		  fitmid->SetParameter(1, parameters[1]);
		  hist->Fit(fitmid, "0R");
		  fitmid->GetParameters(&parameters[3]);

		  TF1 * fit2 = new TF1("fit2", "gaus", center-outer, center+outer );
		  fit2->SetParameter(1, parameters[1]);
		  hist->Fit(fit2, "0R");
		  fit2->GetParameters(&parameters[6]);

		  TF1 * fitproper = new TF1("fitproper", "gaus(0)+gaus(3)+gaus(6)", center-outer, center+outer);
		  fitproper->SetParameters(parameters);
		  fitproper->SetParName(0, "Const.(inner)");
		  fitproper->SetParName(1, "Mean (inner)");
		  fitproper->SetParName(2, "Sigma (inner)");
		  fitproper->SetParName(3, "Const.(middle)");
		  fitproper->SetParName(4, "Mean (middle)");
		  fitproper->SetParName(5, "Sigma (middle)");
		  fitproper->SetParName(6, "Const.(outer)");
		  fitproper->SetParName(7, "Mean (outer)");
		  fitproper->SetParName(8, "Sigma (outer)");

		  hist->Fit(fitproper, "0R");

//		  fitproper->SetLineColor(hist->GetLineColor());
//		  fitproper->SetLineWidth(hist->GetLineWidth());
////		  fitproper->SetLineStyle(2);

		  return fitproper;

	  }

		TF1 * gaussFit(TH1 * hist, double inner) {
			Double_t parameters[6] = {0,0,0,0,0,0};


			double center =0;// hist->GetMean();

			TF1 * fit = new TF1("fit", "gaus", center-inner, center+inner );

			hist->Fit(fit, "Q0R");
			fit->GetParameters(&parameters[0]);
			fit->SetParameters(parameters);
			fit->SetParName(0, "Const.");
			fit->SetParName(1, "Mean");
			fit->SetParName(2, "Sigma");

			return fit;
		}

void plotMttResolution() {

  setTDRStyle();

  TCanvas *c = new TCanvas("c", "c", 800, 800);

  TH1F *mtt = (TH1F*) _file0->Get("mtt_resolution_mva");
  mtt->SetTitle("");
  mtt->SetName("m_{t#bar{t}}");
  mtt->GetXaxis()->SetTitle("m_{t#bar{t}} - m_{t#bar{t}}^{gen} (GeV)");
  mtt->SetLineColor(TColor::GetColor("#542437"));
  mtt->SetLineColor(TColor::GetColor("#8A9B0F"));
  mtt->SetLineColor(TColor::GetColor("#C02942"));
  mtt->SetFillColor(mtt->GetLineColor());
  mtt->SetFillStyle(3004);
  mtt->GetXaxis()->SetTitleOffset(1.2);

  TF1 *voigt = new TF1("voigt", "gaus", -100, 100);
	voigt->SetParName(0, "amp");
	voigt->SetParName(2, "mean");
	voigt->SetParName(1, "sigma");

  voigt->SetParameter(0, mtt->GetMaximum());
	voigt->SetParameter(2, 0);
	voigt->SetParameter(1, 100);

  mtt->Fit(voigt, "RVN");

  voigt->SetLineColor(TColor::GetColor("#8A9B0F"));
  voigt->SetLineWidth(2);

  mtt->Draw();
  voigt->Draw("same");

  TLatex* latex = new TLatex();
  latex->SetTextFont(42);
  latex->SetTextSize(0.033);
  latex->DrawLatexNDC(0.25, 0.84, TString::Format("mean = %.2f", voigt->GetParameter(1)));
  latex->DrawLatexNDC(0.25, 0.8, TString::Format("#sigma = %.2f", voigt->GetParameter(2)));

  gPad->Modified();
  gPad->Update();

  c->Print("mtt_resolution.pdf");
}

Double_t fitfulllang( char* hs ) {

  TH1 *h = (TH1*)gDirectory->Get(hs);

  if( h == NULL ){
    cout << hs << " does not exist\n";
    return 0;
  }

  double aa = h->GetEntries();//normalization

  // find peak:
  int ipk = h->GetMaximumBin();
  double xpk = h->GetBinCenter(ipk);
  double sm = xpk / 9; // sigma
  double ns = sm; // noise

  // fit range:
  int ib0 = ipk/2;
  int ib9 = h->GetNbinsX() - 1;

  double x0 = h->GetBinLowEdge(ib0);
  double x9 = h->GetBinLowEdge(ib9) + h->GetBinWidth(ib9);

  // create a TF1 with the range from x0 to x9 and 4 parameters
  TF1 *fitFcn = new TF1( "fitFcn", fitLandauGauss, x0, x9, 4 );

  fitFcn->SetParName( 0, "peak" );
  fitFcn->SetParName( 1, "sigma" );
  fitFcn->SetParName( 2, "area" );
  fitFcn->SetParName( 3, "smear" );

  fitFcn->SetNpx(500);
  fitFcn->SetLineWidth(4);
  fitFcn->SetLineColor(kMagenta);

  // set start values:
  fitFcn->SetParameter( 0, xpk ); // peak position, defined above
  fitFcn->SetParameter( 1, sm ); // width
  fitFcn->SetParameter( 2, aa ); // area
  fitFcn->SetParameter( 3, ns ); // noise

  h->Fit("fitFcn", "NQR", "ep" );// R = range from fitFcn
  return fitFcn->GetParameter(0);
}

	  TF1 * GetVoigtFit(TH1 * hist, double lower, double upper){
		  /**
		   * @brief performes a Fit with a Breit-Wigner-Distribution.
		  */
		  double parameter[3];

		  double max = hist->GetMaximum();

		  int bin1 = hist->FindFirstBinAbove(max-1);

		  double center =  hist->GetBinCenter(bin1);


		  TString func = TString::Format("TMath::Voigt(x-%.4f, [0], [1])", center);

		  TF1 * fit = new TF1("fit", func, lower, upper);

		  hist->Fit(fit, "0R");

		  fit->GetParameters(&parameter[0]);


		  TF1 * voigt = new TF1("voigt", func+"*[2]", lower, upper);

		  voigt->SetParameter(0, parameter[0]);
		  voigt->SetParameter(1, parameter[1]);

		  voigt->SetParName(0, "#sigma");
		  voigt->SetParName(1, "#Gamma");
		  voigt->SetParName(2, "A");
//		  voigt->SetParName(3, "B");
//		  voigt->SetParName(4, "C");

		  hist->Fit(voigt, "0R");

		  return voigt;



	  }

void countZs(const char* filename, int run1=1, int run2=0, float lumi=0, bool add=false) {
   TFile *file = new TFile(filename);
   if (!file) return;
   
   TCanvas *c1 = new TCanvas();
   c1->cd();

   TH1F *h = (TH1F*) file->Get("Run summary/DiMuonHistograms/GlbGlbMuon_HM");
   h->GetYaxis()->SetRangeUser(0,200);
   h->Draw();

   TF1 *f = new TF1("f","gaus(0)+pol2(3)",70,110);
   f->SetParName(0,"norm");
   f->SetParName(1,"mass");
   f->SetParName(2,"width");
   f->SetParameters(1.16028e+02,9.07785e+01,2.11556e+00,-2.17600e+01,9.67443e-01,-6.98382e-03);
   f->SetParLimits(0,0,1e6);
   f->SetParLimits(1,85,95);
   f->SetParLimits(2,1.5,5);
   h->Fit(f,"","",70,110);
   h->GetYaxis()->SetRangeUser(0,200);
   c1->Update();

   TF1 *f2 = new TF1("f2","gaus(0)",70,110);
   f2->SetParameters(f->GetParameter(0),f->GetParameter(1),f->GetParameter(2));
   int integral = (int) (f2->Integral(70,110) / h->GetBinWidth(2));

   TLatex *tl = new TLatex();
   tl->DrawLatexNDC(0.1,0.95,Form("pPb 8TeV [%i-%i, %.1f nb^{-1}]",run1,run2,lumi));
   tl->DrawLatexNDC(0.18,0.75,Form("%i Z bosons",integral));
   tl->DrawLatexNDC(0.18,0.7,Form("(#sigma = %.1f nb)",integral/lumi));

   if (!add) c1->SaveAs(Form("plotZ_%i_%i.pdf",run1,run2));
   else c1->SaveAs("plotZ.gif+100");

   file->Close();
}

std::vector<double> fitskwpol( char* hs, double xl=-0.1, double xr=0.1) {
  TH1 *h = (TH1*)gDirectory->Get(hs);

  if( h == NULL ) {
    cout << hs << " does not exist\n";
    return std::vector<double>();
  }

  int nb = h->GetNbinsX();
  double x1 = h->GetBinCenter(1); // first
  double x9 = h->GetBinCenter(nb); // last

  if( xl > x1 && xl < x9 ) x1 = xl; // left
  if( xr > x1 && xr < x9 ) x9 = xr; // right

  // create a TF1 with the range from x1 to x9 and 3 parameters
  TF1 *tanhFcn = new TF1( "tanhFcn", fitSkwPol, x1, x9, 4 );
  tanhFcn->SetParName( 0, "dyoff" );
  tanhFcn->SetParName( 1, "dyslp" );
  tanhFcn->SetParName( 2, "dypar" );
  tanhFcn->SetParName( 3, "dyhyp" );

  // set start values:
  tanhFcn->SetParameter( 0, 0 ); // dy off [um]
  tanhFcn->SetParameter( 1, 99 ); // dy slope [um/skw]
  tanhFcn->SetParameter( 2, 0 ); // dy parabola
  tanhFcn->SetParameter( 3, 0 ); // dy hyperbola


  h->Fit( "tanhFcn", "R Q", "p" );// R = range from tanhFcn

  std::vector<double> result;
  for(size_t i = 0; i < 4; i++) {
    result.push_back(tanhFcn->GetParameter(i));
  }
  return result;
}

void fitTProfileCustom(TProfile* p)
{

  // Display fit info on canvas.
  gStyle->SetOptFit(0011);
  // Fit things.
  //        TF1* fit = new TF1("fit", "[0]*exp([4]*x)*sin([1]*x**0.5 + [2]) + [3]", -20, 20);
  TF1* fit = new TF1("fit", "[0]*(x**2)*sin([1]*x**0.5 + [2]) + [3]", -20, 20);

  fit->SetParName(0, "Scale");
  fit->SetParName(1, "Frequency");
  fit->SetParName(2, "Phase");
  fit->SetParName(3, "Shift");
  //        fit->SetParName(4, "Decay Length");

  fit->SetParameter(0, 0.006 );
  //        fit->SetParameter(0, 0.01 );
  fit->SetParameter(1, -7); 
  fit->SetParameter(2, 2.5);
  fit->SetParameter(3, 91);
  //        fit->SetParameter(4, 0.000795);

  bool converged = 0;
  int ntries = 0;

  // Make sure the fit converges.
  while(!converged)
    {   
      if(ntries >= 50) break;
      std::cout << "==== " << p->GetTitle() << " ====" << std::endl;
      p->Fit("fit", "", "", 20, 61);
      TString sconverge = gMinuit->fCstatu.Data();
      converged = sconverge.Contains(TString("CONVERGED"));
      ntries++;
    }   

}

TF1* Mathematics::gamma_func( int OutputLevel )
{
	//      the gamma function in ROOT has issues with verbosity, let's silence it and then return the verbosity back at the end

	streambuf *nullbuf=NULL, *cout_bak=NULL, *cerr_bak=NULL, *clog_bak=NULL;
	ofstream filestr;
	filestr.open ("/dev/null");
	//      If the user wanted silence we point the Std Output Streams to /dev/null
	if( OutputLevel <= -1 )
	{
		cout_bak = cout.rdbuf();
		cerr_bak = cerr.rdbuf();
		clog_bak = clog.rdbuf();
		nullbuf = filestr.rdbuf();
		cout.rdbuf(nullbuf);
		cerr.rdbuf(nullbuf);
		clog.rdbuf(nullbuf);
	}

	TF1* output = new TF1( "gammaf", "[0]*TMath::GammaDist( x, ([1]*[1])/([2]*[2]), 0, ([2]*[2])/[1] )" );

	//  TF1 inherits from TFormula so we can use it's functions to rename the parameters to have consistancy with gaus / landau functions
	output->SetParName( 0, "Constant" );
	output->SetParName( 1, "Mean" );
	output->SetParName( 2, "Sigma" );

	//      Reset Std Output Streams
	if( OutputLevel <= -1 )
	{
		cout.rdbuf(cout_bak);
		cerr.rdbuf(cerr_bak);
		clog.rdbuf(clog_bak);
	}

	return output;
}