StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Public Member Functions | Static Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | List of all members
PeakWindow Class Reference

#include <PeakWindow.h>

Inheritance diagram for PeakWindow:

Public Member Functions

 PeakWindow ()
 Default Constructor. More...
 
 PeakWindow (Double_t start, Double_t end)
 Construct with known start and end points, peak gets set to imposible values.
 
 PeakWindow (const PeakWindow &oldpeak)
 Copy Constructor.
 
PeakWindowoperator= (const PeakWindow &rhs)
 Assignment operator.
 
virtual ~PeakWindow ()
 Destructor.
 
virtual void Copy (TObject &obj) const
 Only copies variables, to copy TLines use Clone()
 
virtual TObject * Clone (const char *newname="") const
 Clone whole object, name is irrelevant.
 
void SetWindow (Double_t s, Double_t e)
 
void GetWindow (Double_t &s, Double_t &e) const
 
void SetPeak (TGraph *gdata)
 sets mPeakX based on mP_Peak using line of slopes from points (mP_Peak-1,mP_Peak) and (mP_Peak,mP_Peak+1) More...
 
virtual void Combine (const PeakWindow &other, bool keepthis=true)
 merges one PeakWindow into another More...
 
virtual UShort_t CompareTo (const PeakWindow &other) const
 compare two PeakWindow objects More...
 
Double_t StartEndLineSlope () const
 Computes the slope of the line formed by the points (mStartX,mStartY) and (mEndX,mEndY)
 
Double_t StartEndSlopeUncertainty (Double_t sigma) const
 Uncertainty int the slope of the line formed by the points (mStartX,mStartY) and (mEndX,mEndY)
 
Double_t StartEndLineYint () const
 Computes the y-intercept of the line formed by the points (mStartX,mStartY) and (mEndX,mEndY)
 
Double_t MidPoint (TGraph *graph=0) const
 Computes the the line formed by the points (mStartX,mStartY) and (mEndX,mEndY) and evaluates that line at mPeakX. More...
 
virtual Double_t SlopeChirality (Double_t scale) const
 
virtual Double_t PeakChirality (Double_t slopescale, Double_t peakscale) const
 
virtual Double_t PeakChiralityProb (Double_t probscale, Double_t chirality) const
 
virtual Double_t PeakChiralityProb (Double_t probscale, Double_t peakscale, Double_t chirscale) const
 
virtual Double_t PeakTunnelProb (TGraph *graph, Double_t scale=1.0, Double_t sigma=1.0) const
 Compute probablity that a given PeakWindow is a real peak. More...
 
virtual void Reset (Double_t start, Double_t end)
 Reset PeakWindow to constructor state.
 
virtual void Print (Option_t *opt="") const
 Prints information about PeakWindow.
 
virtual void Draw (Option_t *opt="")
 Draw the PeakWindow. More...
 
virtual void Paint (Option_t *opt="")
 paint method see Draw() for options
 
TLine * GetStartLine (Double_t ymin=0, Double_t ymax=0)
 Create and return a TLine for the start of the peak window.
 
Color_t GetStartLineColor () const
 
Style_t GetStartLineStyle () const
 
Width_t GetStartLineWidth () const
 
void SetStartLineColor (Color_t color)
 
void SetStartLineColorAlpha (Color_t color, Float_t alpha)
 
void SetStartLineStyle (Style_t style)
 
void SetStartLineWidth (Width_t width)
 
TMarker * GetPeakMarker ()
 Create and return a TMarker to mark the location of the peak.
 
Color_t GetPeakMarkerColor () const
 
Style_t GetPeakMarkerStyle () const
 
Size_t GetPeakMarkerSize () const
 
void SetPeakMarkerColor (Color_t color)
 
void SetPeakMarkerColorAlpha (Color_t color, Float_t alpha)
 
void SetPeakMarkerStyle (Style_t style)
 
void SetPeakMarkerSize (Size_t size)
 
TLine * GetEndLine (Double_t ymin=0, Double_t ymax=0)
 Create and return a TLine for the end of the peak window.
 
Color_t GetEndLineColor () const
 
Style_t GetEndLineStyle () const
 
Width_t GetEndLineWidth () const
 
void SetEndLineColor (Color_t color)
 
void SetEndLineColorAlpha (Color_t color, Float_t alpha)
 
void SetEndLineStyle (Style_t style)
 
void SetEndLineWidth (Width_t width)
 

Static Public Member Functions

static PeakWindow Combine (const PeakWindow &leftpeak, const PeakWindow &rightpeak, bool keepleft=true)
 combine two PeakWindow objects More...
 

Public Attributes

Double_t mStartX
 x value for start of the peak window
 
Double_t mEndX
 x value for end of the peak window
 
Double_t mStartY
 y value associated with mStartX
 
Double_t mEndY
 y value associated with mEndX
 
Int_t mP_Peak
 Point Number of peak in a TGraph object (P for point), point is such that slope with previous point will be positive and next point will be negative.
 
Double_t mPeakX
 x-value of peak position as determined by SetPeak()
 
Double_t mPeakY
 y-value at mP_Peak
 

Protected Member Functions

 ClassDef (PeakWindow, 3)
 

Protected Attributes

TLine * mStartLine
 TLine for drawing the start of the peak window.
 
TMarker * mPeakMarker
 TMarker for drawing the peak location.
 
TLine * mEndLine
 TLine for drawing the end of the peak window.
 

Detailed Description

Data class to hold properties of a found peak like the peak position and its start and end points. Mostly used as a helper class for PeakAna.

A "peak window" consists of 3 points: where a peak starts (first positive slope after negative slopes leading up to peak), where it plateaus (slope is zero i.e. the peak position), and where the peak ends (last negative slope before changing to positive slopes)

Definition at line 55 of file PeakWindow.h.

Constructor & Destructor Documentation

PeakWindow::PeakWindow ( )

Default Constructor.

Start and end points of peak are set to a range of 0 to 1 and mP_Peak, mPeakX, and mPeakY equal to -1 (impossible point)

Definition at line 5 of file PeakWindow.cxx.

Referenced by Clone().

Member Function Documentation

PeakWindow PeakWindow::Combine ( const PeakWindow leftpeak,
const PeakWindow rightpeak,
bool  keepleft = true 
)
static

combine two PeakWindow objects

Parameters
leftpeakone of the peaks to be combined and assumed to have the lower x-value
rightpeakone of the peaks to be combined and assumed to have the higher x-value
keepleftboolean to determine which peak position should be kept in the combined peak, true means "leftpeak", false means "rightpeak"
Returns
combined PeakWindow object is returned

Definition at line 156 of file PeakWindow.cxx.

References mEndX, mEndY, mP_Peak, mPeakX, mPeakY, mStartX, and mStartY.

Referenced by PeakAna::GetPossiblePeaks().

void PeakWindow::Combine ( const PeakWindow other,
bool  keepthis = true 
)
virtual

merges one PeakWindow into another

Changes "this" peak to include the "other" peak.

Parameters
otherPeakWindow to merge with this one
keepthistrue means keep "this" peak's position, false means keep "other" peak's position

Definition at line 137 of file PeakWindow.cxx.

References mEndX, mEndY, mP_Peak, mPeakX, mPeakY, mStartX, and mStartY.

UShort_t PeakWindow::CompareTo ( const PeakWindow other) const
virtual

compare two PeakWindow objects

Compares two PeakWindow's and returns a value based on how different they are

Returns
comparison flag
0 = different PeakWindow's
1 = same mStartX and mStartY only
2 = "1" + same mP_Peak only
3 = "1" + "2" + same mStartY and mEndY only
4 = "1" + "2" + "3" + same mPeakX only
5 = same PeakWindow's ("1"+"2"+"3"+"4"+same mPeakY)

Definition at line 271 of file PeakWindow.cxx.

References mEndX, mEndY, mP_Peak, mPeakX, mPeakY, mStartX, and mStartY.

void PeakWindow::Draw ( Option_t *  opt = "")
virtual

Draw the PeakWindow.

Parameters
optoptions for drawing:
  • "" (empty string) means draw start line, peak marker and end line.
  • "s" means draw start line
  • "p" means draw peak marker
  • "e" means draw end line

Definition at line 288 of file PeakWindow.cxx.

void PeakWindow::GetWindow ( Double_t &  s,
Double_t &  e 
) const
Parameters
sget x-value for start of peak
eget x-value for end of peak

Definition at line 76 of file PeakWindow.cxx.

References mEndX, and mStartX.

Double_t PeakWindow::MidPoint ( TGraph *  graph = 0) const

Computes the the line formed by the points (mStartX,mStartY) and (mEndX,mEndY) and evaluates that line at mPeakX.

The function is mostly needed in computing the probablity for peak tunneling since the difference of MidPoint() and mPeakY is used in the probability formula.

Parameters
graphif graph given use that graph's x-value at mP_Peak rather than mPeakX
Returns
y-value evaluated at mPeakX of line formed by points (mStartX,mStartY) and (mEndX,mEndY)

Definition at line 194 of file PeakWindow.cxx.

References mEndX, mEndY, mP_Peak, mPeakX, mStartX, and mStartY.

Referenced by PeakAna::ConvertPeaksToAna().

Double_t PeakWindow::PeakTunnelProb ( TGraph *  graph,
Double_t  scale = 1.0,
Double_t  sigma = 1.0 
) const
virtual

Compute probablity that a given PeakWindow is a real peak.

\begin{equation} Probability = 1/(scale*StartEndDiff^2+1)*Erfc(PeakHeightDiff/(sqrt(2)*sigma)) \end{equation}

or 1 if PeakHeightDiff<=0
Erfc = complimentary error function
StartEndDiff = mEndX-mStartX
PeakHeightDiff = mPeakY - MidPoint()
This function will work even if SetPeak() is not called since it requires a TGraph and will read the values from there

Parameters
graphTGraph of data points
scaleStartEndDiff scale in formula of probability
sigmasigma of Erfc to use in formula of probability

Definition at line 247 of file PeakWindow.cxx.

References mEndX, mEndY, mP_Peak, mStartX, and mStartY.

Referenced by PeakAna::GetPossiblePeaks(), StFcsPulseAna::MergeByProbability(), PeakAna::PeakProb(), and PeakAna::PeakTunnel().

void PeakWindow::SetPeak ( TGraph *  gdata)

sets mPeakX based on mP_Peak using line of slopes from points (mP_Peak-1,mP_Peak) and (mP_Peak,mP_Peak+1)

This function is used to set mPeakX to a value from the left and right slopes of mP_Peak to correct for any discretization coming from points on a graph. Requires mP_Peak has been set correctly

Parameters
gdataTGraph where data is stored

Definition at line 112 of file PeakWindow.cxx.

References mP_Peak, mPeakX, and mPeakY.

Referenced by PeakAna::GetPossiblePeaks().

void PeakWindow::SetWindow ( Double_t  s,
Double_t  e 
)
Parameters
sset x-value for start of peak
eset x-value for end of peak

Definition at line 70 of file PeakWindow.cxx.

References mEndX, and mStartX.

Referenced by StFcsPulseAna::AnalyzeForPeak(), PeakAna::AnalyzeForPeak(), PeakAna::Init(), PeakAna::SetSearchWindow(), and PeakAna::SetWindow().


The documentation for this class was generated from the following files: