From e780c38282f46b68c29ddc19ad9d03700099df60 Mon Sep 17 00:00:00 2001 From: Sergey Fedorov Date: Thu, 30 Jul 2026 04:18:22 +0000 Subject: [PATCH 29/38] swresample/ppc: add AltiVec resample_linear float + S16P common/linear kernels Extend the PPC resample DSP (previously FLTP resample_common only) with three more kernels, all bit-for-bit faithful to resample_template.c: * resample_common int16 (S16P): 8 s16 taps/iteration via vec_msum (modular, NOT the saturating vec_msums) into 4 s32 lanes, horizontal-summed in int64; FOFFSET added once, av_clip_int16(sum >> 15). The template widens to int64 only at the final combine, so with realistic filter banks no int32 partial ever wraps and the reassociated 8-tap dot is bit-exact vs. the C even/odd split (checkasm's int oracle is exact memcmp). * resample_linear int16 (S16P): two filter rows (phase index and index+1) against one src window; accumulate in uint32 (FELEM2U, mod 2^32, fully associative -> bit-exact unconditionally). Interpolation term copied verbatim from the template, including the signed reinterpretation: val += (int32)(v2 - val) * (int64)frac / c->src_incr. * resample_linear float (FLTP): two vec_madd accumulator sets, then the FILTER_SHIFT==0 branch with inv_src_incr = 1.0/src_incr hoisted out of the loop; interpolation expression copied verbatim (float->double->float). Same alignment discipline as the existing FLTP common kernel: filter rows are 16-byte aligned (plain vec_ld), src is arbitrary so the vec_lvsl splice is recomputed every iteration, and the vector loop stops at filter_length & ~7 (s16) / & ~3 (float) with a scalar tail -- never running into the zero padding, which would read src[] past the C reference's window. S32P (FELEM2 = int64, no 32x32 or s64 in classic AltiVec) and DBLP (no double vector pre-VSX) stay on C on this target and are not attempted. Registered unconditionally (no per-CPU _ARCH_PWR4 gating -- no real-hardware bench data for these yet). checkasm swr_resample: all three new kernels dispatch to AltiVec (s32p/dblp correctly stay on C), pass seeds 1/42/777; full suite 587/587 (was 584, +3). Co-Authored-By: Claude Opus 4.8 (1M context) Claude-Session: https://claude.ai/code/session_017TAQsLz1uCu9jURTT9xdzM --- libswresample/ppc/resample.c | 280 +++++++++++++++++++++++++++++++++++ 1 file changed, 280 insertions(+) diff --git a/libswresample/ppc/resample.c b/libswresample/ppc/resample.c index 8ed4a6f..0724e9a 100644 --- a/libswresample/ppc/resample.c +++ b/libswresample/ppc/resample.c @@ -22,6 +22,7 @@ #include "config.h" #include "libavutil/attributes.h" +#include "libavutil/common.h" #include "libavutil/cpu.h" #include "libavutil/mem_internal.h" #include "libavutil/ppc/cpu.h" @@ -67,6 +68,108 @@ static float apply_filter_x4_float_altivec(const float *src, const float *filter return sums[0] + sums[1] + sums[2] + sums[3]; } +/* + * Same as above but drives TWO filter rows (the interpolation phase and the + * next one, filter and filter + filter_alloc, both 16-aligned) against the one + * unaligned src window in a single pass: resample_linear needs both dot + * products for every output sample. The horizontal sums are returned via *out0 + * / *out1. + */ +static void apply_filter_x4_float_dual_altivec(const float *src, + const float *filter0, + const float *filter1, + int length, + float *out0, float *out1) +{ + DECLARE_ALIGNED(16, float, sums0)[4]; + DECLARE_ALIGNED(16, float, sums1)[4]; + vec_f acc0 = (vec_f)vec_splat_u32(0); + vec_f acc1 = (vec_f)vec_splat_u32(0); + int off; + + for (off = 0; off < length * 4; off += 16) { + vec_f s = vec_perm(vec_ld(off, src), vec_ld(off + 15, src), + vec_lvsl(off, src)); + acc0 = vec_madd(s, vec_ld(off, filter0), acc0); + acc1 = vec_madd(s, vec_ld(off, filter1), acc1); + } + + vec_st(acc0, 0, sums0); + vec_st(acc1, 0, sums1); + *out0 = sums0[0] + sums0[1] + sums0[2] + sums0[3]; + *out1 = sums1[0] + sums1[1] + sums1[2] + sums1[3]; +} + +/* + * int16 leaf: 8 s16 taps per iteration via vec_msum (s16 x s16 -> pairs summed + * into 4 s32 lanes). The C reference keeps the tap products in int32 and only + * widens to int64 at the final val + (int64)val2 combine, so this leaf returns + * the four lanes summed in int64 -- with realistic filter banks (coefficient + * magnitudes summing to ~2^15 times small ripple) no int32 partial ever wraps, + * so the different association vs. the C even/odd split is still bit-exact + * (checkasm's int oracle is exact memcmp, which proves it). + * + * vec_msum is the MODULAR multiply-sum; do NOT use vec_msums, which saturates + * -- the template semantics are plain modular int32 arithmetic. + * + * filter row is 16-aligned (see the float note above; filter_alloc*2 bytes is a + * multiple of 16 for s16 too). src (int16) is arbitrarily aligned, so the same + * per-iteration vec_lvsl splice as the float leaf brings in 8 int16s; off + * advances 16 bytes = 8 samples, changing alignment mod 16 every step, so the + * permute is recomputed each iteration (never reused across steps). + */ +static int64_t apply_filter_x8_s16_altivec(const int16_t *src, + const int16_t *filter, int length) +{ + DECLARE_ALIGNED(16, int32_t, sums)[4]; + vec_s32 acc = vec_splat_s32(0); + int off; + + for (off = 0; off < length * 2; off += 16) { + vec_s16 s = vec_perm(vec_ld(off, src), vec_ld(off + 15, src), + vec_lvsl(off, src)); + vec_s16 f = vec_ld(off, filter); /* filter row is 16-aligned */ + acc = vec_msum(s, f, acc); + } + + vec_st(acc, 0, sums); + return (int64_t)sums[0] + sums[1] + sums[2] + sums[3]; +} + +/* + * int16 linear leaf: two filter rows against one src window, accumulating in + * uint32 lanes. resample_linear's C reference accumulates in FELEM2U (uint32, + * i.e. mod 2^32) which is fully associative, so the 4-lane msum reduction plus + * a wrapping (uint32) horizontal sum is bit-exact UNCONDITIONALLY here -- + * regardless of coefficient magnitudes. vec_msum on s16 lanes produces signed + * int32 products summed with two's-complement wraparound; reinterpreting those + * lanes as uint32 and summing mod 2^32 matches the C exactly. + */ +static void apply_filter_x8_s16_dual_altivec(const int16_t *src, + const int16_t *filter0, + const int16_t *filter1, + int length, + uint32_t *out0, uint32_t *out1) +{ + DECLARE_ALIGNED(16, int32_t, sums0)[4]; + DECLARE_ALIGNED(16, int32_t, sums1)[4]; + vec_s32 acc0 = vec_splat_s32(0); + vec_s32 acc1 = vec_splat_s32(0); + int off; + + for (off = 0; off < length * 2; off += 16) { + vec_s16 s = vec_perm(vec_ld(off, src), vec_ld(off + 15, src), + vec_lvsl(off, src)); + acc0 = vec_msum(s, vec_ld(off, filter0), acc0); + acc1 = vec_msum(s, vec_ld(off, filter1), acc1); + } + + vec_st(acc0, 0, sums0); + vec_st(acc1, 0, sums1); + *out0 = (uint32_t)sums0[0] + sums0[1] + sums0[2] + sums0[3]; + *out1 = (uint32_t)sums1[0] + sums1[1] + sums1[2] + sums1[3]; +} + static int resample_common_float_altivec(ResampleContext *c, void *dest, const void *source, int n, int update_ctx) { @@ -119,6 +222,178 @@ static int resample_common_float_altivec(ResampleContext *c, void *dest, return sample_index; } +static int resample_linear_float_altivec(ResampleContext *c, void *dest, + const void *source, int n, int update_ctx) +{ + float *dst = dest; + const float *src = source; + int dst_index; + int index = c->index; + int frac = c->frac; + int sample_index = 0; + int x4_len = c->filter_length & ~3; + double inv_src_incr = 1.0 / c->src_incr; + + while (index >= c->phase_count) { + sample_index++; + index -= c->phase_count; + } + + for (dst_index = 0; dst_index < n; dst_index++) { + const float *filter = ((const float *)c->filter_bank) + + c->filter_alloc * index; + float val = 0, v2 = 0; + int i = 0; + + if (x4_len >= 4) { + apply_filter_x4_float_dual_altivec(&src[sample_index], filter, + filter + c->filter_alloc, + x4_len, &val, &v2); + i = x4_len; + } + for (; i < c->filter_length; i++) { + val += src[sample_index + i] * filter[i]; + v2 += src[sample_index + i] * filter[i + c->filter_alloc]; + } + + val += (v2 - val) * inv_src_incr * frac; + dst[dst_index] = val; + + frac += c->dst_incr_mod; + index += c->dst_incr_div; + if (frac >= c->src_incr) { + frac -= c->src_incr; + index++; + } + + while (index >= c->phase_count) { + sample_index++; + index -= c->phase_count; + } + } + + if (update_ctx) { + c->frac = frac; + c->index = index; + } + + return sample_index; +} + +static int resample_common_int16_altivec(ResampleContext *c, void *dest, + const void *source, int n, int update_ctx) +{ + int16_t *dst = dest; + const int16_t *src = source; + int dst_index; + int index = c->index; + int frac = c->frac; + int sample_index = 0; + int x8_len = c->filter_length & ~7; + + while (index >= c->phase_count) { + sample_index++; + index -= c->phase_count; + } + + for (dst_index = 0; dst_index < n; dst_index++) { + const int16_t *filter = ((const int16_t *)c->filter_bank) + + c->filter_alloc * index; + int64_t val = 0; + int i = 0; + + if (x8_len >= 8) { + val = apply_filter_x8_s16_altivec(&src[sample_index], filter, x8_len); + i = x8_len; + } + for (; i < c->filter_length; i++) + val += src[sample_index + i] * (int32_t)filter[i]; + + val += 1 << 14; /* FOFFSET */ + dst[dst_index] = av_clip_int16(val >> 15); + + frac += c->dst_incr_mod; + index += c->dst_incr_div; + if (frac >= c->src_incr) { + frac -= c->src_incr; + index++; + } + + while (index >= c->phase_count) { + sample_index++; + index -= c->phase_count; + } + } + + if (update_ctx) { + c->frac = frac; + c->index = index; + } + + return sample_index; +} + +static int resample_linear_int16_altivec(ResampleContext *c, void *dest, + const void *source, int n, int update_ctx) +{ + int16_t *dst = dest; + const int16_t *src = source; + int dst_index; + int index = c->index; + int frac = c->frac; + int sample_index = 0; + int x8_len = c->filter_length & ~7; + + while (index >= c->phase_count) { + sample_index++; + index -= c->phase_count; + } + + for (dst_index = 0; dst_index < n; dst_index++) { + const int16_t *filter = ((const int16_t *)c->filter_bank) + + c->filter_alloc * index; + /* uint32, mod 2^32 -- matches the template's FELEM2U accumulators. */ + uint32_t val = 1 << 14, v2 = 1 << 14; /* FOFFSET */ + int i = 0; + + if (x8_len >= 8) { + uint32_t a0, a1; + apply_filter_x8_s16_dual_altivec(&src[sample_index], filter, + filter + c->filter_alloc, + x8_len, &a0, &a1); + val += a0; + v2 += a1; + i = x8_len; + } + for (; i < c->filter_length; i++) { + val += src[sample_index + i] * (int32_t)filter[i]; + v2 += src[sample_index + i] * (int32_t)filter[i + c->filter_alloc]; + } + + val += (int32_t)(v2 - val) * (int64_t)frac / c->src_incr; + dst[dst_index] = av_clip_int16((int32_t)val >> 15); + + frac += c->dst_incr_mod; + index += c->dst_incr_div; + if (frac >= c->src_incr) { + frac -= c->src_incr; + index++; + } + + while (index >= c->phase_count) { + sample_index++; + index -= c->phase_count; + } + } + + if (update_ctx) { + c->frac = frac; + c->index = index; + } + + return sample_index; +} + #endif /* HAVE_ALTIVEC && HAVE_BIGENDIAN */ av_cold void swri_resample_dsp_ppc_init(ResampleContext *c) @@ -130,6 +405,11 @@ av_cold void swri_resample_dsp_ppc_init(ResampleContext *c) switch (c->format) { case AV_SAMPLE_FMT_FLTP: c->dsp.resample_common = resample_common_float_altivec; + c->dsp.resample_linear = resample_linear_float_altivec; + break; + case AV_SAMPLE_FMT_S16P: + c->dsp.resample_common = resample_common_int16_altivec; + c->dsp.resample_linear = resample_linear_int16_altivec; break; } #endif -- 2.43.0