C++ TH2::Fill方法代码示例

//______________________________________________________________________________
TH2* CreateHisto(const Char_t* name, Int_t nPar, const Double_t* par)
{
    TH2* h;
    if (gIsTAPS)
        h = new TH2F(name, name, 4000, -100, 100, nPar, 0, nPar);
    else
        h = new TH2F(name, name, 2000, -100, 100, nPar, 0, nPar);

    // fill histo
    for (Long64_t i = 0; i < gEntries; i++)
    {
        const Double_t time_1 = gGain[gElem_1[i]] * (gRaw_Time_1[i] - par[gElem_1[i]]);
        const Double_t time_2 = gGain[gElem_2[i]] * (gRaw_Time_2[i] - par[gElem_2[i]]);
        h->Fill(time_1 - time_2, gElem_1[i]);
        h->Fill(time_2 - time_1, gElem_2[i]);
    }

    return h;
}

TH2 *
ReturnCorrFromFit(TH2 *hcorr)
{

  gStyle->SetOptStat(kFALSE);
  gStyle->SetOptFit(kTRUE);

  TObjArray *oa = new TObjArray();
  hcorr->FitSlicesX(0, 0, -1, 0, "QNR", oa);

  TCanvas *cFit = new TCanvas("cFit", "cFit");
  cFit->Divide(2, 2);

  TF1 *fMean = new TF1("fMean", "[0] + [1] * TMath::Power(x, [2])", 1., 100.);
  fMean->SetParameter(0, 0.);
  fMean->SetParameter(1, 1.);
  fMean->SetParameter(2, 1.);
  TH1 *hMean = (TH1 *)oa->At(1);
  hMean->Fit(fMean, "0q", "I", 1., 100.);
  cFit->cd(1)->SetLogx();
  cFit->cd(1)->SetLogy();
  hMean->Draw();
  fMean->Draw("same");

  TF1 *fSigma = new TF1("fSigma", "[0] + [1] * TMath::Power(x, [2])", 1., 100.);
  fSigma->SetParameter(0, 0.);
  fSigma->SetParameter(1, 1.);
  fSigma->SetParameter(2, 0.5);
  TH1 *hSigma = (TH1 *)oa->At(2);
  hSigma->Fit(fSigma, "0q", "", 1., 100.);
  cFit->cd(3)->SetLogx();
  cFit->cd(3)->SetLogy();
  hSigma->Draw();
  fSigma->Draw("same");

  cFit->cd(2)->SetLogx();
  cFit->cd(2)->SetLogy();
  cFit->cd(2)->SetLogz();
  hcorr->Draw("colz");

  TH2 *hcorrfit = (TH2 *)hcorr->Clone("hcorrfit");
  //  hcorrfit->Reset();
  for (Int_t i = 0; i < hcorr->GetNbinsX(); i++) {
    Float_t cent = hcorr->GetXaxis()->GetBinCenter(i + 1);
    Float_t mean = fMean->Eval(cent);
    Float_t sigma = fSigma->Eval(cent);
    if (cent < 25 || cent > 100) continue;
    for (Int_t j = 0; j < 10000; j++) {
      Float_t val = gRandom->Gaus(mean, sigma);
      if (val <= 0.) continue;
      hcorrfit->Fill(val, cent);
    }
  }

  cFit->cd(4)->SetLogx();
  cFit->cd(4)->SetLogy();
  cFit->cd(4)->SetLogz();
  hcorrfit->Draw("colz");

  return hcorrfit;
}

void KVSpIdGUI::SpiderIdentification()
{
   if ((!fHisto) || (!fGrid)) return;

   TVirtualPad* pad = fGrid->GetPad();
   fGrid->UnDraw();

   fZp = fZpEntry->GetIntNumber();
   if (!fUserParameter) fSpFactor = GetFactor();
   else fSpFactor = fSpiderFactorEntry->GetNumber();

   fAnglesUp   = fAngleUpEntry->GetIntNumber();
   fAnglesDown = fAngleDownEntry->GetIntNumber();
   fAlpha      = fApertureUpEntry->GetNumber();
   fPiedType   = fPiedChoice->GetSelected();
   fMatrixType = fTypeChoice->GetSelected();
   Int_t type  = fMatrixType;

   TH2*   tmpHisto = fHisto;
   TList* tmpCut   = 0;
   if (fUseCut) {
      tmpHisto = (TH2*)fHisto->Clone(Form("%s_cut", fHisto->GetName()));
      tmpHisto->Reset();
      for (int i = 1; i <= fHisto->GetNbinsX(); i++) {
         for (int j = 1; j <= fHisto->GetNbinsY(); j++) {
            Stat_t ww = fHisto->GetBinContent(i, j);
            Axis_t x0 = fHisto->GetXaxis()->GetBinCenter(i);
            Axis_t y0 = fHisto->GetYaxis()->GetBinCenter(j);
            if (fGrid->IsIdentifiable(x0, y0)) tmpHisto->Fill(x0, y0, ww);
         }
      }
      tmpCut = (TList*)fGrid->GetCuts()->Clone("tmpCuts");
   }

   fGrid->Clear();

   if (fScaledHisto) delete fScaledHisto;
   KVHistoManipulator hm;

   TF1 RtLt("RtLt", Form("x*%lf", fSfx), 0, tmpHisto->GetXaxis()->GetXmax());
   TF1 RtLty("RtLty", Form("x*%lf", fSfy), 0, tmpHisto->GetXaxis()->GetXmax());
   fScaledHisto = (TH2F*)hm.ScaleHisto(tmpHisto, &RtLt, &RtLty);

   if (fIdentificator) delete fIdentificator;
   fIdentificator = new KVSpiderIdentificator(fScaledHisto, fXm * fSfx, fYm * fSfy);

   switch (fPiedType) {
      case kUser:
         fIdentificator->SetX0(fPdx * fSfx);
         fIdentificator->SetY0(fPdy * fSfy);
         break;
      case kAuto:
         break;
      case kNone:
         fIdentificator->SetX0(0.);
         fIdentificator->SetY0(0.);
   }

   fIdentificator->SetParameters(fSpFactor);
   fIdentificator->SetNangles(fAnglesUp, fAnglesDown);
   fIdentificator->SetAlpha(fAlpha);

   fProgressBar->SetRange(0, fAnglesUp + fAnglesDown + 1);
   fProgressBar->Reset();
   fIdentificator->Connect("Increment(Float_t)", "TGHProgressBar", fProgressBar, "SetPosition(Float_t)");

   fTestButton->SetEnabled(kFALSE);
   fCloseButton->SetEnabled(kFALSE);
   fIdentificator->ProcessIdentification();
   fTestButton->SetEnabled(kTRUE);
   fCloseButton->SetEnabled(kTRUE);

   fIdentificator->Disconnect("Increment(Float_t)", fProgressBar, "SetPosition(Float_t)");
   fProgressBar->Reset();

   if (fDebug) fIdentificator->Draw(fOption.Data());

   TList* ll = (TList*)fIdentificator->GetListOfLines();

   KVIDZALine* TheLine = 0;
   int zmax = 0;

   KVSpiderLine* spline = 0;
   TIter next_line(ll);
   while ((spline = (KVSpiderLine*)next_line())) {
      if ((spline->GetN() > 10)) { //&&(spline->GetX(0)<=fIdentificator->GetX0()+200.))
         TF1* ff1 = 0;
         if (type == kSiCsI) ff1 = spline->GetFunction(fPdx * fSfx, TMath::Max(fScaledHisto->GetXaxis()->GetXmax() * 0.9, spline->GetX(spline->GetN() - 1)));
         else if (type == kSiSi)  ff1 = spline->GetFunction(fPdx * fSfx, TMath::Min(fScaledHisto->GetXaxis()->GetXmax() * 0.9, spline->GetX(spline->GetN() - 1) * 1.5));
         else if (type == kChIoSi)  ff1 = spline->GetFunction(fPdx * fSfx, TMath::Min(fScaledHisto->GetXaxis()->GetXmax() * 0.9, spline->GetX(spline->GetN() - 1) * 1.5));
         else ff1 = spline->GetFunction();
         if ((type == kSiCsI) && (ff1->GetParameter(1) >= 3000. || (ff1->GetParameter(2) <= 0.35) || (ff1->GetParameter(2) >= 1.))) {
            Info("SpiderIdentification", "Z = %d has been rejected (fit parameters)", spline->GetZ());
            continue;
         }
         TheLine = (KVIDZALine*)((KVIDZAGrid*)fGrid)->NewLine("ID");
         TheLine->SetZ(spline->GetZ());
         double min, max;
         ff1->GetRange(min, max);
         double step = TMath::Min((max - min) * 0.05, 20.); //20.;
//.........这里部分代码省略.........

//.........这里部分代码省略.........
     ///////////////////////////////////////////////

     genjet.clear();
     for(int i=0;i < branchGenJet->GetEntriesFast();++i)
     {
       gjet = (Jet*) branchGenJet->At(i);
       genjet.push_back(gjet);
     }

     //////////////////////////////////////////////
     ///colojet bilgilerini alıyor ////////////////
     //////////////////////////////////////////////

     colojet.clear();
     for(int i=0;i < branchJet->GetEntriesFast(); ++i)
     {
      cjet = (Jet*) branchJet->At(i);
      colojet.push_back(cjet);
     }

     /////////////////////////////////////////////
     /////// delta(y) for (60 - (100) - (150) - 250)///////////
     /////////////////////////////////////////////

     double ptbins=0;
     if(genjet.size()==colojet.size())
     {

       for(unsigned i = 0; i< genjet.size();++i)
       {
         ptbins=genjet[i]->PT;
         if( ptbins > 60 && ptbins < 100)
         {
           hist_deltaY_60_100->Fill(genjet[i]->Eta,genjet[i]->Eta - colojet[i]->Eta);
           histevent60100->Fill(genjet[i]->Eta - colojet[i]->Eta);
         }
         if(ptbins > 100 && ptbins < 150)
         {
           hist_deltaY_100_150->Fill(genjet[i]->Eta,genjet[i]->Eta - colojet[i]->Eta);
           histevent100150->Fill(genjet[i]->Eta - colojet[i]->Eta);
         }
         if(ptbins > 150 && ptbins < 200)
         {
           hist_deltaY_150_200->Fill(genjet[i]->Eta,genjet[i]->Eta - colojet[i]->Eta);
           histevent150200->Fill(genjet[i]->Eta - colojet[i]->Eta);
         }
       }
     }



     //////////////////////////////////////////////////////////////////////
     //R parametresi buluyor ve pt resolution yapıyor//////////////////////
     //////////////////////////////////////////////////////////////////////

     double r =0;
     double pt_res_cut=0;
     if(genjet.size()> 1 && colojet.size()>1)
     {
       for(unsigned i = 0; i<3;i++)
       {
         if(abs(genjet[i]->Eta) <= 2.5 && abs(colojet[i]->Eta) <= 2.5)
         {
           double deltaPhi = abs(genjet[i]->Phi - colojet[i]->Phi);
           double deltaEta = abs(genjet[i]->Eta - colojet[i]->Eta);
           r = sqrt( pow(deltaPhi,2) + pow(deltaEta,2));