1 /*
2 * This file is part of Siril, an astronomy image processor.
3 * Copyright (C) 2005-2011 Francois Meyer (dulle at free.fr)
4 * Copyright (C) 2012-2021 team free-astro (see more in AUTHORS file)
5 * Reference site is https://free-astro.org/index.php/Siril
6 *
7 * Siril is free software: you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation, either version 3 of the License, or
10 * (at your option) any later version.
11 *
12 * Siril is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with Siril. If not, see <http://www.gnu.org/licenses/>.
19 */
20 #include <float.h>
21 #include <string.h>
22 #include <gsl/gsl_statistics.h>
23
24 #include "core/siril.h"
25 #include "core/proto.h"
26 #include "core/arithm.h"
27 #include "core/undo.h"
28 #include "core/processing.h"
29 #include "core/OS_utils.h"
30 #include "algos/statistics.h"
31 #include "algos/sorting.h"
32 #include "gui/image_display.h"
33 #include "gui/progress_and_log.h"
34 #include "gui/utils.h"
35 #include "gui/dialogs.h"
36 #include "io/single_image.h"
37 #include "io/image_format_fits.h"
38 #include "io/sequence.h"
39 #include "opencv/opencv.h"
40
41 #include "banding.h"
42
43 /*****************************************************************************
44 * B A N D I N G R E D U C T I O N M A N A G E M E N T *
45 ****************************************************************************/
46
banding_image_hook(struct generic_seq_args * args,int o,int i,fits * fit,rectangle * _)47 int banding_image_hook(struct generic_seq_args *args, int o, int i, fits *fit, rectangle *_) {
48 struct banding_data *banding_args = (struct banding_data *)args->user;
49 return BandingEngine(fit, banding_args->sigma, banding_args->amount,
50 banding_args->protect_highlights, banding_args->applyRotation);
51 }
52
apply_banding_to_sequence(struct banding_data * banding_args)53 void apply_banding_to_sequence(struct banding_data *banding_args) {
54 struct generic_seq_args *args = create_default_seqargs(&com.seq);
55 args->filtering_criterion = seq_filter_included;
56 args->nb_filtered_images = com.seq.selnum;
57 args->prepare_hook = seq_prepare_hook;
58 args->finalize_hook = seq_finalize_hook;
59 args->image_hook = banding_image_hook;
60 args->stop_on_error = FALSE;
61 args->description = _("Banding Reduction");
62 args->has_output = TRUE;
63 args->output_type = get_data_type(args->seq->bitpix);
64 args->new_seq_prefix = banding_args->seqEntry;
65 args->load_new_sequence = TRUE;
66 args->user = banding_args;
67
68 banding_args->fit = NULL; // not used here
69
70 start_in_new_thread(generic_sequence_worker, args);
71 }
72
73 // idle function executed at the end of the BandingEngine processing
end_BandingEngine(gpointer p)74 gboolean end_BandingEngine(gpointer p) {
75 struct banding_data *args = (struct banding_data *) p;
76 stop_processing_thread();// can it be done here in case there is no thread?
77 adjust_cutoff_from_updated_gfit();
78 redraw(com.cvport, REMAP_ALL);
79 redraw_previews();
80 set_cursor_waiting(FALSE);
81
82 free(args);
83 return FALSE;
84 }
85
fmul_layer_ushort(fits * a,int layer,float coeff)86 static int fmul_layer_ushort(fits *a, int layer, float coeff) {
87 WORD *buf;
88 size_t i, n = a->naxes[0] * a->naxes[1];
89
90 if (coeff < 0.0)
91 return 1;
92 buf = a->pdata[layer];
93 for (i = 0; i < n; ++i) {
94 buf[i] = round_to_WORD(buf[i] * coeff);
95 }
96 invalidate_stats_from_fit(a);
97 return 0;
98 }
99
fmul_layer_float(fits * a,int layer,float coeff)100 static int fmul_layer_float(fits *a, int layer, float coeff) {
101 float *buf;
102 size_t i, n = a->naxes[0] * a->naxes[1];
103
104 if (coeff < 0.0)
105 return 1;
106 buf = a->fpdata[layer];
107 for (i = 0; i < n; ++i) {
108 buf[i] = buf[i] * coeff;
109 }
110 invalidate_stats_from_fit(a);
111 return 0;
112 }
113
114 /*** Reduces Banding in Canon DSLR images.
115 * This code come from CanonBandingReduction.js v0.9.1, a script of
116 * PixInsight, originally written by Georg Viehoever and
117 * distributed under the terms of the GNU General Public License ******/
BandingEngineThreaded(gpointer p)118 gpointer BandingEngineThreaded(gpointer p) {
119 struct banding_data *args = (struct banding_data *) p;
120 struct timeval t_start, t_end;
121
122 siril_log_color_message(_("Banding Reducing: processing...\n"), "green");
123 gettimeofday(&t_start, NULL);
124
125 int retval = BandingEngine(args->fit, args->sigma, args->amount, args->protect_highlights, args->applyRotation);
126
127 gettimeofday(&t_end, NULL);
128 show_time(t_start, t_end);
129 siril_add_idle(end_BandingEngine, args);
130
131 return GINT_TO_POINTER(retval);
132 }
133
BandingEngine_ushort(fits * fit,double sigma,double amount,gboolean protect_highlights,gboolean applyRotation)134 static int BandingEngine_ushort(fits *fit, double sigma, double amount, gboolean protect_highlights, gboolean applyRotation) {
135 int chan, row, i, ret = 0;
136 WORD *line, *fixline;
137 double minimum = DBL_MAX, globalsigma = 0.0;
138 fits *fiximage = NULL;
139 double invsigma = 1.0 / sigma;
140
141 if (applyRotation) {
142 point center = {gfit.rx / 2.0, gfit.ry / 2.0};
143 cvRotateImage(fit, center, 90.0, -1, OPENCV_AREA);
144 }
145
146 if (new_fit_image(&fiximage, fit->rx, fit->ry, fit->naxes[2], DATA_USHORT))
147 return 1;
148
149 for (chan = 0; chan < fit->naxes[2]; chan++) {
150 imstats *stat = statistics(NULL, -1, fit, chan, NULL, STATS_BASIC | STATS_MAD, TRUE);
151 if (!stat) {
152 siril_log_message(_("Error: statistics computation failed.\n"));
153 return 1;
154 }
155 double background = stat->median;
156 double *rowvalue = calloc(fit->ry, sizeof(double));
157 if (rowvalue == NULL) {
158 PRINT_ALLOC_ERR;
159 free_stats(stat);
160 return 1;
161 }
162 if (protect_highlights) {
163 globalsigma = stat->mad * MAD_NORM;
164 }
165 free_stats(stat);
166 for (row = 0; row < fit->ry; row++) {
167 line = fit->pdata[chan] + row * fit->rx;
168 WORD *cpyline = calloc(fit->rx, sizeof(WORD));
169 if (cpyline == NULL) {
170 PRINT_ALLOC_ERR;
171 free(rowvalue);
172 return 1;
173 }
174 memcpy(cpyline, line, fit->rx * sizeof(WORD));
175 int n = fit->rx;
176 double median;
177 if (protect_highlights) {
178 quicksort_s(cpyline, n);
179 WORD reject = round_to_WORD(
180 background + invsigma * globalsigma);
181 for (i = fit->rx - 1; i >= 0; i--) {
182 if (cpyline[i] < reject)
183 break;
184 n--;
185 }
186 median = gsl_stats_ushort_median_from_sorted_data(cpyline, 1, n);
187 } else {
188 median = round_to_WORD(quickmedian(cpyline, n));
189 }
190
191 rowvalue[row] = background - median;
192 minimum = min(minimum, rowvalue[row]);
193 free(cpyline);
194 }
195 for (row = 0; row < fit->ry; row++) {
196 fixline = fiximage->pdata[chan] + row * fiximage->rx;
197 for (i = 0; i < fit->rx; i++)
198 fixline[i] = round_to_WORD(rowvalue[row] - minimum);
199 }
200 free(rowvalue);
201 }
202 for (chan = 0; chan < fit->naxes[2]; chan++)
203 fmul_layer_ushort(fiximage, chan, amount);
204 ret = imoper(fit, fiximage, OPER_ADD, FALSE);
205
206 invalidate_stats_from_fit(fit);
207 clearfits(fiximage);
208 if ((!ret) && applyRotation) {
209 point center = {gfit.rx / 2.0, gfit.ry / 2.0};
210 cvRotateImage(fit, center, -90.0, -1, OPENCV_AREA);
211 }
212
213 return ret;
214 }
215
BandingEngine_float(fits * fit,double sigma,double amount,gboolean protect_highlights,gboolean applyRotation)216 static int BandingEngine_float(fits *fit, double sigma, double amount, gboolean protect_highlights, gboolean applyRotation) {
217 int chan, row, i, ret = 0;
218 float *line, *fixline;
219 double minimum = DBL_MAX, globalsigma = 0.0;
220 fits *fiximage = NULL;
221 double invsigma = 1.0 / sigma;
222
223 if (applyRotation) {
224 point center = {gfit.rx / 2.0, gfit.ry / 2.0};
225 cvRotateImage(fit, center, 90.0, -1, OPENCV_AREA);
226 }
227
228 if (new_fit_image(&fiximage, fit->rx, fit->ry, fit->naxes[2], DATA_FLOAT))
229 return 1;
230
231 for (chan = 0; chan < fit->naxes[2]; chan++) {
232 imstats *stat = statistics(NULL, -1, fit, chan, NULL, STATS_BASIC | STATS_MAD, TRUE);
233 if (!stat) {
234 siril_log_message(_("Error: statistics computation failed.\n"));
235 return 1;
236 }
237 double background = stat->median;
238 double *rowvalue = calloc(fit->ry, sizeof(double));
239 if (rowvalue == NULL) {
240 PRINT_ALLOC_ERR;
241 free_stats(stat);
242 return 1;
243 }
244 if (protect_highlights) {
245 globalsigma = stat->mad * MAD_NORM;
246 }
247 free_stats(stat);
248 for (row = 0; row < fit->ry; row++) {
249 line = fit->fpdata[chan] + row * fit->rx;
250 float *cpyline = calloc(fit->rx, sizeof(float));
251 if (cpyline == NULL) {
252 PRINT_ALLOC_ERR;
253 free(rowvalue);
254 return 1;
255 }
256 memcpy(cpyline, line, fit->rx * sizeof(float));
257 int n = fit->rx;
258 double median;
259 if (protect_highlights) {
260 quicksort_f(cpyline, n);
261 float reject = background + invsigma * globalsigma;
262 for (i = fit->rx - 1; i >= 0; i--) {
263 if (cpyline[i] < reject)
264 break;
265 n--;
266 }
267 median = gsl_stats_float_median_from_sorted_data(cpyline, 1, n);
268 } else {
269 median = quickmedian_float(cpyline, n);
270 }
271
272 rowvalue[row] = background - median;
273 minimum = min(minimum, rowvalue[row]);
274 free(cpyline);
275 }
276 for (row = 0; row < fit->ry; row++) {
277 fixline = fiximage->fpdata[chan] + row * fiximage->rx;
278 for (i = 0; i < fit->rx; i++)
279 fixline[i] = rowvalue[row] - minimum;
280 }
281 free(rowvalue);
282 }
283 for (chan = 0; chan < fit->naxes[2]; chan++)
284 fmul_layer_float(fiximage, chan, amount);
285 ret = imoper(fit, fiximage, OPER_ADD, TRUE);
286
287 invalidate_stats_from_fit(fit);
288 clearfits(fiximage);
289 if ((!ret) && applyRotation) {
290 point center = {gfit.rx / 2.0, gfit.ry / 2.0};
291 cvRotateImage(fit, center, -90.0, -1, OPENCV_AREA);
292 }
293
294 return ret;
295 }
296
BandingEngine(fits * fit,double sigma,double amount,gboolean protect_highlights,gboolean applyRotation)297 int BandingEngine(fits *fit, double sigma, double amount, gboolean protect_highlights, gboolean applyRotation) {
298 if (fit->type == DATA_FLOAT)
299 return BandingEngine_float(fit, sigma, amount, protect_highlights, applyRotation);
300 if (fit->type == DATA_USHORT)
301 return BandingEngine_ushort(fit, sigma, amount, protect_highlights, applyRotation);
302 return -1;
303 }
304
305 /***************** GUI for Canon Banding Reduction ********************/
306
on_button_ok_fixbanding_clicked(GtkButton * button,gpointer user_data)307 void on_button_ok_fixbanding_clicked(GtkButton *button, gpointer user_data) {
308 siril_close_dialog("canon_fixbanding_dialog");
309 }
310
on_button_apply_fixbanding_clicked(GtkButton * button,gpointer user_data)311 void on_button_apply_fixbanding_clicked(GtkButton *button, gpointer user_data) {
312 static GtkRange *range_amount = NULL;
313 static GtkRange *range_invsigma = NULL;
314 static GtkToggleButton *toggle_protect_highlights_banding = NULL,
315 *vertical = NULL, *seq = NULL;
316 static GtkEntry *bandingSeqEntry = NULL;
317 double amount, invsigma;
318 gboolean protect_highlights;
319
320 if (get_thread_run()) {
321 PRINT_ANOTHER_THREAD_RUNNING;
322 return;
323 }
324
325 struct banding_data *args = malloc(sizeof(struct banding_data));
326
327 if (range_amount == NULL) {
328 range_amount = GTK_RANGE(lookup_widget("scale_fixbanding_amount"));
329 range_invsigma = GTK_RANGE(lookup_widget("scale_fixbanding_invsigma"));
330 toggle_protect_highlights_banding = GTK_TOGGLE_BUTTON(
331 lookup_widget("checkbutton_fixbanding"));
332 vertical = GTK_TOGGLE_BUTTON(lookup_widget("checkBandingVertical"));
333 seq = GTK_TOGGLE_BUTTON(lookup_widget("checkBandingSeq"));
334 bandingSeqEntry = GTK_ENTRY(lookup_widget("entryBandingSeq"));
335 }
336 amount = gtk_range_get_value(range_amount);
337 invsigma = gtk_range_get_value(range_invsigma);
338 protect_highlights = gtk_toggle_button_get_active(
339 toggle_protect_highlights_banding);
340
341 if (!protect_highlights)
342 undo_save_state(&gfit, _("Canon Banding Reduction (amount=%.2lf)"), amount);
343 else
344 undo_save_state(&gfit, _("Canon Banding Reduction (amount=%.2lf, Protect=TRUE, invsigma=%.2lf)"),
345 amount, invsigma);
346
347 args->fit = &gfit;
348 args->protect_highlights = protect_highlights;
349 args->amount = amount;
350 args->sigma = invsigma;
351 args->applyRotation = gtk_toggle_button_get_active(vertical);
352 args->seqEntry = gtk_entry_get_text(bandingSeqEntry);
353 set_cursor_waiting(TRUE);
354
355 if (gtk_toggle_button_get_active(seq) && sequence_is_loaded()) {
356 if (args->seqEntry && args->seqEntry[0] == '\0')
357 args->seqEntry = "unband_";
358 apply_banding_to_sequence(args);
359 } else {
360 start_in_new_thread(BandingEngineThreaded, args);
361 }
362 }
363
on_checkbutton_fixbanding_toggled(GtkToggleButton * togglebutton,gpointer user_data)364 void on_checkbutton_fixbanding_toggled(GtkToggleButton *togglebutton,
365 gpointer user_data) {
366 static GtkWidget *scalebandingHighlightBox = NULL;
367 static GtkWidget *spinbandingHighlightBox = NULL;
368 gboolean is_active;
369
370 if (scalebandingHighlightBox == NULL) {
371 scalebandingHighlightBox = lookup_widget("scale_fixbanding_invsigma");
372 spinbandingHighlightBox = lookup_widget("spin_fixbanding_invsigma");
373 }
374
375 is_active = gtk_toggle_button_get_active(togglebutton);
376 gtk_widget_set_sensitive(scalebandingHighlightBox, is_active);
377 gtk_widget_set_sensitive(spinbandingHighlightBox, is_active);
378 }
379