the answer may be in the following line in the TSpectrum::Search
source code
for (i=0;i<size;i++) source[i] = hin->GetBinContent(i+1);
which seems to cause a shift by one channel down ( array source is input
of Search1 function ) which is not reverted later.
Btw you can better align your TH1F *hpx vs TGraph *gr
( in your example they are shifted !!! ) by doing something like this
...
ofstream oufile("./peak.1.out") ;
...
TH1F *hpx = new TH1F("hpx","Peaks",npts,0.5,(Float_t) npts + 0.5);
for (Int_t j=0 ; j<npts ; j++) hpx->SetBinContent(j,y[j]) ;
for (Int_t j=0 ; j<npts; j++)
oufile<<j<<" "<<x[j]<<" "<<hpx->GetBinCenter(j)<<" "<<y[j]<<endl;
...
what aligns centers of bins to initial x-values ( see peak.1.out ).
Then your example gives peak positions at 7 and 20 ( actually they are
at 8 and 21 !!! ).
regards
martin
On Tue, 29 May 2001, Antonio Kanaan wrote:
>
>
> Hi All,
>
> I am trying to use TSpectrum to produce a list of peaks in spectra. I
> produced a data file containing two peaks. This is my code's output:
>
> root [0]
> Processing search.cxx...
> Peaks found 2
> Peak at = 6.5
> Peak at = 19.5
> (int)0
> root [1]
>
>
> Almost fine weren't for the peak positions. They are actually at 7.5
> and 20.5. I must be missing something but don't know what.
>
> I copy here the code that does the peak finding using TSpectrum. It
> reads the original x,y table (which is attached as a comment to the
> end of the code), converts it to a histogram (I must confess I am
> still struggling to understand the reason why everything in root uses
> histograms instead of arrays of values, but that's for another SOS
> message) then uses TSpectrum::Search to find the peaks. Histogram,
> x-y table and peaks are shown on the canvas. The red triangles show
> the found peaks which are off the peak by 1 unit.
>
> Any suggestions will be appreciated,
>
> Antonio
>
> ***************************************************************
> Antônio Kanaan
> Departamento de Física - Universidade Federal de Santa Catarina
> CP 476 -- CEP 88040-900 -- Florianópolis -- SC -- Brasil
> http://www.astro.ufsc.br/~kanaan e-mail: kanaan@astro.ufsc.br
> Phone: 55,48,3319069 FAX : 55,48,3319946
> ***************************************************************
>
>
> #include "TCanvas.h"
> #include "TGraph.h"
>
> #include <iostream>
> #include <fstream>
> using namespace std;
>
> const Int_t MAXNPT = 140; /* number of points */
> Float_t x[MAXNPT],y[MAXNPT],z ; /* data arrays */
>
> search()
> {
> Int_t n,npts;
>
> TCanvas *c1 = new TCanvas("c1","the fit canvas",800,600);
>
> // fills up the x and y reading file "
> ifstream infile("./peak.1") ;
> n = 0 ;
> while ( (infile >> x[n] >> y[n++]) )
> npts = n;
>
> TH1F *hpx = new TH1F("hpx","Peaks",npts,0,npts);
> for (Int_t j=0 ; j<npts ; j++) hpx->SetBinContent(j,y[j]) ;
> hpx->Draw() ;
>
> TSpectrum *sp1 = new TSpectrum();
> sp1->Search(hpx,1.0,"") ;
> Int_t npeaks = sp1->GetNPeaks() ;
> cout << "Peaks found " << npeaks << endl ;
>
> Float_t *peaks ;
> peaks = sp1->GetPositionX() ;
> for (j=0 ; j<npeaks ; j++)
> cout << "Peak at = " << peaks[j] << endl ;
>
> // creates a TGraph object to display the computed curves
> TGraph *gr = new TGraph (npts, x, y);
> gr->Draw("C*") ;
>
> return(0) ;
>
> }
>
> /*
>
> File peak.1 to be read by the above code:
>
> 0 10
> 1 10
> 2 10
> 3 10
> 4 10
> 5 20
> 6 30
> 7 50
> 8 100
> 9 50
> 10 30
> 11 20
> 12 10
> 13 10
> 14 10
> 15 10
> 16 10
> 17 10
> 18 20
> 19 30
> 20 50
> 21 100
> 22 50
> 23 30
> 24 20
> 25 10
> 26 10
> 27 10
> 28 10
> 29 10
> 30 10
>
>
> */
>
This archive was generated by hypermail 2b29 : Tue Jan 01 2002 - 17:50:47 MET