CR2RE Pipeline Reference Manual 1.7.0
hdrl_cat_casu_bkg_badpatch-test.c
1/*
2 * This file is part of the HDRL
3 * Copyright (C) 2017 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#ifdef HAVE_CONFIG_H
21#include "config.h"
22#endif
23
24#include <cpl_test.h>
25
26#include "../hdrl_cat_casu.h"
27
28
29#define NTEST 10
30#define COMP_TOLERANCE 1e-2
31#define CORNER_OFFSET 10
32#define CORNER_REL_TOL 1e-2
33#define X_OS_P1 60
34#define X_OS_P2 95
35#define STAR 7
36#define PATCH_SIZE 100 /* even number (patch is a square) */
37#define IMG_XSIZE 1100
38#define IMG_YSIZE 1500
39#define CELL_SIZE 32
40#define MIN_RAMP 10.
41#define MAX_RAMP 100.
42#define APER_FLUX_NUM "Aper_flux_3"
43
44/*----------------------------------------------------------------------------*/
53/*----------------------------------------------------------------------------*/
54cpl_error_code hdrl_casubkg_badpatch_compute(void)
55{
56 double xpos[] = { 100., 200., 300., 400., 500., 600., 700., 800., 900., 1000.};
57 double ypos[] = { 100., 200., 300., 400., 500., 600., 700., 800., 900., 1000.};
58 double norm[] = {1000., 100., 200., 500., 550., 600., 650., 700., 750., 800.};
59
60 /* Check inherited status */
61 hdrl_casu_result *res = cpl_calloc(1, sizeof(hdrl_casu_result));
62 hdrl_casu_result *res_p1 = cpl_calloc(1, sizeof(hdrl_casu_result));
63 hdrl_casu_result *res_p2 = cpl_calloc(1, sizeof(hdrl_casu_result));
64 cpl_test_eq(hdrl_casu_catalogue(NULL, NULL, NULL, 5, 1.5, 1, 3.5, 1, CELL_SIZE,
65 6, 1.0, 2.0, HDRL_SATURATION_INIT, res), CPL_ERROR_NULL_INPUT);
66 cpl_test_error(CPL_ERROR_NULL_INPUT);
67 cpl_test_null(res->catalogue);
68
69 cpl_image_delete(res->segmentation_map);
70 cpl_image_delete(res->background);
71 hdrl_casu_tfits_delete(res->catalogue);
72 cpl_free(res);
73
74 /* Generate a field with some stars and a confidence map */
75 cpl_image *bkg = cpl_image_new(IMG_XSIZE, IMG_YSIZE, CPL_TYPE_DOUBLE);
76 cpl_image *im = cpl_image_new(IMG_XSIZE, IMG_YSIZE, CPL_TYPE_DOUBLE);
77 cpl_image *cnf = cpl_image_new(IMG_XSIZE, IMG_YSIZE, CPL_TYPE_DOUBLE);
78 cpl_image *aux = cpl_image_new(IMG_XSIZE, IMG_YSIZE, CPL_TYPE_DOUBLE);
79
80 double sigma = 2.;
81 double norm2 = 2. * CPL_MATH_PI * sigma * sigma;
82 cpl_image_fill_noise_uniform(bkg, -10., 10.);
83
84 double sky = 500.;
85 cpl_image_add_scalar(bkg, sky);
86
87 cpl_image_fill_noise_uniform(cnf, 99.9, 100.1);
88 cpl_image *cnf_p1 = cpl_image_duplicate(cnf);
89 cpl_image *cnf_p2 = cpl_image_duplicate(cnf);
90
91 double tot[NTEST];
92 for (cpl_size i = 0; i < NTEST; i++) {
93 cpl_image_fill_gaussian(im, xpos[i], ypos[i], norm[i] * norm2, sigma, sigma);
94 tot[i] = cpl_image_get_flux(im);
95 cpl_image_add(bkg, im);
96 }
97 cpl_image_delete(im);
98
99 /* Add gradient */
100 double dlt_ramp = MAX_RAMP - MIN_RAMP;
101 for (cpl_size i = 1; i <= IMG_XSIZE; i++) {
102 for (cpl_size j = 1; j <= IMG_YSIZE; j++) {
103 cpl_image_set(aux, i, j, MIN_RAMP + dlt_ramp * (double)i / IMG_XSIZE);
104 }
105 }
106 cpl_image_add(bkg, aux);
107 cpl_image_delete(aux);
108
109 cpl_size y_os = 0;
110
111 /* Generation of confidence map "Patch Close" */
112 cpl_size x_os = X_OS_P1;
113 for (cpl_size i = 0; i <= PATCH_SIZE; i++) {
114 for (cpl_size j = 0; j <= PATCH_SIZE; j++) {
115 cpl_image_set(cnf_p1, xpos[STAR]-(PATCH_SIZE/2)+i +x_os,
116 ypos[STAR]-(PATCH_SIZE/2)+j +y_os, 0);
117 }
118 }
119
120 /* Generation of confidence map "Patch Away" */
121 x_os = X_OS_P2;
122 for (cpl_size i = 0; i <= PATCH_SIZE; i++) {
123 for (cpl_size j = 0; j <= PATCH_SIZE; j++) {
124 cpl_image_set(cnf_p2, xpos[STAR]-(PATCH_SIZE/2)+i +x_os,
125 ypos[STAR]-(PATCH_SIZE/2)+j +y_os, 0);
126 }
127 }
128
129 hdrl_casu_fits *inf = hdrl_casu_fits_wrap( bkg);
130 hdrl_casu_fits *inconf = hdrl_casu_fits_wrap( cnf);
131 hdrl_casu_fits *inconf_p1 = hdrl_casu_fits_wrap(cnf_p1);
132 hdrl_casu_fits *inconf_p2 = hdrl_casu_fits_wrap(cnf_p2);
133
134 /* Give it a WCS */
135 cpl_propertylist *pl = hdrl_casu_fits_get_ehu(inf);
136 cpl_propertylist_update_string(pl, "CTYPE1", "RA---TAN" );
137 cpl_propertylist_update_string(pl, "CTYPE2", "DEC--TAN" );
138 cpl_propertylist_update_double(pl, "CRVAL1", 30. );
139 cpl_propertylist_update_double(pl, "CRVAL2", 12. );
140 cpl_propertylist_update_double(pl, "CRPIX1", 512. );
141 cpl_propertylist_update_double(pl, "CRPIX2", 512. );
142 cpl_propertylist_update_double(pl, "CD1_1", -1. / 3600.);
143 cpl_propertylist_update_double(pl, "CD1_2", 0. );
144 cpl_propertylist_update_double(pl, "CD2_1", 0. );
145 cpl_propertylist_update_double(pl, "CD2_2", 1. / 3600.);
146 cpl_propertylist_update_int( pl, "NAXIS1", IMG_XSIZE );
147 cpl_propertylist_update_int( pl, "NAXIS2", IMG_YSIZE );
148 cpl_wcs *wcs = cpl_wcs_new_from_propertylist(pl);
149
150 /* Run casu catalogue */
151 res = cpl_calloc(1, sizeof(hdrl_casu_result));
152 cpl_test_eq(hdrl_casu_catalogue(inf, inconf, wcs, 5, 1.5, 0, 5., 1,
153 CELL_SIZE, HDRL_CATALOGUE_ALL, 3., 1., HDRL_SATURATION_INIT, res), CPL_ERROR_NONE);
154 cpl_test_nonnull(res->catalogue);
155 hdrl_casu_fits_delete(inconf);
156
157 /* Run casu catalogue for PC */
158 cpl_test_eq(hdrl_casu_catalogue(inf, inconf_p1, wcs, 5, 1.5, 0, 5., 1,
159 CELL_SIZE, HDRL_CATALOGUE_ALL, 3., 1., HDRL_SATURATION_INIT, res_p1), CPL_ERROR_NONE);
160 cpl_test_nonnull(res_p1->catalogue);
161 hdrl_casu_fits_delete(inconf_p1);
162
163 /* Run casu catalogue for PA*/
164 cpl_test_eq(hdrl_casu_catalogue(inf, inconf_p2, wcs, 5, 1.5, 0, 5., 1,
165 CELL_SIZE, HDRL_CATALOGUE_ALL, 3., 1., HDRL_SATURATION_INIT, res_p2), CPL_ERROR_NONE);
166 cpl_test_nonnull(res_p2->catalogue);
167 hdrl_casu_fits_delete(inconf_p2);
168 hdrl_casu_fits_delete(inf);
169 cpl_wcs_delete(wcs);
170
171 /*** Tests ***/
172
173 cpl_test_image_rel(res->background, res_p1->background, COMP_TOLERANCE);
174 cpl_test_image_rel(res->background, res_p2->background, COMP_TOLERANCE);
175
176 cpl_size corner_x = cpl_image_get_size_x(res->background) -CORNER_OFFSET;
177 cpl_size corner_y = cpl_image_get_size_y(res->background) -CORNER_OFFSET;
178
179 int rej;
180 cpl_test_abs(cpl_image_get(res ->background, corner_x, corner_y, &rej),
181 cpl_image_get(res_p1->background, corner_x, corner_y, &rej), CORNER_REL_TOL);
182 cpl_test_abs(cpl_image_get(res ->background, corner_x, corner_y, &rej),
183 cpl_image_get(res_p2->background, corner_x, corner_y, &rej), CORNER_REL_TOL);
184
185 cpl_test_image_rel(res->segmentation_map, res_p1->segmentation_map, COMP_TOLERANCE);
186 cpl_test_image_rel(res->segmentation_map, res_p2->segmentation_map, COMP_TOLERANCE);
187
188
189 /* Check the results. Start by checking the number of rows and columns. Sort the table by X */
190 cpl_table *tab = hdrl_casu_tfits_get_table(res ->catalogue);
191 cpl_table *tab_p1 = hdrl_casu_tfits_get_table(res_p1->catalogue);
192 cpl_table *tab_p2 = hdrl_casu_tfits_get_table(res_p2->catalogue);
193 cpl_test_nonnull(tab);
194 cpl_test_nonnull(tab_p1);
195 cpl_test_nonnull(tab_p2);
196
197 int nl;
198 for (cpl_size i = 0; i < NTEST; i++) {
199 cpl_test_rel( cpl_table_get_double(tab, APER_FLUX_NUM, i, &nl),
200 cpl_table_get_double(tab_p1, APER_FLUX_NUM, i, &nl), COMP_TOLERANCE );
201 cpl_test_rel( cpl_table_get_double(tab, APER_FLUX_NUM, i, &nl),
202 cpl_table_get_double(tab_p2, APER_FLUX_NUM, i, &nl), COMP_TOLERANCE );
203 }
204
205 cpl_test_eq(cpl_table_get_ncol(tab), 63);
206 cpl_test_eq(cpl_table_get_nrow(tab), NTEST);
207 pl = cpl_propertylist_new();
208 cpl_propertylist_append_bool(pl, "X_coordinate", 0);
209 cpl_table_sort(tab, pl);
210 cpl_propertylist_delete(pl);
211
212 /* Test the column content of the table */
213 for (cpl_size i = 0; i < NTEST; i++) {
214 cpl_test_abs(xpos[i], cpl_table_get_double(tab, "X_coordinate", i, &nl), 0.2);
215 cpl_test_abs(ypos[i], cpl_table_get_double(tab, "Y_coordinate", i, &nl), 0.2);
216 double diff = fabs( cpl_table_get_double(tab, "Aper_flux_5", i, &nl) -tot[i])
217 / cpl_table_get_double(tab, "Aper_flux_5_err", i, &nl);
218 cpl_test_lt(diff, 1.5);
219 cpl_test_eq( cpl_table_get_double(tab, "Classification", i, &nl), -1.);
220 }
221
222 /* Clean up */
223 cpl_image_delete( res->segmentation_map);
224 cpl_image_delete( res->background);
225 hdrl_casu_tfits_delete(res->catalogue);
226 cpl_free(res);
227
228 cpl_image_delete( res_p1->segmentation_map);
229 cpl_image_delete( res_p1->background);
230 hdrl_casu_tfits_delete(res_p1->catalogue);
231 cpl_free(res_p1);
232
233 cpl_image_delete( res_p2->segmentation_map);
234 cpl_image_delete( res_p2->background);
235 hdrl_casu_tfits_delete(res_p2->catalogue);
236 cpl_free(res_p2);
237
238
239 return cpl_error_get_code();
240}
241
242/*----------------------------------------------------------------------------*/
246/*----------------------------------------------------------------------------*/
247int main(void)
248{
249 cpl_test_init(PACKAGE_BUGREPORT, CPL_MSG_WARNING);
250
251 hdrl_casubkg_badpatch_compute();
252
253 cpl_test_error(CPL_ERROR_NONE);
254
255 return cpl_test_end(0);
256}
hdrl_casu_fits_wrap
hdrl_casu_fits * hdrl_casu_fits_wrap(cpl_image *im)
Wrap an image in a hdrl_casu_fits wrapper.
Definition hdrl_cat_casu.c:165
hdrl_casu_tfits_delete
void hdrl_casu_tfits_delete(hdrl_casu_tfits *p)
Free all the workspace associated with a hdrl_casu_fits object.
Definition hdrl_cat_casu.c:139
hdrl_casu_fits_delete
void hdrl_casu_fits_delete(hdrl_casu_fits *p)
Free all the workspace associated with a hdrl_casu_fits object.
Definition hdrl_cat_casu.c:267
hdrl_casu_tfits_get_table
cpl_table * hdrl_casu_tfits_get_table(hdrl_casu_tfits *p)
Get the CPL table from the hdrl_casu_tfits object.
Definition hdrl_cat_casu.c:121
hdrl_casu_fits_get_ehu
cpl_propertylist * hdrl_casu_fits_get_ehu(hdrl_casu_fits *p)
Get the propertylist for the extension header for a given hdrl_casu_fits image.
Definition hdrl_cat_casu.c:226
hdrl_casu_catalogue
cpl_error_code hdrl_casu_catalogue(hdrl_casu_fits *infile, hdrl_casu_fits *conf, const cpl_wcs *wcs, cpl_size ipix, double threshold, cpl_size icrowd, double rcore, cpl_size bkg_subtr, cpl_size nbsize, hdrl_catalogue_options cattype, double filtfwhm, double gainloc, double saturation, hdrl_casu_result *res)
Generate object catalogues from input images.
Definition hdrl_cat_casu.c:321
Generated on Tue Oct 13 2026 12:27:26 for CR2RE Pipeline Reference Manual by doxygen 1.13.2