C++ TH1F::SetName方法代码示例

std::pair<TH1F*,TH1F*> GetPMSyst(TH1F* trueHist, TH1F* intHist, std::string name){
  TH1F *diff = (TH1F*)trueHist->Clone();
  diff->Add(intHist,-1);
  diff->Divide(intHist);
  diff->Scale(3.);
  TH1F *up = (TH1F*)intHist->Clone();
  up->Add(diff,1);
  up->SetName((name+"Up01_sigma").c_str());
  TH1F *down = (TH1F*)intHist->Clone();
  down->Add(diff,-1);
  down->SetName((name+"Down01_sigma").c_str());

  std::pair<TH1F*,TH1F*> result(up,down);
  return result; 
}

TF1 * iterateFitter(TH1F * histo, int i, TF1 * previousResiduals = 0, TCanvas * inputCanvas = 0)
{
  TCanvas * canvas = inputCanvas;
  std::stringstream ss;
  int iCanvas = i;
  ss << i;

  int nBins = histo->GetNbinsX();

  TF1 * functionToFit;
  if( previousResiduals == 0 ) {
    // functionToFit = relativisticBWFit(ss.str());
    // functionToFit = relativisticBWintFit(ss.str());
    // functionToFit = relativisticBWintPhotFit(ss.str());
    functionToFit = expRelativisticBWintPhotFit(ss.str());
    // functionToFit = crystalBallFit();
    // functionToFit = reversedCrystalBallFit();
    // functionToFit = lorentzianFit();
    // functionToFit = relativisticBWFit();
  }
  else {
    functionToFit = combinedFit(previousResiduals->GetParameters(), ss.str());
  }
  histo->Fit(functionToFit, "MN", "", 60, 120);


  double xMin = histo->GetXaxis()->GetXmin();
  double xMax = histo->GetXaxis()->GetXmax();
  double step = (xMax-xMin)/(double(nBins));
  TH1F * functionHisto = new TH1F(("functionHisto"+ss.str()).c_str(), "functionHisto", nBins, xMin, xMax);
  for( int i=0; i<nBins; ++i ) {
    functionHisto->SetBinContent( i+1, functionToFit->Eval(xMin + (i+0.5)*step) );
  }

  if( canvas == 0 ) {
    canvas = new TCanvas(("canvasResiduals"+ss.str()).c_str(), ("canvasResiduals"+ss.str()).c_str(), 1000, 800);
    canvas->Divide(2,1);
    canvas->Draw();
    iCanvas = 0;
  }
  canvas->cd(1+2*iCanvas);
  histo->Draw();
  functionToFit->SetLineColor(kRed);
  functionToFit->Draw("same");
  // functionHisto->Draw("same");
  // functionHisto->SetLineColor(kGreen);
  canvas->cd(2+2*iCanvas);
  TH1F * residuals = (TH1F*)histo->Clone();
  residuals->Add(functionHisto, -1);
  residuals->SetName("Residuals");

  // TF1 * residualFit = new TF1("residualFit", "-[0] + [1]*x+sqrt( ([1]-1)*([1]-1)*x*x + [0]*[0] )", 0., 1000. );
  // TF1 * residualFit = new TF1("residualFit", "[0]*(x - [1])/([2]*x*x + [3]*x + [4])", 0., 1000. );
  TF1 * residualFitFunction = residualFit(ss.str());
  residuals->Fit(residualFitFunction, "ME", "", 90.56, 120);

  residuals->Draw();

  return residualFitFunction;
}

void real_emission_jackknife(TString dir="real_aajj_splitted_for_jackknife", TString nameroot="outphoton_aajj_effunf_mu1", int N=20){

  TFile *outfile = new TFile(Form("%s/%s_jackknifed.root",dir.Data(),nameroot.Data()),"recreate");
  outfile->mkdir("effunf");

  vector<TFile*> files;
  for (int i=1; i<=N; i++){
    files.push_back(new TFile(Form("%s/%s_splitted_%d.root",dir.Data(),nameroot.Data(),i),"read"));
  }

  for (std::vector<TString>::const_iterator it = diffvariables_list.begin(); it!=diffvariables_list.end(); it++){
    for (int reg=0; reg<3; reg++)
      {
	TString diffvar = *it;
	vector<TH1F*> histos;
	TString thisname = Form("htruth_%s_%d",diffvar.Data(),reg);
	TString thisnamedir = Form("effunf/%s",thisname.Data());
	for (uint i=0; i<files.size(); i++){
	  histos.push_back((TH1F*)(files.at(i)->Get(thisnamedir.Data())));
	}
	TH1F *out = JackknifammiStiHistos(histos);
	out->SetName(thisname.Data());
	out->Print();
	outfile->cd("effunf");
	out->Write();
#ifdef DEBUG
	return;
#endif
      } 
  }
  
  outfile->Close();
};

TH1F* fixrange(TH1F* old, int numB) {

  float x1, x2;
  string name = old->GetName();

  if (name.find("Ht")!=string::npos) {
    x1 = 30.;
    x2 = 500.;
    if (numB==2) x2 = 400.;
  } else if (name.find("jet_pt")!=string::npos) {
    x1 = 30.;
    x2 = 300.;
    if (numB==2) {
      if (name.find("first")!=string::npos) x2 = 200.;
      if (name.find("second")!=string::npos) x2 = 120.;
    }
  } else if (name.find("pt_Z")!=string::npos) {
    x1 = 0.;
    x2 = 300.;
    if (numB==2) x2 = 230.;
  } else {
    x1 = old->GetXaxis()->GetBinCenter(1);
    x2 = old->GetXaxis()->GetBinCenter(old->GetNbinsX());
  }

  int nx = old->GetXaxis()->FindBin(x2)-old->GetXaxis()->FindBin(x1)+1;

  x1 = old->GetXaxis()->GetBinLowEdge(old->GetXaxis()->FindBin(x1));
  x2 = old->GetXaxis()->GetBinUpEdge(old->GetXaxis()->FindBin(x2));

  TH1F* tmp = new TH1F("tmp",old->GetTitle(),nx,x1,x2);
  tmp->Sumw2();

  tmp->GetXaxis()->SetTitle(old->GetXaxis()->GetTitle());
  tmp->GetYaxis()->SetTitle(old->GetYaxis()->GetTitle());

  for (int i=0;i<=old->GetNbinsX()+1;i++) {
    int ii = tmp->GetXaxis()->FindBin(old->GetXaxis()->GetBinCenter(i));
    float c1 = tmp->GetBinContent(ii);
    float e1 = tmp->GetBinError(ii);
    float c2 = old->GetBinContent(i);
    float e2 = old->GetBinError(i);

    tmp->SetBinContent(ii,c1+c2);
    tmp->SetBinError(ii,TMath::Sqrt(e1*e1+e2*e2));
  }

  tmp->SetEntries(old->GetEntries());

  old->Delete();
  tmp->SetName(name.c_str());

  return tmp;
}

void makePlot(vector<TTree*> sigTree,vector<double> sigWeight,
	      vector<int> ptHatLo, vector<int> ptHatHi,
              std::string var,TCut cut,TH1F* h,bool norm)
{
   TH1F *hRes = (TH1F*)h->Clone();
   hRes->SetName("hRes");
   hRes->Sumw2();
  
   char tmp[300];
   for (unsigned int i=0; i<sigTree.size(); i++)
     {
       // first determine the pthat cut
       
       sprintf(tmp, "ptHat >= %d && ptHat <= %d",ptHatLo[i], ptHatHi[i]);
       TCut ptHatCut = tmp;
       TCut allCut = cut + ptHatCut;
       cout << "Current cut = " << allCut.GetTitle() << endl;
       TH1F *htmp = (TH1F*)h->Clone();
       htmp->SetName("htmp");
       sigTree[i]->Draw(Form("%s>>htmp",var.data()),allCut);
       htmp->Sumw2();
       htmp->Scale(sigWeight[i]);
       cout << "scale = " << sigWeight[i] << endl;
       cout << "After scaling htmp -> entries() " << htmp->GetEntries() << endl;
       cout << "After scaling htmp -> Integral() " << htmp->Integral() << endl;
       cout << "After scaling htmp -> GetMean()  " << htmp->GetMean() << endl;
       cout << "After scaling htmp -> GetRMS()  " << htmp->GetRMS() << endl;
       hRes->Add(htmp);
       delete htmp;
   }
   h->Sumw2();
   h->Add(hRes);
   if(norm)h->Scale(1.0/(double)h->Integral(0,1000));
   cout << "After scaling h-> entries() " << h->GetEntries() << endl;
   cout << "After scaling h-> Integral() " << h->Integral() << endl;
   cout << "After scaling h -> GetMean()  " << h->GetMean() << endl;
   cout << "After scaling h -> GetRMS()  " << h->GetRMS() << endl;
   
   delete hRes;
}

//Merge Histogram with function histogram
TH1F * merge(std::string iName,double iMergePoint,TH1F *iH,TH1F *iFunc) {
  cout << "====> Name " << iName << " -- " << iFunc << " -- " << iH << endl;
  TH1F *lH = (TH1F*) iH->Clone(iName.c_str());
  lH->SetFillStyle(0);
  int lMergeBin = iH->GetXaxis()->FindBin(iMergePoint);
  double lVal  = iH->GetBinContent(lMergeBin);
  //iFunc->Scale(lVal/iFunc->GetBinContent(lMergeBin));
  iFunc->Scale( (iH->Integral(lMergeBin, iH->GetXaxis()->FindBin(1500))) / (iFunc->Integral(lMergeBin, iFunc->GetXaxis()->FindBin(1500)) )); // felix - last fit bin = 1500; this approach seems to work much better
  for(int i0 = 0;         i0 < lMergeBin;         i0++) lH->SetBinContent(i0,iH->GetBinContent(i0));
  for(int i0 = lMergeBin; i0 < iH->GetNbinsX()+1; i0++) lH->SetBinContent(i0,iFunc->GetBinContent(iFunc->GetXaxis()->FindBin(lH->GetXaxis()->GetBinCenter(i0))));
  lH->SetName(iName.c_str());
  return lH;
}

void storePUHisto(TString inputfilepath, TString outputfilepath, TString histoname = "MCPileUpHisto"){

  TFile *inputFile = TFile::Open(inputfilepath.Data());
  //This should have always the same name
  TH1F *histo = (TH1F*)inputFile->Get("diphotonSignalMCAnalyzer/fpu_n_BeforeCuts");
  //set the name according to the name provided in cfg for crab
  histo->SetName(histoname.Data());
  histo->SetBins(80,0,80);
  //Open a new file to store this histo to
  TFile *outputFile = new TFile(outputfilepath.Data(),"recreate");
  outputFile->cd();
  histo->Write();
  outputFile->cd();
  outputFile->Close();
}

// rebin the spectra
TH1F *rebin_yaxian(TH1F *h, char *histName)
{
  TH1F *hRebin = new TH1F(Form("%s_rebin",h->GetName()),Form("rebin %s",h->GetTitle()),nbins_jetPtBin,boundaries_jetPtBin);
  for (int i=1;i<=h->GetNbinsX();i++)
    {
      double val=h->GetBinContent(i);
      double valErr=h->GetBinError(i);
      int binNum = hRebin->FindBin(h->GetBinCenter(i));
      double val1 = hRebin->GetBinContent(binNum);
      double valErr1 = hRebin->GetBinError(binNum);
      hRebin->SetBinContent(binNum,val+val1);
      hRebin->SetBinError(binNum,sqrt(valErr1*valErr1+valErr*valErr));
    }
  cleanup(hRebin);
  hRebin->SetName(histName);
  return hRebin;
}

std::vector<GaussFitInfo> saveGaussians(TH2F* hp, TList* l, TString savedir)
{
  // Each xbin of the TH2F can be viewed as a 1D histogram of the y info.
  // Fit the histograms in each xbin with a gaussian and save them into a tlist. 
  // Save the histograms as png files as well.

  // Get the title info from the mother histogram.
  TString title = hp->GetTitle();

  // This is only used to align the statboxes.
  DiMuPlottingSystem* dps = new DiMuPlottingSystem();

  // Do the fitting and gather the histos into the tlist.
  std::vector<GaussFitInfo> gfis = std::vector<GaussFitInfo>();
  gfis = dps->plotVarVsGaussFit(hp, l, 1);

  TCanvas* c = new TCanvas();
  c->SetGridx(kTRUE);
  c->SetGridy(kTRUE);

  // Iterate through the list and save the histos as pictures. 
  TIter next(l);
  TObject* obj = 0;
  int i = 0;
  int total = 0;
  while( (obj = next()) )
    {
      TH1F* h = (TH1F*) obj;
      h->SetName(Form("%s_%d_%s", title.Data(), i, hp->GetName()));
      c->cd();
      h->Draw("");
      total+=h->GetEntries();
      std::cout << "    " << i << " num: "<< h->GetEntries() << std::endl;
      dps->arrangeStatBox(c);
      c->Print(savedir+Form("%s_%s_gauss_bins_%d.png", hp->GetName(), title.Data(), i));
      //        c->Print(savedir+Form("%s_%s_gauss_bins.gif+500", hp->GetName(), title.Data()));
      i++;
    }
  c->cd();
  //    c->Print(savedir+Form("%s_%s_gauss_bins.gif++", hp->GetName(), title.Data()));
  std::cout << "nentries total: " << total << std::endl;

  delete c;
  //    delete dps;
  return gfis;
}

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");
}

void CreateFile(const char *fname)
{
   TFile *example = (TFile*)gROOT->ProcessLineFast("hsimple(1)");
   if (!example) return;
   TH1F *hpx = (TH1F*)example->Get("hpx");
   hpx->SetName("hpx1");
   TFile::Cp(example->GetName(), fname);
   TFile *file = TFile::Open(fname, "UPDATE");
   file->mkdir("folder")->cd();
   hpx->Write();
   file->Close();
   example->Close();
   TString sname(fname);
   if (sname.Contains("000")) {
      TFile::Cp(fname, "original.root");
      TFile::Open("original.root");
   }
}

//Rebin the histogram
TH1F* rebin(TH1F* iH,int iNBins,double *iAxis) { 
  std::string lTmp = "tmp"; //Added to avoid Root output errors
  TH1F *lH = new TH1F(lTmp.c_str(),lTmp.c_str(),iNBins,iAxis);
  for(int i0 = 0; i0 < iH->GetNbinsX()+1; i0++) {
    int    lNBin = lH->GetXaxis()->FindBin(iH->GetXaxis()->GetBinCenter(i0));
    double lVal  = iH->GetBinContent(i0);
    double lErr  = iH->GetBinError  (i0);
    double lOldV = lH->GetBinContent(lNBin);
    double lOldE = lH->GetBinError  (lNBin);
    lH->SetBinContent(lNBin,lVal+lOldV);
    lH->SetBinError  (lNBin,sqrt(lOldE*lOldE+lErr*lErr));
  }
  std::string lName2 = iH->GetName();
  std::string fine_binning = "_fine_binning";
  lName2.replace(lName2.find(fine_binning),fine_binning.length(),"");
  lH->SetName (lName2.c_str());
  lH->SetTitle(lName2.c_str());
  delete iH;
  return lH;
}

void billtw(Int_t compress) {
   //write N histograms to a Tree
   timer.Start();
   TFile f("billt.root","recreate","bill benchmark with trees",compress);
   TH1F *h = new TH1F("h","h",1000,-3,3);
   h->FillRandom("gaus",50000);
   TTree *T = new TTree("T","test bill");
   T->Branch("event","TH1F",&h,64000,0);
   for (Int_t i=0;i<N;i++) {
      char name[20];
      sprintf(name,"h%d",i);
      h->SetName(name);
      h->Fill(2*gRandom->Rndm());
      T->Fill();
   }
   T->Write();
   delete T;
   timer.Stop();
   printf("billtw%d : RT=%7.3f s, Cpu=%7.3f s, File size= %9d bytes, CX= %g\n",compress,timer.RealTime(),timer.CpuTime(),
                    (Int_t)f.GetBytesWritten(),f.GetCompressionFactor());
}

TH1F* Interpolate(double massLow, TH1F* low, double massHigh, TH1F* high, double massInt){

  if (low->GetNbinsX()!=high->GetNbinsX()) std::cout << "Cannot interpolate differently binned histograms" << std::endl;
  assert(low->GetNbinsX()==high->GetNbinsX());

  TH1F *interp = (TH1F*)low->Clone();
  for (int i=0; i<interp->GetNbinsX(); i++){
    double OutLow = low->GetBinContent(i+1);
    double OutHigh = high->GetBinContent(i+1);
    double OutInt = (OutHigh*(massInt-massLow)-OutLow*(massInt-massHigh))/(massHigh-massLow);
    interp->SetBinContent(i+1,OutInt);
  }

  std::string name = low->GetName();
  std::string strLow = Form("%3.1f",massLow);
  std::string strInt = Form("%3.1f",massInt);
  int ind = name.rfind(strLow);
  name.replace(ind-6,11,strInt);
  interp->SetName(name.c_str());

  return interp;
   
}

//.........这里部分代码省略.........
      for (int x=1;x<=uhist[i]->hResponse->GetNbinsX();x++) {  	
	double sum=hProj->GetBinContent(y);
	cout <<y<<" "<<x<<" "<<sum<<endl;
	if (sum==0) continue;
	double ratio = uhist[i]->hMeas->GetBinContent(y)/sum;
	if (uhist[i]->hMeas->GetBinContent(y)==0) ratio = 1e-100/sum;
        uhist[i]->hResponse->SetBinContent(x,y,uhist[i]->hResponse->GetBinContent(x,y)*ratio);
	uhist[i]->hResponse->SetBinError(x,y,uhist[i]->hResponse->GetBinError(x,y)*ratio);
      }
    }

    uhist[i]->hResponseNorm = (TH2F*)uhist[i]->hMatrix->Clone(Form("hResponseNorm_cent%d",i));
    for (int x=1;x<=uhist[i]->hResponseNorm->GetNbinsX();x++) {
      double sum=0;
      for (int y=1;y<=uhist[i]->hResponseNorm->GetNbinsY();y++) {
	if (uhist[i]->hResponseNorm->GetBinContent(x,y)<=0*uhist[i]->hResponseNorm->GetBinError(x,y)) {
	  uhist[i]->hResponseNorm->SetBinContent(x,y,0);
	  uhist[i]->hResponseNorm->SetBinError(x,y,0);
	}
	sum+=uhist[i]->hResponseNorm->GetBinContent(x,y);
      }
			
      for (int y=1;y<=uhist[i]->hResponseNorm->GetNbinsY();y++) {  	
	if (sum==0) continue;
	double ratio = 1./sum;
	uhist[i]->hResponseNorm->SetBinContent(x,y,uhist[i]->hResponseNorm->GetBinContent(x,y)*ratio);
	uhist[i]->hResponseNorm->SetBinError(x,y,uhist[i]->hResponseNorm->GetBinError(x,y)*ratio);
      }
    }
		
    uhist[i]->hResponse->Draw("colz");
		
    if (!useMatrixFromFile) uhist[i]->hMatrixFit = uhist[i]->hMatrix;
    uhist[i]->hMatrixFit->SetName(Form("hMatrixFit_cent%d",i));
  }

  pbpb_Unfo->cd();
	
  cout << "==================================== UNFOLD ===================================" << endl;
	
  //char chmet[100]; 
	
  // ======================= Reconstructed pp and PbPb spectra =========================================================
  TCanvas * cPbPb = new TCanvas("cPbPb","Comparison",1200,600);
  cPbPb->Divide(2,1); 
  cPbPb->cd(1);
	
  TH1F *hRecoBW[nbins_cent+1], *hRecoBinByBinBW[nbins_cent+1], *hMeasBW[nbins_cent+1], *hGenBW[nbins_cent+1];


  TCanvas * cPbPbMeas = new TCanvas("cPbPbMeas","Measurement",1200,600);
  cPbPbMeas->Divide(2,1); 
  cPbPbMeas->cd(1);
	

  for (int i=0;i<=nbins_cent;i++) {
    cPbPb->cd(i+1)->SetLogy();   
    // Do Bin-by-bin
    TH1F *hBinByBinCorRaw = (TH1F*)uhist[i]->hResponse->ProjectionY(); 
    TH1F *hMCGen           = (TH1F*)uhist[i]->hResponse->ProjectionX(); // gen
    hBinByBinCorRaw->Divide(hMCGen);
    TF1 *f = new TF1("f","[0]+[1]*x");
    hBinByBinCorRaw->Fit("f","LL ","",90,300);
    TH1F* hBinByBinCor = (TH1F*)hBinByBinCorRaw->Clone();//functionHist(f,hBinByBinCorRaw,Form("hBinByBinCor_cent%d",i));
    uhist[i]->hRecoBinByBin = (TH1F*) uhist[i]->hMeas->Clone(Form("hRecoBinByBin_cent%d",i));
    uhist[i]->hRecoBinByBin->Divide(hBinByBinCor);