Logo ROOT   6.08/07
Reference Guide
TSpectrum3.h
Go to the documentation of this file.
1 // @(#)root/spectrum:$Id$
2 // Author: Miroslav Morhac 25/09/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_TSpectrum3
12 #define ROOT_TSpectrum3
13 
14 #ifndef ROOT_TNamed
15 #include "TNamed.h"
16 #endif
17 
18 class TH1;
19 
20 class TSpectrum3 : public TNamed {
21 protected:
22  Int_t fMaxPeaks; ///< Maximum number of peaks to be found
23  Int_t fNPeaks; ///< number of peaks found
24  Double_t *fPosition; ///< [fNPeaks] array of current peak positions
25  Double_t *fPositionX; ///< [fNPeaks] X positions of peaks
26  Double_t *fPositionY; ///< [fNPeaks] Y positions of peaks
27  Double_t *fPositionZ; ///< [fNPeaks] Z positions of peaks
28  Double_t fResolution; ///< *NOT USED* resolution of the neighboring peaks
29  TH1 *fHistogram; ///< resulting histogram
30 
31 public:
32  enum {
37  };
38 
39  TSpectrum3();
40  TSpectrum3(Int_t maxpositions, Double_t resolution=1); // resolution is *NOT USED*
41  virtual ~TSpectrum3();
42  virtual const char *Background(const TH1 *hist, Int_t niter, Option_t *option="goff");
43  const char *Background(Double_t ***spectrum, Int_t ssizex, Int_t ssizey, Int_t ssizez, Int_t numberIterationsX,Int_t numberIterationsY, Int_t numberIterationsZ, Int_t direction,Int_t filterType);
44  const char *Deconvolution(Double_t ***source, const Double_t ***resp, Int_t ssizex, Int_t ssizey, Int_t ssizez,Int_t numberIterations, Int_t numberRepetitions, Double_t boost);
45  TH1 *GetHistogram() const {return fHistogram;}
46  Int_t GetNPeaks() const {return fNPeaks;}
47  Double_t *GetPositionX() const {return fPositionX;}
48  Double_t *GetPositionY() const {return fPositionY;}
49  Double_t *GetPositionZ() const {return fPositionZ;}
50  virtual void Print(Option_t *option="") const;
51  virtual Int_t Search(const TH1 *hist, Double_t sigma=2, Option_t *option="goff", Double_t threshold=0.05);
52  Int_t SearchFast(const Double_t ***source, Double_t ***dest, Int_t ssizex, Int_t ssizey, Int_t ssizez, Double_t sigma, Double_t threshold, Bool_t markov, Int_t averWindow);
53  Int_t SearchHighRes(const Double_t ***source,Double_t ***dest, Int_t ssizex, Int_t ssizey, Int_t ssizez, Double_t sigma, Double_t threshold, Bool_t backgroundRemove,Int_t deconIterations, Bool_t markov, Int_t averWindow);
54  void SetResolution(Double_t resolution=1); // *NOT USED*
55  const char *SmoothMarkov(Double_t ***source, Int_t ssizex, Int_t ssizey, Int_t ssizez, Int_t averWindow);
56 
57  ClassDef(TSpectrum3,1) //Peak Finder, Background estimator, Markov smoothing and Deconvolution for 3-D histograms
58 };
59 
60 #endif
61 
62 
void SetResolution(Double_t resolution=1)
NOT USED resolution: determines resolution of the neighbouring peaks default value is 1 correspond to...
Definition: TSpectrum3.cxx:227
Double_t * fPosition
[fNPeaks] array of current peak positions
Definition: TSpectrum3.h:24
Int_t fNPeaks
number of peaks found
Definition: TSpectrum3.h:23
const char Option_t
Definition: RtypesCore.h:62
Int_t SearchFast(const Double_t ***source, Double_t ***dest, Int_t ssizex, Int_t ssizey, Int_t ssizez, Double_t sigma, Double_t threshold, Bool_t markov, Int_t averWindow)
THREE-DIMENSIONAL CLASSICAL PEAK SEARCH FUNCTION This function searches for peaks in source spectrum ...
Double_t * GetPositionY() const
Definition: TSpectrum3.h:48
const char * SmoothMarkov(Double_t ***source, Int_t ssizex, Int_t ssizey, Int_t ssizez, Int_t averWindow)
This function calculates smoothed spectrum from source spectrum based on Markov chain method...
Definition: TSpectrum3.cxx:859
const char * Deconvolution(Double_t ***source, const Double_t ***resp, Int_t ssizex, Int_t ssizey, Int_t ssizez, Int_t numberIterations, Int_t numberRepetitions, Double_t boost)
This function calculates deconvolution from source spectrum according to response spectrum The result...
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
#define ClassDef(name, id)
Definition: Rtypes.h:254
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:33
Double_t * fPositionZ
[fNPeaks] Z positions of peaks
Definition: TSpectrum3.h:27
const Double_t sigma
virtual const char * Background(const TH1 *hist, Int_t niter, Option_t *option="goff")
This function calculates background spectrum from source in h.
Definition: TSpectrum3.cxx:115
Double_t * fPositionY
[fNPeaks] Y positions of peaks
Definition: TSpectrum3.h:26
Double_t fResolution
NOT USED resolution of the neighboring peaks
Definition: TSpectrum3.h:28
TH1 * GetHistogram() const
Definition: TSpectrum3.h:45
Double_t * GetPositionZ() const
Definition: TSpectrum3.h:49
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:329
Advanced 3-dimensional spectra processing functions.
Definition: TSpectrum3.h:20
Int_t fMaxPeaks
Maximum number of peaks to be found.
Definition: TSpectrum3.h:22
Int_t SearchHighRes(const Double_t ***source, Double_t ***dest, Int_t ssizex, Int_t ssizey, Int_t ssizez, 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 Double_t
Definition: RtypesCore.h:55
The TH1 histogram class.
Definition: TH1.h:80
TH1 * fHistogram
resulting histogram
Definition: TSpectrum3.h:29
Double_t * fPositionX
[fNPeaks] X positions of peaks
Definition: TSpectrum3.h:25
#define dest(otri, vertexptr)
Definition: triangle.c:1040
virtual Int_t Search(const TH1 *hist, Double_t sigma=2, Option_t *option="goff", Double_t threshold=0.05)
This function searches for peaks in source spectrum in hin The number of found peaks and their positi...
Definition: TSpectrum3.cxx:160
Double_t * GetPositionX() const
Definition: TSpectrum3.h:47
virtual void Print(Option_t *option="") const
Print the array of positions.
Definition: TSpectrum3.cxx:126
Int_t GetNPeaks() const
Definition: TSpectrum3.h:46
TSpectrum3()
Constructor.
Definition: TSpectrum3.cxx:58
virtual ~TSpectrum3()
Destructor.
Definition: TSpectrum3.cxx:97