Linear rebin and rescale of a histogram

ROOT
1

Hello ROOTers,

I’m trying to linear rebin and rescale the yellow histogram in order to match the green spectrum. The code that we wrote is the following:

int Double_TH1F() {

  FILE *file = fopen("./Run0_PHA_HG_0_44.txt", "r");
  FILE *file2 = fopen("./conv_s13_63.txt", "r");

  TCanvas* c1 = new TCanvas("c1","HG Multiphoton Spectrum",10,10,1920,1980);

  c1->SetGridx();
  c1->SetGridy();

  TH1F* htmp = new TH1F("htemp","htemp",2048,0,2048);
  TH1F *hsipm = new TH1F("HG Multiphoton Spectrum","",2048,0,2048);
  //TH1F *hsipm = new TH1F("HG Multiphoton Spectrum","",630,0,21);

  
  Int_t bins = 0;
  Float_t x= 0;
  Int_t i = 0;
  Double_t gain_verde = 33.3085;
  Double_t gain_giallo = 32.5152;
  Double_t ped_verde = 74.9209;
  Double_t ped_giallo = 121.278;
  Double_t spe_verde = 101.996;
  Double_t spe_giallo = 144.391;

  
  while(!feof(file)){
    fscanf(file,"%g\n",&x);
    for(i=0;i<x; ++i){

    htmp->Fill(bins);
    hsipm->Fill(bins);

  }
    ++bins;
  }
  
     
  hsipm->SetFillColorAlpha(kGreen, 0.85);
  hsipm->GetXaxis()->SetRangeUser(0, 550);
  hsipm->GetXaxis()->SetTitle("ADC Channels");
  hsipm->GetYaxis()->SetTitle("Counts");
  
  TH1F* htmp2 = new TH1F("htemp2","htemp2",2048,0,2048);
  TH1F *hsipm2 = new TH1F("HG Multiphoton Spectrum 2","",2048,0,2048);


  Int_t bins2 = 0;
  Float_t y= 0;
  
  while(!feof(file2)){
    fscanf(file2,"%g\n",&y);
    for(i=0;i<y; ++i){

    htmp2->Fill(bins2);
    hsipm2->Fill((bins2-ped_giallo+ped_verde)*(gain_giallo/gain_verde));

  }
    ++bins2;
  }
  
  hsipm2->SetFillColorAlpha(kYellow, 0.85);
  
  Float_t eventi_dt5202 = hsipm->Integral();
  cout << "Eventi DT5202 = " << eventi_dt5202 << endl;
  
  
  Float_t eventi_dt5550w = hsipm2->Integral();
  cout << "Eventi DT5550W = " << eventi_dt5550w << endl;


 // hsipm->SetBins(hsipm->GetNbinsX(),0.,20.);  
  hsipm->Draw();
//  hsipm2->SetBins(hsipm2->GetNbinsX(),0.,1999.22331); 
  hsipm2->Draw("same");
  
  auto legend = new TLegend(0.1, 0.6, 0.25, 0.9);

	legend->AddEntry(hsipm, "CAEN DT5202", "f");
	legend->AddEntry(hsipm2,"CAEN DT5550W", "f");
	legend->Draw();
	
	c1->SaveAs("Confronto_Schede_CAEN.pdf");
  
  return 0;
  
}

We obtain the following output:

c1
c11278×730 96.9 KB

As you can see, the yellow peaks are not superimposed on the green peaks, as you might expect since we scaled the second histogram for the ratio of the two gains. And we notice, also, a problem in the yellow spectra bins (see the extra event in some bins).

The reseult that we want to reproduce is the following that we simply obtained using the same logic with Python.

Data file of the code:
conv_s13_63.txt (9,5 KB)
Run0_PHA_HG_0_44.txt (7,5 KB)

2

I tried an other approach. I am not sure it is correct. It is simply a shift of bins of the 2nd histogram compare to the first one. The shift is equal to the position of the first non empty bin in the first histogram. It is assumed the bins’ spacing is the same in both. The code is:

void Double_TH1F() {

   FILE *file  = fopen("./Run0_PHA_HG_0_44.txt", "r");
   FILE *file2 = fopen("./conv_s13_63.txt", "r");

   TCanvas* c1 = new TCanvas("c1","HG Multiphoton Spectrum",10,10,1920,1980);

   c1->SetGrid();

   TH1F *hsipm = new TH1F("HG Multiphoton Spectrum","",2048,0,2048);

   Int_t bins = 0;
   Float_t x  = 0;
   Int_t i    = 0;

   Int_t fbin = 1;

   while(!feof(file)){
      fscanf(file,"%g\n",&x);
      if(x>0 && fbin==1) fbin=bins;
      for(i=0;i<x; ++i){
         hsipm->Fill(bins);
      }
      ++bins;
   }

   hsipm->SetFillColorAlpha(kGreen, 0.85);
   hsipm->GetXaxis()->SetRangeUser(0, 550);
   hsipm->GetXaxis()->SetTitle("ADC Channels");
   hsipm->GetYaxis()->SetTitle("Counts");

   TH1F *hsipm2 = new TH1F("HG Multiphoton Spectrum 2","",2048,0,2048);

   Int_t bins2 = 0;
   Float_t y= 0;

   while(!feof(file2)){
      fscanf(file2,"%g\n",&y);
      for(i=0;i<y; ++i){
         hsipm2->Fill(bins2-fbin);
      }
      ++bins2;
   }

   hsipm2->SetFillColorAlpha(kYellow, 0.85);

   Float_t eventi_dt5202 = hsipm->Integral();
   cout << "Eventi DT5202 = " << eventi_dt5202 << endl;
   Float_t eventi_dt5550w = hsipm2->Integral();
   cout << "Eventi DT5550W = " << eventi_dt5550w << endl;

   hsipm->Draw();
   hsipm2->Draw("same");

   auto legend = new TLegend(0.1, 0.8, 0.25, 0.9);

   legend->AddEntry(hsipm, "CAEN DT5202", "f");
   legend->AddEntry(hsipm2,"CAEN DT5550W", "f");
   legend->Draw();

   c1->SaveAs("Confronto_Schede_CAEN.pdf");
}

Which gives:

Screenshot 2022-03-16 at 14.02.01
Screenshot 2022-03-16 at 14.02.012866×1934 351 KB

1 Like
3

Thanks for the help and reply. Ok, YOur result is improved but the 1’s and 2nd p.e. are not matched…

We are triying another way:

int Double_TH1F() {

  FILE *file = fopen("./Run0_PHA_HG_0_44.txt", "r");
  FILE *file2 = fopen("./conv_s13_63.txt", "r");

  TCanvas* c1 = new TCanvas("c1","HG Multiphoton Spectrum",10,10,1920,1980);

  c1->SetGridx();
  c1->SetGridy();

  TH1F* htmp = new TH1F("htemp","htemp",2048,0,2048);
  TH1F *hsipm = new TH1F("HG Multiphoton Spectrum","",630,0,12);

  
  Int_t bins = 0;
  Float_t x= 0;
  Int_t i = 0;
  Double_t gain_verde = 33.3085;
  Double_t gain_giallo = 32.5152;
  Double_t ped_verde = 74.9209;
  Double_t ped_giallo = 121.278;
  Double_t spe_verde = 101.996;
  Double_t spe_giallo = 144.391;

  
  while(!feof(file)){
    fscanf(file,"%g\n",&x);
    for(i=0;i<x; ++i){

    htmp->Fill(bins);
    hsipm->Fill((bins-spe_verde)/gain_verde+1);

  }
    ++bins;
  }
  
     
  hsipm->SetFillColorAlpha(kGreen, 0.5);
  hsipm->GetXaxis()->SetRangeUser(0, 12);
  hsipm->GetXaxis()->SetTitle("# p.e.");
  hsipm->GetYaxis()->SetTitle("Counts");
  
  TH1F* htmp2 = new TH1F("htemp2","htemp2",2048,0,2048);
  TH1F *hsipm2 = new TH1F("HG Multiphoton Spectrum 2","",630,0,12);


  Int_t bins2 = 0;
  Float_t y= 0;
  
  while(!feof(file2)){
    fscanf(file2,"%g\n",&y);
    for(i=0;i<y; ++i){

    htmp2->Fill(bins2);
    hsipm2->Fill((bins2-spe_giallo)/gain_giallo+1);

  }
    ++bins2;
  }
  
  hsipm2->SetFillColorAlpha(kYellow, 0.5);
  
  Float_t eventi_dt5202 = hsipm->Integral();
  cout << "Eventi DT5202 = " << eventi_dt5202 << endl;
  
  
  Float_t eventi_dt5550w = hsipm2->Integral();
  cout << "Eventi DT5550W = " << eventi_dt5550w << endl;


 // hsipm->SetBins(hsipm->GetNbinsX(),0.,20.);  
  hsipm->Draw();
//  hsipm2->SetBins(hsipm2->GetNbinsX(),0.,1999.22331); 
  hsipm2->Draw("same");
  
  auto legend = new TLegend(0.1, 0.6, 0.25, 0.9);

	legend->AddEntry(hsipm, "CAEN DT5202", "f");
	legend->AddEntry(hsipm2,"CAEN DT5550W", "f");
	legend->Draw();
	
	c1->SaveAs("Confronto_Schede_CAEN.pdf");
  
  return 0;
  
}

And we get:

immagine
immagine1865×941 17.4 KB

The result is perfect but… not the bin. As you can see there are many empty bin that no have physical sense.
How can i fix this last bug?

4

Well, that’s your logic computing the fill value which provoque this effect. I do not know what is the algorithm behind . For example in the first fill (bins-spe_verde)/gain_verde+1 must be not right. But what to do instead ? you know better than me.Reducing the number of bins helps, but is that correct ? I do not know.

5

With (bins-spe_verde)/gain_verde+1 I want to move the spectrum so that the first photopeak is centred on 1 and rescaled equal to the gain. That way the second photo-peak will be on 2 and so on.

The plus 1 was inserted in order to solve the bin problem, but it doesn’t seem to work.

7

Lets concentrate on the first histogram.

I see that bins varies from 0 to 2048.

spe_verde = 101.996
gain_verde = 33.3085

so (bins-spe_verde)/gain_verde+1 varies from -2.0621613 to 59.423646 … Is it correct ? the histogram range is 0, 12 .

where doe these values for spe_verde and gain_verde come from ?

How is chosen the number of bins 630 ?

I do not understand the logic behind

8

If you see the second code, the bins varies from 0 to 12

The calculus is right.

The values of spe_verde and gain_verde comes from another fit.

The number of bins is choosen in order to plot a significant spectra.

9

In the code given below, you need to improve the two macros which convert “channel -> p.e.” (note: they must be “increasing linear functions”).

int Double_TH1F() {

  FILE *file = fopen("./Run0_PHA_HG_0_44.txt", "r");
  FILE *file2 = fopen("./conv_s13_63.txt", "r");

  TCanvas* c1 = new TCanvas("c1","HG Multiphoton Spectrum",10,10,1920,1980);

  c1->SetGridx();
  c1->SetGridy();

  TH1F* htmp = new TH1F("htemp","htemp",2048,0,2048);
  TH1F *hsipm = new TH1F("HG Multiphoton Spectrum","",2048,0,2048);

  
  Int_t bins = 0;
  Float_t x= 0;
  Int_t i = 0;
  Double_t gain_verde = 33.3085;
  Double_t gain_giallo = 32.5152;
  Double_t ped_verde = 74.9209;
  Double_t ped_giallo = 121.278;
  Double_t spe_verde = 101.996;
  Double_t spe_giallo = 144.391;

  
  while(!feof(file)){
    fscanf(file,"%g\n",&x);
    htmp->SetBinContent(bins + 1, x);
    hsipm->SetBinContent(bins + 1, x);
    ++bins;
  }
  htmp->Sumw2(kFALSE); htmp->ResetStats();
#define _VERDE_(_X_) (((_X_) - spe_verde) / gain_verde + 1.) /* channel -> p.e. */
  hsipm->GetXaxis()->Set(hsipm->GetNbinsX(), // nbins
			 _VERDE_(hsipm->GetXaxis()->GetXmin()), // xmin
			 _VERDE_(hsipm->GetXaxis()->GetXmax())); // xmax
#undef _VERDE_
  hsipm->Sumw2(kFALSE); hsipm->ResetStats();
  
  hsipm->SetLineColor(kGreen);
  //hsipm->SetFillColorAlpha(kGreen, 0.5);
  //hsipm->GetXaxis()->SetRangeUser(0, 12);
  hsipm->GetXaxis()->SetTitle("# p.e.");
  hsipm->GetYaxis()->SetTitle("Counts");
  
  TH1F* htmp2 = new TH1F("htemp2","htemp2",2048,0,2048);
  TH1F *hsipm2 = new TH1F("HG Multiphoton Spectrum 2","",2048,0,2048);


  Int_t bins2 = 0;
  Float_t y= 0;
  
  while(!feof(file2)){
    fscanf(file2,"%g\n",&y);
    htmp2->SetBinContent(bins2 + 1, y);
    hsipm2->SetBinContent(bins2 + 1, y);
    ++bins2;
  }
  htmp2->Sumw2(kFALSE); htmp2->ResetStats();
#define _GIALLO_(_X_) (((_X_) - spe_giallo) / gain_giallo + 1.) /* channel -> p.e. */
  hsipm2->GetXaxis()->Set(hsipm2->GetNbinsX(), // nbins
			  _GIALLO_(hsipm2->GetXaxis()->GetXmin()), // xmin
			  _GIALLO_(hsipm2->GetXaxis()->GetXmax())); // xmax
#undef _GIALLO_
  hsipm2->Sumw2(kFALSE); hsipm2->ResetStats();
  
  hsipm2->SetLineColor(kYellow);
  //hsipm2->SetFillColorAlpha(kYellow, 0.5);
  
  Float_t eventi_dt5202 = hsipm->Integral();
  cout << "Eventi DT5202 = " << eventi_dt5202 << endl;
  
  
  Float_t eventi_dt5550w = hsipm2->Integral();
  cout << "Eventi DT5550W = " << eventi_dt5550w << endl;


 // hsipm->SetBins(hsipm->GetNbinsX(),0.,20.);  
  hsipm->Draw();
//  hsipm2->SetBins(hsipm2->GetNbinsX(),0.,1999.22331); 
  hsipm2->Draw("same");
  
  auto legend = new TLegend(0.1, 0.6, 0.25, 0.9);

	legend->AddEntry(hsipm, "CAEN DT5202", "f");
	legend->AddEntry(hsipm2,"CAEN DT5550W", "f");
	legend->Draw();
	
	c1->SaveAs("Confronto_Schede_CAEN.pdf");
  
  return 0;
  
}
1 Like
10

Your answer is very helpful, thank you very much.

It is a very good starting point that will allow us to improve the macro.

11

I improved the source code in my previous post.

BTW. There was a bug in my original version. I used “Fill(bin+1)” while it should have been “Fill(bin)” (note: “SetBinContent(bin+1, value)” is correct).

1 Like
12

This topic was automatically closed 14 days after the last reply. New replies are no longer allowed.