diff --git a/docs/docs/reference/wasm/modules/image/cielab/api.md b/docs/docs/reference/wasm/modules/image/cielab/api.md index b8b4af65d..e6cc4d475 100644 --- a/docs/docs/reference/wasm/modules/image/cielab/api.md +++ b/docs/docs/reference/wasm/modules/image/cielab/api.md @@ -9,35 +9,56 @@ sidebar_position: 2 ## Overview -The CIELAB module provides bidirectional conversion between 8-bit sRGB and CIELAB (CIE L\*a\*b\*) color spaces. CIELAB is a perceptually uniform color space designed to approximate human vision, where equal Euclidean distances correspond to roughly equal perceived color differences. +The CIELAB module provides bidirectional conversion between sRGB and CIELAB (CIE L\*a\*b\*) color spaces. CIELAB is a perceptually uniform color space designed to approximate human vision, where equal Euclidean distances correspond to roughly equal perceived color differences. + +Functions are templated to accept arbitrary input and output datatypes: + +```cpp +template +``` + +NOTE: Data types for CIELAB color space should support signed values. Suggest using `float` or `double`. ## Functions ### `rgb_to_lab` -Convert 8-bit sRGB to CIELAB color space. +Convert sRGB to CIELAB color space. ```cpp +template void rgb_to_lab( - const uint8_t r_u8, - const uint8_t g_u8, - const uint8_t b_u8, - double& out_l, - double& out_a, - double& out_b + const Tin r_u8, + const Tin g_u8, + const Tin b_u8, + Tout& out_l, + Tout& out_a, + Tout& out_b +); + +template +void rgb_to_lab( + const ImageLib::RGBPixel &rgb, + ImageLib::LABPixel &lab +); + +template +void rgb_to_lab( + const ImageLib::RGBAPixel &rgba, + ImageLib::LABAPixel &laba ); ``` #### Parameters -| Parameter | Type | Range | Description | -| :-------- | :-------- | :----------- | :------------------------------ | -| `r_u8` | `uint8_t` | [0, 255] | Red channel (input) | -| `g_u8` | `uint8_t` | [0, 255] | Green channel (input) | -| `b_u8` | `uint8_t` | [0, 255] | Blue channel (input) | -| `out_l` | `double&` | [0, 100] | L\* lightness (output, clamped) | -| `out_a` | `double&` | ~[-128, 127] | a\* green-red axis (output) | -| `out_b` | `double&` | ~[-128, 127] | b\* blue-yellow axis (output) | +| Parameter | Type | Range | Description | +| :-------- | :------ | :----------- | :------------------------------ | +| `r_u8` | `Tin` | [0, 255] | Red channel (input) | +| `g_u8` | `Tin` | [0, 255] | Green channel (input) | +| `b_u8` | `Tin` | [0, 255] | Blue channel (input) | +| `out_l` | `Tout&` | [0, 100] | L\* lightness (output, clamped) | +| `out_a` | `Tout&` | ~[-128, 127] | a\* green-red axis (output) | +| `out_b` | `Tout&` | ~[-128, 127] | b\* blue-yellow axis (output) | #### Transformation Pipeline @@ -69,29 +90,42 @@ rgb_to_lab(r, g, b, L, a, b_lab); ### `lab_to_rgb` -Convert CIELAB to 8-bit sRGB color space. +Convert CIELAB to sRGB color space. ```cpp +template +void lab_to_rgb( + const Tin L, + const Tin A, + const Tin B, + Tout& r_u8, + Tout& g_u8, + Tout& b_u8 +); + +template +void lab_to_rgb( + const ImageLib::LABPixel &lab, + ImageLib::RGBPixel &rgb +); + +template void lab_to_rgb( - const double L, - const double A, - const double B, - uint8_t& r_u8, - uint8_t& g_u8, - uint8_t& b_u8 + const ImageLib::LABAPixel &laba, + ImageLib::RGBAPixel &rgba ); ``` #### Parameters -| Parameter | Type | Range | Description | -| :-------- | :--------- | :----------- | :------------------------------ | -| `L` | `double` | [0, 100] | L\* lightness (input) | -| `A` | `double` | ~[-128, 127] | a\* green-red axis (input) | -| `B` | `double` | ~[-128, 127] | b\* blue-yellow axis (input) | -| `r_u8` | `uint8_t&` | [0, 255] | Red channel (output, clamped) | -| `g_u8` | `uint8_t&` | [0, 255] | Green channel (output, clamped) | -| `b_u8` | `uint8_t&` | [0, 255] | Blue channel (output, clamped) | +| Parameter | Type | Range | Description | +| :-------- | :------ | :----------- | :------------------------------ | +| `L` | `Tin` | [0, 100] | L\* lightness (input) | +| `A` | `Tin` | ~[-128, 127] | a\* green-red axis (input) | +| `B` | `Tin` | ~[-128, 127] | b\* blue-yellow axis (input) | +| `r_u8` | `Tout&` | [0, 255] | Red channel (output, clamped) | +| `g_u8` | `Tout&` | [0, 255] | Green channel (output, clamped) | +| `b_u8` | `Tout&` | [0, 255] | Blue channel (output, clamped) | #### Transformation Pipeline diff --git a/src/components/WasmImageProcessor.jsx b/src/components/WasmImageProcessor.jsx index 89d9f4efa..fa7f2b5d3 100644 --- a/src/components/WasmImageProcessor.jsx +++ b/src/components/WasmImageProcessor.jsx @@ -13,7 +13,7 @@ const WasmImageProcessor = () => { const inputId = useId(); const inputRef = useRef(null); - const { bilateralFilter, blackThreshold, kmeans, findContours } = useWasmWorker(); + const { bilateralFilter, kmeans, findContours } = useWasmWorker(); const [originalSrc, setOriginalSrc] = useState(null); const [fileData, setFileData] = useState(null); @@ -85,18 +85,19 @@ const WasmImageProcessor = () => { height, }); - step(45); - const thresholded = await blackThreshold({ + // step(45); + // is this needed? num_colors is incorrect should be num_threshold + /*const thresholded = await blackThreshold({ ...fileData, pixels: imgBilateralFiltered, num_colors: 8, - }); + });*/ step(70); const { pixels: kmeansed, labels } = await kmeans({ ...fileData, - pixels: thresholded, - num_colors: 8, + pixels: imgBilateralFiltered, + num_colors: 16, }); const contours = await findContours({ @@ -126,7 +127,7 @@ const WasmImageProcessor = () => { step(0); }, 800); } - }, [fileData, bilateralFilter, blackThreshold, kmeans, findContours, navigate, step]); + }, [fileData, bilateralFilter, kmeans, findContours, navigate, step]); /* Memo'd UI fragments */ const EmptyState = useMemo( diff --git a/src/hooks/useWasmWorker.js b/src/hooks/useWasmWorker.js index b26745e23..5ac1bed48 100644 --- a/src/hooks/useWasmWorker.js +++ b/src/hooks/useWasmWorker.js @@ -42,9 +42,12 @@ export function useWasmWorker() { sigma_spatial = 3.0, sigma_range = 50.0, color_space = 0, + n_threads = 8, }) => { return ( - await call('bilateral_filter', { pixels, width, height, sigma_spatial, sigma_range, color_space }, ['pixels']) + await call('bilateral_filter', { pixels, width, height, sigma_spatial, sigma_range, color_space, n_threads }, [ + 'pixels', + ]) ).output.pixels; }; const blackThreshold = async ({ pixels, width, height, num_colors }) => { @@ -58,13 +61,15 @@ export function useWasmWorker() { out_labels = new Int32Array(width * height), num_colors, max_iter = 100, + color_space = 0, + n_threads = 8, }) => { const result = ( - await call('kmeans', { pixels, out_pixels, out_labels, width, height, num_colors, max_iter }, [ - 'pixels', - 'out_pixels', - 'out_labels', - ]) + await call( + 'kmeans', + { pixels, out_pixels, out_labels, width, height, num_colors, max_iter, color_space, n_threads }, + ['pixels', 'out_pixels', 'out_labels'] + ) ).output; return { pixels: result.out_pixels, diff --git a/src/utils/image-utils.js b/src/utils/image-utils.js index 64c41575c..d01f7045b 100644 --- a/src/utils/image-utils.js +++ b/src/utils/image-utils.js @@ -35,6 +35,7 @@ export const uint8ClampedArrayToSVG = async ({ pixels, width, height }) => { let svgString = ImageTracer.imagedataToSVG(imageData, { ltres: 0, rightangleenhance: false, + // numberofcolors: 128, }); // Inject viewBox diff --git a/src/wasm/modules/image/CMakeLists.txt b/src/wasm/modules/image/CMakeLists.txt index ee65ccd7d..a75eb1914 100644 --- a/src/wasm/modules/image/CMakeLists.txt +++ b/src/wasm/modules/image/CMakeLists.txt @@ -35,21 +35,24 @@ set(COMMON_FLAGS "SHELL:-s MODULARIZE=1" "SHELL:-s EXPORT_ES6=1" "SHELL:-s EXIT_RUNTIME=1" - "SHELL:-s ENVIRONMENT=web" + "SHELL:-s ENVIRONMENT=web,worker" "SHELL:-s EXPORTED_FUNCTIONS=['_malloc','_free']" "SHELL:-s EXPORTED_RUNTIME_METHODS=['ccall','cwrap','getValue','setValue','HEAPU8','HEAP32']" "SHELL:-s INITIAL_MEMORY=1024MB" "SHELL:-s MAXIMUM_MEMORY=2048MB" "SHELL:-s ALLOW_MEMORY_GROWTH=1" "SHELL:-s EXPORT_NAME=create${CAP_MODULE_NAME}Module" + "SHELL:-s USE_PTHREADS=1" + "SHELL:-s PTHREAD_POOL_SIZE=8" ) # Apply common flags -target_link_options(${MODULE_NAME}_wasm PRIVATE ${COMMON_FLAGS}) +target_compile_options(${MODULE_NAME}_wasm PRIVATE -pthread -msimd128) +target_link_options(${MODULE_NAME}_wasm PRIVATE ${COMMON_FLAGS} -pthread -msimd128) # Build-type specific flags if(CMAKE_BUILD_TYPE STREQUAL "Debug") - target_compile_options(${MODULE_NAME}_wasm PRIVATE -O0 -g4) + target_compile_options(${MODULE_NAME}_wasm PRIVATE -O2 -g4) target_link_options(${MODULE_NAME}_wasm PRIVATE "SHELL:-s ASSERTIONS=2" -g4 diff --git a/src/wasm/modules/image/include/LABAPixel.h b/src/wasm/modules/image/include/LABAPixel.h new file mode 100644 index 000000000..80a92e990 --- /dev/null +++ b/src/wasm/modules/image/include/LABAPixel.h @@ -0,0 +1,34 @@ +#ifndef LABAPixel_H +#define LABAPixel_H + +#include "LABPixel.h" + +namespace ImageLib { +template +struct LABAPixel : public ImageLib::LABPixel { + // ----- Members ----- + NumberT alpha; + + // ----- Constructors ----- + constexpr LABAPixel(NumberT l = 0, NumberT a = 0, NumberT b = 0, + NumberT alpha = 255) + : LABPixel(l, a, b), alpha(alpha) {} + + // ----- Modifiers ----- + [[nodiscard]] inline bool operator==(const LABAPixel &other) const { + return LABPixel::operator==(other) && alpha == other.alpha; + } + [[nodiscard]] inline bool operator!=(const LABAPixel &other) const { + return !(*this == other); + } + + // ----- Utilities ----- + inline void setGray(NumberT gray, NumberT alpha = 255) { + LABPixel::setGray(gray); + this->alpha = alpha; + } + +} __attribute__((packed)); +} // namespace ImageLib + +#endif // LABAPixel_H diff --git a/src/wasm/modules/image/include/LABPixel.h b/src/wasm/modules/image/include/LABPixel.h new file mode 100644 index 000000000..b16ac913d --- /dev/null +++ b/src/wasm/modules/image/include/LABPixel.h @@ -0,0 +1,36 @@ +#ifndef LABPIXEL_H +#define LABPIXEL_H + +#include "Pixel.h" + +/* +can support signed data types +preferrably float or double +*/ + +namespace ImageLib { +template struct LABPixel : public Pixel { + // ----- Members ----- + NumberT l, a, b; + + constexpr LABPixel(NumberT l = 0, NumberT a = 0, NumberT b = 0) + : l(l), a(a), b(b) {} + + // ----- Modifiers ----- + [[nodiscard]] inline bool operator==(const LABPixel &other) const { + return l == other.l && a == other.a && b == other.b; + } + [[nodiscard]] inline bool operator!=(const LABPixel &other) const { + return !(*this == other); + } + + // ----- Utilities ----- + inline void setGray(NumberT new_luma) { + l = new_luma; + a = b = 0; + } + +} __attribute__((packed)); +} // namespace ImageLib + +#endif // LABPIXEL_H diff --git a/src/wasm/modules/image/include/bilateral_filter.h b/src/wasm/modules/image/include/bilateral_filter.h index 2b5b938fe..607a3eb24 100644 --- a/src/wasm/modules/image/include/bilateral_filter.h +++ b/src/wasm/modules/image/include/bilateral_filter.h @@ -17,7 +17,7 @@ namespace bilateral { // (radiometric decay) void bilateral_filter(uint8_t *image, size_t width, size_t height, double sigma_spatial, double sigma_range, - uint8_t color_space); + uint8_t color_space, const uint8_t n_threads = 8); } // namespace bilateral diff --git a/src/wasm/modules/image/include/cielab.h b/src/wasm/modules/image/include/cielab.h index 31625f9b4..72e339140 100644 --- a/src/wasm/modules/image/include/cielab.h +++ b/src/wasm/modules/image/include/cielab.h @@ -1,11 +1,39 @@ #ifndef CIELAB_H #define CIELAB_H +#include "LABAPixel.h" +#include "LABPixel.h" +#include "RGBAPixel.h" +#include "RGBPixel.h" #include +#include -void rgb_to_lab(const uint8_t r_u8, const uint8_t g_u8, const uint8_t b_u8, - double &out_l, double &out_a, double &out_b); +// templates must be in headers + +template +void rgb_to_lab(const Tin r_u8, const Tin g_u8, const Tin b_u8, Tout &out_l, + Tout &out_a, Tout &out_b); + +template +void lab_to_rgb(const Tin L, const Tin A, const Tin B, Tout &r_u8, Tout &g_u8, + Tout &b_u8); + +template +void rgb_to_lab(const ImageLib::RGBAPixel &rgba, + ImageLib::LABAPixel &laba); + +template +void rgb_to_lab(const ImageLib::RGBPixel &rgb, + ImageLib::LABPixel &lab); + +template +void lab_to_rgb(const ImageLib::LABAPixel &laba, + ImageLib::RGBAPixel &rgba); + +template +void lab_to_rgb(const ImageLib::LABPixel &lab, + ImageLib::RGBPixel &rgb); + +#include "cielab_impl.h" -void lab_to_rgb(const double L, const double A, const double B, uint8_t &r_u8, - uint8_t &g_u8, uint8_t &b_u8); #endif // CIELAB_H diff --git a/src/wasm/modules/image/src/cielab.cpp b/src/wasm/modules/image/include/cielab_impl.h similarity index 68% rename from src/wasm/modules/image/src/cielab.cpp rename to src/wasm/modules/image/include/cielab_impl.h index 073594f4e..018284680 100644 --- a/src/wasm/modules/image/src/cielab.cpp +++ b/src/wasm/modules/image/include/cielab_impl.h @@ -1,4 +1,10 @@ -#include "cielab.h" +#ifndef CIELAB_IMPL_H +#define CIELAB_IMPL_H + +#ifndef CIELAB_H +#error "cielab_impl.h should not be included directly; include cielab.h instead" +#endif + #include #include @@ -95,12 +101,18 @@ inline double srgb_to_linear(const double c) { SRGB_GAMMA); } -void rgb_to_lab(const uint8_t r_u8, const uint8_t g_u8, const uint8_t b_u8, - double &out_l, double &out_a, double &out_b) { +template +void rgb_to_lab(const Tin r_u8, const Tin g_u8, const Tin b_u8, Tout &out_l, + Tout &out_a, Tout &out_b) { // 1. Convert 8-bit RGB [0, 255] to linear RGB [0.0, 1.0] - double r{srgb_to_linear(r_u8 / 255.0)}; - double g{srgb_to_linear(g_u8 / 255.0)}; - double b{srgb_to_linear(b_u8 / 255.0)}; + + double _r = static_cast(r_u8); + double _g = static_cast(g_u8); + double _b = static_cast(b_u8); + + double r{srgb_to_linear(_r / 255.0)}; + double g{srgb_to_linear(_g / 255.0)}; + double b{srgb_to_linear(_b / 255.0)}; // 2. Convert linear RGB to CIE XYZ (using D65 white point reference) // The matrix below is for sRGB to XYZ (D65) @@ -119,11 +131,11 @@ void rgb_to_lab(const uint8_t r_u8, const uint8_t g_u8, const uint8_t b_u8, const double fz{xyz_to_lab(Zr)}; // 4. Output values - out_l = LAB_L_FACTOR * fy - LAB_L_OFFSET; - out_a = LAB_A_FACTOR * (fx - fy); - out_b = LAB_B_FACTOR * (fy - fz); + out_l = static_cast(LAB_L_FACTOR * fy - LAB_L_OFFSET); + out_a = static_cast(LAB_A_FACTOR * (fx - fy)); + out_b = static_cast(LAB_B_FACTOR * (fy - fz)); - out_l = std::clamp(out_l, 0.0, 100.0); + out_l = std::clamp(out_l, static_cast(0.0), static_cast(100.0)); } constexpr double inverse_xyz_to_lab(double t) { @@ -139,12 +151,18 @@ inline double gamma_encode(double u) { SRGB_GAMMA_OFFSET; } -void lab_to_rgb(const double L, const double A, const double B, - uint8_t &out_r_u8, uint8_t &out_g_u8, uint8_t &out_b_u8) { +template +void lab_to_rgb(const Tin L, const Tin A, const Tin B, Tout &out_r_u8, + Tout &out_g_u8, Tout &out_b_u8) { + + const double _L = static_cast(L); + const double _A = static_cast(A); + const double _B = static_cast(B); + // --- Lab → XYZ (D65 white point) - const double fy{(L + LAB_L_OFFSET) / LAB_L_FACTOR}; - const double fx{fy + A / LAB_A_FACTOR}; - const double fz{fy - B / LAB_B_FACTOR}; + const double fy{(_L + LAB_L_OFFSET) / LAB_L_FACTOR}; + const double fx{fy + _A / LAB_A_FACTOR}; + const double fz{fy - _B / LAB_B_FACTOR}; const double X{D65_Xn * inverse_xyz_to_lab(fx)}; const double Y{D65_Yn * inverse_xyz_to_lab(fy)}; @@ -161,7 +179,41 @@ void lab_to_rgb(const double L, const double A, const double B, b = gamma_encode(std::clamp(b, 0.0, 1.0)); // --- Clamp and convert to 8-bit - out_r_u8 = static_cast(std::round(255.0 * std::clamp(r, 0.0, 1.0))); - out_g_u8 = static_cast(std::round(255.0 * std::clamp(g, 0.0, 1.0))); - out_b_u8 = static_cast(std::round(255.0 * std::clamp(b, 0.0, 1.0))); + out_r_u8 = static_cast(std::round(255.0 * std::clamp(r, 0.0, 1.0))); + out_g_u8 = static_cast(std::round(255.0 * std::clamp(g, 0.0, 1.0))); + out_b_u8 = static_cast(std::round(255.0 * std::clamp(b, 0.0, 1.0))); } + +template +void rgb_to_lab(const ImageLib::RGBAPixel &rgba, + ImageLib::LABAPixel &laba) { + + rgb_to_lab(rgba.red, rgba.green, rgba.blue, laba.l, laba.a, + laba.b); + laba.alpha = static_cast(rgba.alpha); +} + +template +void rgb_to_lab(const ImageLib::RGBPixel &rgb, + ImageLib::LABPixel &lab) { + + rgb_to_lab(rgb.red, rgb.green, rgb.blue, lab.l, lab.a, lab.b); +} + +template +void lab_to_rgb(const ImageLib::LABAPixel &laba, + ImageLib::RGBAPixel &rgba) { + + lab_to_rgb(laba.l, laba.a, laba.b, rgba.red, rgba.green, + rgba.blue); + rgba.alpha = static_cast(laba.alpha); +} + +template +void lab_to_rgb(const ImageLib::LABPixel &lab, + ImageLib::RGBPixel &rgb) { + + lab_to_rgb(lab.l, lab.a, lab.b, rgb.red, rgb.green, rgb.blue); +} + +#endif \ No newline at end of file diff --git a/src/wasm/modules/image/include/image_utils.h b/src/wasm/modules/image/include/image_utils.h index da86ab524..30a286a80 100644 --- a/src/wasm/modules/image/include/image_utils.h +++ b/src/wasm/modules/image/include/image_utils.h @@ -3,7 +3,6 @@ #include "exported.h" // EXPORTED macro -#include #include #include "Image.h" diff --git a/src/wasm/modules/image/include/kmeans.h b/src/wasm/modules/image/include/kmeans.h index f83953f67..e24e5b963 100644 --- a/src/wasm/modules/image/include/kmeans.h +++ b/src/wasm/modules/image/include/kmeans.h @@ -8,11 +8,12 @@ EXPORTED void kmeans(const uint8_t *data, uint8_t *out_data, int32_t *out_labels, const int32_t width, const int32_t height, const int32_t k, - const int32_t max_iter); + const int32_t max_iter, const uint8_t color_space, + const uint8_t n_threads); -EXPORTED void kmeans_clustering_spatial(uint8_t *data, int32_t width, +/*EXPORTED void kmeans_clustering_spatial(uint8_t *data, int32_t width, int32_t height, int32_t k, int32_t max_iter, - float spatial_weight = 1.0); + float spatial_weight = 1.0);*/ #endif diff --git a/src/wasm/modules/image/src/bilateral_filter.cpp b/src/wasm/modules/image/src/bilateral_filter.cpp index dfe807f7f..46de91910 100644 --- a/src/wasm/modules/image/src/bilateral_filter.cpp +++ b/src/wasm/modules/image/src/bilateral_filter.cpp @@ -6,6 +6,9 @@ #include #include #include +#include +#include +#include #include namespace bilateral { @@ -23,6 +26,8 @@ inline double gaussian(double x, double sigma) { return std::exp(-(x * x) / (2.0 * sigma * sigma)); } +std::mutex write_mutex; + /* The Bilateral Filter applies a composite weight based on both spatial distance and radiometric difference (intensity) to return an image that is smoothed while @@ -40,74 +45,18 @@ decay) ├── 0: CIELAB └── 1: RGB */ -void bilateral_filter(uint8_t *image, size_t width, size_t height, - double sigma_spatial, double sigma_range, - uint8_t color_space) { - // bad data -> return - if (sigma_spatial <= 0.0 || sigma_range <= 0.0 || width <= 0 || height <= 0) - return; - if (color_space != COLOR_SPACE_OPTION_CIELAB && - color_space != COLOR_SPACE_OPTION_RGB) - return; - - const int raw_radius{ - static_cast(std::ceil(SIGMA_RADIUS_FACTOR * sigma_spatial))}; - const int radius{std::min(raw_radius, MAX_KERNEL_RADIUS)}; - const int kernel_diameter{2 * radius + 1}; - - std::vector result(width * height * 4); - - std::vector spatial_weights(kernel_diameter * kernel_diameter); - - // Precompute Spatial Weights (Gaussian Kernel) - for (int ky{-radius}; ky <= radius; ++ky) { - for (int kx{-radius}; kx <= radius; ++kx) { - const double dist{static_cast(std::sqrt(kx * kx + ky * ky))}; - spatial_weights[(ky + radius) * kernel_diameter + (kx + radius)] = - gaussian(dist, sigma_spatial); - } - } - - // ========= RGB-only section start ========= - // Precompute Range Weights - std::vector range_lut; - if (color_space == COLOR_SPACE_OPTION_RGB) { - range_lut.resize(MAX_RGB_DIST_SQ + 1); - for (int i{0}; i <= MAX_RGB_DIST_SQ; ++i) { - range_lut[i] = gaussian(static_cast(std::sqrt(i)), sigma_range); - } - } - // ========= RGB-only section end ========= - - // ========= CIELAB section start ========= - // Compute full image RGB - CIELAB conversion - std::vector cie_image; - if (color_space == COLOR_SPACE_OPTION_CIELAB) { - cie_image.resize(width * height * 4); - - for (int y{0}; y < height; y++) { - for (int x{0}; x < width; x++) { - int center_idx{(y * static_cast(width) + x) * 4}; - uint8_t r0{image[center_idx]}; - uint8_t g0{image[center_idx + 1]}; - uint8_t b0{image[center_idx + 2]}; - uint8_t a0{image[center_idx + 3]}; - double L0, A0, B0; - rgb_to_lab(r0, g0, b0, L0, A0, B0); - - cie_image[center_idx] = L0; - cie_image[center_idx + 1] = A0; - cie_image[center_idx + 2] = B0; - cie_image[center_idx + 3] = - 0.0; // unused but keep for indexing purposes - } - } - } - // ========= CIELAB section end ========= +void _process(const uint8_t *image, const std::vector &cie_image, + std::vector &result, + const std::vector &spatial_weights, + const std::vector &range_lut, int radius, + double sigma_range, int start_row, int end_row, size_t height, + size_t width, uint8_t color_space, const uint8_t n_threads) { int h{static_cast(height)}; int w{static_cast(width)}; - for (int y{0}; y < h; ++y) { + const int kernel_diameter{2 * radius + 1}; + + for (int y{start_row}; y < end_row; ++y) { for (int x{0}; x < w; ++x) { size_t center_idx{(y * width + x) * 4}; @@ -184,6 +133,8 @@ void bilateral_filter(uint8_t *image, size_t width, size_t height, weight_acc += w_space * w_range; } } + // writing - must grab mutex + std::unique_lock lock{write_mutex}; switch (color_space) { case COLOR_SPACE_OPTION_RGB: { @@ -201,7 +152,7 @@ void bilateral_filter(uint8_t *image, size_t width, size_t height, double A{weight_acc_channel_1 / weight_acc}; double B{weight_acc_channel_2 / weight_acc}; uint8_t r, g, b; - lab_to_rgb(L, A, B, r, g, b); + lab_to_rgb(L, A, B, r, g, b); result[center_idx] = r; result[center_idx + 1] = g; result[center_idx + 2] = b; @@ -209,8 +160,103 @@ void bilateral_filter(uint8_t *image, size_t width, size_t height, break; } } + // done writing + lock.unlock(); } } +} + +void bilateral_filter(uint8_t *image, size_t width, size_t height, + double sigma_spatial, double sigma_range, + uint8_t color_space, const uint8_t n_threads) { + // bad data -> return + if (sigma_spatial <= 0.0 || sigma_range <= 0.0 || width <= 0 || height <= 0) + return; + if (color_space != COLOR_SPACE_OPTION_CIELAB && + color_space != COLOR_SPACE_OPTION_RGB) + return; + + std::vector threads; + + int rows_per_thread = static_cast(height) / static_cast(n_threads); + + const int raw_radius{ + static_cast(std::ceil(SIGMA_RADIUS_FACTOR * sigma_spatial))}; + const int radius{std::min(raw_radius, MAX_KERNEL_RADIUS)}; + const int kernel_diameter{2 * radius + 1}; + + std::vector result(width * height * 4); + + std::vector spatial_weights(kernel_diameter * kernel_diameter); + + // Precompute Spatial Weights (Gaussian Kernel) + for (int ky{-radius}; ky <= radius; ++ky) { + for (int kx{-radius}; kx <= radius; ++kx) { + const double dist{static_cast(std::sqrt(kx * kx + ky * ky))}; + spatial_weights[(ky + radius) * kernel_diameter + (kx + radius)] = + gaussian(dist, sigma_spatial); + } + } + + // ========= RGB-only section start ========= + // Precompute Range Weights + std::vector range_lut; + if (color_space == COLOR_SPACE_OPTION_RGB) { + range_lut.resize(MAX_RGB_DIST_SQ + 1); + for (int i{0}; i <= MAX_RGB_DIST_SQ; ++i) { + range_lut[i] = gaussian(static_cast(std::sqrt(i)), sigma_range); + } + } + // ========= RGB-only section end ========= + + // ========= CIELAB section start ========= + // Compute full image RGB - CIELAB conversion + std::vector cie_image; + if (color_space == COLOR_SPACE_OPTION_CIELAB) { + cie_image.resize(width * height * 4); + + for (int y{0}; y < height; y++) { + for (int x{0}; x < width; x++) { + int center_idx{(y * static_cast(width) + x) * 4}; + uint8_t r0{image[center_idx]}; + uint8_t g0{image[center_idx + 1]}; + uint8_t b0{image[center_idx + 2]}; + uint8_t a0{image[center_idx + 3]}; + double L0, A0, B0; + rgb_to_lab(r0, g0, b0, L0, A0, B0); + + cie_image[center_idx] = L0; + cie_image[center_idx + 1] = A0; + cie_image[center_idx + 2] = B0; + cie_image[center_idx + 3] = + 0.0; // unused but keep for indexing purposes + } + } + } + // ========= CIELAB section end ========= + if (n_threads > 1) { + for (int i = 0; i < n_threads; ++i) { + int start_row{i * rows_per_thread}; + int end_row{(i == n_threads - 1) ? static_cast(height) + : (i + 1) * rows_per_thread}; + // Launch a thread and add to vector + threads.emplace_back(_process, std::cref(image), std::cref(cie_image), + std::ref(result), std::cref(spatial_weights), + std::cref(range_lut), radius, sigma_range, start_row, + end_row, height, width, color_space, n_threads); + } + + // wait for threads to finish + for (auto &thread : threads) { + if (thread.joinable()) { + thread.join(); + } + } + } else { + _process(image, cie_image, result, spatial_weights, range_lut, radius, + sigma_range, 0, static_cast(height), height, width, + color_space, n_threads); + } std::memcpy(image, result.data(), result.size()); } @@ -220,7 +266,7 @@ void bilateral_filter(uint8_t *image, size_t width, size_t height, // Global wrapper for WASM export EXPORTED void bilateral_filter(uint8_t *image, size_t width, size_t height, double sigma_spatial, double sigma_range, - uint8_t color_space) { + uint8_t color_space, const uint8_t n_threads) { bilateral::bilateral_filter(image, width, height, sigma_spatial, sigma_range, - color_space); -} + color_space, n_threads); +} \ No newline at end of file diff --git a/src/wasm/modules/image/src/graph.cpp b/src/wasm/modules/image/src/graph.cpp index 698ad5934..6c6b6dd0a 100644 --- a/src/wasm/modules/image/src/graph.cpp +++ b/src/wasm/modules/image/src/graph.cpp @@ -1,4 +1,5 @@ #include "graph.h" +#include "Pixel.h" #include #include @@ -6,6 +7,20 @@ *Graph class - manages Node class */ +static inline float colorDistance(const ImageLib::RGBPixel &a, + const ImageLib::RGBPixel &b) { + + ImageLib::RGBPixel af{static_cast(a.red), + static_cast(a.green), + static_cast(a.blue)}; + ImageLib::RGBPixel bf{static_cast(b.red), + static_cast(b.green), + static_cast(b.blue)}; + return std::sqrt((af.red - bf.red) * (af.red - bf.red) + + (af.green - bf.green) * (af.green - bf.green) + + (af.blue - bf.blue) * (af.blue - bf.blue)); +} + /* *To quickly search m_nodes (std::vector) for the index of a node id *create an std::unordered_map of node id - index pairs @@ -128,9 +143,20 @@ void Graph::merge_small_area_nodes(const int32_t min_area) { std::copy(n->edges().begin(), n->edges().end(), std::back_inserter(neighbors)); + ImageLib::RGBPixel col = n->color(); // Sort by size -> a.area < b.area + // std::sort(neighbors.begin(), neighbors.end(), + // [](Node_ptr a, Node_ptr b) { return a->area() < b->area(); + // }); + + // sort by size and color similarity std::sort(neighbors.begin(), neighbors.end(), - [](Node_ptr a, Node_ptr b) { return a->area() < b->area(); }); + [col](Node_ptr a, Node_ptr b) { + float cdista = colorDistance(a->color(), col); + float cdistb = colorDistance(b->color(), col); + return (static_cast(a->area()) + 10.f * cdista) < + (static_cast(b->area()) + 10.f * cdistb); + }); int32_t idx{0}; // find first non-zero area neighbor diff --git a/src/wasm/modules/image/src/kmeans.cpp b/src/wasm/modules/image/src/kmeans.cpp index 66a90ab52..b724fe1a3 100644 --- a/src/wasm/modules/image/src/kmeans.cpp +++ b/src/wasm/modules/image/src/kmeans.cpp @@ -1,23 +1,153 @@ #include "kmeans.h" #include "Image.h" +#include "LABAPixel.h" #include "PixelConverters.h" #include "RGBAPixel.h" +#include "cielab.h" +#include +#include #include #include #include +#include #include +#include +#include +#include #include static inline float colorDistance(const ImageLib::RGBAPixel &a, const ImageLib::RGBAPixel &b) { - return std::sqrt((a.red - b.red) * (a.red - b.red) + - (a.green - b.green) * (a.green - b.green) + - (a.blue - b.blue) * (a.blue - b.blue)); + // sqrt un-necessary + return (a.red - b.red) * (a.red - b.red) + + (a.green - b.green) * (a.green - b.green) + + (a.blue - b.blue) * (a.blue - b.blue); +} + +static inline float colorDistance(const ImageLib::LABAPixel &a, + const ImageLib::LABAPixel &b) { + // sqrt un-necessary + return (a.l - b.l) * (a.l - b.l) + (a.a - b.a) * (a.a - b.a) + + (a.b - b.b) * (a.b - b.b); +} + +static constexpr uint8_t COLOR_SPACE_OPTION_CIELAB{0}; +static constexpr uint8_t COLOR_SPACE_OPTION_RGB{1}; + +template +void _process_dist_per_centroid(const ImageLib::Image &pixels, + const ImageLib::Image ¢roids, + std::vector> &output, + int start_centroid, int end_centroid) { + // threads will not overlap - no need for mutex + std::vector _res(pixels.getPixelCount()); + for (int j{start_centroid}; j < end_centroid; ++j) { + std::transform(pixels.begin(), pixels.end(), _res.begin(), + [¢roids, j](const PixelT &p) { + return colorDistance(p, centroids[j]); + }); + std::copy(_res.begin(), _res.end(), output[j].begin()); + } +} + +template +void _apply_labels(const ImageLib::Image &pixels, + const std::vector> &distances, + std::vector &labels, int start_pixel, int end_pixel, + int k, std::atomic &changed) { + // threads don't overlap + float min_color_dist{std::numeric_limits::max()}; + int32_t best_cluster{0}; + for (int i{start_pixel}; i < end_pixel; ++i) { + min_color_dist = std::numeric_limits::max(); + best_cluster = 0; + for (int j{0}; j < k; ++j) { + if (distances[j][i] < min_color_dist) { + min_color_dist = distances[j][i]; + best_cluster = j; + } + } + if (labels[i] != best_cluster) { + labels[i] = best_cluster; + changed.store(true, std::memory_order_relaxed); + } + } +} + +// The K-Means++ Initialization Function +template +void kMeansPlusPlusInit(const ImageLib::Image &pixels, + ImageLib::Image &out_centroids, int k) { + std::vector centroids; + + int num_pixels = pixels.getSize(); + // Random number generator setup + std::random_device rd; + std::mt19937 gen(rd()); + + // --- Step 1: Choose the first centroid uniformly at random --- + std::uniform_int_distribution<> dis(0, num_pixels - 1); + int first_index = dis(gen); + centroids.push_back(pixels[first_index]); + + // Vector to store the squared distance of each pixel to its NEAREST existing + // centroid. Initialize with max double so the first distance calculation + // always updates it. + std::vector min_dist_sq(num_pixels, + std::numeric_limits::max()); + + // --- Step 2 & 3: Repeat until we have k centroids --- + for (int i = 1; i < k; ++i) { + + double sum_dist_sq = 0.0; + + // Update distances relative to the LAST added centroid (centroids.back()) + // We don't need to recheck previous centroids; min_dist_sq already holds + // the best distance to them. + for (int j = 0; j < num_pixels; ++j) { + double d = colorDistance(pixels[j], centroids.back()); + + // If this new centroid is closer than the previous best, update the min + // distance + if (d < min_dist_sq[j]) { + min_dist_sq[j] = d; + } + sum_dist_sq += min_dist_sq[j]; + } + + // --- Step 3: Choose new center with probability proportional to D(x)^2 --- + // We use a weighted random selection (Roulette Wheel Selection) + std::uniform_real_distribution<> dist_selector(0.0, sum_dist_sq); + double random_value = dist_selector(gen); + + double current_sum = 0.0; + int selected_index = -1; + + // Iterate to find the pixel corresponding to the random_value + for (int j = 0; j < num_pixels; ++j) { + current_sum += min_dist_sq[j]; + if (current_sum >= random_value) { + selected_index = j; + break; + } + } + + // Fallback for floating point rounding errors (pick last one if loop + // finishes) + if (selected_index == -1) { + selected_index = num_pixels - 1; + } + + centroids.push_back(pixels[selected_index]); + } + + std::copy(centroids.begin(), centroids.end(), out_centroids.begin()); } void kmeans(const uint8_t *data, uint8_t *out_data, int32_t *out_labels, const int32_t width, const int32_t height, const int32_t k, - const int32_t max_iter) { + const int32_t max_iter, const uint8_t color_space, + const uint8_t n_threads) { ImageLib::Image> pixels; pixels.loadFromBuffer(data, width, height, ImageLib::RGBA_CONVERTER); const int32_t num_pixels{pixels.getSize()}; @@ -26,41 +156,165 @@ void kmeans(const uint8_t *data, uint8_t *out_data, int32_t *out_labels, // k centroids, initialized to rgba(0,0,0,255) // Init of each pixel is from default in Image constructor ImageLib::Image> centroids{k, 1}; + ImageLib::Image> centroids_lab{k, 1}; std::vector labels(num_pixels, 0); + ImageLib::Image> lab(pixels.getWidth(), + pixels.getHeight()); + if (color_space == COLOR_SPACE_OPTION_CIELAB) { + for (int i{0}; i < pixels.getSize(); ++i) { + rgb_to_lab(pixels[i], lab[i]); + } + } + + std::vector threads; + + int nthreads = std::clamp(static_cast(n_threads), 1, k); + + int pixels_per_thread{num_pixels / nthreads}; + int centroids_per_thread{k / nthreads}; + // Step 2: Initialize centroids randomly - srand(static_cast(time(nullptr))); + // This does not give satisfactory initializations + /*srand(static_cast(time(nullptr))); for (int32_t i{0}; i < k; ++i) { int32_t idx = rand() % num_pixels; - centroids[i] = pixels[idx]; + switch (color_space) { + case COLOR_SPACE_OPTION_RGB : { + centroids[i] = pixels[idx]; + break; + } + case COLOR_SPACE_OPTION_CIELAB : { + centroids_lab[i] = lab[idx]; + break; + } + } + }*/ + + switch (color_space) { + case COLOR_SPACE_OPTION_RGB: { + kMeansPlusPlusInit>(pixels, centroids, k); + break; + } + case COLOR_SPACE_OPTION_CIELAB: { + kMeansPlusPlusInit>(lab, centroids_lab, k); + break; + } } // Step 3: Run k-means iterations + + std::vector> distances; + if (nthreads > 1) { + distances.resize(k); + for (auto &d : distances) { + d.resize(num_pixels, std::numeric_limits::max()); + } + } + + // Assignment step for (int32_t iter{0}; iter < max_iter; ++iter) { bool changed{false}; - // Assignment step - // Iterate over pixels - for (int32_t i{0}; i < num_pixels; ++i) { - float min_color_dist{std::numeric_limits::max()}; - int32_t best_cluster{0}; - - // Iterate over centroids to find centroid with most similar color to - // pixels[i] - for (int32_t j{0}; j < k; ++j) { - float dist{colorDistance(pixels[i], centroids[j])}; - if (dist < min_color_dist) { - min_color_dist = dist; - best_cluster = j; + if (nthreads > 1) { + for (int i = 0; i < nthreads; ++i) { + + int start_c{i * centroids_per_thread}; + int end_c{(i == nthreads - 1) ? k : (i + 1) * centroids_per_thread}; + + switch (color_space) { + case COLOR_SPACE_OPTION_RGB: { + threads.emplace_back( + _process_dist_per_centroid>, + std::cref(pixels), std::cref(centroids), std::ref(distances), + start_c, end_c); + break; + } + case COLOR_SPACE_OPTION_CIELAB: { + threads.emplace_back( + _process_dist_per_centroid>, + std::cref(lab), std::cref(centroids_lab), std::ref(distances), + start_c, end_c); + break; + } } } - if (labels[i] != best_cluster) { - changed = true; - labels[i] = best_cluster; + // wait for threads to finish + for (auto &thread : threads) { + if (thread.joinable()) { + thread.join(); + } + } + // std::cout << "Computed distances" << std::endl; + threads.clear(); + + // assign labels + std::atomic changed_atomic{false}; + for (int i = 0; i < nthreads; ++i) { + int start_pixel{i * pixels_per_thread}; + int end_pixel{(i == nthreads - 1) ? num_pixels + : (i + 1) * pixels_per_thread}; + + switch (color_space) { + case COLOR_SPACE_OPTION_RGB: { + threads.emplace_back(_apply_labels>, + std::cref(pixels), std::cref(distances), + std::ref(labels), start_pixel, end_pixel, k, + std::ref(changed_atomic)); + break; + } + case COLOR_SPACE_OPTION_CIELAB: { + threads.emplace_back(_apply_labels>, + std::cref(lab), std::cref(distances), + std::ref(labels), start_pixel, end_pixel, k, + std::ref(changed_atomic)); + break; + } + } + } + // wait for threads to finish + for (auto &thread : threads) { + if (thread.joinable()) { + thread.join(); + } } - } + changed = changed_atomic.load(); + threads.clear(); + } else { + // Iterate over pixels + for (int32_t i{0}; i < num_pixels; ++i) { + float min_color_dist{std::numeric_limits::max()}; + int32_t best_cluster{0}; + + // Iterate over centroids to find centroid with most similar color to + // pixels[i] + float dist; + for (int32_t j{0}; j < k; ++j) { + // float dist{colorDistance(pixels[i], centroids[j])}; + switch (color_space) { + case COLOR_SPACE_OPTION_RGB: { + dist = colorDistance(pixels[i], centroids[j]); + break; + } + case COLOR_SPACE_OPTION_CIELAB: { + dist = colorDistance(lab[i], centroids_lab[j]); + break; + } + } + if (dist < min_color_dist) { + min_color_dist = dist; + best_cluster = j; + } + } + + if (labels[i] != best_cluster) { + changed = true; + labels[i] = best_cluster; + } + } + } // Stop if no changes if (!changed) { break; @@ -68,13 +322,25 @@ void kmeans(const uint8_t *data, uint8_t *out_data, int32_t *out_labels, // Update step ImageLib::Image> new_centroids(k, 1, 0); + ImageLib::Image> new_centroids_lab(k, 1, 0); std::vector counts(k, 0); for (int32_t i = 0; i < num_pixels; ++i) { int32_t cluster = labels[i]; - new_centroids[cluster].red += pixels[i].red; - new_centroids[cluster].green += pixels[i].green; - new_centroids[cluster].blue += pixels[i].blue; + switch (color_space) { + case COLOR_SPACE_OPTION_RGB: { + new_centroids[cluster].red += pixels[i].red; + new_centroids[cluster].green += pixels[i].green; + new_centroids[cluster].blue += pixels[i].blue; + break; + } + case COLOR_SPACE_OPTION_CIELAB: { + new_centroids_lab[cluster].l += lab[i].l; + new_centroids_lab[cluster].a += lab[i].a; + new_centroids_lab[cluster].b += lab[i].b; + break; + } + } counts[cluster]++; } @@ -84,13 +350,30 @@ void kmeans(const uint8_t *data, uint8_t *out_data, int32_t *out_labels, to it. May be good idea to reinitialize these dead centroids. */ if (counts[j] > 0) { - centroids[j].red = new_centroids[j].red / counts[j]; - centroids[j].green = new_centroids[j].green / counts[j]; - centroids[j].blue = new_centroids[j].blue / counts[j]; + switch (color_space) { + case COLOR_SPACE_OPTION_RGB: { + centroids[j].red = new_centroids[j].red / counts[j]; + centroids[j].green = new_centroids[j].green / counts[j]; + centroids[j].blue = new_centroids[j].blue / counts[j]; + break; + } + case COLOR_SPACE_OPTION_CIELAB: { + centroids_lab[j].l = new_centroids_lab[j].l / counts[j]; + centroids_lab[j].a = new_centroids_lab[j].a / counts[j]; + centroids_lab[j].b = new_centroids_lab[j].b / counts[j]; + break; + } + } } } } + if (color_space == COLOR_SPACE_OPTION_CIELAB) { + for (int32_t i{0}; i < k; ++i) { + lab_to_rgb(centroids_lab[i], centroids[i]); + } + } + // Write the final centroid values to each pixel in the cluster for (int32_t i = 0; i < num_pixels; ++i) { const int32_t cluster = labels[i]; @@ -111,6 +394,7 @@ void kmeans(const uint8_t *data, uint8_t *out_data, int32_t *out_labels, *should actually be taken into account. */ +/* Remove for now struct RGBXY { float r, g, b; float x, y; @@ -211,3 +495,4 @@ void kmeans_clustering_spatial(uint8_t *data, int32_t width, int32_t height, data[i * 4 + 2] = static_cast(centroids[cluster].b * 255); } } +*/ \ No newline at end of file diff --git a/vite.config.js b/vite.config.js index 2fd1937e3..9ee11f951 100644 --- a/vite.config.js +++ b/vite.config.js @@ -46,6 +46,10 @@ export default defineConfig({ watch: { ignored: ['**/docs/**', 'src/wasm/**/*.js', 'src/wasm/**/*.wasm'], }, + headers: { + 'Cross-Origin-Opener-Policy': 'same-origin', + 'Cross-Origin-Embedder-Policy': 'require-corp', + }, }, // Ensure Vite copies .wasm files @@ -109,5 +113,16 @@ export default defineConfig({ files.forEach((file) => server.watcher.add(file)); }, }, + // for multithreading + { + name: 'force-security-headers', + configureServer(server) { + server.middlewares.use((_req, res, next) => { + res.setHeader('Cross-Origin-Opener-Policy', 'same-origin'); + res.setHeader('Cross-Origin-Embedder-Policy', 'require-corp'); + next(); + }); + }, + }, ], });