Logo ROOT   6.16/01
Reference Guide
TSpectrum2.h
Go to the documentation of this file.
1// @(#)root/spectrum:$Id$
2// Author: Miroslav Morhac 17/01/2006
3
4/*************************************************************************
5 * Copyright (C) 1995-2006, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11#ifndef ROOT_TSpectrum2
12#define ROOT_TSpectrum2
13
14#include "TNamed.h"
15
16class TH1;
17
18class TSpectrum2 : public TNamed {
19protected:
20 Int_t fMaxPeaks; ///< Maximum number of peaks to be found
21 Int_t fNPeaks; ///< number of peaks found
22 Double_t *fPosition; ///< [fNPeaks] array of current peak positions
23 Double_t *fPositionX; ///< [fNPeaks] X position of peaks
24 Double_t *fPositionY; ///< [fNPeaks] Y position of peaks
25 Double_t fResolution; ///< *NOT USED* resolution of the neighboring peaks
26 TH1 *fHistogram; ///< resulting histogram
27static Int_t fgAverageWindow; ///< Average window of searched peaks
28static Int_t fgIterations; ///< Maximum number of decon iterations (default=3)
29
30public:
31 enum {
36 };
37
38 TSpectrum2();
39 TSpectrum2(Int_t maxpositions, Double_t resolution=1); // resolution is *NOT USED*
40 virtual ~TSpectrum2();
41 virtual TH1 *Background(const TH1 *hist, Int_t niter=20, Option_t *option="");
42 TH1 *GetHistogram() const {return fHistogram;}
43 Int_t GetNPeaks() const {return fNPeaks;}
44 Double_t *GetPositionX() const {return fPositionX;}
45 Double_t *GetPositionY() const {return fPositionY;}
46 virtual void Print(Option_t *option="") const;
47 virtual Int_t Search(const TH1 *hist, Double_t sigma=2, Option_t *option="", Double_t threshold=0.05);
48 static void SetAverageWindow(Int_t w=3); //set average window
49 static void SetDeconIterations(Int_t n=3); //set max number of decon iterations
50 void SetResolution(Double_t resolution=1); // *NOT USED*
51
52 //new functions January 2006
53 const char *Background(Double_t **spectrum,Int_t ssizex, Int_t ssizey,Int_t numberIterationsX,Int_t numberIterationsY,Int_t direction,Int_t filterType);
54 const char *SmoothMarkov(Double_t **source, Int_t ssizex, Int_t ssizey, Int_t averWindow);
55 const char *Deconvolution(Double_t **source, Double_t **resp, Int_t ssizex, Int_t ssizey,Int_t numberIterations, Int_t numberRepetitions, Double_t boost);
56 Int_t SearchHighRes(Double_t **source,Double_t **dest, Int_t ssizex, Int_t ssizey, Double_t sigma, Double_t threshold, Bool_t backgroundRemove,Int_t deconIterations, Bool_t markov, Int_t averWindow);
57
58 static Int_t StaticSearch(const TH1 *hist, Double_t sigma=2, Option_t *option="goff", Double_t threshold=0.05);
59 static TH1 *StaticBackground(const TH1 *hist,Int_t niter=20, Option_t *option="");
60
61 ClassDef(TSpectrum2,1) //Peak Finder, background estimator, Deconvolution for 2-D histograms
62};
63
64#endif
65
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
double Double_t
Definition: RtypesCore.h:55
const char Option_t
Definition: RtypesCore.h:62
#define ClassDef(name, id)
Definition: Rtypes.h:324
The TH1 histogram class.
Definition: TH1.h:56
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
Advanced 2-dimensional spectra processing.
Definition: TSpectrum2.h:18
Double_t * GetPositionY() const
Definition: TSpectrum2.h:45
Double_t fResolution
NOT USED resolution of the neighboring peaks
Definition: TSpectrum2.h:25
TH1 * fHistogram
resulting histogram
Definition: TSpectrum2.h:26
static Int_t fgIterations
Maximum number of decon iterations (default=3)
Definition: TSpectrum2.h:28
static void SetDeconIterations(Int_t n=3)
static function: Set max number of decon iterations in deconvolution operation see TSpectrum2::Search...
Definition: TSpectrum2.cxx:113
static TH1 * StaticBackground(const TH1 *hist, Int_t niter=20, Option_t *option="")
static function (called by TH1), interface to TSpectrum2::Background
const char * SmoothMarkov(Double_t **source, Int_t ssizex, Int_t ssizey, Int_t averWindow)
This function calculates smoothed spectrum from source spectrum based on Markov chain method.
Definition: TSpectrum2.cxx:736
void SetResolution(Double_t resolution=1)
NOT USED resolution: determines resolution of the neighboring peaks default value is 1 correspond to ...
Definition: TSpectrum2.cxx:287
virtual ~TSpectrum2()
Destructor.
Definition: TSpectrum2.cxx:92
static Int_t fgAverageWindow
Average window of searched peaks.
Definition: TSpectrum2.h:27
Int_t fMaxPeaks
Maximum number of peaks to be found.
Definition: TSpectrum2.h:20
Int_t GetNPeaks() const
Definition: TSpectrum2.h:43
TSpectrum2()
Constructor.
Definition: TSpectrum2.cxx:57
TH1 * GetHistogram() const
Definition: TSpectrum2.h:42
Int_t SearchHighRes(Double_t **source, Double_t **dest, Int_t ssizex, Int_t ssizey, Double_t sigma, Double_t threshold, Bool_t backgroundRemove, Int_t deconIterations, Bool_t markov, Int_t averWindow)
This function searches for peaks in source spectrum It is based on deconvolution method.
Double_t * GetPositionX() const
Definition: TSpectrum2.h:44
@ kBackSuccessiveFiltering
Definition: TSpectrum2.h:34
@ kBackOneStepFiltering
Definition: TSpectrum2.h:35
@ kBackDecreasingWindow
Definition: TSpectrum2.h:33
@ kBackIncreasingWindow
Definition: TSpectrum2.h:32
Int_t fNPeaks
number of peaks found
Definition: TSpectrum2.h:21
virtual TH1 * Background(const TH1 *hist, Int_t niter=20, Option_t *option="")
This function calculates the background spectrum in the input histogram h.
Definition: TSpectrum2.cxx:155
virtual void Print(Option_t *option="") const
Print the array of positions.
Definition: TSpectrum2.cxx:166
static void SetAverageWindow(Int_t w=3)
static function: Set average window of searched peaks see TSpectrum2::SearchHighRes
Definition: TSpectrum2.cxx:104
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...
Definition: TSpectrum2.cxx:206
Double_t * fPositionY
[fNPeaks] Y position of peaks
Definition: TSpectrum2.h:24
static Int_t StaticSearch(const TH1 *hist, Double_t sigma=2, Option_t *option="goff", Double_t threshold=0.05)
static function (called by TH1), interface to TSpectrum2::Search
Double_t * fPosition
[fNPeaks] array of current peak positions
Definition: TSpectrum2.h:22
const char * Deconvolution(Double_t **source, Double_t **resp, Int_t ssizex, Int_t ssizey, Int_t numberIterations, Int_t numberRepetitions, Double_t boost)
This function calculates deconvolution from source spectrum according to response spectrum The result...
Double_t * fPositionX
[fNPeaks] X position of peaks
Definition: TSpectrum2.h:23
const Double_t sigma
const Int_t n
Definition: legend1.C:16
LVector boost(const LVector &v, const BoostVector &b)
Boost a generic Lorentz Vector class using a generic 3D Vector class describing the boost The only re...
Definition: VectorUtil.h:327
#define dest(otri, vertexptr)
Definition: triangle.c:1040