GWOLLUM 4.2.0
Tools for gravitational-wave analyses
Loading...
Searching...
No Matches
Spectrum.h
Go to the documentation of this file.
1
6#ifndef __Spectrum__
7#define __Spectrum__
8
9#include <TMath.h>
10#include <TRandom.h>
11#include <TH1.h>
12#include <TGraph.h>
13#include <TFile.h>
14#include "FFT.h"
15
16using namespace std;
17
58
59 public:
60
76 Spectrum(const unsigned int aSpectrumSize,
77 const unsigned int aDataDuration,
78 const unsigned int aDataSamplingFrequency,
79 const unsigned int aVerbosity=0);
80
84 virtual ~Spectrum(void);
93 void Reset(void);
94
112 bool AddData(const unsigned int aDataSize,
113 const double *aData,
114 const unsigned int aDataStart=0);
115
124 bool LoadData(const unsigned int aDataSize,
125 const double *aData,
126 const unsigned int aDataStart=0);
127
137 double GetPower(const double aFrequency);
138
144 inline double GetPower(const unsigned int aFrequencyIndex){
145 return fft_psd->GetRe_f(aFrequencyIndex);
146 };
147
153 inline double GetAmplitude(const double aFrequency){
154 return TMath::Sqrt(GetPower(aFrequency));
155 };
156
161 TGraph* GetPSD(void);
162
169 TGraph* GetPSD(const double aFrequencyMin, const double aFrequencyMax);
170
175 TGraph* GetASD(void);
176
183 TGraph* GetASD(const double aFrequencyMin, const double aFrequencyMax);
184
191 TGraph* GetPeriodogram(const unsigned int aParity, const unsigned int aIndex);
192
199 TGraph* GetAmplitudePeriodogram(const unsigned int aParity, const unsigned int aIndex);
200
206 bool WritePeriodogram(const string aOutFileName);
207
211 inline unsigned int GetSpectrumSize(void){ return fft_psd->GetSize_f(); };
212
217 inline double GetSpectrumFrequency(const unsigned int aIndex){
218 return (double)aIndex*(double)sampling/(double)fft_psd->GetSize_t();
219 };
220
224 inline double GetSpectrumResolution(void){ return GetSpectrumFrequency(1); };
225
229 inline unsigned int GetSpectrumNyquist(void){
230 return (unsigned int)(GetSpectrumResolution()*(double)fft_psd->GetSize_f()); };
231
236 inline unsigned int GetNSubSegmentsMax(const unsigned int aParity){
237 return nSubSegments[aParity%2];
238 };
239
243 inline unsigned int GetDataBufferLength(void){ return (double)nSubSegments[0]/GetSpectrumResolution(); };
244
245 private:
246
247 unsigned int fVerbosity;
248 unsigned int sampling;
249
250 // derived parameters
251 unsigned int nSubSegments[2];
252
253 // containers
254 double *HannWindow;
255 double **periodogram[2];
256 unsigned int write_index[2];
257
259
260 void ComputeMedianMean(void);
261
262 ClassDef(Spectrum,0)
263};
264
265#endif
266
267
Fast Fourier transform wrapper for FFTW.
Compute the noise power spectral density.
Definition Spectrum.h:57
double GetPower(const double aFrequency)
Returns the current PSD value at a given frequency.
Definition Spectrum.cc:110
unsigned int sampling
Data sampling frequency [Hz].
Definition Spectrum.h:248
TGraph * GetPSD(void)
Returns a copy of the current PSD as a TGraph.
Definition Spectrum.cc:276
void Reset(void)
Resets the current PSD and the periodogram circular buffer.
Definition Spectrum.cc:92
TGraph * GetAmplitudePeriodogram(const unsigned int aParity, const unsigned int aIndex)
Returns a copy of a given (amplitude) periodogram.
Definition Spectrum.cc:348
unsigned int fVerbosity
Verbosity level.
Definition Spectrum.h:247
TGraph * GetASD(void)
Returns a copy of the current ASD as a TGraph.
Definition Spectrum.cc:305
unsigned int GetSpectrumNyquist(void)
Returns the Nyquist frequency of the spectrum [Hz].
Definition Spectrum.h:229
double GetSpectrumResolution(void)
Returns the frequency resolution of the spectrum [Hz]: .
Definition Spectrum.h:224
double ** periodogram[2]
Buffer for periodograms.
Definition Spectrum.h:255
void ComputeMedianMean(void)
Computes median-mean PSD.
Definition Spectrum.cc:204
double * HannWindow
Hann window vector.
Definition Spectrum.h:254
unsigned int GetDataBufferLength(void)
Returns the effective circular buffer length [s] used to compute the PSD.
Definition Spectrum.h:243
Spectrum(const unsigned int aSpectrumSize, const unsigned int aDataDuration, const unsigned int aDataSamplingFrequency, const unsigned int aVerbosity=0)
Constructor of the Spectrum class.
bool WritePeriodogram(const string aOutFileName)
Writes periodograms and the resulting PSD in a ROOT file.
Definition Spectrum.cc:357
TGraph * GetPeriodogram(const unsigned int aParity, const unsigned int aIndex)
Returns a copy of a given periodogram.
Definition Spectrum.cc:334
double GetAmplitude(const double aFrequency)
Returns the current ASD value at a given frequency.
Definition Spectrum.h:153
virtual ~Spectrum(void)
Destructor of the Spectrum class.
Definition Spectrum.cc:80
bool AddData(const unsigned int aDataSize, const double *aData, const unsigned int aDataStart=0)
Loads a data vector and updates the current PSD.
Definition Spectrum.cc:144
bool LoadData(const unsigned int aDataSize, const double *aData, const unsigned int aDataStart=0)
Loads a data vector and computes a new PSD.
Definition Spectrum.cc:195
double GetPower(const unsigned int aFrequencyIndex)
Returns the current PSD value at a given frequency index.
Definition Spectrum.h:144
unsigned int GetNSubSegmentsMax(const unsigned int aParity)
Returns the maximum number of sub-segments in the buffer.
Definition Spectrum.h:236
unsigned int write_index[2]
Write index in the buffer.
Definition Spectrum.h:256
double GetSpectrumFrequency(const unsigned int aIndex)
Returns the frequency of a given spectrum frequency index [Hz].
Definition Spectrum.h:217
unsigned int GetSpectrumSize(void)
Returns the number of points in the spectrum: .
Definition Spectrum.h:211
fft * fft_psd
FFT plan to compute the PSD.
Definition Spectrum.h:258
unsigned int nSubSegments[2]
Number of even/odd sub-subsegments.
Definition Spectrum.h:251
Wrap and optimize FFTW.
Definition FFT.h:19
unsigned int GetSize_f(void)
Returns the size of the frequency-domain vector.
Definition FFT.h:266
double GetRe_f(const unsigned int aIndex)
Returns the real part of a frequency-domain data vector element.
Definition FFT.h:141
unsigned int GetSize_t(void)
Returns the size of the time-domain vector.
Definition FFT.h:261