GRAVI Pipeline Reference Manual 1.8.0
Loading...
Searching...
No Matches
gravi_pca.h
Go to the documentation of this file.
1/*
2 * This file is part of the GRAVI Pipeline
3 * Copyright (C) 2002,2003 European Southern Observatory
4 *
5 * This program is free software; you can redistribute it and/or modify
6 * it under the terms of the GNU General Public License as published by
7 * the Free Software Foundation; either version 2 of the License, or
8 * (at your option) any later version.
9 *
10 * This program is distributed in the hope that it will be useful,
11 * but WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 * GNU General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public License
16 * along with this program; if not, write to the Free Software
17 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
18 */
19
20#ifndef GRAVI_PCA_H
21#define GRAVI_PCA_H
22
23#include <cpl.h>
24
27
28gravi_pca_result *gravi_pca_create_result(const cpl_matrix *data, const cpl_matrix *mask) CPL_ATTR_ALLOC;
30cpl_error_code gravi_pca_set_component_signs(gravi_pca_result *self, const cpl_array *wave, int num_components);
31
32gravi_pca_model *gravi_pca_load_model(const cpl_matrix *components) CPL_ATTR_ALLOC;
33gravi_pca_model *gravi_pca_create_model(const gravi_pca_result **results, int num_results, int num_components) CPL_ATTR_ALLOC;
34
35cpl_error_code gravi_pca_fit_components_polynomial(gravi_pca_model *self, const cpl_array *wave, int degree, int num_components);
36cpl_error_code gravi_pca_fit_components_bspline(gravi_pca_model *self, const cpl_array *wave, int degree, int num_components);
37cpl_error_code gravi_pca_refine_mean(gravi_pca_result *self, const cpl_matrix *residual);
38cpl_error_code gravi_pca_fit_model(gravi_pca_result *self, const gravi_pca_model *model, cpl_boolean fit_mean_subtracted, cpl_boolean verbose);
39
40cpl_vector *gravi_pca_get_component(const gravi_pca_result *self, int component) CPL_ATTR_ALLOC;
41cpl_vector *gravi_pca_get_component_median(const gravi_pca_model *self, int component) CPL_ATTR_ALLOC;
42cpl_vector *gravi_pca_get_component_fit(const gravi_pca_model *self, int component) CPL_ATTR_ALLOC;
43
44cpl_matrix *gravi_pca_get_data_fit(const gravi_pca_result *self) CPL_ATTR_ALLOC;
45cpl_matrix *gravi_pca_get_data_residual(const gravi_pca_result *self) CPL_ATTR_ALLOC;
46
49
50#endif // GRAVI_PCA_H
cpl_error_code gravi_pca_refine_mean(gravi_pca_result *self, const cpl_matrix *residual)
Override the mean component calculated from the PCA decomposition. In conjunction with the fit_mean_s...
Definition: gravi_pca.c:652
cpl_matrix * gravi_pca_get_data_residual(const gravi_pca_result *self) CPL_ATTR_ALLOC
Get residual (data - model).
Definition: gravi_pca.c:956
cpl_error_code gravi_pca_set_component_signs(gravi_pca_result *self, const cpl_array *wave, int num_components)
Impose sign convention on PCA components.
Definition: gravi_pca.c:393
cpl_matrix * gravi_pca_get_data_fit(const gravi_pca_result *self) CPL_ATTR_ALLOC
Get noise-free model.
Definition: gravi_pca.c:935
void gravi_pca_result_delete(gravi_pca_result *self)
Deallocate a gravi_pca_result object.
Definition: gravi_pca.c:99
cpl_error_code gravi_pca_fit_components_polynomial(gravi_pca_model *self, const cpl_array *wave, int degree, int num_components)
Fit polynomial model to each of a set of median-averaged PCA components.
Definition: gravi_pca.c:587
cpl_vector * gravi_pca_get_component(const gravi_pca_result *self, int component) CPL_ATTR_ALLOC
Get components from PCA decomposition.
Definition: gravi_pca.c:907
cpl_error_code gravi_pca_fit_components_bspline(gravi_pca_model *self, const cpl_array *wave, int degree, int num_components)
Fit B-spline model to each of a set of median-averaged PCA components.
Definition: gravi_pca.c:511
gravi_pca_model * gravi_pca_create_model(const gravi_pca_result **results, int num_results, int num_components) CPL_ATTR_ALLOC
Compute median values of PCA components over a set of decomposition results.
Definition: gravi_pca.c:423
cpl_vector * gravi_pca_get_component_fit(const gravi_pca_model *self, int component) CPL_ATTR_ALLOC
Get fit to components from PCA decomposition.
Definition: gravi_pca.c:716
gravi_pca_model * gravi_pca_load_model(const cpl_matrix *components) CPL_ATTR_ALLOC
Create PCA model from existing set of components.
Definition: gravi_pca.c:483
cpl_error_code gravi_pca_decomp_matrix_svd(gravi_pca_result *self)
Perform PCA decomposition by calculating singular value decomposition.
Definition: gravi_pca.c:283
void gravi_pca_model_delete(gravi_pca_model *self)
Deallocate a gravi_pca_model object.
Definition: gravi_pca.c:118
gravi_pca_result * gravi_pca_create_result(const cpl_matrix *data, const cpl_matrix *mask) CPL_ATTR_ALLOC
Construct a new gravi_pca_result object from a matrix of visphi data.
Definition: gravi_pca.c:219
cpl_vector * gravi_pca_get_component_median(const gravi_pca_model *self, int component) CPL_ATTR_ALLOC
Get median-averaged component from a set of PCA decompositions.
Definition: gravi_pca.c:694
cpl_error_code gravi_pca_fit_model(gravi_pca_result *self, const gravi_pca_model *model, cpl_boolean fit_mean_subtracted, cpl_boolean verbose)
Fit model formed from linear combination of PCA components to data.
Definition: gravi_pca.c:792
Type to hold average (median) components obtained from a set of PCA decompositions and/or best-fit mo...
Definition: gravi_pca.c:87
Type to hold results of a PCA decomposition.
Definition: gravi_pca.c:64