Loading [MathJax]/extensions/TeX/AMSsymbols.js
LALPulsar 7.1.1.1-ea7c608
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages

Generate N samples of B-statistic (and F-statistic) values drawn from their respective distributions, assuming Gaussian noise, for given signal parameters. More...

Prototypes

int main (int argc, char *argv[])
 MAIN function Generates samples of B-stat and F-stat according to their pdfs for given signal-params. More...
 
int initUserVars (UserInput_t *uvar)
 Register all our "user-variables" that can be specified from cmd-line and/or config-file. More...
 
int InitCode (ConfigVariables *cfg, const UserInput_t *uvar)
 Initialized Fstat-code: handle user-input and set everything up. More...
 
int XLALsynthesizeSignals (gsl_matrix **A_Mu_i, gsl_matrix **s_mu_i, gsl_matrix **Amp_i, gsl_vector **rho2_i, const gsl_matrix *M_mu_nu, AmpParamsRange_t AmpRange, gsl_rng *rng, UINT4 numDraws)
 Generate random signal draws with uniform priors in given bands [h0, cosi, psi, phi0], and return list of 'numDraws' {s_mu} vectors. More...
 
int XLALsynthesizeNoise (gsl_matrix **n_mu_i, const gsl_matrix *M_mu_nu, gsl_rng *rng, UINT4 numDraws)
 Generate random-noise draws and combine with (FIXME: single!) signal. More...
 
int XLALcomputeLogLikelihood (gsl_vector **lnL_i, const gsl_matrix *A_Mu_i, const gsl_matrix *s_mu_i, const gsl_matrix *x_mu_i)
 Compute log-likelihood function for given input data. More...
 
int XLALcomputeFstatistic (gsl_vector **Fstat_i, gsl_matrix **A_Mu_MLE_i, const gsl_matrix *M_mu_nu, const gsl_matrix *x_mu_i)
 Compute F-statistic for given input data. More...
 
int XLALcomputeBstatisticMC (gsl_vector **Bstat_i, const gsl_matrix *M_mu_nu, const gsl_matrix *x_mu_i, gsl_rng *rng, UINT4 numMCpoints)
 Compute the B-statistic for given input data, using Monte-Carlo integration for the marginalization over {cosi, psi}, while {h0, phi0} have been marginalized analytically. More...
 
int XLALcomputeBstatisticGauss (gsl_vector **Bstat_i, const gsl_matrix *M_mu_nu, const gsl_matrix *x_mu_i)
 Compute the B-statistic for given input data, using standard Gauss-Kronod integration (gsl_integration_qng) for the marginalization over {cosi, psi}, while {h0, phi0} have been marginalized analytically. More...
 
gsl_vector * XLALcomputeBhatStatistic (const gsl_matrix *M_mu_nu, const gsl_matrix *x_mu_i)
 Compute an approximation to the full B-statistic, without using any integrations! Returns Bhat vector of dimensions (numDraws x 1), or NULL on error. More...
 
double BstatIntegrandOuter (double cosi, void *p)
 log likelihood ratio lnL marginalized over {h0, phi0} (analytical) and integrated over psi in [-pi/4,pi/4], for given cosi: BstatIntegrandOuter(cosi) = int lnL dh0 dphi0 dpsi More...
 
double BstatIntegrandInner (double psi, void *p)
 log likelihood ratio lnL marginalized over {h0, phi0} (analytical) for given psi and pars->cosi BstatIntegrandInner(cosi,psi) = int lnL dh0 dphi0 More...
 
double BstatIntegrand (double Amp[], size_t dim, void *p)
 compute log likelihood ratio lnL for given Amp = {h0, cosi, psi, phi0} and M_{mu,nu}. More...
 
REAL8 XLALComputeBhatCorrection (const gsl_vector *A_Mu_in, const gsl_matrix *M_mu_nu)
 Compute 'Bstat' approximate correction terms with respect to F, ie deltaB in Bstat = F + deltaB. More...
 

Detailed Description

Generate N samples of B-statistic (and F-statistic) values drawn from their respective distributions, assuming Gaussian noise, for given signal parameters.

Author
R. Prix

This is mostly meant to be used for Monte-Carlos studies of ROC curves

Definition in file synthesizeBstatMC.c.

Go to the source code of this file.

Data Structures

struct  AmpParamsRange_t
 Signal (amplitude) parameter ranges. More...
 
struct  ConfigVariables
 Configuration settings required for and defining a coherent pulsar search. More...
 
struct  UserInput_t
 User-variables: can be set from config-file or command-line.
 
struct  integrationParams_t
 

Macros

#define SQ(x)   ((x)*(x))
 

Variables

int vrbflg
 defined in lal/lib/std/LALError.c More...
 

Macro Definition Documentation

◆ SQ

#define SQ (   x)    ((x)*(x))

Definition at line 63 of file synthesizeBstatMC.c.

Function Documentation

◆ main()

int main ( int  argc,
char argv[] 
)

MAIN function Generates samples of B-stat and F-stat according to their pdfs for given signal-params.

< various derived configuration settings

< numDraws signal amplitude-params {h0Nat, cosi, psi, phi0}

< list of 'numDraws' signal amplitude vectors {A^mu}

< list of 'numDraws' (covariant) signal amplitude vectors {s_mu = (s|h_mu) = M_mu_nu A^nu}

< list of 'numDraws' (covariant) noise vectors {n_mu = (n|h_mu)}

< list of 'numDraws' (covariant) data vectors {x_mu = n_mu + s_mu}

< list of 'numDraws' (contravariant) data-vectors x^mu = A^mu_MLE = M^{mu nu} x_nu

< list of 'numDraws' optimal SNRs^2

< list of 'numDraws' log-likelihood statistics

< list of 'numDraws' F-statistics

< list of 'numDraws' B-statistics

< list of 'numDraws' approximat B-statistics 'Bhat'

Definition at line 165 of file synthesizeBstatMC.c.

◆ initUserVars()

int initUserVars ( UserInput_t uvar)

Register all our "user-variables" that can be specified from cmd-line and/or config-file.

Here we set defaults for some user-variables and register them with the UserInput module.

Definition at line 321 of file synthesizeBstatMC.c.

◆ InitCode()

int InitCode ( ConfigVariables cfg,
const UserInput_t uvar 
)

Initialized Fstat-code: handle user-input and set everything up.

Definition at line 368 of file synthesizeBstatMC.c.

◆ XLALsynthesizeSignals()

int XLALsynthesizeSignals ( gsl_matrix **  A_Mu_i,
gsl_matrix **  s_mu_i,
gsl_matrix **  Amp_i,
gsl_vector **  rho2_i,
const gsl_matrix *  M_mu_nu,
AmpParamsRange_t  AmpRange,
gsl_rng *  rng,
UINT4  numDraws 
)

Generate random signal draws with uniform priors in given bands [h0, cosi, psi, phi0], and return list of 'numDraws' {s_mu} vectors.

Parameters
A_Mu_i[OUT] list of numDraws 4D line-vectors {A^nu}
s_mu_i[OUT] list of numDraws 4D line-vectors {s_mu = M_mu_nu A^nu}
Amp_i[OUT] list of numDraws 4D amplitude-parameters {h0, cosi, psi, phi}
rho2_i[OUT] optimal SNR^2
M_mu_nuantenna-pattern matrix M_mu_nu
AmpRangesignal amplitude-parameters ranges: lower bound + bands
rnggsl random-number generator
numDrawsnumber of random draws to synthesize

Definition at line 445 of file synthesizeBstatMC.c.

◆ XLALsynthesizeNoise()

int XLALsynthesizeNoise ( gsl_matrix **  n_mu_i,
const gsl_matrix *  M_mu_nu,
gsl_rng *  rng,
UINT4  numDraws 
)

Generate random-noise draws and combine with (FIXME: single!) signal.

Returns a list of numDraws vectors {x_mu}

Parameters
n_mu_i[OUT] list of numDraws 4D line-vectors of noise-components {n_mu}
M_mu_nu4x4 antenna-pattern matrix
rnggsl random-number generator
numDrawsnumber of random draws to synthesize

Definition at line 599 of file synthesizeBstatMC.c.

◆ XLALcomputeLogLikelihood()

int XLALcomputeLogLikelihood ( gsl_vector **  lnL_i,
const gsl_matrix *  A_Mu_i,
const gsl_matrix *  s_mu_i,
const gsl_matrix *  x_mu_i 
)

Compute log-likelihood function for given input data.

Parameters
lnL_i[OUT] log-likelihood vector
A_Mu_i4D amplitude-vector (FIXME: numDraws)
s_mu_i4D signal-component vector s_mu = (s|h_mu) [FIXME]
x_mu_inumDraws x 4D data-vectors x_mu

Definition at line 698 of file synthesizeBstatMC.c.

◆ XLALcomputeFstatistic()

int XLALcomputeFstatistic ( gsl_vector **  Fstat_i,
gsl_matrix **  A_Mu_MLE_i,
const gsl_matrix *  M_mu_nu,
const gsl_matrix *  x_mu_i 
)

Compute F-statistic for given input data.

Parameters
Fstat_i[OUT] F-statistic vector
A_Mu_MLE_i[OUT] vector of {A^mu_MLE} amplitude-vectors
M_mu_nuantenna-pattern matrix M_mu_nu
x_mu_idata-vectors x_mu: numDraws x 4

Definition at line 785 of file synthesizeBstatMC.c.

◆ XLALcomputeBstatisticMC()

int XLALcomputeBstatisticMC ( gsl_vector **  Bstat_i,
const gsl_matrix *  M_mu_nu,
const gsl_matrix *  x_mu_i,
gsl_rng *  rng,
UINT4  numMCpoints 
)

Compute the B-statistic for given input data, using Monte-Carlo integration for the marginalization over {cosi, psi}, while {h0, phi0} have been marginalized analytically.

Currently uses the Vegas Monte-Carlo integrator, which samples more densely where the integrand is larger.

Parameters
Bstat_i[OUT] vector of numDraws B-statistic values
M_mu_nuantenna-pattern matrix M_mu_nu
x_mu_idata-vectors x_mu: numDraws x 4
rnggsl random-number generator
numMCpointsnumber of points to use in Monte-Carlo integration

Definition at line 915 of file synthesizeBstatMC.c.

◆ XLALcomputeBstatisticGauss()

int XLALcomputeBstatisticGauss ( gsl_vector **  Bstat_i,
const gsl_matrix *  M_mu_nu,
const gsl_matrix *  x_mu_i 
)

Compute the B-statistic for given input data, using standard Gauss-Kronod integration (gsl_integration_qng) for the marginalization over {cosi, psi}, while {h0, phi0} have been marginalized analytically.

Parameters
Bstat_i[OUT] vector of numDraws B-statistic values
M_mu_nuantenna-pattern matrix M_mu_nu
x_mu_idata-vectors x_mu: numDraws x 4

Definition at line 1018 of file synthesizeBstatMC.c.

◆ XLALcomputeBhatStatistic()

gsl_vector * XLALcomputeBhatStatistic ( const gsl_matrix *  M_mu_nu,
const gsl_matrix *  x_mu_i 
)

Compute an approximation to the full B-statistic, without using any integrations! Returns Bhat vector of dimensions (numDraws x 1), or NULL on error.

Parameters
M_mu_nuantenna-pattern matrix M_mu_nu
x_mu_idata-vectors x_mu: numDraws x 4

Definition at line 1259 of file synthesizeBstatMC.c.

◆ BstatIntegrandOuter()

double BstatIntegrandOuter ( double  cosi,
void *  p 
)

log likelihood ratio lnL marginalized over {h0, phi0} (analytical) and integrated over psi in [-pi/4,pi/4], for given cosi: BstatIntegrandOuter(cosi) = int lnL dh0 dphi0 dpsi

This function is of type gsl_function for gsl-integration over cosi

Definition at line 1120 of file synthesizeBstatMC.c.

◆ BstatIntegrandInner()

double BstatIntegrandInner ( double  psi,
void *  p 
)

log likelihood ratio lnL marginalized over {h0, phi0} (analytical) for given psi and pars->cosi BstatIntegrandInner(cosi,psi) = int lnL dh0 dphi0

This function is of type gsl_function for gsl-integration over psi at fixed cosi, and represents a simple wrapper around BstatIntegrand() for gsl-integration.

Definition at line 1174 of file synthesizeBstatMC.c.

◆ BstatIntegrand()

double BstatIntegrand ( double  Amp[],
size_t  dim,
void *  p 
)

compute log likelihood ratio lnL for given Amp = {h0, cosi, psi, phi0} and M_{mu,nu}.

computes lnL = A^mu x_mu - 1/2 A^mu M_mu_nu A^nu.

This function is of type gsl_monte_function for gsl Monte-Carlo integration.

Definition at line 1197 of file synthesizeBstatMC.c.

◆ XLALComputeBhatCorrection()

REAL8 XLALComputeBhatCorrection ( const gsl_vector *  A_Mu,
const gsl_matrix *  M_mu_nu 
)

Compute 'Bstat' approximate correction terms with respect to F, ie deltaB in Bstat = F + deltaB.

Definition at line 1383 of file synthesizeBstatMC.c.

Variable Documentation

◆ vrbflg

int vrbflg
extern

defined in lal/lib/std/LALError.c