Spectrum does not find and show right peaks

ROOT
1

Dear All, I am using Root TSpectrum2 to find peaks in 2D histograms. I am attaching here my root macro (taken from Root example file peaks2.C). I am also attaching the plot where peaks markers are not at the right place. The output of the root macro is:

Peak 0 x = 3.63385e+36 y = 3.2518e+39
Peak 1 x = 1.8365e+35 y = 3.2518e+39
Peak 2 x = 5.92494e-315 y = 5.97478e-315
Peak 3 x = 0 y = 0
Peak 4 x = 0 y = 0

Number of positions = 5
x[0] = 46666.7, y[0] = 150676
x[1] = 68541.7, y[1] = 167027
x[2] = 37916.7, y[2] = 141757
x[3] = 49583.3, y[3] = 167027
x[4] = 64166.7, y[4] = 152162
Info in TCanvas::Print: pdf file peaks.pdf has been created

The first list of peak are from s->GetPositionX() and s->GetPositionY() and the second list is printed from TSpectrum2 method s->Print(); .

Questions:

  1. Why TSpectrum2 is not finding the right peaks (peaks are not in the hot region)? I think TSpectrum2 can find bins with maximum contents and set initial values of fit parameters without I set them myself.

  2. Why when I get peaks positions with GetPosition method I do not have the same ones with Print(); method. How ever the number of peaks is good.

Thank you, Rachid
findlabmscanpeaks.C (1.9 KB) peaks.pdf (17.4 KB)

Please read tips for efficient and successful posting and posting code

_ROOT Version: 5.18.00
_Platform: Ubuntu 16.04.5 LTS
_Compiler: gcc/g++ 5.4.0

2

Your macro has input parameters. How should we run it ? it seems also some data are needed,

3

Dear Couet, First of all thank for your fast answer.Yes you are right. I upload here the root file where the histogram is housed. Forgot then about all parameters except use ipm=4.
the file is 2d-scan-24-Nov-2019-21-53-12-imot01-imot02-steps-rate-over-current.root (285.7 KB)

Also I solved the peaks position values I get from GetPosition they are correct now after I went to root version v6-16-00

But peaks still at the wrong positions

Thank you, rachid

4

What should be the value of maxpeaks ?

5

Sorry just found them by hand:

Peak x y
1 64k 128k
2 65k 122k
3 67.5k 150k
4 66k 122k
5 66.5k 126k

You can also look at the attached pdf having all histograms in the root file plotted with “contz” option. You look at the first plot in page 2 with title starting with PM 4: 2d-scan-24-Nov-2019-21-53-12-rate-over-current-steps-cont.pdf (117.3 KB)

Tspectrum finds the right number of peaks !! so I think maybe the problem is related to coordinates types? are there histogram coordinates or in other space?

Thank you, Rachid

6

So, is you problem solved ?

7

No not solved: I sent you the peaks I found by hand on the contour plot using a ruler. But the problem of TSpectrum2 is still because if you compare the peaks positions I found by hand from the plot and the ones given by TSpectrum2 you find they are not the same?? However TSpectrum found the same number of peaks that are the same I found by vision.
As I said above is this related to that TSpectrum2 maybe found the right peaks positions but not giving them in the histogram coordinate system ?

8

To be more clear I give you first the peak I found by hand on the plot:

Peak x y
1 64k 128k
2 65k 122k
3 67.5k 150k
4 66k 122k
5 66.5k 126k

and these are the peaks positions found by TSpectrum2:

Number of positions = 5
x[0] = 46666.7, y[0] = 150676
x[1] = 68541.7, y[1] = 167027
x[2] = 37916.7, y[2] = 141757
x[3] = 49583.3, y[3] = 167027
x[4] = 64166.7, y[4] = 152162

The are not the same ?

Rachid

9

to run your macro I need to know which value you give to maxpeaks

10

Yes this is just a safety number. Me I gave 50

11

I am not completely sure if that will help you but playing with the Search parameters I managed to get a better result:

void findPeak2() {
   auto s = new TSpectrum2(200);
   TH2F *h2 = 0;

   TFile* file = TFile::Open("2d-scan-24-Nov-2019-21-53-12-imot01-imot02-steps-rate-over-current.root");
   file->GetObject("PM 2;2", h2);

   int nfound = s->Search(h2,1,"colz",0.4);
   Double_t *xpeaks = s->GetPositionX();
   Double_t *ypeaks = s->GetPositionY();

   for (int pf=0;pf<nfound;pf++) {
      cout << " Peak " << pf << " x = " << xpeaks[pf] << " y = " << ypeaks[pf] << "\n";
   }
}
12

If you want to make your macro valid for ROOT 5 and ROOT 6, instead of:

   Double_t *xpeaks = s->GetPositionX();
   Double_t *ypeaks = s->GetPositionY();

use:

#if ROOT_VERSION_CODE >= ROOT_VERSION(6,00,00)
   Double_t *xpeaks, *ypeaks; // ROOT 6
#else
   Float_t *xpeaks, *ypeaks; // ROOT 5
#endif
   xpeaks = s->GetPositionX();
   ypeaks = s->GetPositionY();
13

The option noMarkov seems to find the correct peaks, but somehow removes the axis labels.

   Int_t nfound = s->Search(h2,1,"colz,noMarkov");
14

Thank you Olivier, It seems we should work with TSpectrum2 parameters. Now with the input parameters you used, the peaks are in the hot region but will play more with the parameters. Please remind me where I can find Doc about TSpectrum2 input parameter.

Thank you Wile to state about different variables type of TSpectrum peaks positions

Here attached a pdf file showing where are the peaks now with Olivier modification

Rachid

peaks.pdf (18.1 KB)

15

The option noMarkov seems to find the correct peaks, but somehow removes the axis labels.

   Int_t nfound = s->Search(h2,1,"colz,noMarkov");

Thank you Dastudillo indeed it does find the right peaks.

for the axis, the important issue are the peaks and will enter them by hand in the plot

Rachid

peaks-dastudillo.pdf (17.2 KB)

16

https://root.cern.ch/doc/master/classTSpectrum2.html

17

Hi Dastudillo, That Doc is very good. It seems we can improve to find all peaks even with small yield if you choose the right threshold. One last question: Is there a histogram ready after TSpectrum2 call that has only peaks without background as I read in the Doc you sent. In the Doc they are finding it by subtracting background but maybe TSpectrum2 has the histogram without background ready?

For the axies with noMarkov option I have drawn the histogram again with marks like:

auto m = new TMarker();
m->SetMarkerStyle(23);
m->SetMarkerColor(kRed);
gPad->SetLogz();
h2->Draw(“contz”);
c1->Update();
int i;
for (i=0;i<nfound;i++) {
m->DrawMarker(xpeaks[i],ypeaks[i]);
}

The plot is here attached finding small yield peaks using low threshold in TSpectrum2
peaks.pdf (236.9 KB)
Rachid

18

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