Loading [MathJax]/extensions/TeX/AMSsymbols.js
LALPulsar 7.1.1.1-5e288d3
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
MetricUtils.h
Go to the documentation of this file.
1//
2// Copyright (C) 2011--2015 Reinhard Prix, Karl Wette
3//
4// This program is free software; you can redistribute it and/or modify
5// it under the terms of the GNU General Public License as published by
6// the Free Software Foundation; either version 2 of the License, or
7// (at your option) any later version.
8//
9// This program is distributed in the hope that it will be useful,
10// but WITHOUT ANY WARRANTY; without even the implied warranty of
11// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12// GNU General Public License for more details.
13//
14// You should have received a copy of the GNU General Public License
15// along with with program; see the file COPYING. If not, write to the
16// Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
17// MA 02110-1301 USA
18//
19
20#ifndef _METRICUTILS_H
21#define _METRICUTILS_H
22
23#include <gsl/gsl_vector.h>
24#include <gsl/gsl_matrix.h>
25#include <lal/LALStdlib.h>
26#include <lal/UniversalDopplerMetric.h>
27
28#ifdef __cplusplus
29extern "C" {
30#endif
31
32///
33/// \defgroup MetricUtils_h Header MetricUtils.h
34/// \ingroup lalpulsar_metric
35/// \author Reinhard Prix, Karl Wette
36/// \brief Various useful utility functions for working with CW parameter-space metrics.
37///
38/// @{
39///
40
41REAL8 XLALCompareMetrics( const gsl_matrix *g1_ij, const gsl_matrix *g2_ij );
42
43double XLALMetricDeterminant( const gsl_matrix *g_ij );
44gsl_vector *XLALMetricEllipseBoundingBox( const gsl_matrix *g_ij, const double max_mismatch );
45
46int XLALProjectMetric( gsl_matrix **p_gpr_ij, const gsl_matrix *g_ij, const UINT4 c );
47int XLALCholeskyLDLTDecompMetric( gsl_matrix **p_cholesky, const gsl_matrix *g_ij );
48
49int XLALTransformMetric( gsl_matrix **p_gpr_ij, const gsl_matrix *transform, const gsl_matrix *g_ij );
50int XLALInverseTransformMetric( gsl_matrix **p_gpr_ij, const gsl_matrix *transform, const gsl_matrix *g_ij );
51
52int XLALDiagNormalizeMetric( gsl_matrix **p_gpr_ij, gsl_matrix **p_transform, const gsl_matrix *g_ij );
53int XLALNaturalizeMetric( gsl_matrix **p_gpr_ij, gsl_matrix **p_transform, const gsl_matrix *g_ij, const DopplerMetricParams *metricParams );
54int XLALChangeMetricReferenceTime( gsl_matrix **p_gpr_ij, gsl_matrix **p_transform, const gsl_matrix *g_ij, const DopplerCoordinateSystem *coordSys, const double Dtau );
55
56/// @}
57
58#ifdef __cplusplus
59}
60#endif
61
62#endif // _METRICUTILS_H
double REAL8
uint32_t UINT4
int XLALChangeMetricReferenceTime(gsl_matrix **p_gpr_ij, gsl_matrix **p_transform, const gsl_matrix *g_ij, const DopplerCoordinateSystem *coordSys, const double Dtau)
Compute the transform which changes the metric reference time .
Definition: MetricUtils.c:578
int XLALProjectMetric(gsl_matrix **p_gpr_ij, const gsl_matrix *g_ij, const UINT4 c)
Calculate the projected metric onto the subspace orthogonal to coordinate-axis c, namely ,...
Definition: MetricUtils.c:182
REAL8 XLALCompareMetrics(const gsl_matrix *g1_ij, const gsl_matrix *g2_ij)
Flexible comparison function between two metrics and .
Definition: MetricUtils.c:52
int XLALCholeskyLDLTDecompMetric(gsl_matrix **p_cholesky, const gsl_matrix *g_ij)
Decompose a metric as , where is a lower-triangular matrix with unit diagonal, and is a diagonal ...
Definition: MetricUtils.c:231
int XLALInverseTransformMetric(gsl_matrix **p_gpr_ij, const gsl_matrix *transform, const gsl_matrix *g_ij)
Apply the inverse of a transform to a metric such that , or equivalently .
Definition: MetricUtils.c:341
int XLALDiagNormalizeMetric(gsl_matrix **p_gpr_ij, gsl_matrix **p_transform, const gsl_matrix *g_ij)
Diagonally-normalize a metric .
Definition: MetricUtils.c:391
double XLALMetricDeterminant(const gsl_matrix *g_ij)
Compute the determinant of a metric .
Definition: MetricUtils.c:98
int XLALTransformMetric(gsl_matrix **p_gpr_ij, const gsl_matrix *transform, const gsl_matrix *g_ij)
Apply a transform to a metric such that , or equivalently .
Definition: MetricUtils.c:286
gsl_vector * XLALMetricEllipseBoundingBox(const gsl_matrix *g_ij, const double max_mismatch)
Compute the extent of the bounding box of the mismatch ellipse of a metric .
Definition: MetricUtils.c:131
int XLALNaturalizeMetric(gsl_matrix **p_gpr_ij, gsl_matrix **p_transform, const gsl_matrix *g_ij, const DopplerMetricParams *metricParams)
Return a metric in naturalized coordinates.
Definition: MetricUtils.c:456
type describing a Doppler coordinate system: lists the number of dimensions and the symbolic names of...
meta-info specifying a Doppler-metric