From a38137830d2917248784b4e07cd39c13da0f70f9 Mon Sep 17 00:00:00 2001 From: James Jackson-South Date: Wed, 26 Aug 2026 00:40:36 +1000 Subject: [PATCH] Implement SIMD-first HEVC intra prediction --- HEIF_IMPLEMENTATION_PLAN.md | 8 +- .../Hevc/HevcIntraPredictor.Operations.cs | 423 ++++++++++++++++++ ...cIntraPredictor.OperatorImplementations.cs | 280 ++++++++++++ .../Formats/Heif/Hevc/HevcIntraPredictor.cs | 412 +++++++++++++++++ .../Heif/HevcIntraPredictionBenchmarks.cs | 116 +++++ .../Heif/Hevc/HevcIntraPredictorTests.cs | 396 ++++++++++++++++ 6 files changed, 1634 insertions(+), 1 deletion(-) create mode 100644 src/ImageSharp/Formats/Heif/Hevc/HevcIntraPredictor.Operations.cs create mode 100644 src/ImageSharp/Formats/Heif/Hevc/HevcIntraPredictor.OperatorImplementations.cs create mode 100644 src/ImageSharp/Formats/Heif/Hevc/HevcIntraPredictor.cs create mode 100644 tests/ImageSharp.Benchmarks/Codecs/Heif/HevcIntraPredictionBenchmarks.cs create mode 100644 tests/ImageSharp.Tests/Formats/Heif/Hevc/HevcIntraPredictorTests.cs diff --git a/HEIF_IMPLEMENTATION_PLAN.md b/HEIF_IMPLEMENTATION_PLAN.md index ffe49bc14e..78f07240d2 100644 --- a/HEIF_IMPLEMENTATION_PLAN.md +++ b/HEIF_IMPLEMENTATION_PLAN.md @@ -29,7 +29,7 @@ Checkboxes may be marked complete only when the implementation and the verificat ## Active and queued work -- [ ] **Active:** complete the SIMD-first shared HEIF color pipeline, including ICC color management after codec color conversion and image composition. +- [ ] **Active:** complete HEVC still-picture reconstruction and connect its native planes to the verified shared SIMD-first HEIF color pipeline. - [x] Correct AV1 clipped-edge partition entropy handling and verify every block-mode cell from a real libavif AVIF against libaom 3.14.1. - [x] Preserve the exact ICC payload from an independently encoded AVIF primary item. - [x] Prove that a genuine non-sRGB AVIF profile changes decoded pixels and matches the independently converted source image within the documented AV1 tolerance. @@ -42,6 +42,8 @@ Checkboxes may be marked complete only when the implementation and the verificat - [x] Implement SIMD 4:2:0 and 4:2:2 encoder downsampling with odd-width and odd-height tails, and verify the stored 8/12-bit chroma samples against independently encoded full-resolution planes. - [x] Move the H.273 parameter resolver, closed color operators, transfer functions, sample loading/storage, chroma filtering, and RGB packing into one shared HEIF color pipeline used by both AV1 and HEVC rather than maintaining codec-specific arithmetic copies. - [x] Implement sequential pooled HEVC conversion in both directions for monochrome, 4:2:0, 4:2:2, and 4:4:4; independent luma/chroma precision; 8/10/12-bit full/limited ranges; and all six progressive 4:2:0 chroma sample locations. Eight-bit paths use JPEG's optimized RGB plane contracts, while high-bit-depth paths retain 16-bit packed pixels. + - [x] Implement allocation-free SIMD-first HEVC intra prediction for all 35 modes, 4/8/16/32 blocks, and 8/10/12-bit samples, including three-tap and strong-bilinear reference filtering, negative-angle reference extension, luma boundary filters, and SIMD horizontal transposition. Verify the normal and forced-scalar paths against one scalar oracle. + - [ ] Connect prediction to coding-tree and transform-unit traversal, reconstructed-plane reference availability, and the bounded item decoder lifecycle. - [ ] Complete the SIMD YUV/CICP paths for every supported AV1 bit depth, chroma format, range, matrix, transfer function, color primary, and chroma position, with scalar fallback only when hardware vectorization is unavailable or the operation is inherently non-vectorizable. - [x] Apply ICC conversion only after the SIMD YUV/CICP stage, alpha composition, grid assembly, and presentation transforms have produced the presented RGB image; retain ImageSharp's shared ICC converter and optimize reusable bulk kernels rather than creating a HEIF-specific color-management implementation. - [ ] Verify ICC preservation, conversion, compaction, and metadata skipping for grids, alpha-composited images, every presented sequence frame, and the completed HEVC path. @@ -243,6 +245,7 @@ This snapshot pins or classifies the available references and failures; it does | `HevcPlane` and `HevcPictureBuffer` | HEVC sections 6.2 and 6.3 source and decoded picture sample-array dimensions for monochrome, 4:2:0, 4:2:2, 4:4:4, and separate-color-plane coding | HM `source/Lib/TLibCommon/TComPicYuv.cpp` at `9c1f298659ab0cee9dc13d23d0304221575410b9`; Android `libhevc` decoded-picture buffer layout paths at `c83a76b084498d55f252f48b2e3786804cdf24b7`; ImageSharp `MemoryAllocator` and `Buffer2D` | Allocate one allocator-owned native `ushort` plane set for the selected still picture so the same reconstruction path preserves every supported 8-through-16-bit sample without byte-backed reinterpretation. Derive chroma dimensions by ceiling division for each sampling layout and treat separate color planes as full-resolution independently coded arrays. The buffer contains one picture only and introduces no decoded-picture buffer, reference lifetime, frame queue, or playback state. | | `HevcCodingTreeState` | HEVC sections 6.4 coding-tree block and coding-block availability plus section 9.3 split-flag context derivation | HM `source/Lib/TLibCommon/TComDataCU.cpp` function `getCtxSplitFlag` and `source/Lib/TLibDecoder/TDecCu.cpp` coding-tree traversal at `9c1f298659ab0cee9dc13d23d0304221575410b9`; Android `libhevc` `decoder/ihevcd_parse_slice.c` coding-tree paths at `c83a76b084498d55f252f48b2e3786804cdf24b7`; ImageSharp `MemoryAllocator` and `Buffer2D` | Store only leaf depth, effective luma QP, transquant-bypass, and PCM state at minimum-coding-block resolution for the selected picture. Derive split contexts from caller-approved left and above availability so slice and tile boundaries remain owned by traversal rather than hidden in a general block graph. Clip edge writes to the coded still-picture state and add no prediction units, motion fields, references, decoded-picture queue, or sequence lifetime. | | `HevcIntraPredictionState` | HEVC sections 8.4.2 and 9.3 intra luma/chroma prediction-mode derivation and binarization | HM `source/Lib/TLibCommon/TComDataCU.cpp` functions `getIntraDirPredictor` and `getAllowedChromaDir`, plus `source/Lib/TLibDecoder/TDecSbac.cpp` functions `parseIntraDirLumaAng` and `parseIntraDirChroma`, at `9c1f298659ab0cee9dc13d23d0304221575410b9`; Android `libhevc` intra-mode parsing and neighbor derivation paths at `c83a76b084498d55f252f48b2e3786804cdf24b7`; ImageSharp `MemoryAllocator` and `Buffer2D` | Preserve the normative two-pass luma flag/suffix order, spatial most-probable-mode derivation, omitted-mode reinsertion, four-way minimum-CU partition order, explicit chroma candidate substitution, and derived-chroma mode at 4x4 luma resolution. Accept slice/tile availability from traversal, fill only the selected still-picture map, and add no inter prediction unit, motion, reference, or sequence state. | +| `HevcIntraPredictor` | HEVC section 8.4.4.2 intra sample prediction and section 8.4.4.2.3 reference-sample substitution and filtering | HM `source/Lib/TLibCommon/TComPrediction.cpp` functions `predIntraGetPredValDC`, `xPredIntraAng`, `xPredIntraPlanar`, and `xDCPredFiltering`, plus `source/Lib/TLibCommon/TComPattern.cpp` function `fillReferenceSamples` and its reference-filter selection at `9c1f298659ab0cee9dc13d23d0304221575410b9`; Android `libhevc` intra-prediction kernels at `c83a76b084498d55f252f48b2e3786804cdf24b7` | Predict from caller-prepared references with planar, DC, and all 33 angular modes; negative-angle extension; three-tap and strong-bilinear smoothing; and the normative luma boundary filters. Traverse `Vector512`, `Vector256`, and `Vector128` widths before one scalar tail, use caller-owned reusable scratch for horizontal transposition, and retain no generic video, inter-picture, or reference-picture state. | | `HeifContentColorVolume`, `HeifItem.ContentColorVolume`, and `HeifDecoderCore` content color-volume parsing and presentation | HEIF content color-volume item property; AVIF 1.2 content color-volume requirements; ITU-T H.274 (V4) content colour volume syntax and semantics | libavif `src/read.c` function `avifSkipContentColourVolume` at `092276ce89098ead06db80975173191e5fee1826`; official ITU-T H.274 (V4), January 2026 | Decode only the bounded per-image `cclv` property: require zero cancellation, persistence, and reserved bits; preserve optional signed G/B/R primary coordinates and normalized minimum, maximum, and average luminance values; and validate their registered ranges and ordering. Expose the effective grid-or-tile still-image value through `HeifMetadata`. Do not add SEI persistence, retained video state, tracks, samples, timing, or a generic ISO BMFF color-volume box model. | | `HeifAmbientViewingEnvironment`, `HeifReferenceViewingEnvironment`, `HeifNominalDiffuseWhite`, and their per-item presentation metadata | ISOBMFF ambient viewing environment; ITU-T H.274 (V4) section 8.13; HEIF Amendment 1 sections 6.5.44 and 6.5.45; AVIF 1.2 image-item box requirements | libavif `src/read.c` functions `avifSkipAmbientViewingEnvironment`, `avifSkipReferenceViewingEnvironment`, and `avifSkipNominalDiffuseWhite` at `092276ce89098ead06db80975173191e5fee1826`; official ITU-T H.274 (V4), January 2026 | Decode the fixed-size `amve`, version-zero `reve`, and version-zero `ndwt` properties only when associated with the presented still-image item. Preserve physical illuminance and luminance units, distinct surround/periphery chromaticities, and the coded zero that requests the standard nominal diffuse-white default. Validate registered coordinate ranges and duplicate associations. Retain no video-SEI persistence, visual sample entry, display pipeline, track, timing, or generic viewing-environment box model. | | `GridHeifItemDecoder` and `HeifDecoderCore` grid/thumbnail selection | ISO/IEC 23008-12 section 6.6.2.3 image-grid syntax and MIAF grid-cell constraints | libavif `src/read.c` functions `avifParseImageGridBox`, `avifDecoderDataAllocateImagePlanes`, and `avifDecoderDataCopyTileToImage` at `092276ce89098ead06db80975173191e5fee1826` | Parse version-zero 16-bit and 32-bit grid descriptors, preserve row-major `dimg` order, require the declared tile count and one coding format, validate canvas coverage and edge overlap, and crop only the rightmost column and bottom row while copying through ImageSharp pixel buffers. A primary grid whose tile codec is unavailable may use only a decodable thumbnail that explicitly references that grid. | @@ -458,6 +461,8 @@ Implement and verify in dependency order: - [ ] CABAC arithmetic decoding and every required context transition. - [ ] Coding-tree, coding-unit, prediction-unit, and transform-unit traversal across all permitted sizes and partition modes. - [ ] Intra prediction for every luma and chroma mode, including strong intra smoothing and constrained prediction rules. + - [x] Implement and verify the allocation-free SIMD-first predictor primitive for all 35 modes, 4/8/16/32 blocks, and 8/10/12-bit samples. + - [ ] Build reference availability from reconstructed-plane and coding state, and connect the predictor to transform-unit traversal. - [ ] Scaling lists, inverse quantization, transform skip, every required inverse transform, range-extension precision, and lossless reconstruction. - [ ] Deblocking and sample-adaptive offset for every signaled luma/chroma and bit-depth path. - [ ] Tiles, wavefront entry points, dependent slices, and all other parallelization syntax permitted by the exposed still-image profiles. @@ -552,6 +557,7 @@ Tasks: - [ ] Remove known avoidable allocations first: per-transform arrays, the intermediate RGB image, repeated block scratch arrays, and file-sized buffering. - [x] Complete the active AV1 forward and inverse transform-family checklist above. - [ ] Benchmark codec-specific costs for CABAC/range decode, inverse transforms, still-image prediction, deblocking, SAO, CDEF, restoration, chroma upsampling, color conversion, alpha packing, and grid copies. + - [x] Add a permanent frame-wide HEVC intra-prediction benchmark. On .NET 10, SIMD planar, vertical-angular, and horizontal-angular prediction measured 155.2, 150.7, and 255.0 microseconds per padded 1920x1088 frame, compared with forced-scalar timings of 1.890, 1.081, and 0.979 milliseconds: 12.2, 7.2, and 3.8 times faster with zero managed allocations. - [ ] Implement vector paths only for confirmed hot loops, using existing `Vector128`, `Vector256`, and `Vector512` helper and dispatch patterns where supported. - [ ] Prioritize shared color conversion and pixel packing, chroma upsampling, inverse-transform add-and-clip, intra predictors, HEVC deblock/SAO, AV1 loop filter/CDEF/restoration, and contiguous grid copies. - [ ] Benchmark the complete decode color pipeline on representative 8/10/12-bit AVIF and HEIC images with and without embedded ICC profiles. Report absolute end-to-end timings and allocations in addition to the isolated YUV/CICP and ICC stage costs. diff --git a/src/ImageSharp/Formats/Heif/Hevc/HevcIntraPredictor.Operations.cs b/src/ImageSharp/Formats/Heif/Hevc/HevcIntraPredictor.Operations.cs new file mode 100644 index 0000000000..2e4c4c8169 --- /dev/null +++ b/src/ImageSharp/Formats/Heif/Hevc/HevcIntraPredictor.Operations.cs @@ -0,0 +1,423 @@ +// Copyright (c) Six Labors. +// Licensed under the Six Labors Split License. + +using System.Numerics; +using System.Runtime.CompilerServices; +using System.Runtime.InteropServices; +using System.Runtime.Intrinsics; +using SixLabors.ImageSharp.Common.Helpers; + +namespace SixLabors.ImageSharp.Formats.Heif.Hevc; + +/// +/// Provides shared SIMD operations used by the closed prediction operators. +/// +internal static partial class HevcIntraPredictor +{ + /// + /// Calculates one 512-bit half of a planar prediction row. + /// + /// The top reference samples. + /// The zero-based X coordinates. + /// The left reference sample for the row. + /// The top-right reference sample. + /// The bottom-left reference sample. + /// The top-reference weight. + /// The bottom-left-reference weight. + /// The square block side. + /// The division rounding constant. + /// The division shift. + /// The predicted samples as widened lanes. + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private static Vector512 CalculatePlanarVector( + Vector512 top, + Vector512 indices, + uint left, + uint topRight, + uint bottomLeft, + uint topWeight, + uint bottomWeight, + uint size, + uint rounding, + int shift) + { + Vector512 horizontal = ((Vector512.Create(size - 1) - indices) * left) + ((indices + Vector512.One) * topRight); + Vector512 vertical = (top * topWeight) + Vector512.Create(bottomLeft * bottomWeight); + return (horizontal + vertical + Vector512.Create(rounding)) >> shift; + } + + /// + /// Calculates one 256-bit half of a planar prediction row. + /// + /// The top reference samples. + /// The zero-based X coordinates. + /// The left reference sample for the row. + /// The top-right reference sample. + /// The bottom-left reference sample. + /// The top-reference weight. + /// The bottom-left-reference weight. + /// The square block side. + /// The division rounding constant. + /// The division shift. + /// The predicted samples as widened lanes. + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private static Vector256 CalculatePlanarVector( + Vector256 top, + Vector256 indices, + uint left, + uint topRight, + uint bottomLeft, + uint topWeight, + uint bottomWeight, + uint size, + uint rounding, + int shift) + { + Vector256 horizontal = ((Vector256.Create(size - 1) - indices) * left) + ((indices + Vector256.One) * topRight); + Vector256 vertical = (top * topWeight) + Vector256.Create(bottomLeft * bottomWeight); + return (horizontal + vertical + Vector256.Create(rounding)) >> shift; + } + + /// + /// Calculates one 128-bit half of a planar prediction row. + /// + /// The top reference samples. + /// The zero-based X coordinates. + /// The left reference sample for the row. + /// The top-right reference sample. + /// The bottom-left reference sample. + /// The top-reference weight. + /// The bottom-left-reference weight. + /// The square block side. + /// The division rounding constant. + /// The division shift. + /// The predicted samples as widened lanes. + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private static Vector128 CalculatePlanarVector( + Vector128 top, + Vector128 indices, + uint left, + uint topRight, + uint bottomLeft, + uint topWeight, + uint bottomWeight, + uint size, + uint rounding, + int shift) + { + Vector128 horizontal = ((Vector128.Create(size - 1) - indices) * left) + ((indices + Vector128.One) * topRight); + Vector128 vertical = (top * topWeight) + Vector128.Create(bottomLeft * bottomWeight); + return (horizontal + vertical + Vector128.Create(rounding)) >> shift; + } + + /// + /// Sums reconstructed reference samples without overflowing their 16-bit storage. + /// + /// The samples to sum. + /// The exact unsigned sum. + private static uint SumSamples(ReadOnlySpan samples) + { + ref ushort samplesBase = ref MemoryMarshal.GetReference(samples); + uint sum = 0; + int i = 0; + if (Vector512.IsHardwareAccelerated) + { + int oneVectorFromEnd = samples.Length - Vector512.Count; + for (; i <= oneVectorFromEnd; i += Vector512.Count) + { + (Vector512 low, Vector512 high) = Vector512.Widen(Vector512.LoadUnsafe(ref samplesBase, (nuint)i)); + sum += Vector512.Sum(low) + Vector512.Sum(high); + } + } + + if (Vector256.IsHardwareAccelerated) + { + int oneVectorFromEnd = samples.Length - Vector256.Count; + for (; i <= oneVectorFromEnd; i += Vector256.Count) + { + (Vector256 low, Vector256 high) = Vector256.Widen(Vector256.LoadUnsafe(ref samplesBase, (nuint)i)); + sum += Vector256.Sum(low) + Vector256.Sum(high); + } + } + + if (Vector128.IsHardwareAccelerated) + { + int oneVectorFromEnd = samples.Length - Vector128.Count; + for (; i <= oneVectorFromEnd; i += Vector128.Count) + { + (Vector128 low, Vector128 high) = Vector128.Widen(Vector128.LoadUnsafe(ref samplesBase, (nuint)i)); + sum += Vector128.Sum(low) + Vector128.Sum(high); + } + } + + for (; i < samples.Length; i++) + { + sum += Unsafe.Add(ref samplesBase, i); + } + + return sum; + } + + /// + /// Copies the top reference into every row and optionally filters the first column. + /// + /// The top reference samples. + /// The left reference samples. + /// The destination block origin. + /// The destination row stride. + /// The square block side. + /// The reconstructed component precision. + /// Whether the vertical luma edge filter applies. + private static void PredictVertical( + ReadOnlySpan top, + ReadOnlySpan left, + Span destination, + int destinationStride, + int size, + int bitDepth, + bool filterPredictionEdges) + { + ReadOnlySpan row = top.Slice(1, size); + int maximum = (1 << bitDepth) - 1; + for (int y = 0; y < size; y++) + { + row.CopyTo(destination.Slice(y * destinationStride, size)); + if (filterPredictionEdges) + { + int sample = destination[y * destinationStride] + ((left[y + 1] - left[0]) >> 1); + destination[y * destinationStride] = (ushort)Math.Clamp(sample, 0, maximum); + } + } + } + + /// + /// Fills each row from its left reference and optionally filters the first row. + /// + /// The top reference samples. + /// The left reference samples. + /// The destination block origin. + /// The destination row stride. + /// The square block side. + /// The reconstructed component precision. + /// Whether the horizontal luma edge filter applies. + private static void PredictHorizontal( + ReadOnlySpan top, + ReadOnlySpan left, + Span destination, + int destinationStride, + int size, + int bitDepth, + bool filterPredictionEdges) + { + for (int y = 0; y < size; y++) + { + destination.Slice(y * destinationStride, size).Fill(left[y + 1]); + } + + if (!filterPredictionEdges) + { + return; + } + + int maximum = (1 << bitDepth) - 1; + for (int x = 0; x < size; x++) + { + int sample = destination[x] + ((top[x + 1] - top[0]) >> 1); + destination[x] = (ushort)Math.Clamp(sample, 0, maximum); + } + } + + /// + /// Generates a vertical-oriented angular block using contiguous SIMD interpolation within each row. + /// + /// The main reference beginning at logical index zero. + /// The span index corresponding to logical reference index zero. + /// The contiguous destination or transposition scratch block. + /// The destination row stride. + /// The square block side. + /// The signed prediction displacement in thirty-second-sample units. + private static void PredictAngularRows( + ReadOnlySpan main, + int mainOrigin, + Span destination, + int destinationStride, + int size, + int angle) + { + for (int y = 0, deltaPosition = angle; y < size; y++, deltaPosition += angle) + { + int deltaInteger = deltaPosition >> 5; + int deltaFraction = deltaPosition & 31; + int sourceOffset = mainOrigin + deltaInteger + 1; + Span row = destination.Slice(y * destinationStride, size); + if (deltaFraction == 0) + { + main.Slice(sourceOffset, size).CopyTo(row); + } + else + { + InterpolateAngularRow(main[sourceOffset..], row, deltaFraction); + } + } + } + + /// + /// Interpolates one angular prediction row between consecutive main-reference samples. + /// + /// The first main-reference sample for the row. + /// The destination prediction row. + /// The right-hand weight with a denominator of thirty-two. + private static void InterpolateAngularRow(ReadOnlySpan source, Span destination, int fraction) + { + ref ushort sourceBase = ref MemoryMarshal.GetReference(source); + ref ushort destinationBase = ref MemoryMarshal.GetReference(destination); + uint leftWeight = (uint)(32 - fraction); + uint rightWeight = (uint)fraction; + int i = 0; + + if (Vector512.IsHardwareAccelerated) + { + int oneVectorFromEnd = destination.Length - Vector512.Count; + for (; i <= oneVectorFromEnd; i += Vector512.Count) + { + Vector512 left = Vector512.LoadUnsafe(ref sourceBase, (nuint)i); + Vector512 right = Vector512.LoadUnsafe(ref sourceBase, (nuint)(i + 1)); + (Vector512 leftLow, Vector512 leftHigh) = Vector512.Widen(left); + (Vector512 rightLow, Vector512 rightHigh) = Vector512.Widen(right); + Vector512 low = ((leftLow * leftWeight) + (rightLow * rightWeight) + Vector512.Create(16U)) >> 5; + Vector512 high = ((leftHigh * leftWeight) + (rightHigh * rightWeight) + Vector512.Create(16U)) >> 5; + Vector512.Narrow(low, high).StoreUnsafe(ref Unsafe.Add(ref destinationBase, i)); + } + } + + if (Vector256.IsHardwareAccelerated) + { + int oneVectorFromEnd = destination.Length - Vector256.Count; + for (; i <= oneVectorFromEnd; i += Vector256.Count) + { + Vector256 left = Vector256.LoadUnsafe(ref sourceBase, (nuint)i); + Vector256 right = Vector256.LoadUnsafe(ref sourceBase, (nuint)(i + 1)); + (Vector256 leftLow, Vector256 leftHigh) = Vector256.Widen(left); + (Vector256 rightLow, Vector256 rightHigh) = Vector256.Widen(right); + Vector256 low = ((leftLow * leftWeight) + (rightLow * rightWeight) + Vector256.Create(16U)) >> 5; + Vector256 high = ((leftHigh * leftWeight) + (rightHigh * rightWeight) + Vector256.Create(16U)) >> 5; + Vector256.Narrow(low, high).StoreUnsafe(ref Unsafe.Add(ref destinationBase, i)); + } + } + + if (Vector128.IsHardwareAccelerated) + { + int oneVectorFromEnd = destination.Length - Vector128.Count; + for (; i <= oneVectorFromEnd; i += Vector128.Count) + { + Vector128 left = Vector128.LoadUnsafe(ref sourceBase, (nuint)i); + Vector128 right = Vector128.LoadUnsafe(ref sourceBase, (nuint)(i + 1)); + (Vector128 leftLow, Vector128 leftHigh) = Vector128.Widen(left); + (Vector128 rightLow, Vector128 rightHigh) = Vector128.Widen(right); + Vector128 low = ((leftLow * leftWeight) + (rightLow * rightWeight) + Vector128.Create(16U)) >> 5; + Vector128 high = ((leftHigh * leftWeight) + (rightHigh * rightWeight) + Vector128.Create(16U)) >> 5; + Vector128.Narrow(low, high).StoreUnsafe(ref Unsafe.Add(ref destinationBase, i)); + } + } + + for (; i < destination.Length; i++) + { + Unsafe.Add(ref destinationBase, i) = (ushort)(((source[i] * leftWeight) + (source[i + 1] * rightWeight) + 16) >> 5); + } + } + + /// + /// Transposes a square horizontal prediction block into the reconstructed destination. + /// + /// The contiguous transposed prediction block. + /// The destination block origin. + /// The destination row stride. + /// The square block side. + private static void TransposeBlock(ReadOnlySpan source, Span destination, int destinationStride, int size) + { + if (Vector128.IsHardwareAccelerated && size >= Vector128.Count) + { + for (int y = 0; y < size; y += Vector128.Count) + { + for (int x = 0; x < size; x += Vector128.Count) + { + Transpose8x8(source, destination, destinationStride, size, x, y); + } + } + + return; + } + + for (int y = 0; y < size; y++) + { + for (int x = 0; x < size; x++) + { + destination[(x * destinationStride) + y] = source[(y * size) + x]; + } + } + } + + /// + /// Transposes one eight-by-eight tile of 16-bit prediction samples. + /// + /// The contiguous source block. + /// The destination block origin. + /// The destination row stride. + /// The contiguous source row stride. + /// The tile X coordinate in the source block. + /// The tile Y coordinate in the source block. + private static void Transpose8x8( + ReadOnlySpan source, + Span destination, + int destinationStride, + int sourceStride, + int x, + int y) + { + ref ushort sourceBase = ref MemoryMarshal.GetReference(source); + ref ushort destinationBase = ref MemoryMarshal.GetReference(destination); + Vector128 row0 = Vector128.LoadUnsafe(ref sourceBase, (nuint)(((y + 0) * sourceStride) + x)).AsInt16(); + Vector128 row1 = Vector128.LoadUnsafe(ref sourceBase, (nuint)(((y + 1) * sourceStride) + x)).AsInt16(); + Vector128 row2 = Vector128.LoadUnsafe(ref sourceBase, (nuint)(((y + 2) * sourceStride) + x)).AsInt16(); + Vector128 row3 = Vector128.LoadUnsafe(ref sourceBase, (nuint)(((y + 3) * sourceStride) + x)).AsInt16(); + Vector128 row4 = Vector128.LoadUnsafe(ref sourceBase, (nuint)(((y + 4) * sourceStride) + x)).AsInt16(); + Vector128 row5 = Vector128.LoadUnsafe(ref sourceBase, (nuint)(((y + 5) * sourceStride) + x)).AsInt16(); + Vector128 row6 = Vector128.LoadUnsafe(ref sourceBase, (nuint)(((y + 6) * sourceStride) + x)).AsInt16(); + Vector128 row7 = Vector128.LoadUnsafe(ref sourceBase, (nuint)(((y + 7) * sourceStride) + x)).AsInt16(); + + // Three zip stages exchange one, two, then four 16-bit coordinates. The resulting vectors are the eight + // source columns in row order, so each can be stored contiguously into one destination row. + Vector128 pair0 = Vector128_.UnpackLow(row0, row1); + Vector128 pair1 = Vector128_.UnpackHigh(row0, row1); + Vector128 pair2 = Vector128_.UnpackLow(row2, row3); + Vector128 pair3 = Vector128_.UnpackHigh(row2, row3); + Vector128 pair4 = Vector128_.UnpackLow(row4, row5); + Vector128 pair5 = Vector128_.UnpackHigh(row4, row5); + Vector128 pair6 = Vector128_.UnpackLow(row6, row7); + Vector128 pair7 = Vector128_.UnpackHigh(row6, row7); + Vector128 quad0 = Vector128_.UnpackLow(pair0.AsInt32(), pair2.AsInt32()); + Vector128 quad1 = Vector128_.UnpackHigh(pair0.AsInt32(), pair2.AsInt32()); + Vector128 quad2 = Vector128_.UnpackLow(pair1.AsInt32(), pair3.AsInt32()); + Vector128 quad3 = Vector128_.UnpackHigh(pair1.AsInt32(), pair3.AsInt32()); + Vector128 quad4 = Vector128_.UnpackLow(pair4.AsInt32(), pair6.AsInt32()); + Vector128 quad5 = Vector128_.UnpackHigh(pair4.AsInt32(), pair6.AsInt32()); + Vector128 quad6 = Vector128_.UnpackLow(pair5.AsInt32(), pair7.AsInt32()); + Vector128 quad7 = Vector128_.UnpackHigh(pair5.AsInt32(), pair7.AsInt32()); + Vector128 column0 = Vector128_.UnpackLow(quad0.AsInt64(), quad4.AsInt64()).AsUInt16(); + Vector128 column1 = Vector128_.UnpackHigh(quad0.AsInt64(), quad4.AsInt64()).AsUInt16(); + Vector128 column2 = Vector128_.UnpackLow(quad1.AsInt64(), quad5.AsInt64()).AsUInt16(); + Vector128 column3 = Vector128_.UnpackHigh(quad1.AsInt64(), quad5.AsInt64()).AsUInt16(); + Vector128 column4 = Vector128_.UnpackLow(quad2.AsInt64(), quad6.AsInt64()).AsUInt16(); + Vector128 column5 = Vector128_.UnpackHigh(quad2.AsInt64(), quad6.AsInt64()).AsUInt16(); + Vector128 column6 = Vector128_.UnpackLow(quad3.AsInt64(), quad7.AsInt64()).AsUInt16(); + Vector128 column7 = Vector128_.UnpackHigh(quad3.AsInt64(), quad7.AsInt64()).AsUInt16(); + column0.StoreUnsafe(ref destinationBase, (nuint)(((x + 0) * destinationStride) + y)); + column1.StoreUnsafe(ref destinationBase, (nuint)(((x + 1) * destinationStride) + y)); + column2.StoreUnsafe(ref destinationBase, (nuint)(((x + 2) * destinationStride) + y)); + column3.StoreUnsafe(ref destinationBase, (nuint)(((x + 3) * destinationStride) + y)); + column4.StoreUnsafe(ref destinationBase, (nuint)(((x + 4) * destinationStride) + y)); + column5.StoreUnsafe(ref destinationBase, (nuint)(((x + 5) * destinationStride) + y)); + column6.StoreUnsafe(ref destinationBase, (nuint)(((x + 6) * destinationStride) + y)); + column7.StoreUnsafe(ref destinationBase, (nuint)(((x + 7) * destinationStride) + y)); + } +} diff --git a/src/ImageSharp/Formats/Heif/Hevc/HevcIntraPredictor.OperatorImplementations.cs b/src/ImageSharp/Formats/Heif/Hevc/HevcIntraPredictor.OperatorImplementations.cs new file mode 100644 index 0000000000..1dfb66a28b --- /dev/null +++ b/src/ImageSharp/Formats/Heif/Hevc/HevcIntraPredictor.OperatorImplementations.cs @@ -0,0 +1,280 @@ +// Copyright (c) Six Labors. +// Licensed under the Six Labors Split License. + +using System.Numerics; +using System.Runtime.CompilerServices; +using System.Runtime.InteropServices; +using System.Runtime.Intrinsics; + +namespace SixLabors.ImageSharp.Formats.Heif.Hevc; + +/// +/// Provides the closed planar, DC, and angular prediction operators. +/// +internal static partial class HevcIntraPredictor +{ + /// + /// Implements planar interpolation between the top, left, bottom-left, and top-right references. + /// + private readonly struct PlanarPredictionOperator : IHevcIntraPredictionOperator + { + /// + public static void Predict( + ReadOnlySpan top, + ReadOnlySpan left, + Span destination, + int destinationStride, + int size, + int mode, + int bitDepth, + bool filterPredictionEdges, + Span scratch) + { + ref ushort topBase = ref MemoryMarshal.GetReference(top); + ref ushort destinationBase = ref MemoryMarshal.GetReference(destination); + uint bottomLeft = left[size]; + uint topRight = top[size]; + int shift = BitOperations.Log2((uint)size) + 1; + uint rounding = (uint)size; + + for (int y = 0; y < size; y++) + { + uint leftSample = left[y + 1]; + uint topWeight = (uint)(size - y - 1); + uint bottomWeight = (uint)(y + 1); + ref ushort rowBase = ref Unsafe.Add(ref destinationBase, y * destinationStride); + int x = 0; + + // The two widened halves carry consecutive X coordinates. Each lane evaluates the normative + // horizontal and vertical ramps, then narrows after the common rounded power-of-two division. + if (Vector512.IsHardwareAccelerated) + { + Vector512 indices = CreateIndicesVector512(); + int oneVectorFromEnd = size - Vector512.Count; + for (; x <= oneVectorFromEnd; x += Vector512.Count) + { + Vector512 topSamples = Vector512.LoadUnsafe(ref topBase, (nuint)(x + 1)); + (Vector512 topLow, Vector512 topHigh) = Vector512.Widen(topSamples); + Vector512 lowIndices = indices + Vector512.Create((uint)x); + Vector512 highIndices = lowIndices + Vector512.Create((uint)Vector512.Count); + Vector512 low = CalculatePlanarVector( + topLow, + lowIndices, + leftSample, + topRight, + bottomLeft, + topWeight, + bottomWeight, + (uint)size, + rounding, + shift); + + Vector512 high = CalculatePlanarVector( + topHigh, + highIndices, + leftSample, + topRight, + bottomLeft, + topWeight, + bottomWeight, + (uint)size, + rounding, + shift); + + Vector512.Narrow(low, high).StoreUnsafe(ref Unsafe.Add(ref rowBase, x)); + } + } + + if (Vector256.IsHardwareAccelerated) + { + Vector256 indices = CreateIndicesVector256(); + int oneVectorFromEnd = size - Vector256.Count; + for (; x <= oneVectorFromEnd; x += Vector256.Count) + { + Vector256 topSamples = Vector256.LoadUnsafe(ref topBase, (nuint)(x + 1)); + (Vector256 topLow, Vector256 topHigh) = Vector256.Widen(topSamples); + Vector256 lowIndices = indices + Vector256.Create((uint)x); + Vector256 highIndices = lowIndices + Vector256.Create((uint)Vector256.Count); + Vector256 low = CalculatePlanarVector( + topLow, + lowIndices, + leftSample, + topRight, + bottomLeft, + topWeight, + bottomWeight, + (uint)size, + rounding, + shift); + + Vector256 high = CalculatePlanarVector( + topHigh, + highIndices, + leftSample, + topRight, + bottomLeft, + topWeight, + bottomWeight, + (uint)size, + rounding, + shift); + + Vector256.Narrow(low, high).StoreUnsafe(ref Unsafe.Add(ref rowBase, x)); + } + } + + if (Vector128.IsHardwareAccelerated) + { + Vector128 indices = CreateIndicesVector128(); + int oneVectorFromEnd = size - Vector128.Count; + for (; x <= oneVectorFromEnd; x += Vector128.Count) + { + Vector128 topSamples = Vector128.LoadUnsafe(ref topBase, (nuint)(x + 1)); + (Vector128 topLow, Vector128 topHigh) = Vector128.Widen(topSamples); + Vector128 lowIndices = indices + Vector128.Create((uint)x); + Vector128 highIndices = lowIndices + Vector128.Create((uint)Vector128.Count); + Vector128 low = CalculatePlanarVector( + topLow, + lowIndices, + leftSample, + topRight, + bottomLeft, + topWeight, + bottomWeight, + (uint)size, + rounding, + shift); + + Vector128 high = CalculatePlanarVector( + topHigh, + highIndices, + leftSample, + topRight, + bottomLeft, + topWeight, + bottomWeight, + (uint)size, + rounding, + shift); + + Vector128.Narrow(low, high).StoreUnsafe(ref Unsafe.Add(ref rowBase, x)); + } + } + + for (; x < size; x++) + { + uint horizontal = ((uint)(size - x - 1) * leftSample) + ((uint)(x + 1) * topRight); + uint vertical = ((uint)(size - y - 1) * top[x + 1]) + ((uint)(y + 1) * bottomLeft); + Unsafe.Add(ref rowBase, x) = (ushort)((horizontal + vertical + (uint)size) >> shift); + } + } + } + } + + /// + /// Implements DC prediction and its optional luma boundary filter. + /// + private readonly struct DcPredictionOperator : IHevcIntraPredictionOperator + { + /// + public static void Predict( + ReadOnlySpan top, + ReadOnlySpan left, + Span destination, + int destinationStride, + int size, + int mode, + int bitDepth, + bool filterPredictionEdges, + Span scratch) + { + uint sum = SumSamples(top.Slice(1, size)) + SumSamples(left.Slice(1, size)); + ushort dc = (ushort)((sum + (uint)size) >> (BitOperations.Log2((uint)size) + 1)); + for (int y = 0; y < size; y++) + { + destination.Slice(y * destinationStride, size).Fill(dc); + } + + if (!filterPredictionEdges) + { + return; + } + + destination[0] = (ushort)((top[1] + left[1] + (2 * dc) + 2) >> 2); + for (int x = 1; x < size; x++) + { + destination[x] = (ushort)((top[x + 1] + (3 * dc) + 2) >> 2); + } + + for (int y = 1; y < size; y++) + { + destination[y * destinationStride] = (ushort)((left[y + 1] + (3 * dc) + 2) >> 2); + } + } + } + + /// + /// Implements the thirty-three directional intra-prediction modes. + /// + private readonly struct AngularPredictionOperator : IHevcIntraPredictionOperator + { + /// + public static void Predict( + ReadOnlySpan top, + ReadOnlySpan left, + Span destination, + int destinationStride, + int size, + int mode, + int bitDepth, + bool filterPredictionEdges, + Span scratch) + { + if (mode == VerticalMode) + { + PredictVertical(top, left, destination, destinationStride, size, bitDepth, filterPredictionEdges); + return; + } + + if (mode == HorizontalMode) + { + PredictHorizontal(top, left, destination, destinationStride, size, bitDepth, filterPredictionEdges); + return; + } + + bool vertical = mode >= FirstVerticalMode; + int angleMode = vertical ? mode - VerticalMode : HorizontalMode - mode; + int absoluteAngleMode = Math.Abs(angleMode); + int angle = PredictionAngles[absoluteAngleMode] * Math.Sign(angleMode); + ReadOnlySpan main = vertical ? top : left; + ReadOnlySpan side = vertical ? left : top; + Span temporaryBlock = scratch[..(size * size)]; + Span extendedReference = scratch.Slice(size * size, (4 * size) + 1); + int mainOrigin = 0; + + if (angle < 0) + { + mainOrigin = size * 2; + main[..(size + 1)].CopyTo(extendedReference[mainOrigin..]); + int inverseAngle = InversePredictionAngles[absoluteAngleMode]; + int inverseAngleSum = 128; + int minimumIndex = (size * angle) >> 5; + for (int index = -1; index > minimumIndex; index--) + { + inverseAngleSum += inverseAngle; + extendedReference[mainOrigin + index] = side[inverseAngleSum >> 8]; + } + + main = extendedReference; + } + + Span prediction = vertical ? destination : temporaryBlock; + int predictionStride = vertical ? destinationStride : size; + PredictAngularRows(main, mainOrigin, prediction, predictionStride, size, angle); + if (!vertical) + { + TransposeBlock(temporaryBlock, destination, destinationStride, size); + } + } + } +} diff --git a/src/ImageSharp/Formats/Heif/Hevc/HevcIntraPredictor.cs b/src/ImageSharp/Formats/Heif/Hevc/HevcIntraPredictor.cs new file mode 100644 index 0000000000..b79ee9db15 --- /dev/null +++ b/src/ImageSharp/Formats/Heif/Hevc/HevcIntraPredictor.cs @@ -0,0 +1,412 @@ +// Copyright (c) Six Labors. +// Licensed under the Six Labors Split License. + +using System.Numerics; +using System.Runtime.CompilerServices; +using System.Runtime.InteropServices; +using System.Runtime.Intrinsics; +using SixLabors.ImageSharp.Common.Helpers; + +namespace SixLabors.ImageSharp.Formats.Heif.Hevc; + +/// +/// Reconstructs HEVC intra-prediction blocks from prepared neighboring samples. +/// +internal static partial class HevcIntraPredictor +{ + /// + /// The HEVC planar prediction mode. + /// + private const int PlanarMode = 0; + + /// + /// The HEVC DC prediction mode. + /// + private const int DcMode = 1; + + /// + /// The HEVC horizontal prediction mode. + /// + private const int HorizontalMode = 10; + + /// + /// The first prediction mode whose main reference is the top row. + /// + private const int FirstVerticalMode = 18; + + /// + /// The HEVC vertical prediction mode. + /// + private const int VerticalMode = 26; + + /// + /// The largest transform-block side supported by HEVC intra prediction. + /// + private const int MaximumBlockSize = 32; + + /// + /// Defines one closed intra-prediction operation selected by the decoded mode. + /// + /// The implementing operator type. + private interface IHevcIntraPredictionOperator + where TOperator : struct, IHevcIntraPredictionOperator + { + /// + /// Reconstructs one square prediction block. + /// + /// The top-left, top, and top-right reference samples. + /// The top-left, left, and below-left reference samples. + /// The destination buffer beginning at the block origin. + /// The destination row stride in samples. + /// The square block side in samples. + /// The decoded prediction mode. + /// The reconstructed component precision. + /// Whether the luma edge filter applies to the selected block. + /// The caller-owned block and extended-reference scratch space. + public static abstract void Predict( + ReadOnlySpan top, + ReadOnlySpan left, + Span destination, + int destinationStride, + int size, + int mode, + int bitDepth, + bool filterPredictionEdges, + Span scratch); + } + + /// + /// Gets the angle selected by each absolute angular-mode displacement. + /// + private static ReadOnlySpan PredictionAngles => [0, 2, 5, 9, 13, 17, 21, 26, 32]; + + /// + /// Gets the reciprocal angle used to extend the main reference for negative directions. + /// + private static ReadOnlySpan InversePredictionAngles => [0, 4096, 1638, 910, 630, 482, 390, 315, 256]; + + /// + /// Gets the scratch length required to predict a block of the specified size. + /// + /// The base-two logarithm of the square block side. + /// The required number of elements. + public static int GetScratchLength(int log2Size) + { + DebugGuard.MustBeBetweenOrEqualTo(log2Size, 2, 5, nameof(log2Size)); + int size = 1 << log2Size; + return (size * size) + (4 * size) + 1; + } + + /// + /// Reconstructs one square intra-prediction block using a closed operator selected by the decoded mode. + /// + /// The top-left, top, and top-right reference samples. + /// The top-left, left, and below-left reference samples. + /// The destination buffer beginning at the block origin. + /// The destination row stride in samples. + /// The base-two logarithm of the square block side. + /// The decoded prediction mode in the inclusive range zero through thirty-four. + /// The reconstructed component precision. + /// Whether the luma edge filter applies to the selected block. + /// The caller-owned scratch returned by . + public static void Predict( + ReadOnlySpan top, + ReadOnlySpan left, + Span destination, + int destinationStride, + int log2Size, + int mode, + int bitDepth, + bool filterPredictionEdges, + Span scratch) + { + DebugGuard.MustBeBetweenOrEqualTo(log2Size, 2, 5, nameof(log2Size)); + DebugGuard.MustBeBetweenOrEqualTo(mode, PlanarMode, 34, nameof(mode)); + int size = 1 << log2Size; + + switch (mode) + { + case PlanarMode: + Predict( + top, + left, + destination, + destinationStride, + size, + mode, + bitDepth, + filterPredictionEdges, + scratch); + break; + case DcMode: + Predict( + top, + left, + destination, + destinationStride, + size, + mode, + bitDepth, + filterPredictionEdges, + scratch); + break; + default: + Predict( + top, + left, + destination, + destinationStride, + size, + mode, + bitDepth, + filterPredictionEdges, + scratch); + break; + } + } + + /// + /// Filters prepared reference samples using the normative three-tap or strong bilinear filter. + /// + /// The unfiltered top-left, top, and top-right samples. + /// The unfiltered top-left, left, and below-left samples. + /// The destination top reference. + /// The destination left reference. + /// The base-two logarithm of the square prediction-block side. + /// The reconstructed luma precision. + /// Whether the sequence permits strong intra smoothing. + public static void FilterReferenceSamples( + ReadOnlySpan top, + ReadOnlySpan left, + Span filteredTop, + Span filteredLeft, + int log2Size, + int bitDepth, + bool strongIntraSmoothingEnabled) + { + DebugGuard.MustBeBetweenOrEqualTo(log2Size, 2, 5, nameof(log2Size)); + int size = 1 << log2Size; + int referenceLength = (size * 2) + 1; + bool useStrongSmoothing = strongIntraSmoothingEnabled && size == MaximumBlockSize; + if (useStrongSmoothing) + { + int threshold = 1 << (bitDepth - 5); + int last = referenceLength - 1; + bool leftIsBilinear = Math.Abs((left[last] + left[0]) - (2 * left[size])) < threshold; + bool topIsBilinear = Math.Abs((top[0] + top[last]) - (2 * top[size])) < threshold; + useStrongSmoothing = leftIsBilinear && topIsBilinear; + } + + if (useStrongSmoothing) + { + FilterReferenceBilinear(top[..referenceLength], filteredTop, size); + FilterReferenceBilinear(left[..referenceLength], filteredLeft, size); + return; + } + + // The corner belongs to both references. Filtering it once from the first samples on both sides keeps the + // two logical arrays identical at index zero before their independent one-dimensional filters continue. + ushort filteredCorner = (ushort)((left[1] + (2 * top[0]) + top[1] + 2) >> 2); + filteredTop[0] = filteredCorner; + filteredLeft[0] = filteredCorner; + FilterReferenceThreeTap(top[..referenceLength], filteredTop); + FilterReferenceThreeTap(left[..referenceLength], filteredLeft); + } + + /// + /// Invokes one statically selected prediction operator without interface dispatch in the block loop. + /// + /// The selected prediction operator. + /// The prepared top reference. + /// The prepared left reference. + /// The destination block origin. + /// The destination row stride in samples. + /// The square block side in samples. + /// The decoded prediction mode. + /// The reconstructed component precision. + /// Whether the luma edge filter applies. + /// The caller-owned prediction scratch. + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private static void Predict( + ReadOnlySpan top, + ReadOnlySpan left, + Span destination, + int destinationStride, + int size, + int mode, + int bitDepth, + bool filterPredictionEdges, + Span scratch) + where TOperator : struct, IHevcIntraPredictionOperator + => TOperator.Predict( + top, + left, + destination, + destinationStride, + size, + mode, + bitDepth, + filterPredictionEdges, + scratch); + + /// + /// Applies the strong bilinear filter between the reference endpoints. + /// + /// The complete unfiltered reference. + /// The complete filtered reference. + /// The prediction-block side in samples. + private static void FilterReferenceBilinear(ReadOnlySpan source, Span destination, int size) + { + int last = source.Length - 1; + destination[0] = source[0]; + destination[last] = source[last]; + ref ushort sourceBase = ref MemoryMarshal.GetReference(source); + ref ushort destinationBase = ref MemoryMarshal.GetReference(destination); + uint first = source[0]; + uint final = source[last]; + int shift = BitOperations.Log2((uint)(size * 2)); + uint rounding = (uint)size; + int i = 1; + + // Each widened lane represents one reference coordinate. The weights sum to 2N, so narrowing is exact + // after the rounded shift for every supported 8, 10, and 12-bit sample. + if (Vector512.IsHardwareAccelerated) + { + Vector512 indices = CreateIndicesVector512(); + int oneVectorFromEnd = last - Vector512.Count; + for (; i <= oneVectorFromEnd; i += Vector512.Count) + { + Vector512 lowerIndices = indices + Vector512.Create((uint)i); + Vector512 upperIndices = lowerIndices + Vector512.Create((uint)Vector512.Count); + Vector512 lower = (((Vector512.Create((uint)last) - lowerIndices) * first) + (lowerIndices * final) + Vector512.Create(rounding)) >> shift; + Vector512 upper = (((Vector512.Create((uint)last) - upperIndices) * first) + (upperIndices * final) + Vector512.Create(rounding)) >> shift; + Vector512.Narrow(lower, upper).StoreUnsafe(ref Unsafe.Add(ref destinationBase, i)); + } + } + + if (Vector256.IsHardwareAccelerated) + { + Vector256 indices = CreateIndicesVector256(); + int oneVectorFromEnd = last - Vector256.Count; + for (; i <= oneVectorFromEnd; i += Vector256.Count) + { + Vector256 lowerIndices = indices + Vector256.Create((uint)i); + Vector256 upperIndices = lowerIndices + Vector256.Create((uint)Vector256.Count); + Vector256 lower = (((Vector256.Create((uint)last) - lowerIndices) * first) + (lowerIndices * final) + Vector256.Create(rounding)) >> shift; + Vector256 upper = (((Vector256.Create((uint)last) - upperIndices) * first) + (upperIndices * final) + Vector256.Create(rounding)) >> shift; + Vector256.Narrow(lower, upper).StoreUnsafe(ref Unsafe.Add(ref destinationBase, i)); + } + } + + if (Vector128.IsHardwareAccelerated) + { + Vector128 indices = CreateIndicesVector128(); + int oneVectorFromEnd = last - Vector128.Count; + for (; i <= oneVectorFromEnd; i += Vector128.Count) + { + Vector128 lowerIndices = indices + Vector128.Create((uint)i); + Vector128 upperIndices = lowerIndices + Vector128.Create((uint)Vector128.Count); + Vector128 lower = (((Vector128.Create((uint)last) - lowerIndices) * first) + (lowerIndices * final) + Vector128.Create(rounding)) >> shift; + Vector128 upper = (((Vector128.Create((uint)last) - upperIndices) * first) + (upperIndices * final) + Vector128.Create(rounding)) >> shift; + Vector128.Narrow(lower, upper).StoreUnsafe(ref Unsafe.Add(ref destinationBase, i)); + } + } + + for (; i < last; i++) + { + Unsafe.Add(ref destinationBase, i) = (ushort)((((last - i) * first) + (i * final) + rounding) >> shift); + } + } + + /// + /// Applies the normal three-tap reference filter to every non-endpoint sample. + /// + /// The complete unfiltered reference. + /// The complete filtered reference with its corner already initialized. + private static void FilterReferenceThreeTap(ReadOnlySpan source, Span destination) + { + int last = source.Length - 1; + destination[last] = source[last]; + ref ushort sourceBase = ref MemoryMarshal.GetReference(source); + ref ushort destinationBase = ref MemoryMarshal.GetReference(destination); + int i = 1; + + if (Vector512.IsHardwareAccelerated) + { + int oneVectorFromEnd = last - Vector512.Count; + for (; i <= oneVectorFromEnd; i += Vector512.Count) + { + Vector512 previous = Vector512.LoadUnsafe(ref sourceBase, (nuint)(i - 1)); + Vector512 current = Vector512.LoadUnsafe(ref sourceBase, (nuint)i); + Vector512 next = Vector512.LoadUnsafe(ref sourceBase, (nuint)(i + 1)); + (Vector512 previousLow, Vector512 previousHigh) = Vector512.Widen(previous); + (Vector512 currentLow, Vector512 currentHigh) = Vector512.Widen(current); + (Vector512 nextLow, Vector512 nextHigh) = Vector512.Widen(next); + Vector512 low = (previousLow + (currentLow << 1) + nextLow + Vector512.Create(2U)) >> 2; + Vector512 high = (previousHigh + (currentHigh << 1) + nextHigh + Vector512.Create(2U)) >> 2; + Vector512.Narrow(low, high).StoreUnsafe(ref Unsafe.Add(ref destinationBase, i)); + } + } + + if (Vector256.IsHardwareAccelerated) + { + int oneVectorFromEnd = last - Vector256.Count; + for (; i <= oneVectorFromEnd; i += Vector256.Count) + { + Vector256 previous = Vector256.LoadUnsafe(ref sourceBase, (nuint)(i - 1)); + Vector256 current = Vector256.LoadUnsafe(ref sourceBase, (nuint)i); + Vector256 next = Vector256.LoadUnsafe(ref sourceBase, (nuint)(i + 1)); + (Vector256 previousLow, Vector256 previousHigh) = Vector256.Widen(previous); + (Vector256 currentLow, Vector256 currentHigh) = Vector256.Widen(current); + (Vector256 nextLow, Vector256 nextHigh) = Vector256.Widen(next); + Vector256 low = (previousLow + (currentLow << 1) + nextLow + Vector256.Create(2U)) >> 2; + Vector256 high = (previousHigh + (currentHigh << 1) + nextHigh + Vector256.Create(2U)) >> 2; + Vector256.Narrow(low, high).StoreUnsafe(ref Unsafe.Add(ref destinationBase, i)); + } + } + + if (Vector128.IsHardwareAccelerated) + { + int oneVectorFromEnd = last - Vector128.Count; + for (; i <= oneVectorFromEnd; i += Vector128.Count) + { + Vector128 previous = Vector128.LoadUnsafe(ref sourceBase, (nuint)(i - 1)); + Vector128 current = Vector128.LoadUnsafe(ref sourceBase, (nuint)i); + Vector128 next = Vector128.LoadUnsafe(ref sourceBase, (nuint)(i + 1)); + (Vector128 previousLow, Vector128 previousHigh) = Vector128.Widen(previous); + (Vector128 currentLow, Vector128 currentHigh) = Vector128.Widen(current); + (Vector128 nextLow, Vector128 nextHigh) = Vector128.Widen(next); + Vector128 low = (previousLow + (currentLow << 1) + nextLow + Vector128.Create(2U)) >> 2; + Vector128 high = (previousHigh + (currentHigh << 1) + nextHigh + Vector128.Create(2U)) >> 2; + Vector128.Narrow(low, high).StoreUnsafe(ref Unsafe.Add(ref destinationBase, i)); + } + } + + for (; i < last; i++) + { + Unsafe.Add(ref destinationBase, i) = (ushort)((source[i - 1] + (2 * source[i]) + source[i + 1] + 2) >> 2); + } + } + + /// + /// Creates the zero-through-fifteen lane indices used by 512-bit weighted interpolation. + /// + /// The ordered lane indices. + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private static Vector512 CreateIndicesVector512() + => Vector512.Create(0U, 1U, 2U, 3U, 4U, 5U, 6U, 7U, 8U, 9U, 10U, 11U, 12U, 13U, 14U, 15U); + + /// + /// Creates the zero-through-seven lane indices used by 256-bit weighted interpolation. + /// + /// The ordered lane indices. + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private static Vector256 CreateIndicesVector256() => Vector256.Create(0U, 1U, 2U, 3U, 4U, 5U, 6U, 7U); + + /// + /// Creates the zero-through-three lane indices used by 128-bit weighted interpolation. + /// + /// The ordered lane indices. + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private static Vector128 CreateIndicesVector128() => Vector128.Create(0U, 1U, 2U, 3U); +} diff --git a/tests/ImageSharp.Benchmarks/Codecs/Heif/HevcIntraPredictionBenchmarks.cs b/tests/ImageSharp.Benchmarks/Codecs/Heif/HevcIntraPredictionBenchmarks.cs new file mode 100644 index 0000000000..94c885478b --- /dev/null +++ b/tests/ImageSharp.Benchmarks/Codecs/Heif/HevcIntraPredictionBenchmarks.cs @@ -0,0 +1,116 @@ +// Copyright (c) Six Labors. +// Licensed under the Six Labors Split License. + +using BenchmarkDotNet.Attributes; +using SixLabors.ImageSharp.Formats.Heif.Hevc; + +namespace SixLabors.ImageSharp.Benchmarks.Codecs.Heif; + +/// +/// Measures complete coded-frame traversal through representative HEVC intra-prediction modes. +/// +[MemoryDiagnoser(displayGenColumns: false)] +public class HevcIntraPredictionBenchmarks +{ + /// + /// The coded frame width, which is an exact multiple of the maximum transform-block side. + /// + private const int Width = 1920; + + /// + /// The coded frame height including the final padded coding-tree row for a 1080-line presentation. + /// + private const int Height = 1088; + + /// + /// The base-two logarithm of the benchmark prediction-block side. + /// + private const int BlockLog2 = 5; + + /// + /// The prediction-block side in samples. + /// + private const int BlockSize = 1 << BlockLog2; + + /// + /// The prepared top reference shared by deterministic benchmark blocks. + /// + private readonly ushort[] top = new ushort[(BlockSize * 2) + 1]; + + /// + /// The prepared left reference shared by deterministic benchmark blocks. + /// + private readonly ushort[] left = new ushort[(BlockSize * 2) + 1]; + + /// + /// The frame-wide reconstructed prediction samples. + /// + private readonly ushort[] destination = new ushort[Width * Height]; + + /// + /// The maximum-block scratch reused throughout each coded frame. + /// + private readonly ushort[] scratch = new ushort[HevcIntraPredictor.GetScratchLength(BlockLog2)]; + + /// + /// Populates deterministic twelve-bit reference samples outside the measured frame traversal. + /// + [GlobalSetup] + public void Setup() + { + this.top[0] = this.left[0] = 1365; + for (int i = 1; i < this.top.Length; i++) + { + this.top[i] = (ushort)((1365 + (37 * i)) & 4095); + this.left[i] = (ushort)((1365 + (53 * i)) & 4095); + } + } + + /// + /// Measures frame-wide planar prediction with 512-, 256-, and 128-bit row dispatch where available. + /// + /// The final reconstructed sample, keeping the frame output observable. + [Benchmark(Baseline = true)] + public ushort PredictPlanarFrame() => this.PredictFrame(0); + + /// + /// Measures frame-wide fractional vertical prediction using contiguous SIMD interpolation. + /// + /// The final reconstructed sample, keeping the frame output observable. + [Benchmark] + public ushort PredictVerticalAngularFrame() => this.PredictFrame(30); + + /// + /// Measures frame-wide horizontal prediction including the SIMD block transposition stage. + /// + /// The final reconstructed sample, keeping the frame output observable. + [Benchmark] + public ushort PredictHorizontalAngularFrame() => this.PredictFrame(2); + + /// + /// Reconstructs every maximum-size prediction block in the coded benchmark frame. + /// + /// The HEVC intra-prediction mode. + /// The final reconstructed sample. + private ushort PredictFrame(int mode) + { + for (int y = 0; y < Height; y += BlockSize) + { + for (int x = 0; x < Width; x += BlockSize) + { + HevcIntraPredictor.Predict( + this.top, + this.left, + this.destination.AsSpan((y * Width) + x), + Width, + BlockLog2, + mode, + 12, + true, + this.scratch); + } + } + + return this.destination[^1]; + } +} diff --git a/tests/ImageSharp.Tests/Formats/Heif/Hevc/HevcIntraPredictorTests.cs b/tests/ImageSharp.Tests/Formats/Heif/Hevc/HevcIntraPredictorTests.cs new file mode 100644 index 0000000000..8378e77f0f --- /dev/null +++ b/tests/ImageSharp.Tests/Formats/Heif/Hevc/HevcIntraPredictorTests.cs @@ -0,0 +1,396 @@ +// Copyright (c) Six Labors. +// Licensed under the Six Labors Split License. + +using System.Numerics; +using SixLabors.ImageSharp.Formats.Heif.Hevc; + +namespace SixLabors.ImageSharp.Tests.Formats.Heif.Hevc; + +/// +/// Verifies HEVC planar, DC, angular, reference-filter, and SIMD prediction behavior. +/// +[Trait("Format", "Heic")] +public class HevcIntraPredictorTests +{ + /// + /// Verifies fixed four-by-four prediction results derived from the HEVC intra-prediction equations. + /// + [Fact] + public void PredictsKnownFourByFourBlocks() + { + ushort[] top = [64, 80, 96, 112, 128, 144, 160, 176, 192]; + ushort[] left = [64, 70, 76, 82, 88, 94, 100, 106, 112]; + int[] modes = [0, 1, 2, 9, 18, 30, 34]; + ushort[][] expected = + [ + [83, 97, 110, 123, 87, 97, 108, 118, 90, 98, 105, 113, 93, 98, 103, 108], + [84, 93, 97, 101, 88, 92, 92, 92, 90, 92, 92, 92, 91, 92, 92, 92], + [76, 82, 88, 94, 82, 88, 94, 100, 88, 94, 100, 106, 94, 100, 106, 112], + [70, 71, 71, 72, 76, 77, 77, 78, 82, 83, 83, 84, 88, 89, 89, 90], + [64, 80, 96, 112, 70, 64, 80, 96, 76, 70, 64, 80, 82, 76, 70, 64], + [87, 103, 119, 135, 93, 109, 125, 141, 100, 116, 132, 148, 106, 122, 138, 154], + [96, 112, 128, 144, 112, 128, 144, 160, 128, 144, 160, 176, 144, 160, 176, 192] + ]; + + const int size = 4; + const int stride = 7; + ushort[] destination = new ushort[stride * size]; + ushort[] scratch = new ushort[HevcIntraPredictor.GetScratchLength(2)]; + for (int caseIndex = 0; caseIndex < modes.Length; caseIndex++) + { + destination.AsSpan().Fill(ushort.MaxValue); + int mode = modes[caseIndex]; + HevcIntraPredictor.Predict(top, left, destination, stride, 2, mode, 8, mode == 1, scratch); + + for (int y = 0; y < size; y++) + { + ReadOnlySpan expectedRow = expected[caseIndex].AsSpan(y * size, size); + ReadOnlySpan actualRow = destination.AsSpan(y * stride, size); + Assert.True(expectedRow.SequenceEqual(actualRow), $"Mode {mode}, row {y} did not match the fixed HEVC result."); + } + } + } + + /// + /// Verifies the optional luma boundary filter for the pure horizontal and vertical modes. + /// + [Fact] + public void FiltersPureDirectionPredictionEdges() + { + ushort[] top = [64, 80, 96, 112, 128, 144, 160, 176, 192]; + ushort[] left = [64, 70, 76, 82, 88, 94, 100, 106, 112]; + ushort[] scratch = new ushort[HevcIntraPredictor.GetScratchLength(2)]; + ushort[] horizontal = new ushort[16]; + ushort[] vertical = new ushort[16]; + + HevcIntraPredictor.Predict(top, left, horizontal, 4, 2, 10, 8, true, scratch); + HevcIntraPredictor.Predict(top, left, vertical, 4, 2, 26, 8, true, scratch); + + ushort[] expectedHorizontal = [78, 86, 94, 102, 76, 76, 76, 76, 82, 82, 82, 82, 88, 88, 88, 88]; + ushort[] expectedVertical = [83, 96, 112, 128, 86, 96, 112, 128, 89, 96, 112, 128, 92, 96, 112, 128]; + Assert.True(expectedHorizontal.AsSpan().SequenceEqual(horizontal)); + Assert.True(expectedVertical.AsSpan().SequenceEqual(vertical)); + } + + /// + /// Verifies exact three-tap filtering, including the shared top-left sample. + /// + [Fact] + public void FiltersReferenceSamplesWithThreeTapKernel() + { + ushort[] top = [64, 80, 96, 112, 128, 144, 160, 176, 192]; + ushort[] left = [64, 70, 76, 82, 88, 94, 100, 106, 112]; + ushort[] filteredTop = new ushort[top.Length]; + ushort[] filteredLeft = new ushort[left.Length]; + + HevcIntraPredictor.FilterReferenceSamples(top, left, filteredTop, filteredLeft, 2, 8, true); + + ushort[] expectedTop = [70, 80, 96, 112, 128, 144, 160, 176, 192]; + ushort[] expectedLeft = [70, 70, 76, 82, 88, 94, 100, 106, 112]; + Assert.True(expectedTop.AsSpan().SequenceEqual(filteredTop)); + Assert.True(expectedLeft.AsSpan().SequenceEqual(filteredLeft)); + } + + /// + /// Verifies that eligible thirty-two-sample references use strong bilinear smoothing rather than local three-tap filtering. + /// + [Fact] + public void StrongSmoothingReplacesEligibleNonlinearReferences() + { + const int size = 32; + ushort[] top = new ushort[(size * 2) + 1]; + ushort[] left = new ushort[top.Length]; + ushort[] filteredTop = new ushort[top.Length]; + ushort[] filteredLeft = new ushort[left.Length]; + for (int i = 0; i < top.Length; i++) + { + top[i] = (ushort)(100 + i + (i % 3)); + left[i] = (ushort)(100 + (2 * i) + (i % 5)); + } + + // Strong smoothing is selected from the endpoint/midpoint test, so keep those six values exactly bilinear + // while the remaining samples deliberately differ from the expected straight lines. + top[0] = left[0] = 100; + top[size] = 132; + top[size * 2] = 164; + left[size] = 164; + left[size * 2] = 228; + + HevcIntraPredictor.FilterReferenceSamples(top, left, filteredTop, filteredLeft, 5, 10, true); + + for (int i = 0; i < top.Length; i++) + { + Assert.Equal((ushort)(100 + i), filteredTop[i]); + Assert.Equal((ushort)(100 + (2 * i)), filteredLeft[i]); + } + } + + /// + /// Compares every prediction mode and block width with a specification-shaped scalar oracle. + /// + /// The base-two logarithm of the tested block side. + /// The reconstructed component precision. + [Theory] + [InlineData(2, 8)] + [InlineData(3, 10)] + [InlineData(4, 12)] + [InlineData(5, 12)] + public void EveryModeMatchesScalarOracle(int log2Size, int bitDepth) + { + int size = 1 << log2Size; + int maximum = (1 << bitDepth) - 1; + int referenceLength = (size * 2) + 1; + ushort[] top = new ushort[referenceLength]; + ushort[] left = new ushort[referenceLength]; + top[0] = left[0] = (ushort)(maximum / 3); + for (int i = 1; i < referenceLength; i++) + { + top[i] = (ushort)((top[0] + (37 * i) + (3 * size)) & maximum); + left[i] = (ushort)((left[0] + (53 * i) + (5 * size)) & maximum); + } + + int stride = size + 3; + ushort[] expected = new ushort[stride * size]; + ushort[] actual = new ushort[stride * size]; + ushort[] scratch = new ushort[HevcIntraPredictor.GetScratchLength(log2Size)]; + for (int mode = 0; mode <= 34; mode++) + { + expected.AsSpan().Clear(); + actual.AsSpan().Clear(); + PredictScalar(top, left, expected, stride, size, mode, bitDepth, true); + HevcIntraPredictor.Predict(top, left, actual, stride, log2Size, mode, bitDepth, true, scratch); + Assert.True(expected.AsSpan().SequenceEqual(actual), $"Mode {mode}, size {size}, and bit depth {bitDepth} did not match the scalar oracle."); + } + } + + /// + /// Reconstructs one block directly from the HEVC planar, DC, and angular prediction equations. + /// + /// The top reference samples. + /// The left reference samples. + /// The destination block origin. + /// The destination row stride. + /// The square block side. + /// The prediction mode. + /// The reconstructed component precision. + /// Whether the luma edge filter applies. + private static void PredictScalar( + ReadOnlySpan top, + ReadOnlySpan left, + Span destination, + int destinationStride, + int size, + int mode, + int bitDepth, + bool filterPredictionEdges) + { + if (mode == 0) + { + int shift = BitOperations.Log2((uint)size) + 1; + for (int y = 0; y < size; y++) + { + for (int x = 0; x < size; x++) + { + int horizontal = ((size - x - 1) * left[y + 1]) + ((x + 1) * top[size]); + int vertical = ((size - y - 1) * top[x + 1]) + ((y + 1) * left[size]); + destination[(y * destinationStride) + x] = (ushort)((horizontal + vertical + size) >> shift); + } + } + + return; + } + + if (mode == 1) + { + PredictDcScalar(top, left, destination, destinationStride, size, filterPredictionEdges); + return; + } + + if (mode == 10) + { + PredictHorizontalScalar(top, left, destination, destinationStride, size, bitDepth, filterPredictionEdges); + return; + } + + if (mode == 26) + { + PredictVerticalScalar(top, left, destination, destinationStride, size, bitDepth, filterPredictionEdges); + return; + } + + PredictAngularScalar(top, left, destination, destinationStride, size, mode); + } + + /// + /// Reconstructs a scalar DC block and its optional boundary filter. + /// + /// The top reference samples. + /// The left reference samples. + /// The destination block origin. + /// The destination row stride. + /// The square block side. + /// Whether the luma edge filter applies. + private static void PredictDcScalar( + ReadOnlySpan top, + ReadOnlySpan left, + Span destination, + int destinationStride, + int size, + bool filterPredictionEdges) + { + int sum = 0; + for (int i = 1; i <= size; i++) + { + sum += top[i] + left[i]; + } + + ushort dc = (ushort)((sum + size) >> (BitOperations.Log2((uint)size) + 1)); + for (int y = 0; y < size; y++) + { + destination.Slice(y * destinationStride, size).Fill(dc); + } + + if (!filterPredictionEdges) + { + return; + } + + destination[0] = (ushort)((top[1] + left[1] + (2 * dc) + 2) >> 2); + for (int i = 1; i < size; i++) + { + destination[i] = (ushort)((top[i + 1] + (3 * dc) + 2) >> 2); + destination[i * destinationStride] = (ushort)((left[i + 1] + (3 * dc) + 2) >> 2); + } + } + + /// + /// Reconstructs scalar horizontal prediction and its optional boundary filter. + /// + /// The top reference samples. + /// The left reference samples. + /// The destination block origin. + /// The destination row stride. + /// The square block side. + /// The reconstructed component precision. + /// Whether the luma edge filter applies. + private static void PredictHorizontalScalar( + ReadOnlySpan top, + ReadOnlySpan left, + Span destination, + int destinationStride, + int size, + int bitDepth, + bool filterPredictionEdges) + { + for (int y = 0; y < size; y++) + { + destination.Slice(y * destinationStride, size).Fill(left[y + 1]); + } + + if (filterPredictionEdges) + { + int maximum = (1 << bitDepth) - 1; + for (int x = 0; x < size; x++) + { + destination[x] = (ushort)Math.Clamp(destination[x] + ((top[x + 1] - top[0]) >> 1), 0, maximum); + } + } + } + + /// + /// Reconstructs scalar vertical prediction and its optional boundary filter. + /// + /// The top reference samples. + /// The left reference samples. + /// The destination block origin. + /// The destination row stride. + /// The square block side. + /// The reconstructed component precision. + /// Whether the luma edge filter applies. + private static void PredictVerticalScalar( + ReadOnlySpan top, + ReadOnlySpan left, + Span destination, + int destinationStride, + int size, + int bitDepth, + bool filterPredictionEdges) + { + int maximum = (1 << bitDepth) - 1; + for (int y = 0; y < size; y++) + { + top.Slice(1, size).CopyTo(destination[(y * destinationStride)..]); + if (filterPredictionEdges) + { + int offset = y * destinationStride; + destination[offset] = (ushort)Math.Clamp(destination[offset] + ((left[y + 1] - left[0]) >> 1), 0, maximum); + } + } + } + + /// + /// Reconstructs a scalar angular block, including negative-reference extension and horizontal transposition. + /// + /// The top reference samples. + /// The left reference samples. + /// The destination block origin. + /// The destination row stride. + /// The square block side. + /// The angular prediction mode. + private static void PredictAngularScalar( + ReadOnlySpan top, + ReadOnlySpan left, + Span destination, + int destinationStride, + int size, + int mode) + { + ReadOnlySpan angles = [0, 2, 5, 9, 13, 17, 21, 26, 32]; + ReadOnlySpan inverseAngles = [0, 4096, 1638, 910, 630, 482, 390, 315, 256]; + bool vertical = mode >= 18; + int angleMode = vertical ? mode - 26 : 10 - mode; + int absoluteAngleMode = Math.Abs(angleMode); + int angle = angles[absoluteAngleMode] * Math.Sign(angleMode); + ReadOnlySpan main = vertical ? top : left; + ReadOnlySpan side = vertical ? left : top; + int mainOrigin = size * 2; + int[] extendedMain = new int[(4 * size) + 1]; + for (int i = 0; i < main.Length; i++) + { + extendedMain[mainOrigin + i] = main[i]; + } + + if (angle < 0) + { + int inverseAngleSum = 128; + for (int index = -1; index > ((size * angle) >> 5); index--) + { + inverseAngleSum += inverseAngles[absoluteAngleMode]; + extendedMain[mainOrigin + index] = side[inverseAngleSum >> 8]; + } + } + + ushort[] temporary = new ushort[size * size]; + for (int y = 0, deltaPosition = angle; y < size; y++, deltaPosition += angle) + { + int deltaInteger = deltaPosition >> 5; + int deltaFraction = deltaPosition & 31; + for (int x = 0; x < size; x++) + { + int index = mainOrigin + x + deltaInteger + 1; + temporary[(y * size) + x] = deltaFraction == 0 + ? (ushort)extendedMain[index] + : (ushort)(((extendedMain[index] * (32 - deltaFraction)) + (extendedMain[index + 1] * deltaFraction) + 16) >> 5); + } + } + + for (int y = 0; y < size; y++) + { + for (int x = 0; x < size; x++) + { + int sourceIndex = vertical ? (y * size) + x : (x * size) + y; + destination[(y * destinationStride) + x] = temporary[sourceIndex]; + } + } + } +}