From 00fc90b4b9e8e8834a6ca50597dccbc8eb407d62 Mon Sep 17 00:00:00 2001 From: piratenpanda Date: Wed, 28 Oct 2026 09:57:00 +0100 Subject: [PATCH] spektrafilm: tile the Young-van Vliet row pass through local memory The row kernels ran one work-item per row, each walking its own row, so neighbouring items accessed memory one row apart, and the row pass was far slower than the column pass over the same data. A work-group now covers several rows and walks them in chunks of 64 columns: each chunk is loaded and stored by the whole group in contiguous segments through local memory, and one item per row runs the unchanged recursion on the chunk. The work-group is 8 x 16 items, or smaller as dt_opencl_local_buffer_opt() finds the device's work-group and local memory limits require. Output is bit-identical for every tile size. Row pass on 6000x4000, median of 10, built with -cl-fast-relaxed-math, before -> after: device / driver float4 float Apple M1 Pro, rusticl (asahi) 100.3 -> 22.1 ms 13.8 -> 6.0 ms RX 6700 XT, rusticl (radeonsi) 8.8 -> 6.7 ms 5.4 -> 1.5 ms RX 6700 XT, ROCm 9.1 -> 6.3 ms 5.4 -> 1.7 ms --- data/kernels/spektrafilm.cl | 139 ++++++++++++++++++++++++++---------- src/iop/spektrafilm.c | 19 ++++- 2 files changed, 117 insertions(+), 41 deletions(-) diff --git a/data/kernels/spektrafilm.cl b/data/kernels/spektrafilm.cl index 7cbad79ec5..d5f0a80a41 100644 --- a/data/kernels/spektrafilm.cl +++ b/data/kernels/spektrafilm.cl @@ -883,8 +883,8 @@ __kernel void spektrafilm_passthrough(__read_only image2d_t in, __write_only ima /* Young-van Vliet order-3 recursive Gaussian, one work-item per line. Same coefficients and same edge-replicated seeding as sf_gauss_yvv_coeffs / _sf_gauss_iir_1d on the CPU, so both paths deliver the identical blur above - SF_GAUSS_EXACT_MAX_SIGMA. Launch with global size (h, 1) for rows and - (w, 1) for columns. */ + SF_GAUSS_EXACT_MAX_SIGMA. Columns launch with global size (w, 1); rows are + tiled, see spektrafilm_yvv_row_4c. */ /* Incremental form and explicit fma(), in the same nesting as _sf_gauss_iir_1d() in spektra_core.c -- see the comment there: v = w1 + (B*(x - w1) + B2*(w2 - w1) + B3*(w3 - w1)), which uses @@ -893,30 +893,66 @@ __kernel void spektrafilm_passthrough(__read_only image2d_t in, __write_only ima the result back in, so one differing rounding between host and device becomes a different filter rather than a one-ulp difference; FP_CONTRACT pragmas do not settle this on either side, an explicit fma does. */ +/* Rows: one work-group of tw x th items (its local size) covers th rows. The + rows are walked in chunks of YVV_CHUNK columns; each chunk is loaded and + stored by the whole group in contiguous segments through `tile`, local + memory of th * YVV_TILE_PITCH elements, and item 0 of each row runs the + serial recursion over the chunk in local memory. Launch with global size + (tw, h rounded up to th). */ +#define YVV_CHUNK 64 +#define YVV_TILE_PITCH (YVV_CHUNK + 1) + __kernel void spektrafilm_yvv_row_4c(__global const float4 *src, __global float4 *dst, const int w, const int h, const float B, const float B2, - const float B3) + const float B3, __local float4 *tile) { - /* One work-item per line. The launch asks for a second global dimension of - 1, but dt_opencl_enqueue_kernel_2d_args rounds every dimension up to the - device's preferred multiple, so this is entered by a whole row of items - per line -- all of which would otherwise run the same serial recursion - over the same addresses and race. Only item 0 may proceed. */ - const int row = get_global_id(0); - if(row >= h || get_global_id(1) != 0) return; - __global const float4 *s = src + (size_t)row * w; - __global float4 *d = dst + (size_t)row * w; - float4 w1 = s[0], w2 = s[0], w3 = s[0]; - for(int j = 0; j < w; j++) + const int tw = get_local_size(0); + const int lx = get_local_id(0); + const int ly = get_local_id(1); + const int row = get_global_id(1); + const int valid = row < h; + __global const float4 *s = src + (size_t)min(row, h - 1) * w; + __global float4 *d = dst + (size_t)min(row, h - 1) * w; + __local float4 *t = tile + ly * YVV_TILE_PITCH; + const int last = ((w - 1) / YVV_CHUNK) * YVV_CHUNK; + + float4 w1 = s[0], w2 = w1, w3 = w1; + for(int c0 = 0; c0 < w; c0 += YVV_CHUNK) { - const float4 v = w1 + fma(B, s[j] - w1, fma(B2, w2 - w1, B3 * (w3 - w1))); - d[j] = v; w3 = w2; w2 = w1; w1 = v; + const int n = min(YVV_CHUNK, w - c0); + for(int j = lx; valid && j < n; j += tw) t[j] = s[c0 + j]; + barrier(CLK_LOCAL_MEM_FENCE); + if(valid && lx == 0) + { + for(int j = 0; j < n; j++) + { + const float4 v = w1 + fma(B, t[j] - w1, fma(B2, w2 - w1, B3 * (w3 - w1))); + t[j] = v; w3 = w2; w2 = w1; w1 = v; + } + } + barrier(CLK_LOCAL_MEM_FENCE); + for(int j = lx; valid && j < n; j += tw) d[c0 + j] = t[j]; + barrier(CLK_LOCAL_MEM_FENCE); } - float4 v1 = d[w - 1], v2 = v1, v3 = v1; - for(int j = w - 1; j >= 0; j--) + + /* w1 holds d[w - 1] */ + float4 v1 = w1, v2 = v1, v3 = v1; + for(int c0 = last; c0 >= 0; c0 -= YVV_CHUNK) { - const float4 v = v1 + fma(B, d[j] - v1, fma(B2, v2 - v1, B3 * (v3 - v1))); - d[j] = v; v3 = v2; v2 = v1; v1 = v; + const int n = min(YVV_CHUNK, w - c0); + for(int j = lx; valid && j < n; j += tw) t[j] = d[c0 + j]; + barrier(CLK_LOCAL_MEM_FENCE); + if(valid && lx == 0) + { + for(int j = n - 1; j >= 0; j--) + { + const float4 v = v1 + fma(B, t[j] - v1, fma(B2, v2 - v1, B3 * (v3 - v1))); + t[j] = v; v3 = v2; v2 = v1; v1 = v; + } + } + barrier(CLK_LOCAL_MEM_FENCE); + for(int j = lx; valid && j < n; j += tw) d[c0 + j] = t[j]; + barrier(CLK_LOCAL_MEM_FENCE); } } @@ -949,28 +985,55 @@ __kernel void spektrafilm_yvv_col_4c(__global const float4 *src, __global float4 __kernel void spektrafilm_yvv_row_1c(__global const float *src, __global float *dst, const int w, const int h, const float B, const float B2, - const float B3) + const float B3, __local float *tile) { - /* One work-item per line. The launch asks for a second global dimension of - 1, but dt_opencl_enqueue_kernel_2d_args rounds every dimension up to the - device's preferred multiple, so this is entered by a whole row of items - per line -- all of which would otherwise run the same serial recursion - over the same addresses and race. Only item 0 may proceed. */ - const int row = get_global_id(0); - if(row >= h || get_global_id(1) != 0) return; - __global const float *s = src + (size_t)row * w; - __global float *d = dst + (size_t)row * w; - float w1 = s[0], w2 = s[0], w3 = s[0]; - for(int j = 0; j < w; j++) + const int tw = get_local_size(0); + const int lx = get_local_id(0); + const int ly = get_local_id(1); + const int row = get_global_id(1); + const int valid = row < h; + __global const float *s = src + (size_t)min(row, h - 1) * w; + __global float *d = dst + (size_t)min(row, h - 1) * w; + __local float *t = tile + ly * YVV_TILE_PITCH; + const int last = ((w - 1) / YVV_CHUNK) * YVV_CHUNK; + + float w1 = s[0], w2 = w1, w3 = w1; + for(int c0 = 0; c0 < w; c0 += YVV_CHUNK) { - const float v = w1 + fma(B, s[j] - w1, fma(B2, w2 - w1, B3 * (w3 - w1))); - d[j] = v; w3 = w2; w2 = w1; w1 = v; + const int n = min(YVV_CHUNK, w - c0); + for(int j = lx; valid && j < n; j += tw) t[j] = s[c0 + j]; + barrier(CLK_LOCAL_MEM_FENCE); + if(valid && lx == 0) + { + for(int j = 0; j < n; j++) + { + const float v = w1 + fma(B, t[j] - w1, fma(B2, w2 - w1, B3 * (w3 - w1))); + t[j] = v; w3 = w2; w2 = w1; w1 = v; + } + } + barrier(CLK_LOCAL_MEM_FENCE); + for(int j = lx; valid && j < n; j += tw) d[c0 + j] = t[j]; + barrier(CLK_LOCAL_MEM_FENCE); } - float v1 = d[w - 1], v2 = v1, v3 = v1; - for(int j = w - 1; j >= 0; j--) + + /* w1 holds d[w - 1] */ + float v1 = w1, v2 = v1, v3 = v1; + for(int c0 = last; c0 >= 0; c0 -= YVV_CHUNK) { - const float v = v1 + fma(B, d[j] - v1, fma(B2, v2 - v1, B3 * (v3 - v1))); - d[j] = v; v3 = v2; v2 = v1; v1 = v; + const int n = min(YVV_CHUNK, w - c0); + for(int j = lx; valid && j < n; j += tw) t[j] = d[c0 + j]; + barrier(CLK_LOCAL_MEM_FENCE); + if(valid && lx == 0) + { + for(int j = n - 1; j >= 0; j--) + { + const float v = v1 + fma(B, t[j] - v1, fma(B2, v2 - v1, B3 * (v3 - v1))); + t[j] = v; v3 = v2; v2 = v1; v1 = v; + } + } + barrier(CLK_LOCAL_MEM_FENCE); + for(int j = lx; valid && j < n; j += tw) d[c0 + j] = t[j]; + barrier(CLK_LOCAL_MEM_FENCE); } } diff --git a/src/iop/spektrafilm.c b/src/iop/spektrafilm.c index dde517abf8..daeef4c1cf 100644 --- a/src/iop/spektrafilm.c +++ b/src/iop/spektrafilm.c @@ -2590,12 +2590,25 @@ static cl_int _sf_yvv_blur_cl(const int devid, sf_gauss_yvv_coeffs(pass_sigma, b); const int krow = (ch == 4) ? gd->kernel_yvv_row_4c : gd->kernel_yvv_row_1c; const int kcol = (ch == 4) ? gd->kernel_yvv_col_4c : gd->kernel_yvv_col_1c; + /* the row kernels work on tiles of sizex x sizey items, up to 8 x 16, and + keep sizey rows of YVV_TILE_PITCH (65) elements in local memory */ + dt_opencl_local_buffer_t locopt + = (dt_opencl_local_buffer_t){ .xoffset = 65, .xfactor = 0, + .yoffset = 0, .yfactor = 1, + .cellsize = ch * sizeof(float), .overhead = 0, + .sizex = 8, .sizey = 16 }; + if(dt_opencl_local_buffer_opt(devid, krow, &locopt) != CL_SUCCESS) + return DT_OPENCL_DEFAULT_ERROR; + const size_t row_sizes[2] = { locopt.sizex, ROUNDUP(h, locopt.sizey) }; + const size_t row_local[2] = { locopt.sizex, locopt.sizey }; + const size_t tile_bytes = (size_t)locopt.sizey * 65 * ch * sizeof(float); cl_int e = CL_SUCCESS; for(int pass = 0; pass < passes && e == CL_SUCCESS; pass++) { - e = dt_opencl_enqueue_kernel_2d_args(devid, krow, h, 1, CLARG(buf), CLARG(tmp), - CLARG(w), CLARG(h), CLARG(b[0]), CLARG(b[2]), - CLARG(b[3])); + e = dt_opencl_enqueue_kernel_2d_local_args(devid, krow, row_sizes, row_local, + CLARG(buf), CLARG(tmp), + CLARG(w), CLARG(h), CLARG(b[0]), CLARG(b[2]), + CLARG(b[3]), CLLOCAL(tile_bytes)); if(e == CL_SUCCESS) e = dt_opencl_enqueue_kernel_2d_args(devid, kcol, w, 1, CLARG(tmp), CLARG(buf), CLARG(w), CLARG(h), CLARG(b[0]), CLARG(b[2]),