Omicron
4.2.0
An algorithm to detect and characterize transient events in gravitational-wave detectors
|
Ox cross-correlation engine. More...
#include <OxCorr.h>
Public Member Functions | |
double | GetTimeDelay (const unsigned int aTimeDelayBinIndex) |
Returns the time delay of a bin in the cross-correlation map [s]. More... | |
unsigned int | GetTimeDelayBinIndex (const double aTimeDelay) |
Returns the time-delay bin index for a given time delay. More... | |
unsigned int | GetTimeDelayBinsN (void) |
Returns the number of time-delay bins. More... | |
double | GetTimeDelayResolution (void) |
Returns the time-delay resolution [s]. More... | |
double | GetXi (const unsigned int aQindex, const unsigned int aFrequencyIndex, const unsigned int aTimeDelayIndex) |
Returns the cross-correlation coefficient \(\xi(Q,f,dt)\). More... | |
bool | Process (void) |
Runs the cross-correlation analysis. More... | |
bool | Transform (const unsigned int aDetectorIndex) |
Transforms the Omicron spectrogram in the Fourier domain. More... | |
Constructors and destructors | |
OxCorr (const string aOptionFile, const unsigned int aGpsRef) | |
Constructor of the OxCorr class. More... | |
virtual | ~OxCorr (void) |
Destructor of the OxCorr class. More... | |
Public Member Functions inherited from OxOmicron | |
bool | CreateOutputFile (const unsigned int aTime) |
Creates an ouput file. More... | |
unsigned int | GetAnalysisDuration (void) |
Returns the Ox analysis duration [s]. More... | |
Long64_t | GetEntries (const unsigned int aDetectorIndex) |
Returns the number of Omicron results. More... | |
int | GetEntry (const unsigned int aDetectorIndex, Long64_t aEntry) |
Loads an entry in the Omicron results. More... | |
double | GetLightTravelTime (void) |
Returns the light travel time between the two detectors. More... | |
string | GetOmicronChannelName (const unsigned int aDetectorIndex) |
Returns the Omicron channel name. More... | |
unsigned int | GetOmicronChunkDuration (void) |
Returns the Omicron chunk duration [s]. More... | |
unsigned int | GetOmicronEndTime (const unsigned int aDetectorIndex) |
Returns the end time of the current Omicron entry. More... | |
unsigned int | GetOmicronOverlapDuration (void) |
Returns the Omicron overlap duration [s]. More... | |
unsigned int | GetOmicronQ (const unsigned int aQindex) |
Returns Q value of a given Q plane. More... | |
unsigned int | GetOmicronQN (void) |
Returns the number of Omicron Q planes. More... | |
unsigned int | GetOmicronStartTime (const unsigned int aDetectorIndex) |
Returns the start time of the current Omicron entry. More... | |
bool | Process (const unsigned int aDetectorIndex, const unsigned int aTimeStart, const bool aResetPsd) |
Runs the Omicron analysis. More... | |
OxOmicron (const string aOptionFile, const unsigned int aGpsRef) | |
Constructor of the OxOmicron class. More... | |
virtual | ~OxOmicron (void) |
Destructor of the OxOmicron class. More... | |
Public Member Functions inherited from OxInit | |
void | AttachTree (TTree *aTree) |
Attach a TTree to the output file. More... | |
void | CloseOutputFile (void) |
Closes the current ouput file (if any). More... | |
bool | CreateOutputFile (const unsigned int aTime) |
Creates an ouput file. More... | |
string | GetName (void) |
Returns the object name. More... | |
string | GetOutputDirectory (void) |
Returns the path to the output directory. More... | |
bool | GetStatus (void) |
Returns the class status. More... | |
OxInit (const string aOptionFile) | |
Constructor of the OxInit class. More... | |
virtual | ~OxInit (void) |
Destructor of the OxInit class. More... | |
Private Attributes | |
fft *** | map_fft [2] |
Map FFT plans: detector 1, detector 2 (cross-correlation). More... | |
Additional Inherited Members | |
Protected Member Functions inherited from OxInit | |
void | AddOptions (void) |
Adds options from the text file. More... | |
void | DefineOption (const string aTag, const string aKey, const double aDefaultValue, const unsigned int aSize) |
Defines an option (double). More... | |
void | DefineOption (const string aTag, const string aKey, const int aDefaultValue, const unsigned int aSize) |
Defines an option (integer). More... | |
void | DefineOption (const string aTag, const string aKey, const string aDefaultValue, const unsigned int aSize) |
Defines an option (string). More... | |
void | DefineOption (const string aTag, const string aKey, const unsigned int aDefaultValue, const unsigned int aSize) |
Defines an option (unsigned integer). More... | |
void | OverloadOption (const string aTag, const string aKey) |
Overloads an option using the text file as a reference. More... | |
Protected Attributes inherited from OxOmicron | |
Omicron * | omicron [2] |
Omicron objects (detector 1 and 2) - do not delete. More... | |
TTree * | oxo_tree [2] |
TTree: Omicron results. Do not delete as it is owned by the output file. More... | |
TChain * | oxo_tree_read [2] |
TChain: Omicron results (read-only). More... | |
Protected Attributes inherited from OxInit | |
GwollumOptions * | ox_opt |
List of options. More... | |
bool | status |
Class status. More... | |
Ox cross-correlation engine.
This class is designed to compute the cross-correlation of Omicron spectrograms. The one-dimensional cross-correlation is performed in the Fourier domain along the time direction.
The Omicron analysis is managed by the OxOmicron class. Use it to compute the the spectrograms for both detectors. Then, call Transform() to Fourier-transform the spectrograms. Finally, the cross-correlation is computed with Process().
OxCorr::OxCorr | ( | const string | aOptionFile, |
const unsigned int | aGpsRef | ||
) |
Constructor of the OxCorr class.
FFT plans are created. They cover the duration of the Omicron spectrograms excluding half the overlap at both sides.
An option file must be provided to set the cross-correlation parameters:
OXC TIMERES
Time resolution [s] of the cross-correlation. The time resolution cannot be lower than the resolution of the \(Q=0\) Omicron time-frequency plane. Set this parameter to a negative value to use the maximal time resolution.OXC SPECTROPOWER
The amplitude spectrograms from Omicron are applied a power law: \((A^{2}/norm)^{\alpha}\), where \(\alpha\) is the power-law index and \(norm\) is a normalization factor. Provide \(\alpha\) and \(norm\) with this option. See Transform().[in] | aOptionFile | Path to the option file. |
[in] | aGpsRef | Reference GPS time to intiate Omicron objects. |
|
virtual |
Destructor of the OxCorr class.
|
inline |
Returns the time delay of a bin in the cross-correlation map [s].
The requested time delay is given by the bin low edge.
When the two spectrograms are mis-aligned in time, an extra time delay is added. It is derived from the Omicron spectrogram times.
[in] | aTimeDelayBinIndex | Time-delay bin index (starts at 0). |
|
inline |
Returns the time-delay bin index for a given time delay.
When the two spectrograms are mis-aligned in time, the extra time delay is removed. It is derived from the Omicron spectrogram times.
[in] | aTimeDelay | Time delay [s]. |
|
inline |
Returns the number of time-delay bins.
|
inline |
Returns the time-delay resolution [s].
This is the size of one time-delay bin.
|
inline |
Returns the cross-correlation coefficient \(\xi(Q,f,dt)\).
[in] | aQindex | Q plane index (must be valid!). |
[in] | aFrequencyIndex | Frequency row index (must be valid!). |
[in] | aTimeDelayIndex | Time-delay index (must be valid!). |
bool OxCorr::Process | ( | void | ) |
Runs the cross-correlation analysis.
After the Omicron spectrograms have been transformed in the Fourier domain (Transform()), this function runs the cross-correlation analysis. The cross-correlation is performed in the Fourier domain. After Fourier-transforming back the data, the cross-correlation map is obtained.
bool OxCorr::Transform | ( | const unsigned int | aDetectorIndex | ) |
Transforms the Omicron spectrogram in the Fourier domain.
The Omicron time-frequency maps are resampled at the OxCorr time resolution and the Omicron overlap is removed at both sides of the map. The resulting time-frequency maps are Fourier-transform along the time direction. This is done one frequency row at a time. Each frequency row is normalized:
[in] | aDetectorIndex | Detector index: 0 = first detector, 1 = second detector. |
|
private |
Map FFT plans: detector 1, detector 2 (cross-correlation).