C++ TH2D::GetEntries方法代码示例

void DrawCurrHistogram() {


  TObjString* filename = histograms->At(currPos);
  TString f(fFileDir);
  f.Append("/");
  f.Append(filename->GetString());

  TFile* file = new TFile(f, "read");
  // check file
  if (!file) return;
  if (file->IsZombie()) return;

  TH2D* hist = (TH2D*) gFile->Get(fHistName);

  if (!hist) continue;
  if (!hist->GetEntries()) continue;

  // copy the histogram to the root memory
  //gROOT->cd();
  //TH2D* hist_copy = new TH2D(*hist);
  TString t(filename->GetString());
  t.Append(": ");
  t.Append(hist->GetTitle());
  hist->SetTitle(t);
  //histograms->Add(hist_copy);
  //TH1* h = histograms->At(currPos);
  c->cd(1);
  //gPad->cd();

  hist->Draw("colz");
  // file can be closed now...
  //file->Close();

  if(currHist != 0) {
    PreserveAxis(currHist->GetXaxis(), hist->GetXaxis());
    PreserveAxis(currHist->GetYaxis(), hist->GetYaxis());
  }

  if(currFile != 0)
    currFile->Close();
  currFile = file;
  currHist = hist;
  //gPad->Modified();
  //gPad->Update();

  //cout << "Added histogram from " << filename << endl;

}

void HISTndist(){
	TSystemDirectory dir(".",".");
	TList *files = dir.GetListOfFiles();
	TFile *outfile = new TFile("ndist.root","RECREATE");
	TH2D *ndist = new TH2D("ndist", "Neutron Distribution", 301, 0.001, 0.601, 1201, -0.499, 1.901);
	ifstream infile;
	Int_t xi, yi, hits;
	Double_t x, y, prob;
	for (Int_t i = 0; ; i++){
		TObject *f;
		if (f = files->At(i)){
			TString filename = f->GetName();
			if (filename.EndsWith("ndist.out") && TString(filename(0,8)).IsDigit()){
				cout << "Adding " << filename << endl;
				infile.open(filename.Data());
				infile.ignore(9999, '\n');
				while (infile){
					infile >> xi >> yi >> x >> y >> prob >> hits;
					if (infile && (hits != 0)) ndist->Fill(x,y,prob);
				}
				infile.close();
				cout << "Entries: " << ndist->GetEntries() << endl;
			}
		}

void addn2() {

  // open existing f1:

  char* fn[99];
  int ni;

  // SR90 maps, mod D003

  ni = -1;
  ni++; fn[ni] = "SR90-map-00a.root";
  ni++; fn[ni] = "SR90-map-00b.root";
  ni++; fn[ni] = "SR90-map-01a.root";
  ni++; fn[ni] = "SR90-map-01b.root";
  ni++; fn[ni] = "SR90-map-02.root";
  ni++; fn[ni] = "SR90-map-03.root";
  ni++; fn[ni] = "SR90-map-05a.root";
  ni++; fn[ni] = "SR90-map-05b.root";
  ni++; fn[ni] = "SR90-map-07a.root";
  ni++; fn[ni] = "SR90-map-07b.root";
  ni++; fn[ni] = "SR90-map-07c.root";
  ni++; fn[ni] = "SR90-map-07d.root";

  // direct X-rays 9.4.2014 module D0003

  ni = -1;
  ni++; fn[ni] = "X-ray-20-modtd40000.root";
  ni++; fn[ni] = "X-ray-20-vthr60-modtd1000.root";
  ni++; fn[ni] = "X-ray-27-vthr60-modtd1000.root";
  ni++; fn[ni] = "X-ray-35-modtd40000.root";
  ni++; fn[ni] = "X-ray-35-vthr60-close-modtd65000.root";
  ni++; fn[ni] = "X-ray-35-vthr60-midpos-0p1mA-modtd10000.root";
  ni++; fn[ni] = "X-ray-35-vthr60-midpos-halfcur-modtd40000.root";
  ni++; fn[ni] = "X-ray-35-vthr60-modtd40000.root";

  int nmax = ni;

  TFile f1(fn[0]);

  if( f1.IsZombie() ) {
    cout << "Error opening " << fn[0] << endl;
    return;
  }
  cout << "opened " << fn[0] << endl;

  //--------------------------------------------------------------------
  // create f0:

  TFile f0("fileA.root", "RECREATE");
  cout << "created ";
  gDirectory->pwd();

  /*
  TFile options:  
  NEW or CREATE   create a new file and open it for writing,
                  if the file already exists the file is
                  not opened.
  RECREATE        create a new file, if the file already
                  exists it will be overwritten.
  UPDATE          open an existing file for writing.
                  if no file exists, it is created.
  READ            open an existing file for reading (default).
  NET             used by derived remote file access
                  classes, not a user callable option
  WEB             used by derived remote http access
                  class, not a user callable option
  "" (default), READ is assumed.
  */

  //--------------------------------------------------------------------
  // copy f1 to f0:

  f1.cd();

  cout << "keys:\n";
  f1.GetListOfKeys()->Print();

  cout << "pwd: ";
  f1.pwd();

  cout << "ls: \n";
  f1.ls();

  // f1 has sub-dir:

  cout << "First: " << f1.GetListOfKeys()->First()->GetName() << endl;
  cout << "First: " << f1.GetListOfKeys()->First()->ClassName() << endl;
  char* dir1 = f1.GetListOfKeys()->First()->GetName();
  cout << "cd to " << dir1 << endl;
  f1.cd( dir1 );
  cout << "we are in ";
  gDirectory->pwd();

  gDirectory->ReadAll(); // load histos

  TList * lst = gDirectory->GetList();
  cout << lst->GetName() << endl;
  cout << lst->GetTitle() << endl;
  cout << "size    " << lst->GetSize() << endl;
  cout << "entries " << lst->GetEntries() << endl;
//.........这里部分代码省略.........

// infilename - root file with relevant histograms
// system - PP,APAP,PP
// status - Pass,Fail
// rWrite - 0-no,1-png,2-eps
// rPerformance - 0-no,1-yes (ALICE logo etc.)
// bin: 0 - all, 1- 0:5, 2- 5:10, etc
void drawPID(const char* infilename, const char* system, const char* status, Int_t rWrite, Int_t rPerformance, Int_t bin)
{


    TFile *f = new TFile(infilename, "read");

    // TPC dEdx
    TH2D* TPCdEdx =(TH2D*)f->Get(Form("TPCdEdxcut%s1%stpcM%i",status, system,0));

    if (!bin) {
        int minMultBin = 0;
        int maxMultBin = 6; // 8
    }
    else {
        int minMultBin = bin-1;
        int maxMultBin = bin; // 8
    }

    double EvMultall = 0;

    for (int i = minMultBin; i < maxMultBin; i++) {

        TH2D* TPCdEdxN =(TH2D*)f->Get(Form("TPCdEdxcut%s1%stpcM%i",status,system,i));
        TPCdEdx->Add(TPCdEdxN);

        cout << i << " " << TPCdEdxN->GetEntries() << endl;

        //delete hEvMult;
    }

    TCanvas *c2 = new TCanvas("TPC dEdx", "TPC dEdx");
    c2->SetGridx();
    c2->SetGridy();
    c2->SetFillColor(10);
    c2->SetRightMargin(1.9);
    c2->SetLogz();

    TPCdEdx->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
    TPCdEdx->GetXaxis()->SetRangeUser(0.0,6.0);
    TPCdEdx->GetYaxis()->SetTitle("dE/dx");
    TPCdEdx->GetZaxis()->SetLabelSize(0.03);
    TPCdEdx->Draw("colz");

// 	double a1 = -3000.0; double b1 =  1280.0;
// 	double a2 = -312.5;  double b2 =  312.5;
// 	double a3 = -200.0;  double b3 =  240.0;

//    TF1 *fa1 = new TF1("fa1","-1800*x+940",0.3,0.4);
//   fa1->Draw("same");
//   TF1 *fa2 = new TF1("fa2","-500.0*x+420.0",0.4,0.6);
//   fa2->Draw("same");
// TF1 *fa3 = new TF1("fa3","-216.7*x+250.0",0.6,0.9);
//    fa3->Draw("same");
// TF1 *fa4 = new TF1("fa4","-566.7*x+570.0",0.6,0.75);
//    fa4->Draw("same");
// TF1 *fa5 = new TF1("fa5","-2076.92*x+1476.15",0.47,0.6);
//    fa5->Draw("same");
//   cout<<TPCdEdx->GetNbinsX()<<endl;
//   cout<<TPCdEdx->GetNbinsY()<<endl;
//   for (int ii=0;ii<TPCdEdx->GetNbinsX();ii++){

//     for (int jj=0;jj<TPCdEdx->GetNbinsY();jj++){

//       cout<<"binX: "<<ii<<endl;
//       cout<<"binY: "<<jj<<endl;
//       cout<<"val: "<<TPCdEdx->GetBinContent(ii,jj)<<endl;
//     }
//   }

// 	TH1D *py = TPCdEdx->ProjectionY("py", 230, 232); // where firstYbin = 0 and lastYbin = 9
// 	TCanvas *c22 = new TCanvas("TPC2", "TPC2");
// 	py->Draw();

    postprocess(c2,Form("TPCdEdx%s",status),rWrite,rPerformance,system);

    // TPC Nsigma
    TH2D* TPCNsigma =(TH2D*)f->Get(Form("TPCNSigmacut%s1%stpcM%i",status,system,0));

    // int minMultBin = 0;
    // int maxMultBin = 2;
    double EvMultall = 0;

    for(int i = minMultBin; i<maxMultBin; i++) {
        //all
        TH2D* TPCNsigmaN =(TH2D*)f->Get(Form("TPCNSigmacut%s1%stpcM%i",status,system,i));
        TPCNsigma->Add(TPCNsigmaN);
        //delete hEvMult;
    }

    TCanvas *c3 = new TCanvas("TPC Nsigma", "TPC Nsigma");
    c3->SetGridx();
    c3->SetGridy();
    c3->SetFillColor(10);
    c3->SetRightMargin(1.7);
//.........这里部分代码省略.........

/**
 * Loads a histogram from a ROOT-file, scales it, and either returns it as \c
 * returnedHistogram or adds it to \c returnedHistogram . The first if
 * \c returnedHistogram is \c NULL otherwise the latter.
 *
 * @param histogramName
 *                              name of the histogram.
 * @param inputFilename
 *                              name of the ROOT-file.
 * @param scale
 *                              scale factor for histogram.
 * @param returnedHistogram
 *                              the returned histogram (input, pass-by-reference does not work).
 * @param debug
 *                              switch on debug output, defaults to \c false
 *
 * @return
 *                              the returned histogram (output, pass-by-reference does not work).
 */
void
LoadHistogramTH2D(const TString& histogramName, const TString& inputFilename, double scale, TH2D*& returnedHistogram, bool debug = false)
{
// {{{
    TFile inputFile(inputFilename);
    if (!inputFile.IsOpen())
    {
        cerr << "Could not open '" << inputFilename << "' for reading." << endl;
    }
    else
    {
        TH2D* histogram = dynamic_cast<TH2D*>( inputFile.Get(histogramName) );
        if (!histogram)
        {
            cerr << "No histogram named '" << histogramName << "' in file '" << inputFilename << "'" << endl;
        }
        else
        {
            if (debug) cerr << inputFilename << " " << histogramName << " entries=" << histogram->GetEntries() << " integral=" << histogram->Integral() << " scale*integral=" << scale*histogram->Integral() << endl;
            histogram->Scale(scale);
            if (!returnedHistogram)
            {
                returnedHistogram = new TH2D(*histogram);
                returnedHistogram->SetDirectory(0); // otherwise "TFile inputFile" owns this returnedHistogram and deletes it when "TFile inputFile" goes out of scope
            }
            else
            {
                returnedHistogram->Add(histogram);
            }
        }
        inputFile.Close();
    }
// }}}
};

void combine_mistag_rates_2d( )
{

  const char * names [] = {
    "mistag_parameterization_qcd_230.root", 
    "mistag_parameterization_qcd_300.root", 
    "mistag_parameterization_qcd_380.root", 
    "mistag_parameterization_qcd_470.root", 
    "mistag_parameterization_qcd_600.root", 
    "mistag_parameterization_qcd_800.root", 
    "mistag_parameterization_qcd_1000.root", 
    "mistag_parameterization_qcd_1400.root", 
    "mistag_parameterization_qcd_1800.root", 
    "mistag_parameterization_qcd_2200.root", 
    "mistag_parameterization_qcd_2600.root", 
    "mistag_parameterization_qcd_3000.root", 
    "mistag_parameterization_qcd_3500.root"
  };


  double xs[] = {
    10623.2,
    2634.94,
    722.099,
    240.983,
    62.4923,
    9.42062,
    2.34357,
    0.1568550,
    0.013811,
    0.00129608,
    0.00011404,
    0.0000084318,
    0.00000018146
  };


  int nevents[] = {
    54000,
    54000,
    51840,
    27648,
    28620,
    20880,
    24640,
    27744,
    22848,
    22560,
    22800,
    20880,
    34320
  };


  const char * filetitles[] = {
    "QCD Dijets, #hat{pt} = 230-300",
    "QCD Dijets, #hat{pt} = 300-380",
    "QCD Dijets, #hat{pt} = 380-470",
    "QCD Dijets, #hat{pt} = 470-600",
    "QCD Dijets, #hat{pt} = 600-800",
    "QCD Dijets, #hat{pt} = 800-1000",
    "QCD Dijets, #hat{pt} = 1000-1400",
    "QCD Dijets, #hat{pt} = 1400-1800",
    "QCD Dijets, #hat{pt} = 1800-2200",
    "QCD Dijets, #hat{pt} = 2200-2600",
    "QCD Dijets, #hat{pt} = 2600-3000",
    "QCD Dijets, #hat{pt} = 3000-3500",
    "QCD Dijets, #hat{pt} = 3500-Inf"
  };


  static const int N = sizeof( names ) / sizeof ( const char * );

//   palette(N);

  TFile * output = new TFile("mistag_parameterization2d.root", "RECREATE");
  TH2D * numerator_sum = 0;
  TH2D * denominator_sum = 0;

  for ( int i = 0; i < N; ++i ) {
    TFile * f = new TFile(names[i]);
    TH2D * numerator = (TH2D*) f->Get("numerator2d");
    TH2D * denominator = (TH2D*) f->Get("denominator2d");

    char buff[1000];
    sprintf(buff, "%30s & %6.0f & %6.0f ", names[i], denominator->GetEntries(), numerator->GetEntries() );
    cout << buff << endl;

    if ( i == 0 ) {
      output->cd();
      numerator_sum = new TH2D(*numerator);
      denominator_sum = new TH2D(*denominator);
    } else {
      numerator_sum->Add( numerator );
      denominator_sum->Add( denominator );
    }
  }



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

//.........这里部分代码省略.........
                zu[countu] = (*itClus)->z();
                xu[countu] = (*itClus)->position();
                countu++;
            }
        }

        if( (*itClus)->view()== Minerva::IDCluster::V ){
            debug() <<  " StartPoint V view, pe " << (*itClus)->pe() << "; z = " << (*itClus)->z()
                << "; coord " << (*itClus)->position() << endmsg;
            Dx = (*itClus)->position() -  m_mathTool->calcVfromXY(vert.x(),vert.y());
            distance = sqrt( pow(Dx,2) + pow(Dz,2) );
            if ( is_vertex ){
                hV->Fill( (*itClus)->z()-12,(*itClus)->tpos1(), (*itClus)->pe()/distance );
                hV->Fill( (*itClus)->z()+12,(*itClus)->tpos2(), (*itClus)->pe()/distance );
            }
            hV->Fill( (*itClus)->z(),(*itClus)->position(), (*itClus)->pe()/distance );
            if ( countv < 2 ){
                zv[countu] = (*itClus)->z();
                xv[countu] = (*itClus)->position();
                countv++;
            }
        }
    }

    TF1 *fU, *fV;

    double slopeu = -9999;
    double slopev = -9999;
    double bu = -9999;
    double bv = -9999;

    bool goodFit_U = false;
    bool goodFit_V = false;
    if ( hU->GetEntries() > 3 ){
        hU->Fit("pol1","Q0");
        fU = hU->GetFunction("pol1");
        bu = fU->GetParameter(0);
        slopeu = fU->GetParameter(1);
        goodFit_U = true;

        delete fU;
    }
    else if ( hU->GetEntries() == 2 ){ // to deal with 2 clusters on NCPi0
        if ( zu[0] > zu[1] ){
            slopeu = (xu[0] - xu[1]) / (zu[0] - zu[1]);
            bu = xu[1] - zu[1]*slopeu;
            goodFit_U = true;
        }
        else if (zu[0] < zu[1] ) {
            slopeu = (xu[1] - xu[0]) / (zu[1] - zu[0]);
            bu = xu[0] - zu[0]*slopeu;
            goodFit_U = true;
        }
    }

    if ( hV->GetEntries() > 3 ){
        hV->Fit("pol1","Q0");
        fV = hV->GetFunction("pol1");
        bv = fV->GetParameter(0);
        slopev = fV->GetParameter(1);
        goodFit_V = true;

        delete fV;
    }
    else if ( hV->GetEntries() == 2 ){ // to deal with 2 clusters on NCPi0
        if ( zv[0] > zv[1] ){

//.........这里部分代码省略.........
      for(unsigned int kParticle=0; kParticle<particlesPerEvent; kParticle++) {
	particlePhiLabPreviousEvent[kParticle] = particlePhiLab[kParticle];
	particleEtaPreviousEvent[kParticle] = particleEta[kParticle];
      }

    } // loop over events
    cout << "\n Finished NTUPLE scan" << endl;

    //
    // Normalize the delta-phi histograms with the mixed-event yields
    //
    float sumSameEventPhi1D = 0.0;
    float sumMixedEventPhi1D = 0.0;
    for(unsigned int iBinX=1; iBinX<nBinsOneD; iBinX++) {
      sumSameEventPhi1D += pairDeltaPhiHistogram->GetBinContent(iBinX);
      sumMixedEventPhi1D += pairDeltaPhiHistogramMixed->GetBinContent(iBinX);
    }
    if(sumSameEventPhi1D > 0.0 && sumMixedEventPhi1D > 0.0) {
      float normalizationFactor1D = sumMixedEventPhi1D/sumSameEventPhi1D;
      for(unsigned int iBinX=0; iBinX<nBinsOneD; iBinX++) {
	pairDeltaPhiHistogramNormalized->SetBinContent(iBinX,0.0);
	float sameEventContent = pairDeltaPhiHistogram->GetBinContent(iBinX);
	float mixedEventContent = pairDeltaPhiHistogramMixed->GetBinContent(iBinX);
	if(mixedEventContent > 0.0) {
	  float normalizedYield = normalizationFactor1D*sameEventContent/mixedEventContent;
	  pairDeltaPhiHistogramNormalized->SetBinContent(iBinX, normalizedYield);
	} // check for non-zero denominator
      } // loop over phi bins
    } // check for non-zero totals
		
    //
    // Normalize the delta-eta histograms with the mixed-event yields
    //
    float sumSameEvent1D = pairDeltaEtaHistogram->GetEntries();
    float sumMixedEvent1D = pairDeltaEtaHistogramMixed->GetEntries();
    if(sumSameEvent1D > 0.0 && sumMixedEvent1D > 0.0) {
      float normalizationFactor1D = sumMixedEvent1D/sumSameEvent1D;
      for(unsigned int iBinX=0; iBinX<nBinsOneD; iBinX++) {
	pairDeltaEtaHistogramNormalized->SetBinContent(iBinX,0.0);
	float sameEventContent = pairDeltaEtaHistogram->GetBinContent(iBinX);
	float mixedEventContent = pairDeltaEtaHistogramMixed->GetBinContent(iBinX);
	if(mixedEventContent > 0.0) {
	  float normalizedYield = normalizationFactor1D*sameEventContent/mixedEventContent;
	  pairDeltaEtaHistogramNormalized->SetBinContent(iBinX, normalizedYield);
	} // check for non-zero denominator
      } // loop over eta bins
    } // check for non-zero totals
		
    float sumSameEvent2D = pairDeltaEtaDeltaPhiHistogram->GetEntries();
    float sumMixedEvent2D = pairDeltaEtaDeltaPhiHistogramMixed->GetEntries();

    if(sumSameEvent2D > 0.0 && sumMixedEvent2D > 0.0) {
      float normalizationFactor2D = sumMixedEvent2D/sumSameEvent2D;
      for(unsigned int iBinX=0; iBinX<nBinsX; iBinX++) {
	for(unsigned int iBinY=0; iBinY<nBinsY; iBinY++) {
	  pairDeltaEtaDeltaPhiHistogramNormalized->SetBinContent(iBinX, iBinY, 0.0);
	  if(iBinX<iBinXMinPlot || iBinX>iBinXMaxPlot ||
	     iBinY<iBinYMinPlot || iBinY>iBinYMaxPlot)
	    continue;  // limit the fiducial area of the plot because of possible low statistics

	  float sameEventContent = pairDeltaEtaDeltaPhiHistogram->GetBinContent(iBinX, iBinY);
	  float mixedEventContent = pairDeltaEtaDeltaPhiHistogramMixed->GetBinContent(iBinX, iBinY);
	  if(mixedEventContent > 0.0) {
	    float normalizedYield = normalizationFactor2D*sameEventContent/mixedEventContent;
	    pairDeltaEtaDeltaPhiHistogramNormalized->SetBinContent(iBinX, iBinY, normalizedYield);
	  } // check for non-zero mixed event content

// infilename - root file with relevant histograms
// system - PP,APAP,PP
// status - Pass,Fail
// rWrite - 0-no,1-png,2-eps
// rPerformance - 0-no,1-yes (ALICE logo etc.)
// bin: 0 - all, 1- 0:5, 2- 5:10, etc
void drawDCA(const char* infilename, const char* system, const char* status, Int_t rWrite, Int_t rPerformance, int isMC, Int_t bin, Int_t ptrange)
{

    myOptions(0);

    gROOT->ForceStyle();
    gStyle->SetPalette(1.0);

    TDatime now;
    int iDate = now.GetDate();
    int iYear=iDate/10000;
    int iMonth=(iDate%10000)/100;
    int iDay=iDate%100;
    char* cMonth[12]={"Jan","Feb","Mar","Apr","May","Jun",
                      "Jul","Aug","Sep","Oct","Nov","Dec"};
    char cStamp1[25],cStamp2[25];
    sprintf(cStamp1,"%i %s %i",iDay, cMonth[iMonth-1], iYear);
    sprintf(cStamp2,"%i/%.2d/%i",iDay, iMonth, iYear);

    TFile *f = new TFile(infilename, "read");

    // DCA xy
    TH2D* DCAxy =(TH2D*)f->Get(Form("DCARPtcut%s1%stpcM%d","Pass", system,0));


    if (!bin) {
        int minMultBin = 0;
        int maxMultBin = 6; // 8
    }
    else {
        int minMultBin = bin-1;
        int maxMultBin = bin; // 8
    }

    // int minMultBin = 0;
    // int maxMultBin = 0; // 8

    double EvMultall = 0;

    for(int i = minMultBin; i<maxMultBin; i++) {

        TH2D* DCAxyN = (TH2D*)f->Get(Form("DCARPtcut%s1%stpcM%d",status, system,i));
        DCAxy->Add(DCAxyN);

        cout<<i<<" "<<DCAxyN->GetEntries()<<endl;

        //delete hEvMult;
    }


    if (!isMC) {
        TCanvas *c2 = new TCanvas("DCA xy prim", "DCA xy prim");
        c2->SetGridx();
        c2->SetGridy();
        c2->SetFillColor(10);
        c2->SetRightMargin(1.9);
        c2->SetLogz();

        DCAxy->GetXaxis()->SetTitle("DCA_{XY} (cm)");
        DCAxy->GetXaxis()->SetRangeUser(-5.0,5.0);
        DCAxy->GetYaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
        // DCAxy->GetZaxis()->SetLabelSize(0.05);
        DCAxy->Draw("colz");
        postprocess(c2,Form("DCAxy%s",status),rWrite,rPerformance,system);

        // TCanvas *c4 = new TCanvas("DCA xy Projection X", "DCA xy Projection X");
        // c4->SetGridx();
        // c4->SetGridy();
        // c4->SetFillColor(10);
        // c4->SetRightMargin(1.9);
        // c4->SetLogy();

        gStyle->SetOptTitle(1);
        TCanvas *myCan = new TCanvas("myCan",cStamp1,600,400);
        myCan->Draw();
        myCan->cd();

        TPad *myPad = new TPad("myPad", "The pad",0,0,1,1);
        myPadSetUp(myPad,0.15,0.04,0.04,0.15);
        myPad->Draw();
        myPad->SetLogy();

        myPad->cd();

        if (ptrange == 0) {
            TH1D* pripp = (TH1D*)DCAxy->ProjectionX("zxc1",1,100);
            pripp->SetTitle("0.5 < #it{p}_{T} < 3 GeV/#it{c}");
        }
        else if (ptrange == 1) {
            TH1D* pripp = (TH1D*)DCAxy->ProjectionX("zxc1",15,33);
            pripp->SetTitle("0.5 < #it{p}_{T} < 1 GeV/#it{c}");
        }
        else if (ptrange == 2) {
            TH1D* pripp = (TH1D*)DCAxy->ProjectionX("zxc1",33,100);
//.........这里部分代码省略.........

void plotXY(){
  gROOT->Reset();
  gROOT->SetStyle("Plain");
  gStyle->SetHistMinimumZero(kFALSE);

  Int_t status = gSystem->Load("../rootils_C.so");
  std::cout << "status: " << status << std::endl;

  // TString dir("allStat_6hits/");
  // TString dir("allStat_noWeiBug/");
  TString dir("./");

  //Agguanta gli istogrammi
  TFile * SimF = new TFile(dir+"Sim_XY.root");
  TH2D * SimH = (TH2D*) SimF->Get("Sim_XY");;
  //
  TFile * MCF = new TFile(dir+"MC_XY.root");
  TH2D * MCH = (TH2D*) MCF->Get("MC_XY");;
  //
  //  TFile * MCFsF = new TFile(dir+"rMCFs_XY.root");
  //  TH2D * MCFsH = (TH2D*) MCFsF->Get("rMCFs_XY");;
  //
  TFile * FakeF = new TFile(dir+"MCFake_XY.root");
  TH2D * FakeH = (TH2D*) FakeF->Get("MCFake_XY");;
  //
  TFile * DataF = new TFile(dir+"Data_XY.root");
  TH2D * DataH = (TH2D*) DataF->Get("Data_XY");;
  //
  //  TFile * DataFsF = new TFile(dir+"rDataFs_XY.root");
  //  TH2D * DataFsH = (TH2D*) DataFsF->Get("rDataFs_XY");;
  //
  
  cout << MCH->GetEntries() << " " << MCH->GetEffectiveEntries() << " " << MCH->Integral() << endl;
  cout << DataH->GetEntries() << " " << DataH->GetEffectiveEntries() << " " << DataH->Integral() << endl;

  //Normalize to plot entries
  SimH->Scale(1./SimH->Integral()); 
  MCH->Scale(1./MCH->Integral());
  FakeH->Scale(1./FakeH->Integral());
  DataH->Scale(1./DataH->Integral());

  cout << MCH->GetEntries() << " " << MCH->GetEffectiveEntries() << " " << MCH->Integral() << endl;
  cout << DataH->GetEntries() << " " << DataH->GetEffectiveEntries() << " " << DataH->Integral() << endl;

  //Max
  Double_t simMax=SimH->GetMaximum();
  Double_t simPxlMax=GetMaxWithinRadius(SimH, 15.);
  Double_t simPxlInt=GetIntWithinRZ(SimH, 15., 0.);
  Double_t MCMax=MCH->GetMaximum();
  Double_t MCPxlMax=GetMaxWithinRadius(MCH, 15.);
  Double_t MCPxlInt=GetIntWithinRZ(MCH, 15., 0.);
  Double_t DataMax=DataH->GetMaximum();
  Double_t DataPxlMax=GetMaxWithinRadius(DataH, 15.);
  Double_t DataPxlInt=GetIntWithinRZ(DataH, 15., 0.);
  Double_t mAx=max(MCMax,DataMax);


  std::cout << " max     >>>> mc " << MCMax << " data " << DataMax << " Sim " << simMax << std::endl;
  std::cout << " pxl max >>>> mc " << MCPxlMax << " data " << DataPxlMax << " Sim " << simPxlMax << std::endl;
  std::cout << " pxl int >>>> mc " << MCPxlInt << " data " << DataPxlInt << " Sim " << simPxlInt << std::endl;

  // This would be need to normalize SimH to the pixel region material (taking into account that SimH and DataH/MCH have different binning!)
  //
  //  SimH->Scale(DataPxlInt*(DataH->GetXaxis()->GetBinWidth(1)*DataH->GetYaxis()->GetBinWidth(1))/simPxlInt/(SimH->GetXaxis()->GetBinWidth(1)*SimH->GetYaxis()->GetBinWidth(1)));
  // I prefer to have Sim plots and Data/MC plot to look similar by appropriately scaling the SimH range, not by scaling the histo

  mAx=0.001;
  //  mAx=0.011;
  //  mAx=0.15*0.5*(MCPxlMax+DataPxlMax);

  Double_t simRangeScale = 1./(DataPxlInt*(DataH->GetXaxis()->GetBinWidth(1)*DataH->GetYaxis()->GetBinWidth(1))/simPxlInt/(SimH->GetXaxis()->GetBinWidth(1)*SimH->GetYaxis()->GetBinWidth(1)));

  std::cout << simRangeScale << endl;
  std::cout << (DataH->GetXaxis()->GetBinWidth(1)*DataH->GetYaxis()->GetBinWidth(1)) << endl;
  std::cout << (SimH->GetXaxis()->GetBinWidth(1)*SimH->GetYaxis()->GetBinWidth(1)) << endl;

  SimH->SetMaximum(mAx*simRangeScale); 
  MCH->SetMaximum(mAx);  
  DataH->SetMaximum(mAx);

  //

  makeColorTable();

  Plot2D(SimH, "#gamma conv., MC Truth", dir);
  Plot2D(MCH, "#gamma conv., MC Reco #sqrt{s}=8TeV", dir);
  Plot2D(DataH, "#gamma conv., Data #sqrt{s}=8TeV", dir);

  //

  /*
  TH2D *MCCoarse = MCH->Rebin2D(2,2,"MCCoarse");
  TH2D *DataCoarse = DataH->Rebin2D(2,2,"DataCoarse");
  */

  /*
  cout << MCCoarse->GetEntries() << " " << MCCoarse->GetEffectiveEntries() << endl;
  cout << DataCoarse->GetEntries() << " " << DataCoarse->GetEffectiveEntries() << endl;
  */

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