1 /*
2 This file is part of darktable,
3 Copyright (C) 2010-2020 darktable developers.
4
5 darktable 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 3 of the License, or
8 (at your option) any later version.
9
10 darktable 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 darktable. If not, see <http://www.gnu.org/licenses/>.
17 */
18 #ifdef HAVE_CONFIG_H
19 #include "config.h"
20 #endif
21 #include <assert.h>
22 #include <math.h>
23 #include <stdlib.h>
24 #include <string.h>
25 #include "bauhaus/bauhaus.h"
26 #include "common/opencl.h"
27 #include "control/control.h"
28 #include "develop/develop.h"
29 #include "develop/imageop.h"
30 #include "develop/imageop_math.h"
31 #include "develop/imageop_gui.h"
32 #include "develop/tiling.h"
33 #include "gui/accelerators.h"
34 #include "gui/gtk.h"
35 #include "iop/iop_api.h"
36
37 #include <gtk/gtk.h>
38 #include <inttypes.h>
39
40
41 DT_MODULE_INTROSPECTION(2, dt_iop_highlights_params_t)
42
43 typedef enum dt_iop_highlights_mode_t
44 {
45 DT_IOP_HIGHLIGHTS_CLIP = 0, // $DESCRIPTION: "clip highlights"
46 DT_IOP_HIGHLIGHTS_LCH = 1, // $DESCRIPTION: "reconstruct in LCh"
47 DT_IOP_HIGHLIGHTS_INPAINT = 2, // $DESCRIPTION: "reconstruct color"
48 } dt_iop_highlights_mode_t;
49
50 typedef struct dt_iop_highlights_params_t
51 {
52 dt_iop_highlights_mode_t mode; // $DEFAULT: DT_IOP_HIGHLIGHTS_CLIP $DESCRIPTION: "method"
53 float blendL; // unused $DEFAULT: 1.0
54 float blendC; // unused $DEFAULT: 0.0
55 float blendh; // unused $DEFAULT: 0.0
56 float clip; // $MIN: 0.0 $MAX: 2.0 $DEFAULT: 1.0 $DESCRIPTION: "clipping threshold"
57 } dt_iop_highlights_params_t;
58
59 typedef struct dt_iop_highlights_gui_data_t
60 {
61 GtkWidget *clip;
62 GtkWidget *mode;
63 } dt_iop_highlights_gui_data_t;
64
65 typedef dt_iop_highlights_params_t dt_iop_highlights_data_t;
66
67 typedef struct dt_iop_highlights_global_data_t
68 {
69 int kernel_highlights_1f_clip;
70 int kernel_highlights_1f_lch_bayer;
71 int kernel_highlights_1f_lch_xtrans;
72 int kernel_highlights_4f_clip;
73 } dt_iop_highlights_global_data_t;
74
75
name()76 const char *name()
77 {
78 return _("highlight reconstruction");
79 }
80
description(struct dt_iop_module_t * self)81 const char *description(struct dt_iop_module_t *self)
82 {
83 return dt_iop_set_description(self, _("avoid magenta highlights and try to recover highlights colors"),
84 _("corrective"),
85 _("linear, raw, scene-referred"),
86 _("reconstruction, raw"),
87 _("linear, raw, scene-referred"));
88 }
89
default_group()90 int default_group()
91 {
92 return IOP_GROUP_BASIC | IOP_GROUP_TECHNICAL;
93 }
94
flags()95 int flags()
96 {
97 return IOP_FLAGS_SUPPORTS_BLENDING | IOP_FLAGS_ALLOW_TILING | IOP_FLAGS_ONE_INSTANCE;
98 }
99
default_colorspace(dt_iop_module_t * self,dt_dev_pixelpipe_t * pipe,dt_dev_pixelpipe_iop_t * piece)100 int default_colorspace(dt_iop_module_t *self, dt_dev_pixelpipe_t *pipe, dt_dev_pixelpipe_iop_t *piece)
101 {
102 return iop_cs_RAW;
103 }
104
legacy_params(dt_iop_module_t * self,const void * const old_params,const int old_version,void * new_params,const int new_version)105 int legacy_params(dt_iop_module_t *self, const void *const old_params, const int old_version,
106 void *new_params, const int new_version)
107 {
108 if(old_version == 1 && new_version == 2)
109 {
110 memcpy(new_params, old_params, sizeof(dt_iop_highlights_params_t) - sizeof(float));
111 dt_iop_highlights_params_t *n = (dt_iop_highlights_params_t *)new_params;
112 n->clip = 1.0f;
113 return 0;
114 }
115 return 1;
116 }
117
118 #ifdef HAVE_OPENCL
process_cl(struct dt_iop_module_t * self,dt_dev_pixelpipe_iop_t * piece,cl_mem dev_in,cl_mem dev_out,const dt_iop_roi_t * const roi_in,const dt_iop_roi_t * const roi_out)119 int process_cl(struct dt_iop_module_t *self, dt_dev_pixelpipe_iop_t *piece, cl_mem dev_in, cl_mem dev_out,
120 const dt_iop_roi_t *const roi_in, const dt_iop_roi_t *const roi_out)
121 {
122 dt_iop_highlights_data_t *d = (dt_iop_highlights_data_t *)piece->data;
123 dt_iop_highlights_global_data_t *gd = (dt_iop_highlights_global_data_t *)self->global_data;
124
125 cl_int err = -999;
126 cl_mem dev_xtrans = NULL;
127
128 const int devid = piece->pipe->devid;
129 const int width = roi_in->width;
130 const int height = roi_in->height;
131
132 const float clip = d->clip
133 * fminf(piece->pipe->dsc.processed_maximum[0],
134 fminf(piece->pipe->dsc.processed_maximum[1], piece->pipe->dsc.processed_maximum[2]));
135
136 const uint32_t filters = piece->pipe->dsc.filters;
137
138 if(!filters)
139 {
140 // non-raw images use dedicated kernel which just clips
141 size_t sizes[] = { ROUNDUPWD(width), ROUNDUPHT(height), 1 };
142 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_4f_clip, 0, sizeof(cl_mem), (void *)&dev_in);
143 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_4f_clip, 1, sizeof(cl_mem), (void *)&dev_out);
144 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_4f_clip, 2, sizeof(int), (void *)&width);
145 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_4f_clip, 3, sizeof(int), (void *)&height);
146 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_4f_clip, 4, sizeof(int), (void *)&d->mode);
147 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_4f_clip, 5, sizeof(float), (void *)&clip);
148 err = dt_opencl_enqueue_kernel_2d(devid, gd->kernel_highlights_4f_clip, sizes);
149 if(err != CL_SUCCESS) goto error;
150 }
151 else if(d->mode == DT_IOP_HIGHLIGHTS_CLIP)
152 {
153 // raw images with clip mode (both bayer and xtrans)
154 size_t sizes[] = { ROUNDUPWD(width), ROUNDUPHT(height), 1 };
155 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_clip, 0, sizeof(cl_mem), (void *)&dev_in);
156 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_clip, 1, sizeof(cl_mem), (void *)&dev_out);
157 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_clip, 2, sizeof(int), (void *)&width);
158 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_clip, 3, sizeof(int), (void *)&height);
159 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_clip, 4, sizeof(float), (void *)&clip);
160 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_clip, 5, sizeof(int), (void *)&roi_out->x);
161 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_clip, 6, sizeof(int), (void *)&roi_out->y);
162 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_clip, 7, sizeof(int), (void *)&filters);
163 err = dt_opencl_enqueue_kernel_2d(devid, gd->kernel_highlights_1f_clip, sizes);
164 if(err != CL_SUCCESS) goto error;
165 }
166 else if(d->mode == DT_IOP_HIGHLIGHTS_LCH && filters != 9u)
167 {
168 // bayer sensor raws with LCH mode
169 size_t sizes[] = { ROUNDUPWD(width), ROUNDUPHT(height), 1 };
170 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_bayer, 0, sizeof(cl_mem), (void *)&dev_in);
171 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_bayer, 1, sizeof(cl_mem), (void *)&dev_out);
172 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_bayer, 2, sizeof(int), (void *)&width);
173 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_bayer, 3, sizeof(int), (void *)&height);
174 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_bayer, 4, sizeof(float), (void *)&clip);
175 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_bayer, 5, sizeof(int), (void *)&roi_out->x);
176 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_bayer, 6, sizeof(int), (void *)&roi_out->y);
177 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_bayer, 7, sizeof(int), (void *)&filters);
178 err = dt_opencl_enqueue_kernel_2d(devid, gd->kernel_highlights_1f_lch_bayer, sizes);
179 if(err != CL_SUCCESS) goto error;
180 }
181 else if(d->mode == DT_IOP_HIGHLIGHTS_LCH && filters == 9u)
182 {
183 // xtrans sensor raws with LCH mode
184 int blocksizex, blocksizey;
185
186 dt_opencl_local_buffer_t locopt
187 = (dt_opencl_local_buffer_t){ .xoffset = 2 * 2, .xfactor = 1, .yoffset = 2 * 2, .yfactor = 1,
188 .cellsize = sizeof(float), .overhead = 0,
189 .sizex = 1 << 8, .sizey = 1 << 8 };
190
191 if(dt_opencl_local_buffer_opt(devid, gd->kernel_highlights_1f_lch_xtrans, &locopt))
192 {
193 blocksizex = locopt.sizex;
194 blocksizey = locopt.sizey;
195 }
196 else
197 blocksizex = blocksizey = 1;
198
199 dev_xtrans
200 = dt_opencl_copy_host_to_device_constant(devid, sizeof(piece->pipe->dsc.xtrans), piece->pipe->dsc.xtrans);
201 if(dev_xtrans == NULL) goto error;
202
203 size_t sizes[] = { ROUNDUP(width, blocksizex), ROUNDUP(height, blocksizey), 1 };
204 size_t local[] = { blocksizex, blocksizey, 1 };
205 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_xtrans, 0, sizeof(cl_mem), (void *)&dev_in);
206 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_xtrans, 1, sizeof(cl_mem), (void *)&dev_out);
207 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_xtrans, 2, sizeof(int), (void *)&width);
208 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_xtrans, 3, sizeof(int), (void *)&height);
209 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_xtrans, 4, sizeof(float), (void *)&clip);
210 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_xtrans, 5, sizeof(int), (void *)&roi_out->x);
211 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_xtrans, 6, sizeof(int), (void *)&roi_out->y);
212 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_xtrans, 7, sizeof(cl_mem), (void *)&dev_xtrans);
213 dt_opencl_set_kernel_arg(devid, gd->kernel_highlights_1f_lch_xtrans, 8,
214 sizeof(float) * (blocksizex + 4) * (blocksizey + 4), NULL);
215
216 err = dt_opencl_enqueue_kernel_2d_with_local(devid, gd->kernel_highlights_1f_lch_xtrans, sizes, local);
217 if(err != CL_SUCCESS) goto error;
218 }
219
220 // update processed maximum
221 const float m = fmaxf(fmaxf(piece->pipe->dsc.processed_maximum[0], piece->pipe->dsc.processed_maximum[1]),
222 piece->pipe->dsc.processed_maximum[2]);
223 for(int k = 0; k < 3; k++) piece->pipe->dsc.processed_maximum[k] = m;
224
225 dt_opencl_release_mem_object(dev_xtrans);
226 return TRUE;
227
228 error:
229 dt_opencl_release_mem_object(dev_xtrans);
230 dt_print(DT_DEBUG_OPENCL, "[opencl_highlights] couldn't enqueue kernel! %d\n", err);
231 return FALSE;
232 }
233 #endif
234
tiling_callback(struct dt_iop_module_t * self,struct dt_dev_pixelpipe_iop_t * piece,const dt_iop_roi_t * roi_in,const dt_iop_roi_t * roi_out,struct dt_develop_tiling_t * tiling)235 void tiling_callback(struct dt_iop_module_t *self, struct dt_dev_pixelpipe_iop_t *piece,
236 const dt_iop_roi_t *roi_in, const dt_iop_roi_t *roi_out,
237 struct dt_develop_tiling_t *tiling)
238 {
239 dt_iop_highlights_data_t *d = (dt_iop_highlights_data_t *)piece->data;
240 const uint32_t filters = piece->pipe->dsc.filters;
241
242 tiling->factor = 2.0f; // in + out
243 tiling->maxbuf = 1.0f;
244 tiling->overhead = 0;
245
246 if(filters == 9u)
247 {
248 // xtrans
249 tiling->xalign = 6;
250 tiling->yalign = 6;
251 tiling->overlap = (d->mode == DT_IOP_HIGHLIGHTS_LCH) ? 2 : 0;
252 }
253 else if(filters)
254 {
255 // bayer
256 tiling->xalign = 2;
257 tiling->yalign = 2;
258 tiling->overlap = (d->mode == DT_IOP_HIGHLIGHTS_LCH) ? 1 : 0;
259 }
260 else
261 {
262 // non-raw
263 tiling->xalign = 1;
264 tiling->yalign = 1;
265 tiling->overlap = 0;
266 }
267 }
268
269 /* interpolate value for a pixel, ideal via ratio to nearby pixel */
interp_pix_xtrans(const int ratio_next,const ssize_t offset_next,const float clip0,const float clip_next,const float * const in,const float * const ratios)270 static inline float interp_pix_xtrans(const int ratio_next,
271 const ssize_t offset_next,
272 const float clip0, const float clip_next,
273 const float *const in,
274 const float *const ratios)
275 {
276 assert(ratio_next != 0);
277 // it's OK to exceed clipping of current pixel's color based on a
278 // neighbor -- that is the purpose of interpolating highlight
279 // colors
280 const float clip_val = fmaxf(clip0, clip_next);
281 if(in[offset_next] >= clip_next - 1e-5f)
282 {
283 // next pixel is also clipped
284 return clip_val;
285 }
286 else
287 {
288 // set this pixel in ratio to the next
289 assert(ratio_next != 0);
290 if (ratio_next > 0)
291 return fminf(in[offset_next] / ratios[ratio_next], clip_val);
292 else
293 return fminf(in[offset_next] * ratios[-ratio_next], clip_val);
294 }
295 }
296
interpolate_color_xtrans(const void * const ivoid,void * const ovoid,const dt_iop_roi_t * const roi_in,const dt_iop_roi_t * const roi_out,int dim,int dir,int other,const float * const clip,const uint8_t (* const xtrans)[6],const int pass)297 static inline void interpolate_color_xtrans(const void *const ivoid, void *const ovoid,
298 const dt_iop_roi_t *const roi_in,
299 const dt_iop_roi_t *const roi_out,
300 int dim, int dir, int other,
301 const float *const clip,
302 const uint8_t (*const xtrans)[6],
303 const int pass)
304 {
305 // In Bayer each row/col has only green/red or green/blue
306 // transitions, hence can reconstruct color by single ratio per
307 // row. In x-trans there can be transitions between arbitrary colors
308 // in a row/col (and 2x2 green blocks which provide no color
309 // transition information). Hence calculate multiple color ratios
310 // for each row/col.
311
312 // Lookup for color ratios, e.g. red -> blue is roff[0][2] and blue
313 // -> red is roff[2][0]. Returned value is an index into ratios. If
314 // negative, then need to invert the ratio. Identity color
315 // transitions aren't used.
316 const int roff[3][3] = {{ 0, -1, -2},
317 { 1, 0, -3},
318 { 2, 3, 0}};
319 // record ratios of color transitions 0:unused, 1:RG, 2:RB, and 3:GB
320 dt_aligned_pixel_t ratios = {1.0f, 1.0f, 1.0f, 1.0f};
321
322 // passes are 0:+x, 1:-x, 2:+y, 3:-y
323 // dims are 0:traverse a row, 1:traverse a column
324 // dir is 1:left to right, -1: right to left
325 int i = (dim == 0) ? 0 : other;
326 int j = (dim == 0) ? other : 0;
327 const ssize_t offs = (ssize_t)(dim ? roi_out->width : 1) * ((dir < 0) ? -1 : 1);
328 const ssize_t offl = offs - (dim ? 1 : roi_out->width);
329 const ssize_t offr = offs + (dim ? 1 : roi_out->width);
330 int beg, end;
331 if(dir == 1)
332 {
333 beg = 0;
334 end = (dim == 0) ? roi_out->width : roi_out->height;
335 }
336 else
337 {
338 beg = ((dim == 0) ? roi_out->width : roi_out->height) - 1;
339 end = -1;
340 }
341
342 float *in, *out;
343 if(dim == 1)
344 {
345 out = (float *)ovoid + (size_t)i + (size_t)beg * roi_out->width;
346 in = (float *)ivoid + (size_t)i + (size_t)beg * roi_in->width;
347 }
348 else
349 {
350 out = (float *)ovoid + (size_t)beg + (size_t)j * roi_out->width;
351 in = (float *)ivoid + (size_t)beg + (size_t)j * roi_in->width;
352 }
353
354 for(int k = beg; k != end; k += dir)
355 {
356 if(dim == 1)
357 j = k;
358 else
359 i = k;
360
361 const uint8_t f0 = FCxtrans(j, i, roi_in, xtrans);
362 const uint8_t f1 = FCxtrans(dim ? (j + dir) : j, dim ? i : (i + dir), roi_in, xtrans);
363 const uint8_t fl = FCxtrans(dim ? (j + dir) : (j - 1), dim ? (i - 1) : (i + dir), roi_in, xtrans);
364 const uint8_t fr = FCxtrans(dim ? (j + dir) : (j + 1), dim ? (i + 1) : (i + dir), roi_in, xtrans);
365 const float clip0 = clip[f0];
366 const float clip1 = clip[f1];
367 const float clipl = clip[fl];
368 const float clipr = clip[fr];
369 const float clip_max = fmaxf(fmaxf(clip[0], clip[1]), clip[2]);
370
371 if(i == 0 || i == roi_out->width - 1 || j == 0 || j == roi_out->height - 1)
372 {
373 if(pass == 3) out[0] = fminf(clip_max, in[0]);
374 }
375 else
376 {
377 // ratio to next pixel if this & next are unclamped and not in
378 // 2x2 green block
379 if ((f0 != f1) &&
380 (in[0] < clip0 && in[0] > 1e-5f) &&
381 (in[offs] < clip1 && in[offs] > 1e-5f))
382 {
383 const int r = roff[f0][f1];
384 assert(r != 0);
385 if (r > 0)
386 ratios[r] = (3.f * ratios[r] + (in[offs] / in[0])) / 4.f;
387 else
388 ratios[-r] = (3.f * ratios[-r] + (in[0] / in[offs])) / 4.f;
389 }
390
391 if(in[0] >= clip0 - 1e-5f)
392 {
393 // interplate color for clipped pixel
394 float add;
395 if(f0 != f1)
396 // next pixel is different color
397 add =
398 interp_pix_xtrans(roff[f0][f1], offs, clip0, clip1, in, ratios);
399 else
400 // at start of 2x2 green block, look diagonally
401 add = (fl != f0) ?
402 interp_pix_xtrans(roff[f0][fl], offl, clip0, clipl, in, ratios) :
403 interp_pix_xtrans(roff[f0][fr], offr, clip0, clipr, in, ratios);
404
405 if(pass == 0)
406 out[0] = add;
407 else if(pass == 3)
408 out[0] = fminf(clip_max, (out[0] + add) / 4.0f);
409 else
410 out[0] += add;
411 }
412 else
413 {
414 // pixel is not clipped
415 if(pass == 3) out[0] = in[0];
416 }
417 }
418 out += offs;
419 in += offs;
420 }
421 }
422
interpolate_color(const void * const ivoid,void * const ovoid,const dt_iop_roi_t * const roi_out,int dim,int dir,int other,const float * clip,const uint32_t filters,const int pass)423 static inline void interpolate_color(const void *const ivoid, void *const ovoid,
424 const dt_iop_roi_t *const roi_out, int dim, int dir, int other,
425 const float *clip, const uint32_t filters, const int pass)
426 {
427 float ratio = 1.0f;
428 float *in, *out;
429
430 int i = 0, j = 0;
431 if(dim == 0)
432 j = other;
433 else
434 i = other;
435 ssize_t offs = dim ? roi_out->width : 1;
436 if(dir < 0) offs = -offs;
437 int beg, end;
438 if(dim == 0 && dir == 1)
439 {
440 beg = 0;
441 end = roi_out->width;
442 }
443 else if(dim == 0 && dir == -1)
444 {
445 beg = roi_out->width - 1;
446 end = -1;
447 }
448 else if(dim == 1 && dir == 1)
449 {
450 beg = 0;
451 end = roi_out->height;
452 }
453 else if(dim == 1 && dir == -1)
454 {
455 beg = roi_out->height - 1;
456 end = -1;
457 }
458 else
459 return;
460
461 if(dim == 1)
462 {
463 out = (float *)ovoid + i + (size_t)beg * roi_out->width;
464 in = (float *)ivoid + i + (size_t)beg * roi_out->width;
465 }
466 else
467 {
468 out = (float *)ovoid + beg + (size_t)j * roi_out->width;
469 in = (float *)ivoid + beg + (size_t)j * roi_out->width;
470 }
471 for(int k = beg; k != end; k += dir)
472 {
473 if(dim == 1)
474 j = k;
475 else
476 i = k;
477 const float clip0 = clip[FC(j, i, filters)];
478 const float clip1 = clip[FC(dim ? (j + 1) : j, dim ? i : (i + 1), filters)];
479 if(i == 0 || i == roi_out->width - 1 || j == 0 || j == roi_out->height - 1)
480 {
481 if(pass == 3) out[0] = in[0];
482 }
483 else
484 {
485 if(in[0] < clip0 && in[0] > 1e-5f)
486 { // both are not clipped
487 if(in[offs] < clip1 && in[offs] > 1e-5f)
488 { // update ratio, exponential decay. ratio = in[odd]/in[even]
489 if(k & 1)
490 ratio = (3.0f * ratio + in[0] / in[offs]) / 4.0f;
491 else
492 ratio = (3.0f * ratio + in[offs] / in[0]) / 4.0f;
493 }
494 }
495
496 if(in[0] >= clip0 - 1e-5f)
497 { // in[0] is clipped, restore it as in[1] adjusted according to ratio
498 float add = 0.0f;
499 if(in[offs] >= clip1 - 1e-5f)
500 add = fmaxf(clip0, clip1);
501 else if(k & 1)
502 add = in[offs] * ratio;
503 else
504 add = in[offs] / ratio;
505
506 if(pass == 0)
507 out[0] = add;
508 else if(pass == 3)
509 out[0] = (out[0] + add) / 4.0f;
510 else
511 out[0] += add;
512 }
513 else
514 {
515 if(pass == 3) out[0] = in[0];
516 }
517 }
518 out += offs;
519 in += offs;
520 }
521 }
522
523 /*
524 * these 2 constants were computed using following Sage code:
525 *
526 * sqrt3 = sqrt(3)
527 * sqrt12 = sqrt(12) # 2*sqrt(3)
528 *
529 * print 'sqrt3 = ', sqrt3, ' ~= ', RealField(128)(sqrt3)
530 * print 'sqrt12 = ', sqrt12, ' ~= ', RealField(128)(sqrt12)
531 */
532 #define SQRT3 1.7320508075688772935274463415058723669L
533 #define SQRT12 3.4641016151377545870548926830117447339L // 2*SQRT3
534
process_lch_bayer(dt_iop_module_t * self,dt_dev_pixelpipe_iop_t * piece,const void * const ivoid,void * const ovoid,const dt_iop_roi_t * const roi_in,const dt_iop_roi_t * const roi_out,const float clip)535 static void process_lch_bayer(dt_iop_module_t *self, dt_dev_pixelpipe_iop_t *piece, const void *const ivoid,
536 void *const ovoid, const dt_iop_roi_t *const roi_in,
537 const dt_iop_roi_t *const roi_out, const float clip)
538 {
539 const uint32_t filters = piece->pipe->dsc.filters;
540
541 #ifdef _OPENMP
542 #pragma omp parallel for default(none) \
543 dt_omp_firstprivate(clip, filters, ivoid, ovoid, roi_out) \
544 schedule(static) collapse(2)
545 #endif
546 for(int j = 0; j < roi_out->height; j++)
547 {
548 for(int i = 0; i < roi_out->width; i++)
549 {
550 float *const out = (float *)ovoid + (size_t)roi_out->width * j + i;
551 const float *const in = (float *)ivoid + (size_t)roi_out->width * j + i;
552
553 if(i == roi_out->width - 1 || j == roi_out->height - 1)
554 {
555 // fast path for border
556 out[0] = MIN(clip, in[0]);
557 }
558 else
559 {
560 int clipped = 0;
561
562 // sample 1 bayer block. thus we will have 2 green values.
563 float R = 0.0f, Gmin = FLT_MAX, Gmax = -FLT_MAX, B = 0.0f;
564 for(int jj = 0; jj <= 1; jj++)
565 {
566 for(int ii = 0; ii <= 1; ii++)
567 {
568 const float val = in[(size_t)jj * roi_out->width + ii];
569
570 clipped = (clipped || (val > clip));
571
572 const int c = FC(j + jj + roi_out->y, i + ii + roi_out->x, filters);
573 switch(c)
574 {
575 case 0:
576 R = val;
577 break;
578 case 1:
579 Gmin = MIN(Gmin, val);
580 Gmax = MAX(Gmax, val);
581 break;
582 case 2:
583 B = val;
584 break;
585 }
586 }
587 }
588
589 if(clipped)
590 {
591 const float Ro = MIN(R, clip);
592 const float Go = MIN(Gmin, clip);
593 const float Bo = MIN(B, clip);
594
595 const float L = (R + Gmax + B) / 3.0f;
596
597 float C = SQRT3 * (R - Gmax);
598 float H = 2.0f * B - Gmax - R;
599
600 const float Co = SQRT3 * (Ro - Go);
601 const float Ho = 2.0f * Bo - Go - Ro;
602
603 if(R != Gmax && Gmax != B)
604 {
605 const float ratio = sqrtf((Co * Co + Ho * Ho) / (C * C + H * H));
606 C *= ratio;
607 H *= ratio;
608 }
609
610 dt_aligned_pixel_t RGB = { 0.0f, 0.0f, 0.0f };
611
612 /*
613 * backtransform proof, sage:
614 *
615 * R,G,B,L,C,H = var('R,G,B,L,C,H')
616 * solve([L==(R+G+B)/3, C==sqrt(3)*(R-G), H==2*B-G-R], R, G, B)
617 *
618 * result:
619 * [[R == 1/6*sqrt(3)*C - 1/6*H + L, G == -1/6*sqrt(3)*C - 1/6*H + L, B == 1/3*H + L]]
620 */
621 RGB[0] = L - H / 6.0f + C / SQRT12;
622 RGB[1] = L - H / 6.0f - C / SQRT12;
623 RGB[2] = L + H / 3.0f;
624
625 out[0] = RGB[FC(j + roi_out->y, i + roi_out->x, filters)];
626 }
627 else
628 {
629 out[0] = in[0];
630 }
631 }
632 }
633 }
634 }
635
process_lch_xtrans(dt_iop_module_t * self,dt_dev_pixelpipe_iop_t * piece,const void * const ivoid,void * const ovoid,const dt_iop_roi_t * const roi_in,const dt_iop_roi_t * const roi_out,const float clip)636 static void process_lch_xtrans(dt_iop_module_t *self, dt_dev_pixelpipe_iop_t *piece, const void *const ivoid,
637 void *const ovoid, const dt_iop_roi_t *const roi_in,
638 const dt_iop_roi_t *const roi_out, const float clip)
639 {
640 const uint8_t(*const xtrans)[6] = (const uint8_t(*const)[6])piece->pipe->dsc.xtrans;
641
642 #ifdef _OPENMP
643 #pragma omp parallel for default(none) \
644 dt_omp_firstprivate(clip, ivoid, ovoid, roi_in, roi_out, xtrans) \
645 schedule(static)
646 #endif
647 for(int j = 0; j < roi_out->height; j++)
648 {
649 float *out = (float *)ovoid + (size_t)roi_out->width * j;
650 float *in = (float *)ivoid + (size_t)roi_in->width * j;
651
652 // bit vector used as ring buffer to remember clipping of current
653 // and last two columns, checking current pixel and its vertical
654 // neighbors
655 int cl = 0;
656
657 for(int i = 0; i < roi_out->width; i++)
658 {
659 // update clipping ring buffer
660 cl = (cl << 1) & 6;
661 if(j >= 2 && j <= roi_out->height - 3)
662 {
663 cl |= (in[-roi_in->width] > clip) | (in[0] > clip) | (in[roi_in->width] > clip);
664 }
665
666 if(i < 2 || i > roi_out->width - 3 || j < 2 || j > roi_out->height - 3)
667 {
668 // fast path for border
669 out[0] = MIN(clip, in[0]);
670 }
671 else
672 {
673 // if current pixel is clipped, always reconstruct
674 int clipped = (in[0] > clip);
675 if(!clipped)
676 {
677 clipped = cl;
678 if(clipped)
679 {
680 // If the ring buffer can't show we are in an obviously
681 // unclipped region, this is the slow case: check if there
682 // is any 3x3 block touching the current pixel which has
683 // no clipping, as then don't need to reconstruct the
684 // current pixel. This avoids zippering in edge
685 // transitions from clipped to unclipped areas. The
686 // X-Trans sensor seems prone to this, unlike Bayer, due
687 // to its irregular pattern.
688 for(int offset_j = -2; offset_j <= 0; offset_j++)
689 {
690 for(int offset_i = -2; offset_i <= 0; offset_i++)
691 {
692 if(clipped)
693 {
694 clipped = 0;
695 for(int jj = offset_j; jj <= offset_j + 2; jj++)
696 {
697 for(int ii = offset_i; ii <= offset_i + 2; ii++)
698 {
699 const float val = in[(ssize_t)jj * roi_in->width + ii];
700 clipped = (clipped || (val > clip));
701 }
702 }
703 }
704 }
705 }
706 }
707 }
708
709 if(clipped)
710 {
711 dt_aligned_pixel_t mean = { 0.0f, 0.0f, 0.0f };
712 dt_aligned_pixel_t RGBmax = { -FLT_MAX, -FLT_MAX, -FLT_MAX };
713 int cnt[3] = { 0, 0, 0 };
714
715 for(int jj = -1; jj <= 1; jj++)
716 {
717 for(int ii = -1; ii <= 1; ii++)
718 {
719 const float val = in[(ssize_t)jj * roi_in->width + ii];
720 const int c = FCxtrans(j+jj, i+ii, roi_in, xtrans);
721 mean[c] += val;
722 cnt[c]++;
723 RGBmax[c] = MAX(RGBmax[c], val);
724 }
725 }
726
727 const float Ro = MIN(mean[0]/cnt[0], clip);
728 const float Go = MIN(mean[1]/cnt[1], clip);
729 const float Bo = MIN(mean[2]/cnt[2], clip);
730
731 const float R = RGBmax[0];
732 const float G = RGBmax[1];
733 const float B = RGBmax[2];
734
735 const float L = (R + G + B) / 3.0f;
736
737 float C = SQRT3 * (R - G);
738 float H = 2.0f * B - G - R;
739
740 const float Co = SQRT3 * (Ro - Go);
741 const float Ho = 2.0f * Bo - Go - Ro;
742
743 if(R != G && G != B)
744 {
745 const float ratio = sqrtf((Co * Co + Ho * Ho) / (C * C + H * H));
746 C *= ratio;
747 H *= ratio;
748 }
749
750 dt_aligned_pixel_t RGB = { 0.0f, 0.0f, 0.0f };
751
752 RGB[0] = L - H / 6.0f + C / SQRT12;
753 RGB[1] = L - H / 6.0f - C / SQRT12;
754 RGB[2] = L + H / 3.0f;
755
756 out[0] = RGB[FCxtrans(j, i, roi_out, xtrans)];
757 }
758 else
759 out[0] = in[0];
760 }
761 out++;
762 in++;
763 }
764 }
765 }
766
767 #undef SQRT3
768 #undef SQRT12
769
process_clip(dt_dev_pixelpipe_iop_t * piece,const void * const ivoid,void * const ovoid,const dt_iop_roi_t * const roi_in,const dt_iop_roi_t * const roi_out,const float clip)770 static void process_clip(dt_dev_pixelpipe_iop_t *piece, const void *const ivoid, void *const ovoid,
771 const dt_iop_roi_t *const roi_in, const dt_iop_roi_t *const roi_out,
772 const float clip)
773 {
774 const float *const in = (const float *const)ivoid;
775 float *const out = (float *const)ovoid;
776
777 if(piece->pipe->dsc.filters)
778 { // raw mosaic
779 #ifdef _OPENMP
780 #pragma omp parallel for SIMD() default(none) \
781 dt_omp_firstprivate(clip, in, out, roi_out) \
782 schedule(static)
783 #endif
784 for(size_t k = 0; k < (size_t)roi_out->width * roi_out->height; k++)
785 {
786 out[k] = MIN(clip, in[k]);
787 }
788 }
789 else
790 {
791 const int ch = piece->colors;
792
793 #ifdef _OPENMP
794 #pragma omp parallel for SIMD() default(none) \
795 dt_omp_firstprivate(ch, clip, in, out, roi_out) \
796 schedule(static)
797 #endif
798 for(size_t k = 0; k < (size_t)ch * roi_out->width * roi_out->height; k++)
799 {
800 out[k] = MIN(clip, in[k]);
801 }
802 }
803 }
804
process(struct dt_iop_module_t * self,dt_dev_pixelpipe_iop_t * piece,const void * const ivoid,void * const ovoid,const dt_iop_roi_t * const roi_in,const dt_iop_roi_t * const roi_out)805 void process(struct dt_iop_module_t *self, dt_dev_pixelpipe_iop_t *piece, const void *const ivoid,
806 void *const ovoid, const dt_iop_roi_t *const roi_in, const dt_iop_roi_t *const roi_out)
807 {
808 const uint32_t filters = piece->pipe->dsc.filters;
809 dt_iop_highlights_data_t *data = (dt_iop_highlights_data_t *)piece->data;
810
811 const float clip
812 = data->clip * fminf(piece->pipe->dsc.processed_maximum[0],
813 fminf(piece->pipe->dsc.processed_maximum[1], piece->pipe->dsc.processed_maximum[2]));
814 // const int ch = piece->colors;
815 if(!filters)
816 {
817 process_clip(piece, ivoid, ovoid, roi_in, roi_out, clip);
818 for(int k=0;k<3;k++)
819 piece->pipe->dsc.processed_maximum[k]
820 = fminf(piece->pipe->dsc.processed_maximum[0],
821 fminf(piece->pipe->dsc.processed_maximum[1], piece->pipe->dsc.processed_maximum[2]));
822 return;
823 }
824
825 switch(data->mode)
826 {
827 case DT_IOP_HIGHLIGHTS_INPAINT: // a1ex's (magiclantern) idea of color inpainting:
828 {
829 const float clips[4] = { 0.987 * data->clip * piece->pipe->dsc.processed_maximum[0],
830 0.987 * data->clip * piece->pipe->dsc.processed_maximum[1],
831 0.987 * data->clip * piece->pipe->dsc.processed_maximum[2], clip };
832
833 if(filters == 9u)
834 {
835 const uint8_t(*const xtrans)[6] = (const uint8_t(*const)[6])piece->pipe->dsc.xtrans;
836 #ifdef _OPENMP
837 #pragma omp parallel for default(none) \
838 dt_omp_firstprivate(clips, filters, ivoid, ovoid, roi_in, roi_out, \
839 xtrans) \
840 schedule(static)
841 #endif
842 for(int j = 0; j < roi_out->height; j++)
843 {
844 interpolate_color_xtrans(ivoid, ovoid, roi_in, roi_out, 0, 1, j, clips, xtrans, 0);
845 interpolate_color_xtrans(ivoid, ovoid, roi_in, roi_out, 0, -1, j, clips, xtrans, 1);
846 }
847 #ifdef _OPENMP
848 #pragma omp parallel for default(none) \
849 dt_omp_firstprivate(clips, filters, ivoid, ovoid, roi_in, roi_out, \
850 xtrans) \
851 schedule(static)
852 #endif
853 for(int i = 0; i < roi_out->width; i++)
854 {
855 interpolate_color_xtrans(ivoid, ovoid, roi_in, roi_out, 1, 1, i, clips, xtrans, 2);
856 interpolate_color_xtrans(ivoid, ovoid, roi_in, roi_out, 1, -1, i, clips, xtrans, 3);
857 }
858 }
859 else
860 {
861 #ifdef _OPENMP
862 #pragma omp parallel for default(none) \
863 dt_omp_firstprivate(clips, filters, ivoid, ovoid, roi_out) \
864 shared(data, piece) \
865 schedule(static)
866 #endif
867 for(int j = 0; j < roi_out->height; j++)
868 {
869 interpolate_color(ivoid, ovoid, roi_out, 0, 1, j, clips, filters, 0);
870 interpolate_color(ivoid, ovoid, roi_out, 0, -1, j, clips, filters, 1);
871 }
872
873 // up/down directions
874 #ifdef _OPENMP
875 #pragma omp parallel for default(none) \
876 dt_omp_firstprivate(clips, filters, ivoid, ovoid, roi_out) \
877 shared(data, piece) \
878 schedule(static)
879 #endif
880 for(int i = 0; i < roi_out->width; i++)
881 {
882 interpolate_color(ivoid, ovoid, roi_out, 1, 1, i, clips, filters, 2);
883 interpolate_color(ivoid, ovoid, roi_out, 1, -1, i, clips, filters, 3);
884 }
885 }
886 break;
887 }
888 case DT_IOP_HIGHLIGHTS_LCH:
889 if(filters == 9u)
890 process_lch_xtrans(self, piece, ivoid, ovoid, roi_in, roi_out, clip);
891 else
892 process_lch_bayer(self, piece, ivoid, ovoid, roi_in, roi_out, clip);
893 break;
894 default:
895 case DT_IOP_HIGHLIGHTS_CLIP:
896 process_clip(piece, ivoid, ovoid, roi_in, roi_out, clip);
897 break;
898 }
899
900 // update processed maximum
901 const float m = fmaxf(fmaxf(piece->pipe->dsc.processed_maximum[0], piece->pipe->dsc.processed_maximum[1]),
902 piece->pipe->dsc.processed_maximum[2]);
903 for(int k = 0; k < 3; k++) piece->pipe->dsc.processed_maximum[k] = m;
904
905 if(piece->pipe->mask_display & DT_DEV_PIXELPIPE_DISPLAY_MASK) dt_iop_alpha_copy(ivoid, ovoid, roi_out->width, roi_out->height);
906 }
907
commit_params(struct dt_iop_module_t * self,dt_iop_params_t * p1,dt_dev_pixelpipe_t * pipe,dt_dev_pixelpipe_iop_t * piece)908 void commit_params(struct dt_iop_module_t *self, dt_iop_params_t *p1, dt_dev_pixelpipe_t *pipe,
909 dt_dev_pixelpipe_iop_t *piece)
910 {
911 dt_iop_highlights_params_t *p = (dt_iop_highlights_params_t *)p1;
912 dt_iop_highlights_data_t *d = (dt_iop_highlights_data_t *)piece->data;
913 memcpy(d, p, sizeof(*p));
914
915 piece->process_cl_ready = 1;
916
917 // no OpenCL for DT_IOP_HIGHLIGHTS_INPAINT yet.
918 if(d->mode == DT_IOP_HIGHLIGHTS_INPAINT) piece->process_cl_ready = 0;
919 }
920
init_global(dt_iop_module_so_t * module)921 void init_global(dt_iop_module_so_t *module)
922 {
923 const int program = 2; // basic.cl, from programs.conf
924 dt_iop_highlights_global_data_t *gd
925 = (dt_iop_highlights_global_data_t *)malloc(sizeof(dt_iop_highlights_global_data_t));
926 module->data = gd;
927 gd->kernel_highlights_1f_clip = dt_opencl_create_kernel(program, "highlights_1f_clip");
928 gd->kernel_highlights_1f_lch_bayer = dt_opencl_create_kernel(program, "highlights_1f_lch_bayer");
929 gd->kernel_highlights_1f_lch_xtrans = dt_opencl_create_kernel(program, "highlights_1f_lch_xtrans");
930 gd->kernel_highlights_4f_clip = dt_opencl_create_kernel(program, "highlights_4f_clip");
931 }
932
cleanup_global(dt_iop_module_so_t * module)933 void cleanup_global(dt_iop_module_so_t *module)
934 {
935 dt_iop_highlights_global_data_t *gd = (dt_iop_highlights_global_data_t *)module->data;
936 dt_opencl_free_kernel(gd->kernel_highlights_4f_clip);
937 dt_opencl_free_kernel(gd->kernel_highlights_1f_lch_bayer);
938 dt_opencl_free_kernel(gd->kernel_highlights_1f_lch_xtrans);
939 dt_opencl_free_kernel(gd->kernel_highlights_1f_clip);
940 free(module->data);
941 module->data = NULL;
942 }
943
init_pipe(struct dt_iop_module_t * self,dt_dev_pixelpipe_t * pipe,dt_dev_pixelpipe_iop_t * piece)944 void init_pipe(struct dt_iop_module_t *self, dt_dev_pixelpipe_t *pipe, dt_dev_pixelpipe_iop_t *piece)
945 {
946 piece->data = malloc(sizeof(dt_iop_highlights_data_t));
947 }
948
cleanup_pipe(struct dt_iop_module_t * self,dt_dev_pixelpipe_t * pipe,dt_dev_pixelpipe_iop_t * piece)949 void cleanup_pipe(struct dt_iop_module_t *self, dt_dev_pixelpipe_t *pipe, dt_dev_pixelpipe_iop_t *piece)
950 {
951 free(piece->data);
952 piece->data = NULL;
953 }
954
gui_update(struct dt_iop_module_t * self)955 void gui_update(struct dt_iop_module_t *self)
956 {
957 dt_iop_highlights_gui_data_t *g = (dt_iop_highlights_gui_data_t *)self->gui_data;
958 dt_iop_highlights_params_t *p = (dt_iop_highlights_params_t *)self->params;
959 dt_bauhaus_slider_set(g->clip, p->clip);
960 dt_bauhaus_combobox_set(g->mode, p->mode);
961 }
962
reload_defaults(dt_iop_module_t * module)963 void reload_defaults(dt_iop_module_t *module)
964 {
965 // enable this per default if raw or sraw,
966 module->default_enabled = dt_image_is_rawprepare_supported(&(module->dev->image_storage));
967 }
968
gui_init(struct dt_iop_module_t * self)969 void gui_init(struct dt_iop_module_t *self)
970 {
971 dt_iop_highlights_gui_data_t *g = IOP_GUI_ALLOC(highlights);
972
973 g->mode = dt_bauhaus_combobox_from_params(self, "mode");
974 gtk_widget_set_tooltip_text(g->mode, _("highlight reconstruction method"));
975
976 g->clip = dt_bauhaus_slider_from_params(self, "clip");
977 dt_bauhaus_slider_set_digits(g->clip, 3);
978 gtk_widget_set_tooltip_text(g->clip, _("manually adjust the clipping threshold against "
979 "magenta highlights (you shouldn't ever need to touch this)"));
980 }
981
982 // modelines: These editor modelines have been set for all relevant files by tools/update_modelines.sh
983 // vim: shiftwidth=2 expandtab tabstop=2 cindent
984 // kate: tab-indents: off; indent-width 2; replace-tabs on; indent-mode cstyle; remove-trailing-spaces modified;
985