Loading [MathJax]/extensions/TeX/AMSsymbols.js
LALInference 4.1.9.1-8a6b96f
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
LALInferencePriorTest.c
Go to the documentation of this file.
1/*
2 * LALInferencePriorTest.c: Testing the Nested Sampling routines in LALInferencePrior.c
3 *
4 * Copyright (C) 2011 Ben Aylott, Ilya Mandel, Chiara Mingarelli, Vivien Raymond, Christian Roever, Marc van der Sluys, John Veitch, Will Vousden
5 *
6 *
7 * This program is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * This program is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with with program; see the file COPYING. If not, write to the
19 * Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
20 * MA 02110-1301 USA
21 */
22
23#include <stdio.h>
24#include <lal/LALInferencePrior.h>
25#include <lal/LALInference.h>
26#include <lal/XLALError.h>
27#include <math.h>
28#include "LALInferenceTest.h"
29
30#define EPSILON 1e-13
31
37
38//Old tests
41
42int main(void)
43{
44 int failureCount = 0;
45
46 failureCount += computePriorMassNormTest();
47 printf("\n");
48 failureCount += LALInferenceRotateInitialPhaseTest();
49 printf("\n");
51 printf("\n");
52 failureCount += LALInferenceDrawFromPriorTest();
53 printf("\n");
54 failureCount += LALInferenceInspiralPriorTest();
55 printf("\n");
56
57 printf("Test results: %i failure(s).\n", failureCount);
58 return failureCount;
59}
60
62{
64
65 int errnum;
66 double result;
67
68 char massRatioName[VARNAME_MAX];
69 REAL8 MMin;
70 REAL8 MMax;
71 REAL8 MTotMax;
72 REAL8 McMin;
73 REAL8 McMax;
74 REAL8 massRatioMin;
75 REAL8 massRatioMax;
76
77 MMin = 1;
78 MMax = 10;
79 MTotMax = 100;
80 McMin = 1;
81 McMax = 2;
82 massRatioMin = 1;
83 massRatioMax = 10;
84 XLAL_TRY(result = LALInferenceComputePriorMassNorm(MMin, MMax, MTotMax, McMin, McMax, massRatioMin, massRatioMax, NULL), errnum);
85 if (!XLAL_IS_REAL8_FAIL_NAN(result) || errnum != XLAL_EFAULT)
86 TEST_FAIL("Null reference check failed.");
87
88 strcpy(massRatioName, "foo");
89 XLAL_TRY(result = LALInferenceComputePriorMassNorm(MMin, MMax, MTotMax, McMin, McMax, massRatioMin, massRatioMax, massRatioName), errnum);
90 if (!XLAL_IS_REAL8_FAIL_NAN(result) || errnum != XLAL_ENAME)
91 TEST_FAIL("Invalid mass ratio name specified but appropriate error not generated.");
92
93 strcpy(massRatioName, "q");
94 MMin = 1;
95 MMax = -10;
96 MTotMax = -100;
97 McMin = -1;
98 McMax = -2;
99 massRatioMin = 1;
100 massRatioMax = 10;
101 XLAL_TRY(result = LALInferenceComputePriorMassNorm(MMin, MMax, MTotMax, McMin, McMax, massRatioMin, massRatioMax, massRatioName), errnum);
102 if (!XLAL_IS_REAL8_FAIL_NAN(result) || errnum == XLAL_SUCCESS)
103 TEST_FAIL("Unphysical masses given but appropriate error not generated.");
104
105 TEST_FOOTER();
106}
107
109{
110 TEST_HEADER();
111
112 int errnum;
113
114 // A basic null reference check.
116 if (errnum != XLAL_EFAULT)
117 TEST_FAIL("Null reference check failed.");
118
119 // Construct a variable list containing phi0 and psi.
120 REAL8 psi = 0;
121 REAL8 phi0 = 0;
122 REAL8 phi0_2 = 0;
126
127 // Check that if psi is in [0,2pi], phi0 remains unchanged.
128 psi = LAL_PI;
129 phi0 = 0;
131 LALInferenceSetVariable(variables, "phi0", &phi0);
133 phi0_2 = *(REAL8 *)LALInferenceGetVariable(variables, "phi0");
134 if (!compareFloats(phi0_2, phi0, EPSILON))
135 TEST_FAIL("Psi in [0,2pi] but phi0 has changed!");
136
137 // Check that if psi is outside [0,2pi], phi0 is rotated.
138 psi = -2 * LAL_TWOPI;
139 phi0 = 0;
141 LALInferenceSetVariable(variables, "phi0", &phi0);
143 phi0_2 = *(REAL8 *)LALInferenceGetVariable(variables, "phi0");
144 if (!compareFloats(phi0_2, phi0 + LAL_PI, EPSILON))
145 TEST_FAIL("Psi outside [0,2pi] but phi0 not rotated by 2pi!");
146 psi = 2 * LAL_TWOPI;
147 phi0 = 0;
149 LALInferenceSetVariable(variables, "phi0", &phi0);
151 phi0_2 = *(REAL8 *)LALInferenceGetVariable(variables, "phi0");
152 if (!compareFloats(phi0_2, phi0 + LAL_PI, EPSILON))
153 TEST_FAIL("Psi outside [0,2pi] but phi0 not rotated by 2pi!");
154
155 // Check boundary cases: if psi=0 or psi=2pi, phi0 shouldn't be rotated.
156 psi = 0;
157 phi0 = 0;
159 LALInferenceSetVariable(variables, "phi0", &phi0);
161 phi0_2 = *(REAL8 *)LALInferenceGetVariable(variables, "phi0");
162 if (phi0 != phi0_2) // Should be exact.)
163 TEST_FAIL("Psi on boundary of [0,2pi] but phi0 has changed!");
164 psi = LAL_TWOPI;
165 phi0 = 0;
167 LALInferenceSetVariable(variables, "phi0", &phi0);
169 phi0_2 = *(REAL8 *)LALInferenceGetVariable(variables, "phi0");
170 if (phi0 != phi0_2) // Should be exact.)
171 TEST_FAIL("Psi on boundary of [0,2pi] but phi0 has changed!");
172
174 TEST_FOOTER();
175}
176
178{
179 TEST_HEADER();
180
181 int errnum;
182 int outcome;
183
185 REAL8 a;
186 REAL8 b;
187 REAL8 a_2;
188 REAL8 b_2;
189 REAL8 a_min;
190 REAL8 a_max;
191 REAL8 a_delta;
192 REAL8 b_min;
193 REAL8 b_max;
194 REAL8 b_delta;
197
198 // A basic null reference check.
199 outcome = 1;
200 XLAL_TRY(LALInferenceCyclicReflectiveBound(NULL, priorArgs), errnum);
201 outcome &= errnum == XLAL_EFAULT;
203 outcome &= errnum == XLAL_EFAULT;
204 if (!outcome)
205 TEST_FAIL("Null reference check failed.");
206
207 // Check some (meaningful) minima/maxima.
208 a_min = -LAL_PI;
209 a_max = LAL_PI;
210 a_delta = a_max - a_min;
211 b_min = -1;
212 b_max = 1;
213 b_delta = b_max - b_min;
214 LALInferenceRemoveMinMaxPrior(priorArgs, "a");
215 LALInferenceRemoveMinMaxPrior(priorArgs, "b");
216 LALInferenceAddMinMaxPrior(priorArgs, "a", &a_min, &a_max, type);
217 LALInferenceAddMinMaxPrior(priorArgs, "b", &b_min, &b_max, type);
218
219 // Variables within [min,max]: should remain unchanged.
220 a = a_min + (a_max - a_min) / 2;
221 b = b_min + (b_max - b_min) / 2;
227 if (!compareFloats(a, a_2, EPSILON) || !compareFloats(b, b_2, EPSILON))
228 TEST_FAIL("Values within bounds should remain unchanged.");
229
230 // Boundary cases (circular): variables on [min, max] boundaries should be equal modulo period.
231 outcome = 1;
232 a = a_min;
236 outcome &= compareFloats(fmod(a - a_2, a_delta), 0, EPSILON);
237 a = a_max;
241 outcome &= compareFloats(fmod(a - a_2, a_delta), 0, EPSILON);
242 if (!outcome)
243 TEST_FAIL("Circular boundary values should remain equal modulo their period.");
244
245 // Boundary cases (linear): variables on [min, max] boundaries should be equal modulo period.
246 outcome = 1;
247 b = b_min;
251 outcome &= compareFloats(b, b_2, EPSILON);
252 b = b_max;
256 outcome &= compareFloats(b, b_2, EPSILON);
257 if (!outcome)
258 TEST_FAIL("Linear boundary values should remain unchanged.");
259
260 // Outside range (circular).
261 outcome = 1;
262 a = a_min - a_delta / 3;
266 outcome &= compareFloats(a_2, a_max - a_delta / 3, EPSILON);
267 a = a_max + a_delta / 3;
271 outcome &= compareFloats(a_2, a_min + a_delta / 3, EPSILON);
272 if (!outcome)
273 TEST_FAIL("Circular values outside range should be correctly modded into range.");
274
275 // Outside range (linear).
276 outcome = 1;
277 b = b_min - 10 * b_delta / 3;
281 outcome &= compareFloats(b_2, b_max - b_delta / 3, EPSILON);
282 b = b_max + 7 * b_delta / 5;
286 outcome &= compareFloats(b_2, b_min + 2 * b_delta / 5, EPSILON);
287 if (!outcome)
288 TEST_FAIL("Linear values outside range should be correctly reflected into range.");
289
291 XLALFree(priorArgs);
292 TEST_FOOTER();
293}
294
296{
297 TEST_HEADER();
298
299 int errnum;
300 const char *name;
301 gsl_rng *rng = gsl_rng_alloc(gsl_rng_default);
302 gsl_rng_set(rng, 0);
303
306
307 // Null reference checks.
308 int outcome = 1;
309 XLAL_TRY(LALInferenceDrawFromPrior(NULL, priorArgs, rng), errnum);
310 outcome &= errnum == XLAL_EFAULT;
311 XLAL_TRY(LALInferenceDrawFromPrior(output, NULL, rng), errnum);
312 outcome &= errnum == XLAL_EFAULT;
313 XLAL_TRY(LALInferenceDrawFromPrior(output, priorArgs, NULL), errnum);
314 outcome &= errnum == XLAL_EFAULT;
315 if (!outcome)
316 TEST_FAIL("Null reference check failed.");
317
318 int i;
319 const char *varyName=NULL;
320 char caseTag[VARNAME_MAX];
323 for (i = 0; i < 2; i++)
324 {
325 switch (i)
326 {
327 case 0:
329 varyName = "linear";
330 break;
331 case 1:
333 varyName = "circular";
334 break;
335 }
336 sprintf(caseTag, "[%s] ", varyName);
337
338 // Try and generate some normally distributed variables for various mu and sigma.
339 REAL8 gaussian = 0;
340 REAL8 mu;
341 REAL8 sigma;
342 name = "gaussian";
343 LALInferenceAddVariable(output, name, &gaussian, type, vary);
344
345 // Zero standard deviation; should always equal mean.
346 mu = -50;
347 sigma = 0;
349 LALInferenceAddGaussianPrior(priorArgs, name, &mu, &sigma, type);
350 XLAL_TRY(LALInferenceDrawFromPrior(output, priorArgs, rng), errnum);
352 if (errnum != XLAL_SUCCESS)
353 {
354 TEST_FAIL("%sFailed to generate Gaussian variable; XLAL error: %s.", caseTag, XLALErrorString(errnum));
355 }
356 else if (!compareFloats(gaussian, mu, EPSILON))
357 {
358 TEST_FAIL("%sGaussian variable with zero standard deviation did not match the mean; X = %f, mu = %f.", caseTag, gaussian, mu);
359 }
360
363
364 // Try a uniform variable!
365 REAL8 uniform = 0;
366 REAL8 min;
367 REAL8 max;
368 name = "uniform";
369 LALInferenceAddVariable(output, name, &uniform, type, vary);
370
371 min = -1;
372 max = 1;
374 LALInferenceAddMinMaxPrior(priorArgs, name, &min, &max, type);
375 XLAL_TRY(LALInferenceDrawFromPrior(output, priorArgs, rng), errnum);
376 if (errnum != XLAL_SUCCESS)
377 TEST_FAIL("%sFailed to generate uniform variable; XLAL error: %s.", caseTag, XLALErrorString(errnum));
378
381
382 // Try a correlated variable!
383 REAL8 correlated = 0;
384 UINT4 idx = 0;
385 name = "correlated";
386 gsl_matrix *covariance = gsl_matrix_calloc(3, 3);
387 LALInferenceAddVariable(output, name, &correlated, type, vary);
389
390 // See what happens when we try to add a non-positive-definite covariance matrix
391 gsl_matrix_set(covariance, 0, 0, -1);
392 XLAL_TRY(LALInferenceAddCorrelatedPrior(priorArgs, name, &covariance, &mu, &sigma, &idx), errnum);
393 if (errnum == XLAL_SUCCESS)
394 TEST_FAIL("%sNon-positive-definite covariance matrix was not rejected.", caseTag);
396
397 // Now try a positive-semi-definite matrix; this should be accepted
398 covariance = gsl_matrix_calloc(3, 3);
399 gsl_matrix_set(covariance, 0, 0, 1);
400 XLAL_TRY(LALInferenceAddCorrelatedPrior(priorArgs, name, &covariance, &mu, &sigma, &idx), errnum);
401 if (errnum != XLAL_SUCCESS)
402 TEST_FAIL("%sCould not add semi-positive-definite covariance matrix.", caseTag);
404
405 // Try a legitimate positive-definite covariance matrix.
406 covariance = gsl_matrix_calloc(3, 3);
407 gsl_matrix_set(covariance, 0, 0, 2);
408 gsl_matrix_set(covariance, 0, 1, 1);
409 gsl_matrix_set(covariance, 0, 2, 0);
410 gsl_matrix_set(covariance, 1, 0, 1);
411 gsl_matrix_set(covariance, 1, 1, 5);
412 gsl_matrix_set(covariance, 1, 2, 1);
413 gsl_matrix_set(covariance, 2, 0, 0);
414 gsl_matrix_set(covariance, 2, 1, 1);
415 gsl_matrix_set(covariance, 2, 2, 1);
416 XLAL_TRY(LALInferenceAddCorrelatedPrior(priorArgs, name, &covariance, &mu, &sigma, &idx), errnum);
417 if (errnum != XLAL_SUCCESS)
418 TEST_FAIL("%sCould not add correlated prior.", caseTag);
419 XLAL_TRY(LALInferenceDrawFromPrior(output, priorArgs, rng), errnum);
420 if (errnum != XLAL_SUCCESS)
421 TEST_FAIL("%sCould not generate correlated variable from positive-definite matrix; XLAL error: %s.", caseTag, XLALErrorString(errnum));
422
425
426 //gsl_matrix_free(covariance);
429 LALInferenceRemoveVariable(output, "correlated");
430 }
431
433 XLALFree(priorArgs);
434 TEST_FOOTER();
435}
436
438{
439 TEST_HEADER();
440
441 int errnum;
442
443 REAL8 result;
447 runState->threads = XLALCalloc(1,sizeof(LALInferenceThreadState));
448 LALInferenceInitThread(&runState->threads[0]);
449
450 // Standard null reference check.
451 int failed = 1;
452 runState->priorArgs = NULL;
453 LALInferenceThreadState *thread = &runState->threads[0];
454 XLAL_TRY(result = LALInferenceInspiralPrior(runState, params, thread->model), errnum);
455 failed &= !XLAL_IS_REAL8_FAIL_NAN(result) || errnum != XLAL_EFAULT;
456 runState->priorArgs = priorArgs;
457 XLAL_TRY(result = LALInferenceInspiralPrior(NULL, params, thread->model), errnum);
458 failed &= !XLAL_IS_REAL8_FAIL_NAN(result) || errnum != XLAL_EFAULT;
459 XLAL_TRY(result = LALInferenceInspiralPrior(runState, NULL, thread->model), errnum);
460 failed &= !XLAL_IS_REAL8_FAIL_NAN(result) || errnum != XLAL_EFAULT;
461 if (failed)
462 TEST_FAIL("Null reference check failed.");
463
464 // Set up parameters.
465 REAL8 value = 0;
475 /*LALInferenceAddVariable(params, "q", &value, LALINFERENCE_REAL8_t, LALINFERENCE_PARAM_LINEAR);*/
477
478 REAL8 min, max;
479 min = 1.0; max = 100.0;
480 LALInferenceAddMinMaxPrior(priorArgs, "distance", &min, &max, LALINFERENCE_REAL8_t);
481 min = log(min); max = log(max);
482 LALInferenceAddMinMaxPrior(priorArgs, "logdistance", &min, &max, LALINFERENCE_REAL8_t);
483 min = 1.0; max = 20.5;
484 LALInferenceAddMinMaxPrior(priorArgs, "chirpmass", &min, &max, LALINFERENCE_REAL8_t);
485 min = log(min); max = log(max);
486 LALInferenceAddMinMaxPrior(priorArgs, "logmc", &min, &max, LALINFERENCE_REAL8_t);
487 min = 1.0; max = 30.0;
488 LALInferenceAddMinMaxPrior(priorArgs, "mass1", &min, &max, LALINFERENCE_REAL8_t);
489 LALInferenceAddMinMaxPrior(priorArgs, "mass2", &min, &max, LALINFERENCE_REAL8_t);
490 /*max *= 2;*/
491 /*LALInferenceAddVariable(priorArgs, "MTotMax", &max, LALINFERENCE_REAL8_t, LALINFERENCE_PARAM_FIXED);*/
492 min = -LAL_PI; max = LAL_PI;
493 LALInferenceAddMinMaxPrior(priorArgs, "inclination", &min, &max, LALINFERENCE_REAL8_t);
494 min = 0; max = LAL_TWOPI;
495 LALInferenceAddMinMaxPrior(priorArgs, "rightascension", &min, &max, LALINFERENCE_REAL8_t);
496 min = -LAL_PI / 2.0; max = LAL_PI / 2.0;
497 LALInferenceAddMinMaxPrior(priorArgs, "declination", &min, &max, LALINFERENCE_REAL8_t);
498 min = -LAL_PI / 2.0; max = LAL_PI / 2.0;
499 LALInferenceAddMinMaxPrior(priorArgs, "theta_spin1", &min, &max, LALINFERENCE_REAL8_t);
500 min = -LAL_PI / 2.0; max = LAL_PI / 2.0;
501 LALInferenceAddMinMaxPrior(priorArgs, "theta_spin2", &min, &max, LALINFERENCE_REAL8_t);
502 min = 0.01; max = 0.25;
503 LALInferenceAddMinMaxPrior(priorArgs, "eta", &min, &max, LALINFERENCE_REAL8_t);
504
505 // Pick a random point in the non-zero region of the parameter space.
506 gsl_rng *rng = gsl_rng_alloc(gsl_rng_default);
507 gsl_rng_set(rng, 0);
508 LALInferenceDrawFromPrior(params, priorArgs, rng);
509
510 // Check that we get a finite log prior.
511 XLAL_TRY(result = LALInferenceInspiralPrior(runState, params, thread->model), errnum);
512
513 if (XLAL_IS_REAL8_FAIL_NAN(result) || errnum != XLAL_SUCCESS)
514 {
515 TEST_FAIL("Could not generate inspiral prior; XLAL error: %s", XLALErrorString(errnum));
516 }
517 else if (isinf(result))
518 {
519 TEST_FAIL("Parameter configuration within specified min/max bounds for each parameter gave zero prior.");
520 }
521
522 // Now set a parameter outside its bounds and see what happens.
523 LALInferenceGetMinMaxPrior(priorArgs, "distance", &min, &max);
524 value = max + (max - min) / 2;
525 LALInferenceSetVariable(params, "distance", &value);
526 XLAL_TRY(result = LALInferenceInspiralPrior(runState, params, thread->model), errnum);
527 if (XLAL_IS_REAL8_FAIL_NAN(result) || errnum != XLAL_SUCCESS)
528 {
529 TEST_FAIL("Could not generate inspiral prior; XLAL error: %s", XLALErrorString(errnum));
530 }
531 else if (isfinite(result))
532 {
533 TEST_FAIL("Distance %f is outside [%f,%f] but prior is non-zero.", value, min, max);
534 }
535
536 // Try another configuration; this time set m1 and m2 such that one is *outside* its bounds,
537 // but the chirp mass and symmetric mass ratio are still OK; this should be picked up and a
538 // zero prior returned.
539 LALInferenceDrawFromPrior(params, priorArgs, rng);
540 LALInferenceGetMinMaxPrior(priorArgs, "mass1", &min, &max);
541 REAL8 m2 = 0.5;
542 REAL8 m1 = 3.82;
543 REAL8 eta = m1 * m2 / pow(m1 + m2, 2);
545 REAL8 Mc = pow(m1 * m2, 3.0 / 5.0) / pow(m1 + m2, 1.0 / 5.0);
546 LALInferenceSetVariable(params, "chirpmass", &Mc);
547 REAL8 logMc = log(Mc);
548 LALInferenceSetVariable(params, "logmc", &logMc);
549 XLAL_TRY(result = LALInferenceInspiralPrior(runState, params, thread->model), errnum);
550 if (XLAL_IS_REAL8_FAIL_NAN(result) || errnum != XLAL_SUCCESS)
551 {
552 TEST_FAIL("Could not generate inspiral prior; XLAL error: %s", XLALErrorString(errnum));
553 }
554 else if (isfinite(result))
555 {
556 TEST_FAIL("Mass ratio %f and chirp mass %f define masses outside bounds [%f,%f], but prior is non-zero.", eta, Mc, min, max);
557 }
558
559 TEST_FOOTER();
560}
561
562/******************************************
563 *
564 * Old tests
565 *
566 ******************************************/
567
568//Test LALInferencePriorFunction
570/****************************************/
571/* Returns unnormalized (!), */
572/* logarithmic (!) prior density. */
573/****************************************/
574/* Assumes the following parameters */
575/* exist (e.g., for TaylorT1): */
576/* chirpmass, massratio, inclination, */
577/* phase, time, rightascension, */
578/* desclination, polarisation, distance.*/
579/* Prior is flat if within range */
580/****************************************/
581{
582 (void) runState; /* avoid warning about unused parameter */
583 REAL8 eta, iota, phi, ra, dec, psi;
584 REAL8 logdensity;
585
586 // UNUSED!!: REAL8 mc = *(REAL8*) LALInferenceGetVariable(params, "chirpmass"); /* solar masses*/
587 eta = *(REAL8*) LALInferenceGetVariable(params, "eta"); /* dim-less */
588 iota = *(REAL8*) LALInferenceGetVariable(params, "inclination"); /* radian */
589 // UNUSED!!: REAL8 tc = *(REAL8*) LALInferenceGetVariable(params, "time"); /* GPS seconds */
590 phi = *(REAL8*) LALInferenceGetVariable(params, "phase"); /* radian */
591 ra = *(REAL8*) LALInferenceGetVariable(params, "rightascension"); /* radian */
592 dec = *(REAL8*) LALInferenceGetVariable(params, "declination"); /* radian */
593 psi = *(REAL8*) LALInferenceGetVariable(params, "polarisation"); /* radian */
594 // UNUSED!!: REAL8 dist = *(REAL8*) LALInferenceGetVariable(params, "distance"); /* Mpc */
595
596 if(eta>0.0 && eta<=0.25 && iota>=0.0 && iota<=LAL_PI && phi>=0.0 && phi<=LAL_TWOPI
597 && ra>=0.0 && ra<=LAL_TWOPI && dec>=-LAL_PI_2 && dec<=LAL_PI_2 && psi>=0.0 && psi<=LAL_PI)
598 logdensity = 0.0;
599 else
600 logdensity = -HUGE_VAL;
601 //TODO: should be properly normalized; pass in range via priorArgs?
602
603 return(logdensity);
604}
605
606//Test LALPriorFunction
608/****************************************/
609/* Prior for two-parameter */
610/* waveform family ASinOmegaT */
611/* Assumes the following parameters */
612/* exist: A, Omega */
613/* Prior is flat if within range */
614/****************************************/
615{
616 (void) runState; /* avoid warning about unused parameter */
617 REAL8 A, Omega;
618 REAL8 logdensity;
619
620 A = *(REAL8*) LALInferenceGetVariable(params, "A"); /* dim-less */
621 Omega = *(REAL8*) LALInferenceGetVariable(params, "Omega"); /* rad/sec */
622
623 if ((A>0.0) & (Omega>0))
624 logdensity = 0.0;
625 else
626 logdensity = -HUGE_VAL;
627
628 return logdensity;
629}
#define EPSILON
int LALInferenceInspiralPriorTest(void)
REAL8 BasicUniformLALPrior(LALInferenceRunState *runState, LALInferenceVariables *params)
int main(void)
int LALInferenceRotateInitialPhaseTest(void)
int LALInferenceCyclicReflectiveBoundTest(void)
int computePriorMassNormTest(void)
REAL8 ASinOmegaTPrior(LALInferenceRunState *runState, LALInferenceVariables *params)
int LALInferenceDrawFromPriorTest(void)
#define max(a, b)
LALInferenceVariables variables
int compareFloats(REAL8 x, REAL8 y, REAL8 epsilon)
#define TEST_FAIL(...)
#define TEST_FOOTER()
#define TEST_HEADER()
#define MMax
const char *const name
#define LAL_PI_2
#define LAL_PI
#define LAL_TWOPI
double REAL8
uint32_t UINT4
#define VARNAME_MAX
Definition: LALInference.h:50
void LALInferenceAddVariable(LALInferenceVariables *vars, const char *name, const void *value, LALInferenceVariableType type, LALInferenceParamVaryType vary)
Add a variable named name to vars with initial value referenced by value.
Definition: LALInference.c:395
void LALInferenceRemoveVariable(LALInferenceVariables *vars, const char *name)
Remove name from vars Frees the memory for the name structure and its contents.
Definition: LALInference.c:450
LALInferenceParamVaryType
An enumerated type for denoting the topology of a parameter.
Definition: LALInference.h:127
void * LALInferenceGetVariable(const LALInferenceVariables *vars, const char *name)
Return a pointer to the memory the variable vars is stored in specified by name User must cast this p...
Definition: LALInference.c:238
LALInferenceThreadState * LALInferenceInitThread(LALInferenceThreadState *thread)
Structure to contain data-related Reduced Order Quadrature quantities.
Definition: LALInference.c:105
LALInferenceVariableType
An enumerated type for denoting the type of a variable.
Definition: LALInference.h:104
void LALInferenceSetVariable(LALInferenceVariables *vars, const char *name, const void *value)
Set a variable named name in vars with a value.
Definition: LALInference.c:352
@ LALINFERENCE_PARAM_CIRCULAR
A parameter that simply has a maximum and a minimum.
Definition: LALInference.h:129
@ LALINFERENCE_PARAM_LINEAR
Definition: LALInference.h:128
@ LALINFERENCE_REAL8_t
Definition: LALInference.h:109
void LALInferenceRemoveCorrelatedPrior(LALInferenceVariables *priorArgs)
Remove the correlation coefficient matrix and index for a parameter from the priorArgs list.
void LALInferenceAddGaussianPrior(LALInferenceVariables *priorArgs, const char *name, REAL8 *mu, REAL8 *sigma, LALInferenceVariableType type)
Function to add the mu and sigma values for the Gaussian prior onto the priorArgs.
REAL8 LALInferenceComputePriorMassNorm(const double MMin, const double MMax, const double MTotMax, const double McMin, const double McMax, const double massRatioMin, const double massRatioMax, const char *massRatioName)
Computes the numerical normalization of the mass prior applying all cuts in the mass plane implied b...
REAL8 LALInferenceInspiralPrior(LALInferenceRunState *runState, LALInferenceVariables *params, LALInferenceModel *model)
Return the logarithmic prior density of the variables specified, for the non-spinning/spinning inspir...
void LALInferenceAddMinMaxPrior(LALInferenceVariables *priorArgs, const char *name, REAL8 *min, REAL8 *max, LALInferenceVariableType type)
Function to add the minimum and maximum values for the uniform prior onto the priorArgs.
void LALInferenceRemoveMinMaxPrior(LALInferenceVariables *priorArgs, const char *name)
Function to remove the minimum and maximum values for the uniform prior onto the priorArgs.
void LALInferenceRemoveGaussianPrior(LALInferenceVariables *priorArgs, const char *name)
Function to remove the mu and sigma values for the Gaussian prior onto the priorArgs.
void LALInferenceCyclicReflectiveBound(LALInferenceVariables *parameter, LALInferenceVariables *priorArgs)
Apply cyclic and reflective boundaries to parameter to bring it back within the allowed prior ranges ...
void LALInferenceRotateInitialPhase(LALInferenceVariables *parameter)
Rotate initial phase if polarisation angle is cyclic around ranges.
void LALInferenceDrawFromPrior(LALInferenceVariables *output, LALInferenceVariables *priorArgs, gsl_rng *rdm)
Draw variables from the prior ranges.
void LALInferenceGetMinMaxPrior(LALInferenceVariables *priorArgs, const char *name, REAL8 *min, REAL8 *max)
Get the minimum and maximum values of the uniform prior from the priorArgs list, given a name.
void LALInferenceAddCorrelatedPrior(LALInferenceVariables *priorArgs, const char *name, gsl_matrix **cor, REAL8 *mu, REAL8 *sigma, UINT4 *idx)
Function to add a correlation matrix and parameter index for a prior defined as part of a multivariat...
void * XLALCalloc(size_t m, size_t n)
void XLALFree(void *p)
static const INT4 a
const char * XLALErrorString(int errnum)
#define XLAL_TRY(statement, errnum)
#define XLAL_IS_REAL8_FAIL_NAN(val)
XLAL_SUCCESS
XLAL_EFAULT
XLAL_ENAME
list mu
type
def Omega(v, m1, m2, S1, S2, Ln)
Preccesion frequency spins Eqs.
Definition: pneqns.py:17
Structure containing inference run state This includes pointers to the function types required to run...
Definition: LALInference.h:592
LALInferenceVariables * priorArgs
Common arguments needed by proposals, to be copied to thread->cycle.
Definition: LALInference.h:603
LALInferenceThreadState * threads
Definition: LALInference.h:612
Structure containing chain-specific variables.
Definition: LALInference.h:541
LALInferenceModel * model
Cycle of proposals to call.
Definition: LALInference.h:548
The LALInferenceVariables structure to contain a set of parameters Implemented as a linked list of LA...
Definition: LALInference.h:170