CR2RE Pipeline Reference Manual 1.6.7
hdrl_resample.h
1/*
2 * This file is part of the HDRL
3 * Copyright (C) 2020 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., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
18 */
19
20#ifndef HDRL_RESAMPLE_H
21#define HDRL_RESAMPLE_H
22
23/*-----------------------------------------------------------------------------
24 Includes
25 -----------------------------------------------------------------------------*/
26
27#include <cpl.h>
28#include <hdrl_parameter.h>
29CPL_BEGIN_DECLS
30
31/*-----------------------------------------------------------------------------
32 Defines
33 -----------------------------------------------------------------------------*/
34
35 /* Table columns of the resampling table. */
36
37#define HDRL_RESAMPLE_TABLE_RA "ra" /* Right ascension column of the table*/
38#define HDRL_RESAMPLE_TABLE_DEC "dec" /* Declination column of the table*/
39#define HDRL_RESAMPLE_TABLE_LAMBDA "lambda" /* wavelength column of the table */
40#define HDRL_RESAMPLE_TABLE_DATA "data" /* data column of the table */
41#define HDRL_RESAMPLE_TABLE_BPM "bpm" /* bpm column of the table */
42#define HDRL_RESAMPLE_TABLE_ERRORS "errors" /* error column of the table */
43
44/* Right ascension column type of the table*/
45#define HDRL_RESAMPLE_TABLE_RA_TYPE CPL_TYPE_DOUBLE
46/* Declination column type of the table*/
47#define HDRL_RESAMPLE_TABLE_DEC_TYPE CPL_TYPE_DOUBLE
48/* wavelength column type of the table */
49#define HDRL_RESAMPLE_TABLE_LAMBDA_TYPE CPL_TYPE_DOUBLE
50/* data column type of the table */
51#define HDRL_RESAMPLE_TABLE_DATA_TYPE CPL_TYPE_DOUBLE
52/* bpm column type of the table */
53#define HDRL_RESAMPLE_TABLE_BPM_TYPE CPL_TYPE_INT
54/* error type column of the table */
55#define HDRL_RESAMPLE_TABLE_ERRORS_TYPE CPL_TYPE_DOUBLE
56
57typedef enum {
58 HDRL_RESAMPLE_METHOD_NEAREST = 0,
59 HDRL_RESAMPLE_METHOD_RENKA = 1,
60 HDRL_RESAMPLE_METHOD_LINEAR,
61 HDRL_RESAMPLE_METHOD_QUADRATIC,
62 HDRL_RESAMPLE_METHOD_DRIZZLE,
63 HDRL_RESAMPLE_METHOD_LANCZOS,
64 HDRL_RESAMPLE_METHOD_NONE /* Can be used for range checking */
65} hdrl_resample_method ;
66
67typedef enum {
68 HDRL_RESAMPLE_OUTGRID_2D,
69 HDRL_RESAMPLE_OUTGRID_3D
70} hdrl_resample_outgrid ;
71
72typedef struct {
73 cpl_propertylist * header; /* cpl propertylist containing the wcs informations */
74 hdrl_imagelist * himlist; /* hdrl imagelist containing the data/errors/bpm */
75} hdrl_resample_result;
76
77/*----------------------------------------------------------------------------
78 RESAMPLING Computation
79 ----------------------------------------------------------------------------*/
80
81cpl_table*
82hdrl_resample_imagelist_to_table(const hdrl_imagelist * himlist,
83 const cpl_wcs* wcs);
84cpl_table*
85hdrl_resample_image_to_table(const hdrl_image * hima, const cpl_wcs* wcs);
86
87hdrl_resample_result *
88hdrl_resample_compute(const cpl_table *restable,
89 hdrl_parameter *method,
90 hdrl_parameter *outputgrid,
91 const cpl_wcs* wcs);
92
93void
94hdrl_resample_result_delete(hdrl_resample_result *resdata);
95
96/*----------------------------------------------------------------------------
97 RESAMPLING Parameters for output definition
98 ----------------------------------------------------------------------------*/
99
100hdrl_parameter*
101hdrl_resample_parameter_create_outgrid2D(const double delta_ra,
102 const double delta_dec);
103
104hdrl_parameter*
105hdrl_resample_parameter_create_outgrid3D(const double delta_ra,
106 const double delta_dec,
107 const double delta_lambda);
108
109hdrl_parameter*
111 const double delta_dec,
112 const double ra_min,
113 const double ra_max,
114 const double dec_min,
115 const double dec_max,
116 const double fieldmargin);
117
118hdrl_parameter*
120 const double delta_dec,
121 const double delta_lambda,
122 const double ra_min,
123 const double ra_max,
124 const double dec_min,
125 const double dec_max,
126 const double lambda_min,
127 const double lambda_max,
128 const double fieldmargin);
129
130/*----------------------------------------------------------------------------
131 RESAMPLING Parameters for method definition
132 ----------------------------------------------------------------------------*/
133
134hdrl_parameter*
135hdrl_resample_parameter_create_renka(const int loop_distance,
136 cpl_boolean use_errorweights,
137 const double critical_radius);
138hdrl_parameter*
139hdrl_resample_parameter_create_drizzle(const int loop_distance,
140 cpl_boolean use_errorweights,
141 const double pix_frac_x,
142 const double pix_frac_y,
143 const double pix_frac_lambda);
144
145hdrl_parameter*
147
148hdrl_parameter*
149hdrl_resample_parameter_create_linear(const int loop_distance,
150 cpl_boolean use_errorweights);
151
152hdrl_parameter*
153hdrl_resample_parameter_create_quadratic(const int loop_distance,
154 cpl_boolean use_errorweights);
155
156hdrl_parameter*
157hdrl_resample_parameter_create_lanczos(const int loop_distance,
158 cpl_boolean use_errorweights,
159 const int kernel_size);
160
161/*----------------------------------------------------------------------------
162 Checks and verifications for the RESAMPLING Parameters
163 ----------------------------------------------------------------------------*/
164
165cpl_error_code
166hdrl_resample_parameter_outgrid_verify(const hdrl_parameter * hp);
167cpl_error_code
168hdrl_resample_parameter_method_verify(const hdrl_parameter * hp);
169int
170hdrl_resample_parameter_outgrid_check(const hdrl_parameter * self);
171int
172hdrl_resample_parameter_method_check(const hdrl_parameter * hp);
173
174/*-----------------------------------------------------------------------------
175 Private declarations - must not be used outside of hdrl
176 -----------------------------------------------------------------------------*/
177
178#ifdef HDRL_USE_PRIVATE
179
180/* Helper Functions */
181
182cpl_error_code
183hdrl_wcs_xy_to_radec(const cpl_wcs *wcs, double x, double y,
184 double *ra, double *dec);
185
186cpl_error_code
187hdrl_wcs_to_propertylist(const cpl_wcs * wcs, cpl_propertylist * header,
188 cpl_boolean only2d);
189
190#endif
191CPL_END_DECLS
192
193#endif /* RECIPES_HDRL_RESAMPLE_H_ */
hdrl_parameter * hdrl_resample_parameter_create_nearest(void)
Creates a resample nearest neighbor hdrl parameter object. The algorithm does not use any weighting f...
cpl_error_code hdrl_resample_parameter_method_verify(const hdrl_parameter *hp)
verify parameters have proper values
int hdrl_resample_parameter_method_check(const hdrl_parameter *hp)
check method is of proper type
hdrl_parameter * hdrl_resample_parameter_create_renka(const int loop_distance, cpl_boolean use_errorweights, const double critical_radius)
Creates a resample renka hdrl parameter object. The algorithm uses a modified Shepard-like distance w...
hdrl_parameter * hdrl_resample_parameter_create_lanczos(const int loop_distance, cpl_boolean use_errorweights, const int kernel_size)
Creates a resample Lanczos hdrl parameter object. The algorithm uses a restricted SINC distance weigh...
hdrl_parameter * hdrl_resample_parameter_create_drizzle(const int loop_distance, cpl_boolean use_errorweights, const double pix_frac_x, const double pix_frac_y, const double pix_frac_lambda)
Creates a resample drizzle hdrl parameter object. The algorithm uses a drizzle-like distance weightin...
hdrl_parameter * hdrl_resample_parameter_create_quadratic(const int loop_distance, cpl_boolean use_errorweights)
Creates a resample quadratic hdrl parameter object. The algorithm uses a quadratic inverse distance w...
hdrl_parameter * hdrl_resample_parameter_create_linear(const int loop_distance, cpl_boolean use_errorweights)
Creates a resample linear hdrl parameter object. The algorithm uses a linear inverse distance weighti...
hdrl_parameter * hdrl_resample_parameter_create_outgrid3D(const double delta_ra, const double delta_dec, const double delta_lambda)
Creates a resample_outgrid hdrl parameter object for a 3 dimensional interpolation,...
hdrl_parameter * hdrl_resample_parameter_create_outgrid3D_userdef(const double delta_ra, const double delta_dec, const double delta_lambda, const double ra_min, const double ra_max, const double dec_min, const double dec_max, const double lambda_min, const double lambda_max, const double fieldmargin)
Creates a resample_outgrid hdrl parameter object for a 3 dimensional interpolation,...
int hdrl_resample_parameter_outgrid_check(const hdrl_parameter *self)
check method is of proper type
hdrl_parameter * hdrl_resample_parameter_create_outgrid2D(const double delta_ra, const double delta_dec)
Creates a resample_outgrid hdrl parameter object for a 2 dimensional interpolation,...
hdrl_parameter * hdrl_resample_parameter_create_outgrid2D_userdef(const double delta_ra, const double delta_dec, const double ra_min, const double ra_max, const double dec_min, const double dec_max, const double fieldmargin)
Creates a resample_outgrid hdrl parameter object for a 2 dimensional interpolation,...
cpl_error_code hdrl_resample_parameter_outgrid_verify(const hdrl_parameter *hp)
verify parameters have proper values
hdrl_resample_result * hdrl_resample_compute(const cpl_table *restable, hdrl_parameter *method, hdrl_parameter *outputgrid, const cpl_wcs *wcs)
High level resampling function.
void hdrl_resample_result_delete(hdrl_resample_result *resdata)
Deallocates the memory associated to a hdrl_resample_result object.
cpl_table * hdrl_resample_image_to_table(const hdrl_image *hima, const cpl_wcs *wcs)
Convert a hdrl image into a cpl table that can be given as input to hdrl_resample_compute()
cpl_table * hdrl_resample_imagelist_to_table(const hdrl_imagelist *himlist, const cpl_wcs *wcs)
Convert a hdrl imagelist into a cpl table that can be given as input to hdrl_resample_compute().