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
LALInferenceHDF5Test.c
Go to the documentation of this file.
1#include <lal/XLALError.h>
2#include <lal/LALInferenceHDF5.h>
3#include <gsl/gsl_test.h>
4
5int main(int argc, char **argv)
6{
7 /* Not used */
8 (void)argc;
9 (void)argv;
11
12 /* Populate example variables array. */
13 UINT4 N = 64;
14 LALInferenceVariables **vars_array = XLALCalloc(
15 N, sizeof(LALInferenceVariables *));
16 for (UINT4 i = 0; i < N; i ++)
17 {
19 vars_array[i] = vars;
28 }
29
30 /* Open HDF5 file for writing and create group structure. */
31 LALH5File *file = XLALH5FileOpen("test.hdf5", "w");
33 file, "lalinference", "lalinference_mcmc");
34
35 /* Write variables array to dataset. */
38
39 /* Free variables array. */
40 for (UINT4 i = 0; i < N; i ++)
41 {
42 LALInferenceClearVariables(vars_array[i]);
43 XLALFree(vars_array[i]);
44 }
45 XLALFree(vars_array);
46
47 /* Close file. */
50
51 /* Open file for reading and find dataset. */
52 file = XLALH5FileOpen("test.hdf5", "r");
54 file, "lalinference/lalinference_mcmc/posterior_samples");
55
56 /* Read dataset back to variables array. */
57 N = 0;
58 vars_array = NULL;
59 LALInferenceH5DatasetToVariablesArray(dataset, &vars_array, &N);
60
61 gsl_test_int(N, 64, "number of rows read back");
62 for (UINT4 i = 0; i < N; i ++)
63 {
64 LALInferenceVariables *vars = vars_array[i];
65 gsl_test_int(LALInferenceGetVariableDimension(vars), 8,
66 "number of columns read back");
67 gsl_test_abs(LALInferenceGetREAL8Variable(vars, "abc"), i, 0,
68 "value of column abc");
69 gsl_test_abs(LALInferenceGetREAL8Variable(vars, "def"), i, 0,
70 "value of column def");
71 gsl_test_abs(LALInferenceGetREAL8Variable(vars, "ghi"), 5, 0,
72 "value of column ghi");
73 gsl_test_abs(LALInferenceGetREAL8Variable(vars, "ijk"), i, 0,
74 "value of column ijk");
75 gsl_test_int(LALInferenceGetINT4Variable (vars, "lmn"), i,
76 "value of column lmn");
77 gsl_test_int(LALInferenceGetINT4Variable (vars, "opq"), i,
78 "value of column opq");
79 gsl_test_int(LALInferenceGetINT4Variable (vars, "rst"), 5,
80 "value of column rst");
81 gsl_test_int(LALInferenceGetINT4Variable (vars, "uvw"), i,
82 "value of column uvw");
83
84 gsl_test_int(LALInferenceGetVariableVaryType(vars, "abc"),
85 LALINFERENCE_PARAM_LINEAR, "vary type of column abc");
86 gsl_test_int(LALInferenceGetVariableVaryType(vars, "def"),
87 LALINFERENCE_PARAM_CIRCULAR, "vary type of column def");
88 gsl_test_int(LALInferenceGetVariableVaryType(vars, "ghi"),
89 LALINFERENCE_PARAM_FIXED, "vary type of column ghi");
90 gsl_test_int(LALInferenceGetVariableVaryType(vars, "ijk"),
91 LALINFERENCE_PARAM_OUTPUT, "vary type of column ijk");
92 gsl_test_int(LALInferenceGetVariableVaryType(vars, "lmn"),
93 LALINFERENCE_PARAM_LINEAR, "vary type of column lmn");
94 gsl_test_int(LALInferenceGetVariableVaryType(vars, "opq"),
95 LALINFERENCE_PARAM_CIRCULAR, "vary type of column opq");
96 gsl_test_int(LALInferenceGetVariableVaryType(vars, "rst"),
97 LALINFERENCE_PARAM_FIXED, "vary type of column rst");
98 gsl_test_int(LALInferenceGetVariableVaryType(vars, "uvw"),
99 LALINFERENCE_PARAM_OUTPUT, "vary type of column uvw");
100 }
101
102 /* Free variables array. */
103 for (UINT4 i = 0; i < N; i ++)
104 {
105 LALInferenceClearVariables(vars_array[i]);
106 XLALFree(vars_array[i]);
107 }
108 XLALFree(vars_array);
109
110 /* Close dataset. */
111 XLALH5DatasetFree(dataset);
112
113 /* Close file. */
115
116 /* Check for memory leaks. */
118
119 /* Done! */
120 return gsl_test_summary();
121}
struct tagLALH5File LALH5File
struct tagLALH5Dataset LALH5Dataset
int LALInferenceH5DatasetToVariablesArray(LALH5Dataset *dataset, LALInferenceVariables ***varsArray, UINT4 *N)
const char LALInferenceHDF5PosteriorSamplesDatasetName[]
int LALInferenceH5VariablesArrayToDataset(LALH5File *h5file, LALInferenceVariables *const *const varsArray, UINT4 N, const char *TableName)
LALH5File * LALInferenceH5CreateGroupStructure(LALH5File *h5file, const char *codename, const char *runID)
Create a HDF5 heirarchy in the given LALH5File reference /codename/runID/ Returns a LALH5File pointer...
int main(int argc, char **argv)
void LALCheckMemoryLeaks(void)
void XLALH5FileClose(LALH5File UNUSED *file)
void XLALH5DatasetFree(LALH5Dataset UNUSED *dset)
LALH5File * XLALH5FileOpen(const char UNUSED *path, const char UNUSED *mode)
LALH5Dataset * XLALH5DatasetRead(LALH5File UNUSED *file, const char UNUSED *name)
uint32_t UINT4
REAL8 LALInferenceGetREAL8Variable(LALInferenceVariables *vars, const char *name)
void LALInferenceAddREAL8Variable(LALInferenceVariables *vars, const char *name, REAL8 value, LALInferenceParamVaryType vary)
INT4 LALInferenceGetINT4Variable(LALInferenceVariables *vars, const char *name)
void LALInferenceAddINT4Variable(LALInferenceVariables *vars, const char *name, INT4 value, LALInferenceParamVaryType vary)
INT4 LALInferenceGetVariableDimension(LALInferenceVariables *vars)
Get number of dimensions in variable vars.
Definition: LALInference.c:250
void LALInferenceClearVariables(LALInferenceVariables *vars)
Delete the variables in this structure.
Definition: LALInference.c:532
LALInferenceParamVaryType LALInferenceGetVariableVaryType(LALInferenceVariables *vars, const char *name)
Get the LALInferenceParamVaryType of the parameter named name in vars see the declaration of LALInfer...
Definition: LALInference.c:226
@ LALINFERENCE_PARAM_OUTPUT
A parameter that never changes, functions should respect this.
Definition: LALInference.h:131
@ LALINFERENCE_PARAM_FIXED
A parameter that is cyclic, such as an angle between 0 and 2pi.
Definition: LALInference.h:130
@ LALINFERENCE_PARAM_CIRCULAR
A parameter that simply has a maximum and a minimum.
Definition: LALInference.h:129
@ LALINFERENCE_PARAM_LINEAR
Definition: LALInference.h:128
void * XLALCalloc(size_t m, size_t n)
void XLALFree(void *p)
XLALErrorHandlerType * XLALSetErrorHandler(XLALErrorHandlerType *newHandler)
void XLALExitErrorHandler(const char *func, const char *file, int line, int errnum)
The LALInferenceVariables structure to contain a set of parameters Implemented as a linked list of LA...
Definition: LALInference.h:170