From e2ff2b7499b77c9c7a310814f6265cbe7159f1bf Mon Sep 17 00:00:00 2001 From: johnnynunez Date: Sun, 28 Jun 2026 18:14:27 +0200 Subject: [PATCH] feat(image-processing): rebuild box filter as FKL DPP --- .../image_processing/box_filter_fast.h | 293 ++++++++++++++++++ .../image_processing/image_processing.h | 1 + .../test_box_filter_quad_dpp.h | 159 ++++++++++ 3 files changed, 453 insertions(+) create mode 100644 include/fused_kernel/algorithms/image_processing/box_filter_fast.h create mode 100644 tests/image_processing/test_box_filter_quad_dpp.h diff --git a/include/fused_kernel/algorithms/image_processing/box_filter_fast.h b/include/fused_kernel/algorithms/image_processing/box_filter_fast.h new file mode 100644 index 00000000..cd433458 --- /dev/null +++ b/include/fused_kernel/algorithms/image_processing/box_filter_fast.h @@ -0,0 +1,293 @@ +/* Copyright 2025-2026 Oscar Amoros Huguet, Johnny Nunez + + Licensed under the Apache License, Version 2.0 (the "License"); + you may not use this file except in compliance with the License. + You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + + Unless required by applicable law or agreed to in writing, software + distributed under the License is distributed on an "AS IS" BASIS, + WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + See the License for the specific language governing permissions and + limitations under the License. */ + +#ifndef FK_BOX_FILTER_FAST_H +#define FK_BOX_FILTER_FAST_H + +#include +#include +#include +#include +#include + +#include + +namespace fk { + +struct BoxFilterQuadDetails { + int width; + int height; + int kernelWidth; + int kernelHeight; + int anchorX; + int anchorY; + + template + FK_HOST_DEVICE_CNST bool validFor( + const int maxRuntimeKernelWidth) const { + const int effectiveWidth = + STATIC_KW > 0 ? STATIC_KW : kernelWidth; + const int effectiveHeight = + STATIC_KH > 0 ? STATIC_KH : kernelHeight; + return width > 0 && height > 0 && + effectiveWidth > 0 && effectiveHeight > 0 && + (STATIC_KW > 0 || effectiveWidth <= maxRuntimeKernelWidth) && + anchorX >= 0 && anchorX < effectiveWidth && + anchorY >= 0 && anchorY < effectiveHeight; + } +}; + +template +struct BoxFilterQuadDPP; + +template +struct BoxFilterQuadDPP { +private: + using SelfType = BoxFilterQuadDPP; + + FK_HOST_FUSE int clamp(const int value, const int upper) { + return value < 0 ? 0 : (value >= upper ? upper - 1 : value); + } + +public: + FK_STATIC_STRUCT(BoxFilterQuadDPP, SelfType) + static constexpr ParArch PAR_ARCH = ParArch::CPU; + static constexpr int MAX_RUNTIME_KERNEL_WIDTH = 31; + + FK_HOST_DEVICE_FUSE bool accepts(const BoxFilterQuadDetails& details) { + return details.template validFor( + MAX_RUNTIME_KERNEL_WIDTH); + } + + template + FK_HOST_FUSE void exec(const BoxFilterQuadDetails& details, + const InIOp& input, + const ComputeIOps& compute, + const OutIOp& output) { + static_assert(isAnyCompleteReadType, + "BoxFilterQuadDPP requires a complete Read IOp"); + static_assert(isAnyWriteType, + "BoxFilterQuadDPP requires a Write IOp"); + static_assert(std::is_same_v< + typename InIOp::Operation::OutputType, float>, + "BoxFilterQuadDPP Read IOp must produce float"); + static_assert(std::is_same_v< + typename OutIOp::Operation::InputType, float>, + "BoxFilterQuadDPP Write IOp must consume float"); + static_assert(EX > 0 && EY > 0, + "BoxFilterQuadDPP output tile must be positive"); + static_assert(KW >= 0 && KH >= 0, + "BoxFilterQuadDPP static kernel sizes cannot be negative"); + static_assert((KW == 0) == (KH == 0), + "BoxFilterQuadDPP kernel dimensions are both static or both runtime"); + static_assert(cn == 1, + "BoxFilterQuadDPP currently supports scalar pixel types"); + if (!accepts(details)) return; + + const int kernelWidth = KW > 0 ? KW : details.kernelWidth; + const int kernelHeight = KH > 0 ? KH : details.kernelHeight; + const auto& combine = get<0>(compute); + const auto& normalize = get<2>(compute); + + for (int y = 0; y < details.height; ++y) { + for (int x = 0; x < details.width; ++x) { + float sum = 0.f; + for (int ky = 0; ky < kernelHeight; ++ky) { + for (int kx = 0; kx < kernelWidth; ++kx) { + const Point source{ + clamp(x + kx - details.anchorX, details.width), + clamp(y + ky - details.anchorY, details.height), 0}; + const float value = + InIOp::Operation::exec(source, input); + sum = make_tuple(sum, value) | combine; + } + } + const float mean = sum | normalize; + OutIOp::Operation::exec(Point{x, y, 0}, mean, output); + } + } + } +}; + +#if defined(__NVCC__) +template +struct BoxFilterQuadDPP { +private: + using SelfType = BoxFilterQuadDPP< + ParArch::GPU_NVIDIA, T, EX, EY, KW, KH>; + static constexpr int MAX_RUNTIME_KERNEL_WIDTH = 31; + static constexpr int MAX_SPAN = + EX + (KW > 0 ? KW : MAX_RUNTIME_KERNEL_WIDTH) - 1; + +public: + FK_STATIC_STRUCT(BoxFilterQuadDPP, SelfType) + static constexpr ParArch PAR_ARCH = ParArch::GPU_NVIDIA; + static constexpr int BLOCK_THREADS = 256; + + FK_HOST_DEVICE_FUSE bool accepts(const BoxFilterQuadDetails& details) { + return details.template validFor( + MAX_RUNTIME_KERNEL_WIDTH); + } + + struct LaunchConfig { + unsigned int blocks; + unsigned int threads; + }; + + FK_HOST_FUSE LaunchConfig launchConfig( + const BoxFilterQuadDetails& details) { + const int columns = (details.width + EX - 1) / EX; + const int rows = (details.height + EY - 1) / EY; + const int workItems = columns * rows; + return {static_cast( + (workItems + BLOCK_THREADS - 1) / BLOCK_THREADS), + static_cast(BLOCK_THREADS)}; + } + + template + FK_DEVICE_FUSE void exec(const BoxFilterQuadDetails& details, + const InIOp& input, + const ComputeIOps& compute, + const OutIOp& output) { + static_assert(isAnyCompleteReadType, + "BoxFilterQuadDPP requires a complete Read IOp"); + static_assert(isAnyWriteType, + "BoxFilterQuadDPP requires a Write IOp"); + static_assert(std::is_same_v< + typename InIOp::Operation::OutputType, float>, + "BoxFilterQuadDPP Read IOp must produce float"); + static_assert(std::is_same_v< + typename OutIOp::Operation::InputType, float>, + "BoxFilterQuadDPP Write IOp must consume float"); + static_assert(EX > 0 && EY > 0, + "BoxFilterQuadDPP output tile must be positive"); + static_assert(KW >= 0 && KH >= 0, + "BoxFilterQuadDPP static kernel sizes cannot be negative"); + static_assert((KW == 0) == (KH == 0), + "BoxFilterQuadDPP kernel dimensions are both static or both runtime"); + static_assert(cn == 1, + "BoxFilterQuadDPP currently supports scalar pixel types"); +#if defined(__CUDA_ARCH__) + if (!accepts(details)) return; + + const int kernelWidth = KW > 0 ? KW : details.kernelWidth; + const int kernelHeight = KH > 0 ? KH : details.kernelHeight; + const auto& combine = get<0>(compute); + const auto& subtract = get<1>(compute); + const auto& normalize = get<2>(compute); + + const int columns = (details.width + EX - 1) / EX; + const int workItem = static_cast( + blockIdx.x * blockDim.x + threadIdx.x); + const int x0 = (workItem % columns) * EX; + const int y0 = (workItem / columns) * EY; + if (x0 >= details.width || y0 >= details.height) return; + + const int span = EX + kernelWidth - 1; + if (span > MAX_SPAN) return; + const int firstColumn = x0 - details.anchorX; + + auto clamp = [](const int value, const int upper) { + return value < 0 ? 0 : (value >= upper ? upper - 1 : value); + }; + auto readSource = [&](const int x, const int y) { + return InIOp::Operation::exec( + Point{clamp(x, details.width), + clamp(y, details.height), 0}, input); + }; + + float columnsSum[MAX_SPAN]; + #pragma unroll + for (int index = 0; index < span; ++index) { + float sum = 0.f; + for (int ky = 0; ky < kernelHeight; ++ky) { + sum = make_tuple( + sum, readSource(firstColumn + index, + y0 + ky - details.anchorY)) | combine; + } + columnsSum[index] = sum; + } + + #pragma unroll + for (int localY = 0; localY < EY; ++localY) { + const int y = y0 + localY; + if (y >= details.height) break; + #pragma unroll + for (int localX = 0; localX < EX; ++localX) { + const int x = x0 + localX; + if (x >= details.width) break; + float sum = 0.f; + for (int kx = 0; kx < kernelWidth; ++kx) { + sum = make_tuple( + sum, columnsSum[localX + kx]) | combine; + } + const float mean = sum | normalize; + OutIOp::Operation::exec(Point{x, y, 0}, mean, output); + } + + if (localY + 1 < EY && y + 1 < details.height) { + const int top = y - details.anchorY; + const int bottom = top + kernelHeight; + #pragma unroll + for (int index = 0; index < span; ++index) { + const float withBottom = make_tuple( + columnsSum[index], + readSource(firstColumn + index, bottom)) | combine; + columnsSum[index] = make_tuple( + withBottom, + readSource(firstColumn + index, top)) | subtract; + } + } + } +#endif // defined(__CUDA_ARCH__) + } +}; + +template +__global__ void launchBoxFilterQuadDPP_Kernel( + const __grid_constant__ BoxFilterQuadDetails details, + const __grid_constant__ IOps... iOps) { + DPP::exec(details, iOps...); +} + +template +FK_HOST_FUSE void executeBoxFilterQuad( + Stream_& stream, + const BoxFilterQuadDetails& details, + const IOps&... iOps) { + static_assert(DPP::PAR_ARCH == ParArch::GPU_NVIDIA, + "GPU stream requires the NVIDIA BoxFilterQuadDPP specialization"); + if (!DPP::accepts(details)) return; + const auto launch = DPP::launchConfig(details); + launchBoxFilterQuadDPP_Kernel + <<>>( + details, iOps...); + gpuErrchk(cudaGetLastError()); +} +#endif // defined(__NVCC__) + +template +FK_HOST_FUSE void executeBoxFilterQuad( + Stream_&, + const BoxFilterQuadDetails& details, + const IOps&... iOps) { + static_assert(DPP::PAR_ARCH == ParArch::CPU, + "CPU stream requires the CPU BoxFilterQuadDPP specialization"); + DPP::exec(details, iOps...); +} + +} // namespace fk + +#endif // FK_BOX_FILTER_FAST_H diff --git a/include/fused_kernel/algorithms/image_processing/image_processing.h b/include/fused_kernel/algorithms/image_processing/image_processing.h index bd7b7d64..42d1a4d7 100644 --- a/include/fused_kernel/algorithms/image_processing/image_processing.h +++ b/include/fused_kernel/algorithms/image_processing/image_processing.h @@ -23,5 +23,6 @@ #include #include #include +#include #endif // FK_IMAGE_PROCESSING diff --git a/tests/image_processing/test_box_filter_quad_dpp.h b/tests/image_processing/test_box_filter_quad_dpp.h new file mode 100644 index 00000000..c6099746 --- /dev/null +++ b/tests/image_processing/test_box_filter_quad_dpp.h @@ -0,0 +1,159 @@ +/* Copyright 2026 Johnny Nunez + + Licensed under the Apache License, Version 2.0 (the "License"); + you may not use this file except in compliance with the License. + You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + + Unless required by applicable law or agreed to in writing, software + distributed under the License is distributed on an "AS IS" BASIS, + WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + See the License for the specific language governing permissions and + limitations under the License. */ + +#include + +#include +#include +#include +#include +#include + +#include +#include + +using namespace fk; + +namespace { + +constexpr unsigned char READ_BIAS = 1; +constexpr unsigned char WRITE_BIAS = 2; + +unsigned char sourceValue(const int x, const int y) { + return static_cast((x * 17 + y * 29 + 3) % 251); +} + +unsigned char oracleAt(const int width, const int height, + const int kernelWidth, const int kernelHeight, + const int anchorX, const int anchorY, + const int x, const int y) { + float sum = 0.f; + for (int ky = 0; ky < kernelHeight; ++ky) { + for (int kx = 0; kx < kernelWidth; ++kx) { + int sx = x + kx - anchorX; + int sy = y + ky - anchorY; + sx = sx < 0 ? 0 : (sx >= width ? width - 1 : sx); + sy = sy < 0 ? 0 : (sy >= height ? height - 1 : sy); + sum += static_cast(sourceValue(sx, sy) + READ_BIAS); + } + } + const auto filtered = static_cast( + sum / static_cast(kernelWidth * kernelHeight)); + return static_cast(filtered + WRITE_BIAS); +} + +struct CaseResult { + bool passed; + std::vector output; +}; + +template +CaseResult runCase(const int width, const int height, + const int runtimeKW, const int runtimeKH, + const int anchorX, const int anchorY) { + using DPP = BoxFilterQuadDPP; + constexpr bool GPU = PA == ParArch::GPU_NVIDIA; + const MemType memoryType = GPU ? MemType::DeviceAndPinned : MemType::Host; + + Ptr2D input(width, height, 0, memoryType); + Ptr2D output(width, height, 0, memoryType); + for (int y = 0; y < height; ++y) { + for (int x = 0; x < width; ++x) { + input.at(Point{x, y, 0}) = sourceValue(x, y); + output.at(Point{x, y, 0}) = 0xA5; + } + } + + Stream_ stream; +#if defined(__NVCC__) + if constexpr (GPU) { + input.upload(stream); + output.upload(stream); + } +#endif + const int kernelWidth = KW > 0 ? KW : runtimeKW; + const int kernelHeight = KH > 0 ? KH : runtimeKH; + const auto read = PerThreadRead::build(input) + .then(Add::build(READ_BIAS)) + .then(Cast::build()); + const auto compute = make_tuple( + Add::build(), + Sub::build(), + Div::build( + static_cast(kernelWidth * kernelHeight))); + const auto write = Cast::build() + .then(Add::build(WRITE_BIAS)) + .then(PerThreadWrite::build(output)); + const BoxFilterQuadDetails details{ + width, height, runtimeKW, runtimeKH, anchorX, anchorY}; + executeBoxFilterQuad(stream, details, read, compute, write); +#if defined(__NVCC__) + if constexpr (GPU) output.download(stream); +#endif + stream.sync(); + + bool passed = true; + std::vector actual( + static_cast(width * height)); + for (int y = 0; y < height; ++y) { + for (int x = 0; x < width; ++x) { + const auto value = output.at(Point{x, y, 0}); + actual[static_cast(y * width + x)] = value; + passed = value == oracleAt( + width, height, kernelWidth, kernelHeight, + anchorX, anchorY, x, y) && passed; + } + } + return {passed, std::move(actual)}; +} + +template +bool verifyCase(const int width, const int height, + const int runtimeKW, const int runtimeKH, + const int anchorX, const int anchorY, + const char* name) { + const auto cpu = runCase( + width, height, runtimeKW, runtimeKH, anchorX, anchorY); + bool ok = cpu.passed; +#if defined(__NVCC__) + const auto gpu = runCase( + width, height, runtimeKW, runtimeKH, anchorX, anchorY); + ok = gpu.passed && gpu.output == cpu.output && ok; + std::printf("BoxFilterQuad %-14s %dx%d k%dx%d CPU/GPU %s\n", + name, width, height, + KW > 0 ? KW : runtimeKW, + KH > 0 ? KH : runtimeKH, + ok ? "PASS" : "FAIL"); +#else + std::printf("BoxFilterQuad CPU %-10s %dx%d k%dx%d %s\n", + name, width, height, + KW > 0 ? KW : runtimeKW, + KH > 0 ? KH : runtimeKH, + ok ? "PASS" : "FAIL"); +#endif + return ok; +} + +} // namespace + +int launch() { + bool ok = true; + ok = verifyCase<4, 4, 3, 3>(257, 129, 3, 3, 1, 1, "3x3-odd") && ok; + ok = verifyCase<4, 4, 5, 5>(192, 108, 5, 5, 2, 2, "5x5") && ok; + ok = verifyCase<2, 3, 7, 7>(65, 37, 7, 7, 3, 3, "7x7-ragged") && ok; + ok = verifyCase<4, 4, 7, 7>(39, 27, 1, 1, 3, 3, "static-details") && ok; + ok = verifyCase<4, 4, 9, 9>(120, 80, 9, 9, 4, 4, "9x9") && ok; + ok = verifyCase<4, 4, 0, 0>(73, 41, 11, 5, 2, 3, "runtime-11x5") && ok; + return ok ? 0 : -1; +}