|
GRAVI Pipeline Reference Manual 1.8.0
|
#include <cpl.h>#include <math.h>#include <time.h>#include <gsl/gsl_matrix.h>#include <gsl/gsl_vector.h>#include <gsl/gsl_multimin.h>#include <gsl/gsl_statistics.h>#include "gravi_demodulate.h"#include "gravi_pfits.h"#include "gravi_utils.h"Go to the source code of this file.
Data Structures | |
| struct | _demodulation_model_params_ |
Macros | |
| #define | FT 0 |
| #define | SC 1 |
| #define | STEPS_PER_SECOND 500 |
| #define | MAX_SECONDS_PER_CHUNK 100 |
| #define | PI 3.14159265359 |
| #define | TWOPI 6.28318530718 |
Typedefs | |
| typedef struct _demodulation_model_params_ | model_params |
Functions | |
| static cpl_size | diode_column_index (int tele, int diode, int side) |
| Return column index in VOLT table for specific diode. | |
| static double | model_x (const gsl_vector *X, int step) |
| static double | model_y (const gsl_vector *X, int step) |
| static double | modulation_model_chi2 (const gsl_vector *X, void *params) |
| Calculate chi^2 statistic for modulation model. | |
| static cpl_error_code | fit_model_modulation (const gsl_vector *vx, const gsl_vector *vy, gsl_vector *Xsolve) |
| Fit the modulation model to voltage data for a single diode. | |
| cpl_error_code | gravi_metrology_demodulate (gravi_data *met_data, cpl_boolean zero_subtracted) |
| Demodulate the metrology. | |
Variables | |
| const cpl_size | ntel = 4 |
| const cpl_size | ndiode = 4 |
| const cpl_size | nside = 2 |
| const double | diode_zeros [] |
| #define FT 0 |
Definition at line 54 of file gravi_demodulate.c.
| #define MAX_SECONDS_PER_CHUNK 100 |
Definition at line 58 of file gravi_demodulate.c.
| #define PI 3.14159265359 |
Definition at line 60 of file gravi_demodulate.c.
| #define SC 1 |
Definition at line 55 of file gravi_demodulate.c.
| #define STEPS_PER_SECOND 500 |
Definition at line 57 of file gravi_demodulate.c.
| #define TWOPI 6.28318530718 |
Definition at line 61 of file gravi_demodulate.c.
| typedef struct _demodulation_model_params_ model_params |
|
static |
Return column index in VOLT table for specific diode.
| tele | Telescope in {0..3} |
| diode | Diode in {0..3} |
| side | 0 (FT) or 1 (SC) |
Definition at line 76 of file gravi_demodulate.c.
References FT.
Referenced by gravi_metrology_demodulate().
|
static |
Fit the modulation model to voltage data for a single diode.
| vx | x voltages | |
| vy | y voltages | |
| zx | x zero offset | |
| zy | y zero offset | |
| [out] | Xsolve | returned best-fit parameters, which must already be allocated. |
Definition at line 404 of file gravi_demodulate.c.
References cpl_msg_debug(), FREE, modulation_model_chi2(), PI, SIGN, TWOPI, and _demodulation_model_params_::volts_x.
Referenced by gravi_metrology_demodulate().
| cpl_error_code gravi_metrology_demodulate | ( | gravi_data * | met_data, |
| cpl_boolean | zero_subtracted | ||
| ) |
Demodulate the metrology.
| [in,out] | met_data | input table of metrology data |
| [in] | zero_subtracted | flag indicating metrology has been zero-subtracted using DARK |
Definition at line 161 of file gravi_demodulate.c.
References cpl_msg_debug(), cpl_msg_info(), CPLCHECK_INT, diode_column_index(), diode_zeros, fit_model_modulation(), FREE, FT, gravi_data_get_header, gravi_data_get_table(), GRAVI_METROLOGY_EXT, gravi_msg_function_exit, gravi_msg_function_start, MAX_SECONDS_PER_CHUNK, ndiode, nside, ntel, SC, and STEPS_PER_SECOND.
Referenced by gravity_vis().
|
static |
Definition at line 117 of file gravi_demodulate.c.
References STEPS_PER_SECOND, and TWOPI.
Referenced by modulation_model_chi2().
|
static |
Definition at line 128 of file gravi_demodulate.c.
References STEPS_PER_SECOND, and TWOPI.
Referenced by modulation_model_chi2().
|
static |
Calculate chi^2 statistic for modulation model.
| X | vector of model parameters [a, b, pha1, pha2] |
| params | pointer to model_params struct |
Definition at line 358 of file gravi_demodulate.c.
References model_x(), model_y(), STEPS_PER_SECOND, _demodulation_model_params_::volts_x, and _demodulation_model_params_::volts_y.
Referenced by fit_model_modulation().
| const double diode_zeros[] |
Diode zero offsets from Stefan Gillesen dated TODO To be updated with new values dated TODO Ordering is such that it matches the VOLT column in the metrology data Used iff no DARK was provided
Definition at line 94 of file gravi_demodulate.c.
Referenced by gravi_metrology_demodulate().
| const cpl_size ndiode = 4 |
Definition at line 85 of file gravi_demodulate.c.
Referenced by gravi_flux_create_met_sc(), gravi_metrology_demodulate(), gravi_metrology_drs(), gravi_metrology_tac(), gravi_metrology_telfc(), and gravi_vis_create_met_sc().
| const cpl_size nside = 2 |
Definition at line 86 of file gravi_demodulate.c.
Referenced by gravi_metrology_demodulate().
| const cpl_size ntel = 4 |
Definition at line 84 of file gravi_demodulate.c.
Referenced by gravi_acqcam_field(), gravi_compute_disp(), gravi_compute_opdc_state(), gravi_compute_outliers(), gravi_compute_p2vmred(), gravi_compute_piezotf(), gravi_compute_qc_injection(), gravi_compute_signals(), gravi_compute_vis(), gravi_compute_vis_qc(), gravi_create_outlier_flag_sc(), gravi_create_p2vm_table(), gravi_disp_cleanup(), gravi_fit_dispersion(), gravi_fit_fddl_lin(), gravi_flux_average_bootstrap(), gravi_flux_create_acq_sc(), gravi_flux_create_average(), gravi_flux_create_fddllin_sc(), gravi_flux_create_fddlpos_sc(), gravi_flux_create_met_sc(), gravi_flux_create_totalflux_sc(), gravi_lkdt_get_sequence(), gravi_metrology_acq(), gravi_metrology_compute_p2vm(), gravi_metrology_create(), gravi_metrology_demodulate(), gravi_metrology_drs(), gravi_metrology_reform(), gravi_metrology_tac(), gravi_metrology_telfc(), gravi_opds_compute_guess(), gravi_p2vm_normalisation(), gravi_phase_correct_closures(), gravi_phase_fit_opdmet(), gravi_reduce_acqcam(), gravi_t3_average_bootstrap(), gravi_vis_compute_column_mean(), gravi_vis_create_acq_sc(), gravi_vis_create_f1f2_ft(), gravi_vis_create_f1f2_sc(), gravi_vis_create_met_ft(), gravi_vis_create_met_sc(), gravi_vis_create_opddisp_sc(), gravi_vis_create_opdguess_sc(), gravi_vis_create_pfactor_ft(), gravi_vis_create_pfactor_sc(), and gravi_vis_erase_obs().