Logo ROOT  
Reference Guide
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
Loading...
Searching...
No Matches
peaks2.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_spectrum
3/// Example to illustrate the 2-d peak finder (class TSpectrum2).
4///
5/// This script generates a random number of 2-d gaussian peaks
6/// The position of the peaks is found via TSpectrum2
7/// To execute this example, do:
8///
9/// ~~~{.cpp}
10/// root > .x peaks2.C (generate up to 50 peaks by default)
11/// root > .x peaks2.C(10) (generate up to 10 peaks)
12/// root > .x peaks2.C+(200) (generate up to 200 peaks via ACLIC)
13/// ~~~
14///
15/// The script will iterate generating a new histogram having
16/// between 5 and the maximun number of peaks specified.
17/// Double Click on the bottom right corner of the pad to go to a new spectrum
18/// To Quit, select the "quit" item in the canvas "File" menu
19///
20/// \macro_image
21/// \macro_code
22///
23/// \author Rene Brun
24
25#include "TSpectrum2.h"
26#include "TCanvas.h"
27#include "TRandom.h"
28#include "TH2.h"
29#include "TF2.h"
30#include "TMath.h"
31#include "TROOT.h"
32
33TSpectrum2 *s;
35Int_t npeaks = 30;
37{
38 Double_t result = 0.1;
39 for (Int_t p = 0; p < npeaks; p++) {
40 Double_t norm = par[5 * p + 0];
41 Double_t mean1 = par[5 * p + 1];
42 Double_t sigma1 = par[5 * p + 2];
43 Double_t mean2 = par[5 * p + 3];
44 Double_t sigma2 = par[5 * p + 4];
46 }
47 return result;
48}
49void findPeak2()
50{
51 printf("Generating histogram with %d peaks\n", npeaks);
52 Int_t nbinsx = 200;
53 Int_t nbinsy = 200;
54 Double_t xmin = 0;
56 Double_t ymin = 0;
58 Double_t dx = (xmax - xmin) / nbinsx;
59 Double_t dy = (ymax - ymin) / nbinsy;
63 // generate n peaks at random
64 Double_t par[3000];
65 Int_t p;
66 for (p = 0; p < npeaks; p++) {
67 par[5 * p + 0] = gRandom->Uniform(0.2, 1);
68 par[5 * p + 1] = gRandom->Uniform(xmin, xmax);
69 par[5 * p + 2] = gRandom->Uniform(dx, 5 * dx);
70 par[5 * p + 3] = gRandom->Uniform(ymin, ymax);
71 par[5 * p + 4] = gRandom->Uniform(dy, 5 * dy);
72 }
73 TF2 *f2 = new TF2("f2", fpeaks2, xmin, xmax, ymin, ymax, 5 * npeaks);
74 f2->SetNpx(100);
75 f2->SetNpy(100);
76 f2->SetParameters(par);
77 TCanvas *c1 = (TCanvas *)gROOT->GetListOfCanvases()->FindObject("c1");
78 if (!c1)
79 c1 = new TCanvas("c1", "c1", 10, 10, 1000, 700);
81
82 // now the real stuff: Finding the peaks
83 Int_t nfound = s->Search(h2, 2, "col");
84
85 // searching good and ghost peaks (approximation)
86 Int_t pf, ngood = 0;
89 for (p = 0; p < npeaks; p++) {
90 for (pf = 0; pf < nfound; pf++) {
91 Double_t diffx = TMath::Abs(xpeaks[pf] - par[5 * p + 1]);
92 Double_t diffy = TMath::Abs(ypeaks[pf] - par[5 * p + 3]);
93 if (diffx < 2 * dx && diffy < 2 * dy)
94 ngood++;
95 }
96 }
97 if (ngood > nfound)
98 ngood = nfound;
99 // Search ghost peaks (approximation)
100 Int_t nghost = 0;
101 for (pf = 0; pf < nfound; pf++) {
102 Int_t nf = 0;
103 for (p = 0; p < npeaks; p++) {
104 Double_t diffx = TMath::Abs(xpeaks[pf] - par[5 * p + 1]);
105 Double_t diffy = TMath::Abs(ypeaks[pf] - par[5 * p + 3]);
106 if (diffx < 2 * dx && diffy < 2 * dy)
107 nf++;
108 }
109 if (nf == 0)
110 nghost++;
111 }
112 c1->Update();
113
114 s->Print();
115 printf("Gener=%d, Found=%d, Good=%d, Ghost=%d\n", npeaks, nfound, ngood, nghost);
116 if (!gROOT->IsBatch()) {
117 printf("\nDouble click in the bottom right corner of the pad to continue\n");
118 c1->WaitPrimitive();
119 }
120}
121void peaks2(Int_t maxpeaks = 50)
122{
123 s = new TSpectrum2(2 * maxpeaks);
124 for (int i = 0; i < 10; ++i) {
126 findPeak2();
127 }
128}
int Int_t
Definition RtypesCore.h:45
double Double_t
Definition RtypesCore.h:59
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
float xmin
float ymin
float xmax
float ymax
#define gROOT
Definition TROOT.h:406
R__EXTERN TRandom * gRandom
Definition TRandom.h:62
The Canvas class.
Definition TCanvas.h:23
virtual void SetNpx(Int_t npx=100)
Set the number of points used to draw the function.
Definition TF1.cxx:3434
virtual void SetParameters(const Double_t *params)
Definition TF1.h:677
A 2-Dim function with parameters.
Definition TF2.h:29
virtual void SetNpy(Int_t npy=100)
Set the number of points used to draw the function.
Definition TF2.cxx:928
2-D histogram with a float per channel (see TH1 documentation)
Definition TH2.h:307
virtual Double_t Uniform(Double_t x1=1)
Returns a uniform deviate on the interval (0, x1).
Definition TRandom.cxx:682
Advanced 2-dimensional spectra processing.
Definition TSpectrum2.h:18
Double_t * GetPositionY() const
Definition TSpectrum2.h:45
Double_t * GetPositionX() const
Definition TSpectrum2.h:44
void Print(Option_t *option="") const override
Print the array of positions.
virtual Int_t Search(const TH1 *hist, Double_t sigma=2, Option_t *option="", Double_t threshold=0.05)
This function searches for peaks in source spectrum in hin The number of found peaks and their positi...
return c1
Definition legend1.C:41
Double_t x[n]
Definition legend1.C:17
Double_t Gaus(Double_t x, Double_t mean=0, Double_t sigma=1, Bool_t norm=kFALSE)
Calculates a gaussian function with mean and sigma.
Definition TMath.cxx:471
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition TMathBase.h:123