data/kernels/capture.cl
Rama master · réplica de lectura del repositorio original.
/*
This file is part of darktable,
copyright (c) 2026 darktable developer.
darktable is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
darktable is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with darktable. If not, see <http://www.gnu.org/licenses/>.
*/
#include "common.h"
#define CAPTURE_KERNEL_ALIGN 32
#define CAPTURE_YMIN 0.001f
#define CAPTURE_SMALL_ULIM 66
__kernel void kernel_9x9_mul(global float *in,
global float *out,
global float *blend,
global float *kernels,
global unsigned char *table,
const int w1,
const int height)
{
const int col = get_global_id(0);
const int row = get_global_id(1);
if(col >= w1 || row >= height) return;
const int i = mad24(row, w1, col);
const int w2 = 2 * w1;
const int w3 = 3 * w1;
const int w4 = 4 * w1;
if(blend[i] <= 0.0f) return;
global const float *kern = kernels + CAPTURE_KERNEL_ALIGN * table[i];
global float *d = in + i;
const bool small = table[i] < CAPTURE_SMALL_ULIM;
const int bd = small ? 2 : 4;
float val = 0.0f;
if(col >= bd && row >= bd && (col < w1 - bd) && (row < height - bd))
{
if(small)
{
val =
kern[ 5+2] * (d[-w2-1] + d[-w2+1] + d[-w1-2] + d[-w1+2] + d[w1-2] + d[w1+2] + d[w2-1] + d[w2+1]) +
kern[ 2] * (d[-w2 ] + d[ -2] + d[ 2] + d[ w2 ]) +
kern[ 5+1] * (d[-w1-1] + d[-w1+1] + d[ w1-1] + d[ w1+1]) +
kern[ 1] * (d[-w1 ] + d[ -1] + d[ 1] + d[ w1 ]) +
kern[ 0] * (d[0]);
}
else
{
val =
kern[10+4] * (d[-w4-2] + d[-w4+2] + d[-w2-4] + d[-w2+4] + d[w2-4] + d[w2+4] + d[w4-2] + d[w4+2]) +
kern[5 +4] * (d[-w4-1] + d[-w4+1] + d[-w1-4] + d[-w1+4] + d[w1-4] + d[w1+4] + d[w4-1] + d[w4+1]) +
kern[4] * (d[-w4 ] + d[ -4] + d[ 4] + d[ w4 ]) +
kern[15+3] * (d[-w3-3] + d[-w3+3] + d[ w3-3] + d[ w3+3]) +
kern[10+3] * (d[-w3-2] + d[-w3+2] + d[-w2-3] + d[-w2+3] + d[w2-3] + d[w2+3] + d[w3-2] + d[w3+2]) +
kern[ 5+3] * (d[-w3-1] + d[-w3+1] + d[-w1-3] + d[-w1+3] + d[w1-3] + d[w1+3] + d[w3-1] + d[w3+1]) +
kern[ 3] * (d[-w3 ] + d[ -3] + d[ 3] + d[ w3 ]) +
kern[10+2] * (d[-w2-2] + d[-w2+2] + d[ w2-2] + d[ w2+2]) +
kern[ 5+2] * (d[-w2-1] + d[-w2+1] + d[-w1-2] + d[-w1+2] + d[w1-2] + d[w1+2] + d[w2-1] + d[w2+1]) +
kern[ 2] * (d[-w2 ] + d[ -2] + d[ 2] + d[ w2 ]) +
kern[ 5+1] * (d[-w1-1] + d[-w1+1] + d[ w1-1] + d[ w1+1]) +
kern[ 1] * (d[-w1 ] + d[ -1] + d[ 1] + d[ w1 ]) +
kern[ 0] * (d[0]);
}
}
else
{
for(int ir = -bd; ir <= bd; ir++)
{
const int irow = row+ir;
if(irow >= 0 && irow < height)
{
for(int ic = -bd; ic <= bd; ic++)
{
const int icol = col+ic;
if(icol >=0 && icol < w1)
val += kern[5 * abs(ir) + abs(ic)] * in[mad24(irow, w1, icol)];
}
}
}
}
out[i] *= val;
}
__kernel void kernel_9x9_div(global float *in,
global float *out,
global float *luminance,
global float *blend,
global float *kernels,
global unsigned char *table,
const int w1,
const int height)
{
const int col = get_global_id(0);
const int row = get_global_id(1);
if(col >= w1 || row >= height) return;
const int i = mad24(row, w1, col);
const int w2 = 2 * w1;
const int w3 = 3 * w1;
const int w4 = 4 * w1;
if(blend[i] <= 0.0f) return;
global const float *kern = kernels + CAPTURE_KERNEL_ALIGN * table[i];
global float *d = in + i;
const bool small = table[i] < CAPTURE_SMALL_ULIM;
const int bd = small ? 2 : 4;
float val = 0.0f;
if(col >= bd && row >= bd && (col < w1 - bd) && (row < height - bd))
{
if(small)
{
val =
kern[ 5+2] * (d[-w2-1] + d[-w2+1] + d[-w1-2] + d[-w1+2] + d[w1-2] + d[w1+2] + d[w2-1] + d[w2+1]) +
kern[ 2] * (d[-w2 ] + d[ -2] + d[ 2] + d[ w2 ]) +
kern[ 5+1] * (d[-w1-1] + d[-w1+1] + d[ w1-1] + d[ w1+1]) +
kern[ 1] * (d[-w1 ] + d[ -1] + d[ 1] + d[ w1 ]) +
kern[ 0] * (d[0]);
}
else
{
val =
kern[10+4] * (d[-w4-2] + d[-w4+2] + d[-w2-4] + d[-w2+4] + d[w2-4] + d[w2+4] + d[w4-2] + d[w4+2]) +
kern[5 +4] * (d[-w4-1] + d[-w4+1] + d[-w1-4] + d[-w1+4] + d[w1-4] + d[w1+4] + d[w4-1] + d[w4+1]) +
kern[4] * (d[-w4 ] + d[ -4] + d[ 4] + d[ w4 ]) +
kern[15+3] * (d[-w3-3] + d[-w3+3] + d[ w3-3] + d[ w3+3]) +
kern[10+3] * (d[-w3-2] + d[-w3+2] + d[-w2-3] + d[-w2+3] + d[w2-3] + d[w2+3] + d[w3-2] + d[w3+2]) +
kern[ 5+3] * (d[-w3-1] + d[-w3+1] + d[-w1-3] + d[-w1+3] + d[w1-3] + d[w1+3] + d[w3-1] + d[w3+1]) +
kern[ 3] * (d[-w3 ] + d[ -3] + d[ 3] + d[ w3 ]) +
kern[10+2] * (d[-w2-2] + d[-w2+2] + d[ w2-2] + d[ w2+2]) +
kern[ 5+2] * (d[-w2-1] + d[-w2+1] + d[-w1-2] + d[-w1+2] + d[w1-2] + d[w1+2] + d[w2-1] + d[w2+1]) +
kern[ 2] * (d[-w2 ] + d[ -2] + d[ 2] + d[ w2 ]) +
kern[ 5+1] * (d[-w1-1] + d[-w1+1] + d[ w1-1] + d[ w1+1]) +
kern[ 1] * (d[-w1 ] + d[ -1] + d[ 1] + d[ w1 ]) +
kern[ 0] * (d[0]);
}
}
else
{
for(int ir = -bd; ir <= bd; ir++)
{
const int irow = row+ir;
if(irow >= 0 && irow < height)
{
for(int ic = -bd; ic <= bd; ic++)
{
const int icol = col+ic;
if(icol >=0 && icol < w1)
val += kern[5 * abs(ir) + abs(ic)] * in[mad24(irow, w1, icol)];
}
}
}
}
out[i] = fp32_div(luminance[i], fmax(val, CAPTURE_YMIN));
}
static float _mask_dilated(global float *in, const int w1)
{
const int w2 = 2*w1;
return ( in[-w2-1] || in[-w2] || in[-w2+1]
|| in[-w1-2] || in[-w1-1] || in[-w1] || in[-w1+1] || in[-w1+2]
|| in[-2] || in[-1] || in[0] || in[1] || in[2]
|| in[w1-2] || in[w1-1] || in[w1] || in[w1+1] || in[w1+2]
|| in[w2-1] || in[w2] || in[w2+1] )
? 1.0f : 0.0f;
}
kernel void capture_dilatemask(global float *in,
global float *out,
const int width,
const int height)
{
const int col = get_global_id(0);
const int row = get_global_id(1);
if(col >= width || row >= height)
return;
int mx = mad24(row, width, col);
bool safe = col >= 2 && row >= 2 && col < width - 2 && row < height - 2;
out[mx] = safe ? _mask_dilated(in + mx, width) : in[mx];
}
__kernel void prepare_blend(__read_only image2d_t cfa,
__read_only image2d_t dev_out,
const int filters,
global const unsigned char (*const xtrans)[6],
global float *mask,
global float *Yold,
const float4 wb,
const int w,
const int height)
{
const int col = get_global_id(0);
const int row = get_global_id(1);
if(col >= w || row >= height) return;
float whites[4] = { wb.x, wb.y, wb.z, wb.w };
float4 rgb = Areadpixel(dev_out, col, row);
// Photometric/digital ITU BT.709
rgb *= (float4)( 0.212671f, 0.715160f, 0.072169f, 0.0f );
const float Y = clipf(rgb.x + rgb.y + rgb.z);
const int k = mad24(row, w, col);
Yold[k] = Y;
const int color = fcol(row, col, filters, xtrans);
const float val = Areadsingle(cfa, col, row);
mask[k] = val > whites[color] || Y < CAPTURE_YMIN ? 0.0f : 1.0f;
}
__kernel void modify_blend(global float *blend,
global float *Yold,
global float *luminance,
const float dthresh,
const int width,
const int height)
{
const int icol = get_global_id(0);
const int irow = get_global_id(1);
if(icol >= width || irow >= height) return;
const int row = clamp(irow, 2, height-3);
const int col = clamp(icol, 2, width-3);
const float threshold = 0.6f * fsquare(dthresh);
const float tscale = 200.0f;
const float offset = -2.5f + tscale * threshold / 2.0f;
float sum = 0.0f;
float sum_sq = 0.0f;
for(int y = row-1; y < row+2; y++)
{
for(int x = col-2; x < col+3; x++)
{
sum += Yold[mad24(y, width, x)];
sum_sq += fsquare(Yold[mad24(y, width, x)]);
}
}
for(int x = col-1; x < col+2; x++)
{
sum += Yold[mad24(row-2, width, x)];
sum_sq += fsquare(Yold[mad24(row-2, width, x)]);
sum += Yold[mad24(row+2, width, x)];
sum_sq += fsquare(Yold[mad24(row+2, width, x)]);
}
const int k = mad24(irow, width, icol);
const float sum_of_squares = fmax(0.0f, sum_sq - fp32_div(fsquare(sum), 21.0f));
const float std_deviation = sqrt(fp32_div(sum_of_squares, 21.0f));
const float modified_coef_variation = fp32_div(std_deviation, sqrt(fmax(NORM_MIN, fp32_div(sum, 21.0f))));
const float t = dtcl_log(1.0f + modified_coef_variation);
const float weight = fp32_div(1.0f, (1.0f + exp(offset - tscale * t)));
blend[k] = clipf(blend[k] * 1.01011f * (weight - 0.01f));
luminance[k] = Yold[k];
}
__kernel void final_blend(global float *blendmask,
global float *unblurred,
int pixels)
{
const int k = get_global_id(0);
if(k >= pixels) return;
const float diff = unblurred[k] - blendmask[k];
const float w_tmp2 = fp32_div(1.0f, (1.0f + exp(5.0f - 10.0f * diff)));
blendmask[k] = clipf(w_tmp2 * unblurred[k] + (1.0f - w_tmp2) * blendmask[k]);
}
__kernel void show_blend_mask(__read_only image2d_t in,
__write_only image2d_t out,
global float *blend_mask,
global unsigned char *sigma_mask,
const int width,
const int height,
const int blender)
{
const int col = get_global_id(0);
const int row = get_global_id(1);
if(col >= width || row >= height) return;
float4 pix = fmax(0.0f, Areadpixel(in, col, row));
const float blend = blender ? blend_mask[mad24(row, width, col)]
: (float)sigma_mask[mad24(row, width, col)] / 255.0f;
pix.w = blend;
write_imagef(out, (int2)(col, row), pix);
}
__kernel void capture_result( __read_only image2d_t in,
__write_only image2d_t out,
global float *blendmask,
global float *luminance,
global float *tmp,
const int width,
const int height)
{
const int col = get_global_id(0);
const int row = get_global_id(1);
if(col >= width || row >= height) return;
float4 pix = Areadpixel(in, col, row);
const int k = mad24(row, width, col);
if(blendmask[k] > 0.0f)
{
const float mixer = clipf(blendmask[k]);
const float luminance_new = mix(luminance[k], tmp[k], mixer);
const float lumcorr = fp32_div(luminance_new, fmax(luminance[k], CAPTURE_YMIN));
pix *= (float4)lumcorr;
}
write_imagef(out, (int2)(col, row), pix);
}
// auto radius: sensor value at (row, col) of the tested area, white balance undone
static inline float _radius_val(read_only image2d_t in,
global const unsigned char (*const xtrans)[6],
const unsigned int filters,
const float4 coeff,
const int dx,
const int dy,
const int row,
const int col)
{
const int c = fcol(row, col, filters, xtrans);
const float f = c == 0 ? coeff.x : (c == 1 ? coeff.y : (c == 2 ? coeff.z : coeff.w));
return f * readsingle(in, col + dx, row + dy);
}
#define RAWEPS 0.005f
#define lowerLimit 0.01f
#define upperLimit 0.9f
#define V(r, c) _radius_val(in, xtrans, filters, coeff, dx, dy, row + (r), col + (c))
// largest ratio of neighbouring greens at one site, same tests as _calcRadiusBayer()
static inline float _radius_bayer(read_only image2d_t in,
global const unsigned char (*const xtrans)[6],
const unsigned int filters,
const float4 coeff,
const int dx,
const int dy,
const int row,
const int col)
{
float maxRatio = 1.0f;
const float val00 = V(0, 0);
if(val00 <= RAWEPS) return maxRatio;
const float val1m1 = V(1, -1);
const float val1p1 = V(1, 1);
const float maxVal0 = fmax(val00, val1m1);
if(val1m1 > RAWEPS && maxVal0 > lowerLimit)
{
const float minVal = fmin(val00, val1m1);
if(maxVal0 > maxRatio * minVal)
{
const float clip = maxVal0 == val00
? fmax(fmax(V(-1, -1), V(-1, 1)), val1p1)
: fmax(fmax(fmax(V(0, -2), val00), V(2, -2)), V(2, 0));
if(clip < upperLimit)
maxRatio = maxVal0 / minVal;
}
}
const float maxVal1 = fmax(val00, val1p1);
if(val1p1 > RAWEPS && maxVal1 > lowerLimit)
{
const float minVal = fmin(val00, val1p1);
if(maxVal1 > maxRatio * minVal)
{
const float clip = maxVal1 == val00
? fmax(fmax(V(-1, -1), V(-1, 1)), val1p1)
: fmax(fmax(fmax(val00, V(0, 2)), V(2, 0)), V(2, 2));
if(clip < upperLimit)
maxRatio = maxVal1 / minVal;
}
}
return maxRatio;
}
static inline float _radius_pair(const float a,
const float b,
const float maxRatio)
{
const float maxVal = fmax(a, b);
const float minVal = fmin(a, b);
return maxVal > lowerLimit && maxVal > maxRatio * minVal ? maxVal / minVal : maxRatio;
}
// same tests as _calcRadiusXtrans()
static inline float _radius_xtrans(read_only image2d_t in,
global const unsigned char (*const xtrans)[6],
const unsigned int filters,
const float4 coeff,
const int dx,
const int dy,
const int row,
const int col)
{
float maxRatio = 1.0f;
const float valp1p1 = V(1, 1);
const bool squareClipped = fmax(fmax(fmax(valp1p1, V(1, 2)), V(2, 1)), V(2, 2)) >= upperLimit;
const float greenSolitary = V(0, 0);
if(greenSolitary > RAWEPS && fmax(V(-1, -1), V(-1, 1)) < upperLimit)
{
if(greenSolitary < upperLimit)
{
const float valp1m1 = V(1, -1);
if(valp1m1 > RAWEPS && fmax(fmax(fmax(V(1, -2), valp1m1), V(2, -2)), V(1, -1)) < upperLimit)
maxRatio = _radius_pair(greenSolitary, valp1m1, maxRatio);
if(valp1p1 > RAWEPS && !squareClipped)
maxRatio = _radius_pair(greenSolitary, valp1p1, maxRatio);
}
}
if(!squareClipped)
{
const float valp2p2 = V(2, 2);
if(valp2p2 > RAWEPS)
{
if(valp1p1 > RAWEPS)
maxRatio = _radius_pair(valp1p1, valp2p2, maxRatio);
const float greenSolitaryRight = V(3, 3);
if(fmax(fmax(greenSolitaryRight, V(4, 2)), V(4, 4)) < upperLimit && greenSolitaryRight > RAWEPS)
maxRatio = _radius_pair(greenSolitaryRight, valp2p2, maxRatio);
}
const float valp1p2 = V(1, 2);
const float valp2p1 = V(2, 1);
if(valp2p1 > RAWEPS)
{
if(valp1p2 > RAWEPS)
maxRatio = _radius_pair(valp1p2, valp2p1, maxRatio);
const float greenSolitaryLeft = V(3, 0);
if(fmax(fmax(greenSolitaryLeft, V(4, -1)), V(4, 1)) < upperLimit && greenSolitaryLeft > RAWEPS)
maxRatio = _radius_pair(greenSolitaryLeft, valp2p1, maxRatio);
}
}
return maxRatio;
}
#undef V
// one work item per tested site, each work group writes the largest ratio it found
__kernel void capture_radius(read_only image2d_t in,
global const unsigned char (*const xtrans)[6],
const unsigned int filters,
const float4 coeff,
const int dx,
const int dy,
const int width,
const int height,
const int startx,
const int starty,
global float *partial,
local float *buffer)
{
const int x = get_global_id(0);
const int y = get_global_id(1);
const int xlsz = get_local_size(0);
const int ylsz = get_local_size(1);
const int xlid = get_local_id(0);
const int ylid = get_local_id(1);
const int l = mad24(ylid, xlsz, xlid);
float maxRatio = 1.0f;
if(filters == 9u)
{
const int row = starty + 2 + 3 * y;
const int col = startx + 2 + 3 * x;
if(row < height - 4 && col < width - 4)
maxRatio = _radius_xtrans(in, xtrans, filters, coeff, dx, dy, row, col);
}
else
{
const int row = 4 + y;
const int col = (filters ? 5 + (FC(row & 1, 0, filters) & 1) : 5) + 2 * x;
if(row < height - 4 && col < width - 4)
maxRatio = _radius_bayer(in, xtrans, filters, coeff, dx, dy, row, col);
}
buffer[l] = maxRatio;
barrier(CLK_LOCAL_MEM_FENCE);
for(int offset = mul24(xlsz, ylsz) / 2; offset > 0; offset = offset / 2)
{
if(l < offset)
buffer[l] = fmax(buffer[l], buffer[l + offset]);
barrier(CLK_LOCAL_MEM_FENCE);
}
if(l == 0)
{
const int xgid = get_group_id(0);
const int ygid = get_group_id(1);
const int xgsz = get_num_groups(0);
partial[mad24(ygid, xgsz, xgid)] = buffer[0];
}
}
#undef RAWEPS
#undef lowerLimit
#undef upperLimit