Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
SearchHR3.C File Reference

Detailed Description

View in nbviewer Open in SWAN Example to illustrate the influence of number of iterations in deconvolution in high resolution peak searching function (class TSpectrum).

Found 11 candidate peaks
void SearchHR3() {
Double_t fPositionX[100];
Double_t fPositionY[100];
Int_t fNPeaks = 0;
Int_t i,nfound,bin;
const Int_t nbins = 1024;
Double_t xmax = nbins;
Double_t source[nbins], dest[nbins];
gROOT->ForceStyle();
TString dir = gROOT->GetTutorialDir();
TString file = dir+"/spectrum/TSpectrum.root";
TFile *f = new TFile(file.Data());
TH1F *h = (TH1F*) f->Get("back2");
h->SetTitle("Influence of # of iterations in deconvolution in peak searching");
h->GetXaxis()->SetRange(1,nbins);
TH1F *d1 = new TH1F("d1","",nbins,xmin,xmax);
TH1F *d2 = new TH1F("d2","",nbins,xmin,xmax);
TH1F *d3 = new TH1F("d3","",nbins,xmin,xmax);
TH1F *d4 = new TH1F("d4","",nbins,xmin,xmax);
TSpectrum *s = new TSpectrum();
for (i = 0; i < nbins; i++) source[i]=h->GetBinContent(i + 1);
nfound = s->SearchHighRes(source, dest, nbins, 8, 2, kTRUE, 3, kTRUE, 3);
Double_t *xpeaks = s->GetPositionX();
for (i = 0; i < nfound; i++) {
a=xpeaks[i];
bin = 1 + Int_t(a + 0.5);
fPositionX[i] = h->GetBinCenter(bin);
fPositionY[i] = h->GetBinContent(bin);
}
TPolyMarker * pm = (TPolyMarker*)h->GetListOfFunctions()->FindObject("TPolyMarker");
if (pm) {
h->GetListOfFunctions()->Remove(pm);
delete pm;
}
pm = new TPolyMarker(nfound, fPositionX, fPositionY);
h->GetListOfFunctions()->Add(pm);
pm->SetMarkerStyle(23);
pm->SetMarkerSize(1.3);
h->Draw("L");
for (i = 0; i < nbins; i++) d1->SetBinContent(i + 1,dest[i]);
d1->Draw("SAME");
for (i = 0; i < nbins; i++) source[i]=h->GetBinContent(i + 1);
s->SearchHighRes(source, dest, nbins, 8, 2, kTRUE, 10, kTRUE, 3);
for (i = 0; i < nbins; i++) d2->SetBinContent(i + 1,dest[i]);
d2->Draw("SAME");
for (i = 0; i < nbins; i++) source[i]=h->GetBinContent(i + 1);
s->SearchHighRes(source, dest, nbins, 8, 2, kTRUE, 100, kTRUE, 3);
for (i = 0; i < nbins; i++) d3->SetBinContent(i + 1,dest[i]);
d3->Draw("SAME");
for (i = 0; i < nbins; i++) source[i]=h->GetBinContent(i + 1);
s->SearchHighRes(source, dest, nbins, 8, 2, kTRUE, 1000, kTRUE, 3);
for (i = 0; i < nbins; i++) d4->SetBinContent(i + 1,dest[i]);
d4->Draw("SAME");
printf("Found %d candidate peaks\n",nfound);
}
#define f(i)
Definition RSha256.hxx:104
#define a(i)
Definition RSha256.hxx:99
#define h(i)
Definition RSha256.hxx:106
int Int_t
Definition RtypesCore.h:45
double Double_t
Definition RtypesCore.h:59
const Bool_t kTRUE
Definition RtypesCore.h:91
@ kRed
Definition Rtypes.h:66
@ kGreen
Definition Rtypes.h:66
@ kMagenta
Definition Rtypes.h:66
@ kBlue
Definition Rtypes.h:66
float xmin
float xmax
#define gROOT
Definition TROOT.h:406
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition TAttLine.h:40
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition TAttMarker.h:38
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition TAttMarker.h:40
virtual void SetMarkerSize(Size_t msize=1)
Set the marker size.
Definition TAttMarker.h:41
A ROOT file is a suite of consecutive data records (TKey instances) with a well defined format.
Definition TFile.h:54
1-D histogram with a float per channel (see TH1 documentation)}
Definition TH1.h:575
virtual void SetTitle(const char *title)
See GetStatOverflows for more information.
Definition TH1.cxx:6678
virtual void SetBinContent(Int_t bin, Double_t content)
Set bin content see convention for numbering bins in TH1::GetBin In case the bin number is greater th...
Definition TH1.cxx:9062
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition TH1.cxx:3073
virtual TObject * FindObject(const char *name) const
Must be redefined in derived classes.
Definition TObject.cxx:323
A PolyMarker is defined by an array on N points in a 2-D space.
Definition TPolyMarker.h:31
Advanced Spectra Processing.
Definition TSpectrum.h:18
Int_t SearchHighRes(Double_t *source, Double_t *destVector, Int_t ssize, Double_t sigma, Double_t threshold, bool backgroundRemove, Int_t deconIterations, bool markov, Int_t averWindow)
One-dimensional high-resolution peak search function.
Double_t * GetPositionX() const
Definition TSpectrum.h:58
Basic string class.
Definition TString.h:136
Definition file.py:1
#define dest(otri, vertexptr)
Definition triangle.c:1040
Authors
Miroslav Morhac, Olivier Couet

Definition in file SearchHR3.C.