20 #ifndef _SUPERSKYMETRICS_H
21 #define _SUPERSKYMETRICS_H
24 #include <gsl/gsl_matrix.h>
25 #include <lal/LALStdlib.h>
26 #include <lal/UniversalDopplerMetric.h>
27 #include <lal/LatticeTiling.h>
28 #include <lal/FITSFileIO.h>
51 typedef enum tagSuperskyMetricType {
60 SWIGLAL( ARRAY_MULTIPLE_LENGTHS( tagSuperskyMetrics, num_segments ) );
62 typedef struct tagSuperskyMetrics {
66 SWIGLAL( ARRAY_1D(
SuperskyMetrics, gsl_matrix *, coh_rssky_metric,
size_t, num_segments ) );
70 SWIGLAL( ARRAY_1D(
SuperskyMetrics, SuperskyTransformData *, coh_rssky_transf,
size_t, num_segments ) );
84 const size_t spindowns,
87 const double fiducial_freq,
137 const double new_fiducial_freq
146 const double coh_max_mismatch,
147 const double semi_max_mismatch
155 const SuperskyTransformData *rssky_transf
162 gsl_vector *out_rssky,
164 const SuperskyTransformData *rssky_transf
172 const gsl_vector *in_rssky,
173 const gsl_vector *ref_rssky,
174 const SuperskyTransformData *rssky_transf
181 gsl_vector *out_rssky,
182 const SuperskyTransformData *out_rssky_transf,
183 const gsl_vector *in_rssky,
184 const gsl_vector *ref_rssky,
185 const SuperskyTransformData *in_rssky_transf
192 SWIGLAL( INOUT_STRUCTS( gsl_matrix **, out_rssky ) );
195 gsl_matrix **out_rssky,
196 const gsl_matrix *in_phys,
197 const SuperskyTransformData *rssky_transf
204 SWIGLAL( INOUT_STRUCTS( gsl_matrix **, out_phys ) );
207 gsl_matrix **out_phys,
208 const gsl_matrix *in_rssky,
209 const SuperskyTransformData *rssky_transf
213 SWIGLAL( COPYINOUT_ARRAYS( gsl_matrix, rssky_metric, rssky_transf ) );
223 LatticeTiling *tiling,
224 gsl_matrix *rssky_metric,
225 SuperskyTransformData *rssky_transf,
238 LatticeTiling *tiling,
239 const gsl_matrix *rssky_metric,
240 const double max_mismatch,
241 const UINT4 patch_count,
242 const UINT4 patch_index
246 SWIGLAL_CLEAR( COPYINOUT_ARRAYS( gsl_matrix, rssky_metric, rssky_transf ) );
254 LatticeTiling *tiling,
255 const SuperskyTransformData *rssky_transf,
266 LatticeTiling *tiling,
267 const SuperskyTransformData *rssky_transf,
281 LatticeTiling *tiling,
282 const SuperskyTransformData *rssky_transf,
292 SWIGLAL( OUTPUT_OWNED_BY_1ST_ARG( gsl_vector **, min_rssky2 ) );
293 SWIGLAL( OUTPUT_OWNED_BY_1ST_ARG( gsl_vector **, max_rssky2 ) );
296 LatticeTiling *tiling,
297 const SuperskyTransformData *rssky_transf,
298 const SuperskyTransformData *rssky2_transf,
299 const gsl_vector **min_rssky2,
300 const gsl_vector **max_rssky2
308 LatticeTiling *tiling,
309 const gsl_vector *min_rssky,
310 const gsl_vector *max_rssky
struct tagFITSFile FITSFile
Representation of a FITS file.
int XLALConvertSuperskyToPhysicalPoints(gsl_matrix **out_phys, const gsl_matrix *in_rssky, const SuperskyTransformData *rssky_transf)
Convert a set of points from supersky to physical coordinates.
int XLALSetSuperskyPhysicalSkyBounds(LatticeTiling *tiling, gsl_matrix *rssky_metric, SuperskyTransformData *rssky_transf, const double alpha1, const double alpha2, const double delta1, const double delta2)
Set parameter-space bounds on the physical sky position for a lattice tiling using the reduced super...
void XLALDestroySuperskyMetrics(SuperskyMetrics *metrics)
Destroy a SuperskyMetrics struct.
int XLALSetSuperskyEqualAreaSkyBounds(LatticeTiling *tiling, const gsl_matrix *rssky_metric, const double max_mismatch, const UINT4 patch_count, const UINT4 patch_index)
Divide the reduced supersky parameter space into patch_count equal-area patches, and set parameter-sp...
int XLALConvertPhysicalToSuperskyPoints(gsl_matrix **out_rssky, const gsl_matrix *in_phys, const SuperskyTransformData *rssky_transf)
Convert a set of points from physical to supersky coordinates.
int XLALSuperskyMetricsDimensions(const SuperskyMetrics *metrics, size_t *spindowns)
Return dimensions of the supersky metrics.
int XLALEqualizeReducedSuperskyMetricsFreqSpacing(SuperskyMetrics *metrics, const double coh_max_mismatch, const double semi_max_mismatch)
Project and rescale the reduced supersky metrics in the frequency dimension, such that all reduced su...
SuperskyMetrics * XLALCopySuperskyMetrics(const SuperskyMetrics *metrics)
Copy a SuperskyMetrics struct.
SuperskyMetricType
Type of supersky metric to compute.
int XLALScaleSuperskyMetricsFiducialFreq(SuperskyMetrics *metrics, const double new_fiducial_freq)
Scale all supersky metrics and their coordinate transform data to a new fiducial frequency.
int XLALSetPhysicalPointSuperskyRefTime(PulsarDopplerParams *out_phys, const SuperskyTransformData *rssky_transf)
Set the reference time of a physical point to that of the reduced supersky coordinates.
int XLALConvertPhysicalToSuperskyPoint(gsl_vector *out_rssky, const PulsarDopplerParams *in_phys, const SuperskyTransformData *rssky_transf)
Convert a point from physical to supersky coordinates.
int XLALConvertSuperskyToPhysicalPoint(PulsarDopplerParams *out_phys, const gsl_vector *in_rssky, const gsl_vector *ref_rssky, const SuperskyTransformData *rssky_transf)
Convert a point from supersky to physical coordinates.
int XLALFITSWriteSuperskyMetrics(FITSFile *file, const SuperskyMetrics *metrics)
Write a SuperskyMetrics struct to a FITS file.
int XLALSetSuperskyRangeBounds(LatticeTiling *tiling, const gsl_vector *min_rssky, const gsl_vector *max_rssky)
Set parameter-space bounds on an entire lattice tiling given minimum and maximum ranges in reduced su...
int XLALRegisterSuperskyLatticePhysicalRangeCallback(LatticeTiling *tiling, const SuperskyTransformData *rssky_transf, const PulsarDopplerParams **min_phys, const PulsarDopplerParams **max_phys)
Register a lattice tiling callback function which computes the physical range covered by a reduced su...
int XLALConvertSuperskyToSuperskyPoint(gsl_vector *out_rssky, const SuperskyTransformData *out_rssky_transf, const gsl_vector *in_rssky, const gsl_vector *ref_rssky, const SuperskyTransformData *in_rssky_transf)
Convert a point between supersky coordinates.
int XLALSetSuperskyPhysicalSpinBoundPadding(LatticeTiling *tiling, const SuperskyTransformData *rssky_transf, const size_t s, const bool padding)
Set parameter-space bound padding on the physical frequency/spindowns for a lattice tiling using the...
int XLALSetSuperskyPhysicalSpinBound(LatticeTiling *tiling, const SuperskyTransformData *rssky_transf, const size_t s, const double bound1, const double bound2)
Set parameter-space bounds on the physical frequency/spindowns for a lattice tiling using the reduce...
SuperskyMetrics * XLALComputeSuperskyMetrics(const SuperskyMetricType type, const size_t spindowns, const LIGOTimeGPS *ref_time, const LALSegList *segments, const double fiducial_freq, const MultiLALDetector *detectors, const MultiNoiseFloor *detector_weights, const DetectorMotionType detector_motion, const EphemerisData *ephemerides)
Compute the supersky metrics, which are returned in a SuperskyMetrics struct.
int XLALRegisterSuperskyLatticeSuperskyRangeCallback(LatticeTiling *tiling, const SuperskyTransformData *rssky_transf, const SuperskyTransformData *rssky2_transf, const gsl_vector **min_rssky2, const gsl_vector **max_rssky2)
Register a lattice tiling callback function which computes the range covered by a reduced supersky la...
int XLALFITSReadSuperskyMetrics(FITSFile *file, SuperskyMetrics **metrics)
Read a SuperskyMetrics struct from a FITS file.
@ SUPERSKY_METRIC_TYPE
Metric for all-sky searches.
DetectorMotionType
Bitfield of different types of detector-motion to use in order to compute the Doppler-metric.
This structure contains all information about the center-of-mass positions of the Earth and Sun,...
array of detectors definitions 'LALDetector'
array of detector-specific 'noise floors' (ie PSD values), assumed constant over the frequency-band o...
Type containing the 'Doppler-parameters' affecting the time-evolution of the phase.
Computed supersky metrics, returned by XLALComputeSuperskyMetrics().
gsl_matrix * semi_rssky_metric
Semicoherent reduced supersky metric (2-dimensional sky)
gsl_matrix ** coh_rssky_metric
Coherent reduced supersky metric (2-dimensional sky) for each segment.
size_t num_segments
Number of segments.
SuperskyTransformData * semi_rssky_transf
Semicoherent reduced supersky metric coordinate transform data.
SuperskyTransformData ** coh_rssky_transf
Coherent reduced supersky metric coordinate transform data for each segment.