fix chroma resampling

NEWS.md

@@ -1,3 +1,12 @@
+unreleased
+----------
+
+- the encoder averages subsampled chroma (was point sampling) and the
+decoder interpolates 4:2:0, 4:2:2 and 4:4:0 chroma like libjpeg for
+full-resolution output (was replicated; also for chroma <= 2 samples wide,
+unlike libjpeg) - output changes, the previous behavior can be selected by
+GPUJPEG_ENC_OPT_CHROMA_FILTER_BOOL and GPUJPEG_DEC_OPT_CHROMA_FILTER_BOOL
+
 2026-03-19 - 0.27.13
 ----------
 

libgpujpeg/gpujpeg_decoder.h

@@ -297,6 +297,10 @@ gpujpeg_decoder_get_image_info(uint8_t *image, size_t image_size, struct gpujpeg
 /// output image should be vertically flipped (bottom-up): values @ref GPUJPEG_VAL_TRUE or @ref GPUJPEG_VAL_FALSE
 #define GPUJPEG_DEC_OPT_FLIPPED_BOOL "dec_opt_flipped"
 
+/// interpolate 4:2:0, 4:2:2 and 4:4:0 chroma for full-resolution output (triangle filter as libjpeg "fancy
+/// upsampling") instead of replicating it; values @ref GPUJPEG_VAL_TRUE (default) or @ref GPUJPEG_VAL_FALSE
+#define GPUJPEG_DEC_OPT_CHROMA_FILTER_BOOL "dec_opt_chroma_filter"
+
 /// @copydoc GPUJPEG_ENC_OPT_CHANNEL_REMAP
 #define GPUJPEG_DEC_OPT_CHANNEL_REMAP "dec_opt_channel_remap"
 

libgpujpeg/gpujpeg_encoder.h

@@ -225,6 +225,9 @@ gpujpeg_encoder_suggest_restart_interval(const struct gpujpeg_image_parameters*
 
 /// input image is vertically flipped (bottom-up): values @ref GPUJPEG_VAL_TRUE or @ref GPUJPEG_VAL_FALSE
 #define GPUJPEG_ENC_OPT_FLIPPED_BOOL "enc_opt_flipped"
+/// average subsampled chroma over its sampling window (centered as JFIF assumes) instead of point sampling; values
+/// @ref GPUJPEG_VAL_TRUE (default) or @ref GPUJPEG_VAL_FALSE
+#define GPUJPEG_ENC_OPT_CHROMA_FILTER_BOOL "enc_opt_chroma_filter"
 /// custom exif tag in format <key>:TYPE=<value>
 #define GPUJPEG_ENC_OPT_EXIF_TAG "enc_exif_tag"
 /// set image orientation - syntax "<name>=<deg>[-]" or "help"; only if header supports (Exif, SPIFF)

src/gpujpeg_common.c

@@ -582,6 +582,7 @@ gpujpeg_coder_init(struct gpujpeg_coder * coder)
     coder->d_data_compressed = NULL;
     coder->d_temp_huffman = NULL;
     coder->data_compressed_allocated_size = 0;
+    coder->preprocessor.chroma_filter = true;
 
     GPUJPEG_CUSTOM_TIMER_CREATE(coder->duration_memory_to, return -1);
     GPUJPEG_CUSTOM_TIMER_CREATE(coder->duration_memory_from, return -1);

src/gpujpeg_decoder.c

@@ -504,6 +504,10 @@ gpujpeg_decoder_set_option(struct gpujpeg_decoder* decoder, const char *opt, con
     if ( strcmp(opt, GPUJPEG_DEC_OPT_CHANNEL_REMAP) == 0 ) {
         return gpujpeg_opt_set_channel_remap(&decoder->coder, val, GPUJPEG_DEC_OPT_CHANNEL_REMAP);
     }
+    if ( strcmp(opt, GPUJPEG_DEC_OPT_CHROMA_FILTER_BOOL) == 0 ) {
+        return gpujpeg_parse_bool_opt(&decoder->coder.preprocessor.chroma_filter, val,
+                                      GPUJPEG_DEC_OPT_CHROMA_FILTER_BOOL);
+    }
     if ( strcmp(opt, GPUJPEG_DEC_OPT_ALIGNMENT_BYTES_INT) == 0 ) {
         decoder->req_alignment = atoi(val);
         return GPUJPEG_NOERR;
@@ -521,6 +525,8 @@ gpujpeg_decoder_print_options()
            "] - vertically flip the output image\n");
     printf("\t" GPUJPEG_DEC_OPT_CHANNEL_REMAP "=XYZ[W] - ouput channel mapping, eg. '210F' for GBRX,\n"
         "\t\t'210' for GBR; special placeholders 'F' and 'Z' to set a channel to all-ones or all-zeros\n");
+    printf("\t" GPUJPEG_DEC_OPT_CHROMA_FILTER_BOOL "=[" GPUJPEG_VAL_FALSE "|" GPUJPEG_VAL_TRUE
+           "] - interpolate subsampled chroma (default) instead of replicating it\n");
     printf("\t" GPUJPEG_DEC_OPT_ALIGNMENT_BYTES_INT "=<num> - output alignment in bytes,\n"
         "\t\tshould not be used for gpujpegtool (cannot process)\n");
 }

src/gpujpeg_encoder.c

@@ -770,6 +770,10 @@ gpujpeg_encoder_set_option(struct gpujpeg_encoder* encoder, const char *opt, con
     if ( strcmp(opt, GPUJPEG_ENC_OPT_CHANNEL_REMAP) == 0 ) {
         return gpujpeg_opt_set_channel_remap(&encoder->coder, val, GPUJPEG_ENC_OPT_CHANNEL_REMAP);
     }
+    if ( strcmp(opt, GPUJPEG_ENC_OPT_CHROMA_FILTER_BOOL) == 0 ) {
+        return gpujpeg_parse_bool_opt(&encoder->coder.preprocessor.chroma_filter, val,
+                                      GPUJPEG_ENC_OPT_CHROMA_FILTER_BOOL);
+    }
     if ( strcmp(opt, GPUJPEG_ENC_OPT_EXIF_TAG) == 0 ) {
         encoder->header_type = GPUJPEG_HEADER_EXIF;
         return gpujpeg_exif_add_tag(&encoder->writer->exif_tags, val) ? GPUJPEG_NOERR : GPUJPEG_ERROR;
@@ -790,6 +794,8 @@ gpujpeg_encoder_print_options() {
            "] - whether is the input image should be vertically flipped (prior encode)\n");
     printf("\t" GPUJPEG_ENC_OPT_CHANNEL_REMAP "=XYZ[W] - input channel mapping, eg. '210F' for GBRX,\n"
         "\t\t'210' for GBR; special placeholders 'F' and 'Z' to set a channel to all-ones or all-zeros\n");
+    printf("\t" GPUJPEG_ENC_OPT_CHROMA_FILTER_BOOL "=[" GPUJPEG_VAL_FALSE "|" GPUJPEG_VAL_TRUE
+           "] - average subsampled chroma (default) instead of point sampling\n");
     printf("\t" GPUJPEG_ENC_OPT_EXIF_TAG "=<key>=<value>|help - custom EXIF tag (use help for syntax)\n");
     printf("\t" GPUJPEG_ENC_OPT_METADATA "=<key>=<value>|help - set image metadata\n");
 }

src/gpujpeg_postprocessor.cu

@@ -112,6 +112,31 @@ struct gpujpeg_preprocessor_comp_to_raw_load
     }
 };
 
+/**
+ * Interpolate a 2x subsampled component at full-resolution pixel (x, y) like libjpeg "fancy" upsampling (triangle
+ * filter, weights 3:1 per subsampled axis); neighbours are clamped to the component edges.
+ */
+template <int fh, int fv>
+static __device__ uint8_t
+gpujpeg_postprocessor_chroma_interpolate(const struct gpujpeg_preprocessor_data_component& comp, int x, int y)
+{
+    const uint8_t* d = comp.d_data;
+    const int pitch = comp.data_width;
+    const int cx = x / fh;
+    const int cy = y / fv;
+    const int nx = fh == 2 ? ((x & 1) ? min(cx + 1, comp.width - 1) : max(cx - 1, 0)) : cx;
+    const int ny = fv == 2 ? ((y & 1) ? min(cy + 1, comp.height - 1) : max(cy - 1, 0)) : cy;
+    if ( fh == 2 && fv == 2 ) {
+        const int this_col = 3 * d[cy * pitch + cx] + d[ny * pitch + cx];
+        const int next_col = 3 * d[cy * pitch + nx] + d[ny * pitch + nx];
+        return (3 * this_col + next_col + ((x & 1) ? 7 : 8)) >> 4;
+    }
+    if ( fh == 2 ) {
+        return (3 * d[cy * pitch + cx] + d[cy * pitch + nx] + ((x & 1) ? 2 : 1)) >> 2;
+    }
+    return (3 * d[cy * pitch + cx] + d[ny * pitch + cx] + ((y & 1) ? 2 : 1)) >> 2;
+}
+
 /// set alpha that may not be included in input data but may be read by the CS conv
 template <enum gpujpeg_pixel_format pixel_format>
 struct pre_load {
@@ -207,6 +232,20 @@ gpujpeg_preprocessor_comp_to_raw_kernel(struct gpujpeg_preprocessor_data data, u
     gpujpeg_preprocessor_comp_to_raw_load<s_comp1_samp_factor_h, s_comp1_samp_factor_v, s_comp2_samp_factor_h, s_comp2_samp_factor_v, s_comp3_samp_factor_h, s_comp3_samp_factor_v, s_comp4_samp_factor_h, s_comp4_samp_factor_v>::perform(r, image_position_x, image_position_y, data);
     post_load<color_space_internal == GPUJPEG_RGB, s_comp2_samp_factor_h>::perform(r, data);
 
+    // Interpolate 2x subsampled chroma (Cb and Cr with equal sampling) for full-resolution output, precompiled
+    // kernels only
+    constexpr bool full_chroma_output = pixel_format == GPUJPEG_444_U8_P012 || pixel_format == GPUJPEG_4444_U8_P0123 ||
+                                        pixel_format == GPUJPEG_444_U8_P0P1P2;
+    constexpr bool interpolate =
+        s_comp2_samp_factor_h <= 2 && s_comp2_samp_factor_v <= 2 && s_comp2_samp_factor_h * s_comp2_samp_factor_v > 1 &&
+        s_comp2_samp_factor_h == s_comp3_samp_factor_h && s_comp2_samp_factor_v == s_comp3_samp_factor_v;
+    if ( full_chroma_output && interpolate && data.chroma_filter ) {
+        constexpr int fh = interpolate ? s_comp2_samp_factor_h : 1;
+        constexpr int fv = interpolate ? s_comp2_samp_factor_v : 1;
+        r.y = gpujpeg_postprocessor_chroma_interpolate<fh, fv>(data.comp[1], image_position_x, image_position_y);
+        r.z = gpujpeg_postprocessor_chroma_interpolate<fh, fv>(data.comp[2], image_position_x, image_position_y);
+    }
+
     // Color transform
     gpujpeg_color_transform<color_space_internal, color_space>::perform(r);
 
@@ -359,6 +398,8 @@ gpujpeg_postprocessor_decoder_init(struct gpujpeg_coder* coder)
         data->comp[comp].sampling_factor.horizontal = coder->sampling_factor.horizontal / coder->component[comp].sampling_factor.horizontal;
         data->comp[comp].sampling_factor.vertical = coder->sampling_factor.vertical / coder->component[comp].sampling_factor.vertical;
         data->comp[comp].data_width = coder->component[comp].data_width;
+        data->comp[comp].width = coder->component[comp].width;
+        data->comp[comp].height = coder->component[comp].height;
     }
 
     if (!gpujpeg_pixel_format_is_interleaved(coder->param_image.pixel_format) &&
@@ -481,6 +522,7 @@ gpujpeg_postprocessor_decode(struct gpujpeg_coder* coder, cudaStream_t stream)
     }
 
     // Run kernel
+    coder->preprocessor.data.chroma_filter = coder->preprocessor.chroma_filter;
     kernel<<<grid, threads, 0, stream>>>(
         coder->preprocessor.data,
         coder->d_data_raw,

src/gpujpeg_preprocessor.cu

@@ -149,6 +149,38 @@ inline __device__ void raw_to_comp_load<GPUJPEG_422_U8_P1020>(const uint8_t* d_d
     }
 }
 
+/**
+ * Set chroma of the top-left pixel of each subsampling window to the window average, so that the stored sample is
+ * centered as JFIF assumes (instead of point sampling). The window is clamped at the right and bottom image edges.
+ */
+template <enum gpujpeg_color_space color_space_internal, enum gpujpeg_color_space color_space,
+          enum gpujpeg_pixel_format pixel_format>
+static __device__ void
+gpujpeg_preprocessor_chroma_box_filter(uchar4& r, const uint8_t* d_data_raw, int image_width_padding, int image_width,
+                                       int image_height, int x, int y, int fh, int fv)
+{
+    if ( x % fh != 0 || y % fv != 0 ) {
+        return;
+    }
+    int sum_cb = 0;
+    int sum_cr = 0;
+    for ( int dy = 0; dy < fv; ++dy ) {
+        for ( int dx = 0; dx < fh; ++dx ) {
+            int xx = min(x + dx, image_width - 1);
+            int yy = min(y + dy, image_height - 1);
+            int offset = (yy * image_width + xx) * unit_size<pixel_format>() + image_width_padding * yy;
+            uchar4 q;
+            raw_to_comp_load<pixel_format>(d_data_raw, image_width, image_height, offset, xx, yy, q);
+            gpujpeg_color_transform<color_space, color_space_internal>::perform(q);
+            sum_cb += q.y;
+            sum_cr += q.z;
+        }
+    }
+    const int n = fh * fv;
+    r.y = (sum_cb + n / 2) / n;
+    r.z = (sum_cr + n / 2) / n;
+}
+
 /**
  * Kernel - Copy raw image source data into three separated component buffers
  */
@@ -192,6 +224,25 @@ gpujpeg_preprocessor_raw_to_comp_kernel(struct gpujpeg_preprocessor_data data, c
     // Color transform
     gpujpeg_color_transform<color_space, color_space_internal>::perform(r);
 
+    // Average chroma if subsampled from full resolution (Cb and Cr with equal sampling)
+    constexpr bool full_chroma_input = pixel_format == GPUJPEG_444_U8_P012 || pixel_format == GPUJPEG_4444_U8_P0123 ||
+                                       pixel_format == GPUJPEG_444_U8_P0P1P2;
+    if ( full_chroma_input && data.chroma_filter ) {
+        const int fh =
+            s_comp2_samp_factor_h == GPUJPEG_DYNAMIC ? data.comp[1].sampling_factor.horizontal : s_comp2_samp_factor_h;
+        const int fv =
+            s_comp2_samp_factor_v == GPUJPEG_DYNAMIC ? data.comp[1].sampling_factor.vertical : s_comp2_samp_factor_v;
+        const int fh_cr =
+            s_comp3_samp_factor_h == GPUJPEG_DYNAMIC ? data.comp[2].sampling_factor.horizontal : s_comp3_samp_factor_h;
+        const int fv_cr =
+            s_comp3_samp_factor_v == GPUJPEG_DYNAMIC ? data.comp[2].sampling_factor.vertical : s_comp3_samp_factor_v;
+        if ( fh * fv > 1 && fh == fh_cr && fv == fv_cr ) {
+            gpujpeg_preprocessor_chroma_box_filter<color_space_internal, color_space, pixel_format>(
+                r, d_data_raw, image_width_padding, image_width, image_height, image_position_x, image_position_y, fh,
+                fv);
+        }
+    }
+
     // Store
     gpujpeg_preprocessor_raw_to_comp_store<s_comp1_samp_factor_h, s_comp1_samp_factor_v, s_comp2_samp_factor_h,
                                            s_comp2_samp_factor_v, s_comp3_samp_factor_h, s_comp3_samp_factor_v,
@@ -394,6 +445,7 @@ gpujpeg_preprocessor_encode_interlaced(struct gpujpeg_encoder * encoder)
     gpujpeg_const_div_prepare(image_width, width_div_mul, width_div_shift);
 
     // Run kernel
+    coder->preprocessor.data.chroma_filter = coder->preprocessor.chroma_filter;
     kernel<<<grid, threads, 0, coder->stream>>>(
         coder->preprocessor.data,
         coder->d_data_raw,

src/gpujpeg_preprocessor.h

@@ -52,6 +52,8 @@ struct gpujpeg_preprocessor_data_component
 {
     uint8_t* d_data;
     int data_width;
+    int width;  ///< component width without MCU padding (set for decoder)
+    int height; ///< component height without MCU padding (set for decoder)
     struct gpujpeg_component_sampling_factor sampling_factor;
 };
 
@@ -61,6 +63,7 @@ struct gpujpeg_preprocessor_data_component
 struct gpujpeg_preprocessor_data
 {
     struct gpujpeg_preprocessor_data_component comp[GPUJPEG_MAX_COMPONENT_COUNT];
+    bool chroma_filter; ///< copy of gpujpeg_preprocessor::chroma_filter set before kernel launch
 };
 
 /**
@@ -70,6 +73,7 @@ struct gpujpeg_preprocessor
 {
     void* kernel;       // function poitner
     bool flipped;               ///< flip image before encode or after decode
+    bool chroma_filter;         ///< box-filter subsampled chroma (encode), interpolate it (decode); default true
     unsigned int channel_remap; ///< remap channels if != 0
                                 ///< format: count_8b | 00000000 | idx0_4b | idx1_4b | idx2_4b | idx3_4b
     struct gpujpeg_preprocessor_data data;

test/unit/run_tests.c

@@ -84,6 +84,10 @@ static void encode_gpu_mem_as_cpu() {
 void
 test_gh_95();
 
+// defined in test_chroma_filter.c
+void
+test_chroma_filter();
+
 static char*
 test_format_number_with_delim() {
         struct {
@@ -127,6 +131,7 @@ int main() {
         encode_gpu_mem_as_cpu();
         test_gh_95();
         test_format_number_with_delim();
+        test_chroma_filter();
         printf("PASSED\n");
 }