Add CMS* prototype, XYB decoding, TOC, convolution and coefficient ordering

* CMS - Color Management System

src/ImageSharp/Formats/Jxl/Cms/JxlCieXyPrimaries.cs

@@ -0,0 +1,27 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+using SixLabors.ImageSharp.ColorProfiles;
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Cms;
+
+/// <summary>
+/// RGB primaries for CIEXY
+/// </summary>
+internal struct JxlCieXyPrimaries
+{
+    /// <summary>
+    /// Gets or sets the R component
+    /// </summary>
+    public CieXyChromaticityCoordinates R { get; set; }
+
+    /// <summary>
+    /// Gets or sets the G component
+    /// </summary>
+    public CieXyChromaticityCoordinates G { get; set; }
+
+    /// <summary>
+    /// Gets or sets the B component
+    /// </summary>
+    public CieXyChromaticityCoordinates B { get; set; }
+}

src/ImageSharp/Formats/Jxl/Cms/JxlColorEncoding.cs

@@ -0,0 +1,77 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Cms;
+
+internal sealed class JxlColorEncoding
+{
+    public JxlWhitePoint WhitePoint { get; set; } = JxlWhitePoint.D65;
+
+    public JxlPrimaries Primaries { get; set; } = JxlPrimaries.SRgb;
+
+    public JxlRenderingIntent RenderingIntent { get; set; } = JxlRenderingIntent.Relative;
+
+    public bool HaveFields { get; set; } = true;
+
+    public JxlIccBytes? Icc { get; set; }
+
+    public JxlColorSpace ColorSpace { get; set; } = JxlColorSpace.Rgb;
+
+    public bool Cmyk { get; set; }
+
+    public JxlCustomTransferFunction TransferFunction { get; set; }
+
+    public JxlCustomXy White { get; set; }
+
+    public JxlCustomXy Red { get; set; }
+
+    public JxlCustomXy Green { get; set; }
+
+    public JxlCustomXy Blue { get; set; }
+
+    public bool HasPrimaries => this.ColorSpace is not (JxlColorSpace.Gray or JxlColorSpace.Xyb);
+
+    public int Channels => (this.ColorSpace == JxlColorSpace.Gray) ? 1 : 3;
+
+    public bool TryGetPrimaries(out JxlCieXyPrimaries xy)
+    {
+        xy = default;
+
+        if (!this.HasPrimaries || !this.HasPrimaries)
+        {
+            return false;
+        }
+
+        switch (this.Primaries)
+        {
+            case JxlPrimaries.Custom:
+                xy.R = this.Red.GetValue();
+                xy.G = this.Green.GetValue();
+                xy.B = this.Blue.GetValue();
+                break;
+
+            case JxlPrimaries.SRgb:
+                xy.R = new(0.639998686f, 0.330010138f);
+                xy.G = new(0.300003784f, 0.600003357f);
+                xy.B = new(0.150002046f, 0.059997204f);
+                break;
+
+            case JxlPrimaries.Bt2020:
+                xy.R = new(0.708f, 0.292f);
+                xy.G = new(0.170f, 0.797f);
+                xy.B = new(0.131f, 0.046f);
+                break;
+
+            case JxlPrimaries.P3:
+                xy.R = new(0.680f, 0.320f);
+                xy.G = new(0.265f, 0.690f);
+                xy.B = new(0.150f, 0.060f);
+                break;
+
+            default:
+                throw new InvalidOperationException("Invalid primaries: " + this.Primaries);
+        }
+
+        return true;
+    }
+}

src/ImageSharp/Formats/Jxl/Cms/JxlColorSpace.cs

@@ -0,0 +1,31 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Cms;
+
+/// <summary>
+/// Supported, JPEG XL-specific color space types.
+/// </summary>
+internal enum JxlColorSpace : byte
+{
+    /// <summary>
+    /// Trichromatic color data. This also includes CMYK if Black
+    /// ExtraChannelInfo is present.
+    /// </summary>
+    Rgb,
+
+    /// <summary>
+    /// Single-channel data.
+    /// </summary>
+    Gray,
+
+    /// <summary>
+    /// Like Rgb but fixed values for primaries.
+    /// </summary>
+    Xyb,
+
+    /// <summary>
+    /// Unknown color space
+    /// </summary>
+    Unknown
+}

src/ImageSharp/Formats/Jxl/Cms/JxlCustomTransferFunction.cs

@@ -0,0 +1,109 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Cms;
+
+internal struct JxlCustomTransferFunction
+{
+    private const uint MaxGamma = 8192;
+    private const uint GammaMultiplier = 10000000;
+
+    public JxlCustomTransferFunction()
+    {
+    }
+
+    public bool HaveGamma { get; set; }
+
+    public uint Gamma { get; set; }
+
+    public JxlTransferFunction TransferFunction { get; set; } = JxlTransferFunction.SRgb;
+
+    public readonly bool IsUnknown => !this.HaveGamma && this.TransferFunction == JxlTransferFunction.Unknown;
+
+    public readonly bool IsSrgb => !this.HaveGamma && this.TransferFunction == JxlTransferFunction.SRgb;
+
+    public readonly bool IsLinear => !this.HaveGamma && this.TransferFunction == JxlTransferFunction.Linear;
+
+    public readonly bool IsPq => !this.HaveGamma && this.TransferFunction == JxlTransferFunction.Pq;
+
+    public readonly bool IsHlg => !this.HaveGamma && this.TransferFunction == JxlTransferFunction.Hlg;
+
+    public readonly bool Is709 => !this.HaveGamma && this.TransferFunction == JxlTransferFunction.Bt709;
+
+    public readonly bool IsDci => !this.HaveGamma && this.TransferFunction == JxlTransferFunction.Dci;
+
+    public readonly JxlTransferFunction GetTransferFunction()
+    {
+        if (this.HaveGamma)
+        {
+            return JxlTransferFunction.Unknown;
+        }
+
+        return this.TransferFunction;
+    }
+
+    public void SetTransferFunction(JxlTransferFunction tf)
+    {
+        this.HaveGamma = false;
+        this.TransferFunction = tf;
+    }
+
+    public readonly float GetGamma()
+    {
+        if (!this.HaveGamma)
+        {
+            return 0.0f;
+        }
+
+        return this.Gamma * (1.0f / GammaMultiplier);
+    }
+
+    public void SetGamma(float newGamma)
+    {
+        if (newGamma is < 1.0f / MaxGamma or > 1.0f)
+        {
+            throw new InvalidOperationException($"Invalid gamma {newGamma}");
+        }
+
+        this.HaveGamma = false;
+
+        if (IsAlmostEqual(newGamma, 1.0f))
+        {
+            this.TransferFunction = JxlTransferFunction.Linear;
+            return;
+        }
+
+        if (IsAlmostEqual(newGamma, 1.0f / 2.6f))
+        {
+            this.TransferFunction = JxlTransferFunction.Dci;
+            return;
+        }
+
+        // Don't translate 0.45.. to kSRGB nor k709 - that might change pixel
+        // values because those curves also have a linear part.
+        this.HaveGamma = true;
+        this.Gamma = (uint)MathF.Round((float)(newGamma * GammaMultiplier));
+        this.TransferFunction = JxlTransferFunction.Unknown;
+    }
+
+    public readonly bool IsSame(JxlCustomTransferFunction other)
+    {
+        if (this.HaveGamma != other.HaveGamma)
+        {
+            return false;
+        }
+
+        if (this.HaveGamma)
+        {
+            return this.Gamma == other.Gamma;
+        }
+
+        return this.TransferFunction == other.TransferFunction;
+    }
+
+    private static bool IsAlmostEqual(float a, float b)
+    {
+        const float dist = 1e-3f;
+        return MathF.Abs(a - b) < dist;
+    }
+}

src/ImageSharp/Formats/Jxl/Cms/JxlCustomXy.cs

@@ -0,0 +1,53 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+using SixLabors.ImageSharp.ColorProfiles;
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Cms;
+
+/// <summary>
+/// A serializable form of CieXyChromaticityCoordinates
+/// </summary>
+internal struct JxlCustomXy
+{
+    private const uint Multiplier = 1000000;
+    private const float RoughLimit = 4.0f;
+    private const int Min = -0x200000;
+    private const int Max = 0x1FFFFF;
+
+    public int X { get; set; }
+
+    public int Y { get; set; }
+
+    public readonly CieXyChromaticityCoordinates GetValue() => new(
+        x: this.X * (1.0f / Multiplier),
+        y: this.Y * (1.0f / Multiplier));
+
+    public bool SetValue(CieXyChromaticityCoordinates xy)
+    {
+        bool ok = (Math.Abs(xy.X) < RoughLimit) && (Math.Abs(xy.Y) < RoughLimit);
+
+        if (!ok)
+        {
+            throw new InvalidOperationException("X or Y is out of bounds");
+        }
+
+        this.X = (int)MathF.Round((float)(xy.X * Multiplier));
+
+        if (this.X is < Min or > Max)
+        {
+            throw new InvalidOperationException("X is out of bounds");
+        }
+
+        this.Y = (int)MathF.Round((float)(xy.Y * Multiplier));
+
+        if (this.Y is < Min or > Max)
+        {
+            throw new InvalidOperationException("Y is out of bounds");
+        }
+
+        return true;
+    }
+
+    public readonly bool IsSame(CieXyChromaticityCoordinates other) => this.X == other.X && this.Y == other.Y;
+}

src/ImageSharp/Formats/Jxl/Cms/JxlPrimaries.cs

@@ -0,0 +1,27 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Cms;
+
+/// <summary>
+/// JPEG XL primaries
+/// </summary>
+internal enum JxlPrimaries : byte
+{
+    /// <summary>
+    /// Same as ITU-R BT.709
+    /// </summary>
+    SRgb = 1,
+
+    /// <summary>
+    /// Values encoded in separate fields
+    /// </summary>
+    Custom = 2,
+
+    /// <summary>
+    /// ITU-R BT.2020
+    /// </summary>
+    Bt2020 = 9,
+
+    P3 = 11,
+}

src/ImageSharp/Formats/Jxl/Cms/JxlRenderingIntent.cs

@@ -0,0 +1,13 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Cms;
+
+internal enum JxlRenderingIntent : byte
+{
+    // Values match ICC sRGB encodings
+    Perceptual, // Good for photos, requires a profile with LUT
+    Relative,   // Good for logos
+    Saturation, // Perhaps useful for CG with fully saturated colors
+    Absolute,   // Leaves white point unchanged; good for proofing
+}

src/ImageSharp/Formats/Jxl/Cms/JxlTransferFunction.cs

@@ -0,0 +1,45 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Cms;
+
+/// <summary>
+/// JPEG XL transfer function type
+/// </summary>
+internal enum JxlTransferFunction : byte
+{
+    /// <summary>
+    /// ITU-R BT.709
+    /// </summary>
+    Bt709 = 1,
+
+    /// <summary>
+    /// Unknown transfer function
+    /// </summary>
+    Unknown = 2,
+
+    /// <summary>
+    /// Linear transfer function
+    /// </summary>
+    Linear = 8,
+
+    /// <summary>
+    /// sRGB
+    /// </summary>
+    SRgb = 13,
+
+    /// <summary>
+    /// From ITU-R BT.2100
+    /// </summary>
+    Pq = 16,
+
+    /// <summary>
+    /// From SMPTE RP 431-2 reference projector
+    /// </summary>
+    Dci = 17,
+
+    /// <summary>
+    /// From ITU-R BT.2100
+    /// </summary>
+    Hlg = 18,
+}

src/ImageSharp/Formats/Jxl/Cms/JxlWhitePoint.cs

@@ -0,0 +1,33 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Cms;
+
+/// <summary>
+/// White point from CICP Color Primaries.
+/// </summary>
+// Note that we define a separate enum instead of using the ColorPrimaries
+// enum from CICP code because JPEG XL doesn't support all color primaries defined
+// by CICP.
+internal enum JxlWhitePoint : byte
+{
+    /// <summary>
+    /// sRGB/ITU-R BT.709/Display P3/ITU-R BT.2020
+    /// </summary>
+    D65 = 1,
+
+    /// <summary>
+    /// Actual values encoded in separate fields
+    /// </summary>
+    Custom = 2,
+
+    /// <summary>
+    /// XYZ
+    /// </summary>
+    E = 10,
+
+    /// <summary>
+    /// DCI-P3
+    /// </summary>
+    Dci = 11,
+}

src/ImageSharp/Formats/Jxl/Cms/README.md

@@ -0,0 +1,4 @@
+# CMS
+This is the JPEG XL Color Management System component.
+
+Not to be confused with Content Management System.

src/ImageSharp/Formats/Jxl/IO/Metadata/JxlImageMetadata.cs

@@ -2,6 +2,7 @@
 // Licensed under the Six Labors Split License.
 
 using System.Diagnostics;
+using SixLabors.ImageSharp.Formats.Jxl.Cms;
 using SixLabors.ImageSharp.Formats.Jxl.Fields;
 
 namespace SixLabors.ImageSharp.Formats.Jxl.IO.Metadata;
@@ -124,5 +125,26 @@ internal sealed class JxlImageMetadata : IJxlFields
         this.Modular16BitBufferSufficient = bits <= 12;
     }
 
+    public void SetIntensityTarget()
+    {
+        JxlCustomTransferFunction? tf = this.ColorEncoding?.TransferFunction;
+
+        if (tf is not null)
+        {
+            if (tf.Value.IsPq)
+            {
+                this.SetIntensityTarget(10000);
+            }
+            else if (tf.Value.IsHlg)
+            {
+                this.SetIntensityTarget(1000);
+            }
+            else
+            {
+                this.SetIntensityTarget(DefaultIntensityTarget);
+            }
+        }
+    }
+
     public bool Visit(JxlVisitor visitor) => throw new NotImplementedException();
 }

src/ImageSharp/Formats/Jxl/IO/Metadata/JxlOpsinInverseMatrix.cs

@@ -12,9 +12,13 @@ internal sealed class JxlOpsinInverseMatrix : IJxlFields
 
     public JxlMatrix3x3F InverseMatrix { get; set; }
 
-    public InlineArray3<float> OpsinBiases { get; set; }
+    // Prefer arrays so we can set values like this:
+    //      JxlOpsinInverseMatrix m = ...;
+    //      m.OpsinBiases[0] = 1f;
+    // An InlineArray can't do that.
+    public float[] OpsinBiases { get; set; } = new float[3];
 
-    public InlineArray4<float> QuantBiases { get; set; }
+    public float[] QuantBiases { get; set; } = new float[4];
 
     public bool Visit(JxlVisitor visitor) => throw new NotImplementedException();
 }

src/ImageSharp/Formats/Jxl/InlineArrays.cs

@@ -13,6 +13,15 @@ internal struct InlineArray3<T>
     private T first;
 }
 
+/// <summary>
+/// Used by JxlOpsinParameters
+/// </summary>
+[InlineArray(36)]
+internal struct InlineArray36<T>
+{
+    private T first;
+}
+
 /// <summary>
 /// Used by JxlCustomTransformData
 /// </summary>

src/ImageSharp/Formats/Jxl/Processing/Decoder/JxlDecoderCore.cs

@@ -0,0 +1,157 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+using SixLabors.ImageSharp.IO;
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Decoder;
+
+internal sealed class JxlDecoderCore : ImageDecoderCore
+{
+    /// <summary>
+    /// Identifies the signature of the JPEG XL file.
+    /// </summary>
+    private enum JxlSignature : byte
+    {
+        /// <summary>
+        /// Error status indicating not enough bytes to detect the signature.
+        /// </summary>
+        NotEnoughBytes,
+
+        /// <summary>
+        /// A JPEG XL code stream.
+        /// </summary>
+        CodeStream,
+
+        /// <summary>
+        /// The signature is invalid.
+        /// </summary>
+        Invalid,
+
+        /// <summary>
+        /// Container format.
+        /// </summary>
+        Container
+    }
+
+    /// <summary>
+    /// Represents a data type.
+    /// </summary>
+    private enum JxlDataType : byte
+    {
+        /// <summary>
+        /// <see cref="byte"/>
+        /// </summary>
+        UInt8,
+
+        /// <summary>
+        /// <see cref="ushort"/>
+        /// </summary>
+        UInt16,
+
+        /// <summary>
+        /// <see cref="float"/>
+        /// </summary>
+        Float,
+
+        /// <summary>
+        /// <see cref="Half"/>
+        /// </summary>
+        Float16
+    }
+
+    public JxlDecoderCore(DecoderOptions options)
+        : base(options)
+    {
+    }
+
+    /// <summary>
+    /// Ensures that the coordinates are not out of bounds.
+    /// </summary>
+    /// <param name="a">First coordinate</param>
+    /// <param name="b">Second coordinate</param>
+    /// <param name="size">Image width</param>
+    /// <returns>Boolean indicating whether the coordinates are out of bounds</returns>
+    private static bool IsOutOfBounds(int a, int b, int size)
+    {
+        int position = a + b;
+
+        return position > size || position < a;
+    }
+
+    private static int InitialBasicInfoSizeHint()
+    {
+        const int containerHeaderSize = 48;
+        const int maxCodestreamBasicInfoSize = 50;
+        return containerHeaderSize + maxCodestreamBasicInfoSize;
+    }
+
+    private static JxlSignature DetectSignature(ReadOnlySpan<byte> buffer, int length, ref int position)
+    {
+        if (position >= length)
+        {
+            return JxlSignature.NotEnoughBytes;
+        }
+
+        buffer = buffer[position..];
+        length -= position;
+
+        // 0xFF 0x0A represents a codestream
+        if (length >= 1 && buffer[0] == 0xFF)
+        {
+            if (length < 2)
+            {
+                // We need at least two bytes for a valid codestream signature
+                return JxlSignature.NotEnoughBytes;
+            }
+            else if (buffer[1] == CodestreamMarker)
+            {
+                position += 2;
+                return JxlSignature.CodeStream;
+            }
+            else
+            {
+                return JxlSignature.Invalid;
+            }
+        }
+
+        // Container?
+        if (length >= 1 && buffer[0] == 0)
+        {
+            if (length < SignatureBox.Length)
+            {
+                return JxlSignature.NotEnoughBytes;
+            }
+            else if (buffer[SignatureBox.Length..].SequenceEqual(SignatureBox))
+            {
+                position += SignatureBox.Length;
+                return JxlSignature.Container;
+            }
+            else
+            {
+                return JxlSignature.Invalid;
+            }
+        }
+
+        // Signature is invalid
+        return JxlSignature.Invalid;
+    }
+
+    private static JxlSignature DetectSignature(ReadOnlySpan<byte> buffer, int length)
+    {
+        int position = 0;
+        return DetectSignature(buffer, length, ref position);
+    }
+
+    private static int BitsPerChannel(JxlDataType dataType)
+        => dataType switch
+        {
+            JxlDataType.UInt8 => 8,
+            JxlDataType.UInt16 or JxlDataType.Float16 => 16,
+            JxlDataType.Float => 32,
+            _ => 0
+        };
+
+    protected override Image<TPixel> Decode<TPixel>(BufferedReadStream stream, CancellationToken cancellationToken) => throw new NotImplementedException();
+
+    protected override ImageInfo Identify(BufferedReadStream stream, CancellationToken cancellationToken) => throw new NotImplementedException();
+}

src/ImageSharp/Formats/Jxl/Processing/Decoder/JxlNoiseDecoder.cs

@@ -0,0 +1,63 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+using System.Numerics;
+using System.Runtime.CompilerServices;
+using System.Runtime.InteropServices;
+using SixLabors.ImageSharp.Formats.Jxl.Memory.ImageTypes;
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Decoder;
+
+internal static class JxlNoiseDecoder
+{
+    [MethodImpl(MethodImplOptions.AggressiveInlining)]
+    private static void BitsToFloatingPoint(ReadOnlySpan<uint> randomBits, Span<float> floats)
+    {
+        Vector<uint> bits = new(randomBits);
+        Vector<float> rand12 = ((bits >> 9) | new Vector<uint>(0x3F800000u)).As<uint, float>();
+        rand12.StoreUnsafe(ref MemoryMarshal.GetReference(floats));
+    }
+
+    public static void GenerateRandomImage(JxlXorShift rng, Rectangle rectangle, JxlImageF noise)
+    {
+        const int floatsPerBatch = JxlXorShift.Generators * sizeof(ulong) / sizeof(float);
+
+        int xSize = rectangle.Width;
+        int ySize = rectangle.Height;
+
+        Span<ulong> batch64 = stackalloc ulong[JxlXorShift.Generators];
+        Span<uint> batch32 = stackalloc uint[JxlXorShift.Generators * 2];
+
+        // stackalloc doesn't zero-initialize, so clear values
+        batch64.Clear();
+        batch32.Clear();
+
+        int n = Vector<float>.Count;
+
+        for (int y = 0; y < ySize; y++)
+        {
+            Span<float> row = noise.GetRow(rectangle, y);
+            int x = 0;
+            for (; x + floatsPerBatch < xSize; x += floatsPerBatch)
+            {
+                rng.Fill(batch64);
+                MemoryMarshal.Cast<uint, ulong>(batch32).CopyTo(batch64);
+                for (int i = 0; i < floatsPerBatch; i += n)
+                {
+                    BitsToFloatingPoint(batch32[i..], row[(x + i)..]);
+                }
+            }
+
+            rng.Fill(batch64);
+            MemoryMarshal.Cast<uint, ulong>(batch32).CopyTo(batch64);
+
+            int batchPos = 0;
+
+            for (; x < xSize; x += n)
+            {
+                BitsToFloatingPoint(batch32[batchPos..], row[x..]);
+                batchPos += n;
+            }
+        }
+    }
+}

src/ImageSharp/Formats/Jxl/Processing/Decoder/JxlOpsinParameters.cs

@@ -0,0 +1,18 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Decoder;
+
+internal sealed class JxlOpsinParameters
+{
+    // Use arrays instead of InlineArrays because, with inline arrays we can't do:
+    //      JxlOpsinParameters parameters = ...;
+    //      parameters.OpsinBiasesCbrt[0] /* <-- error */ = 1.25f;
+    public float[] InverseOpsinMatrix { get; set; } = new float[36];
+
+    public float[] OpsinBiases { get; set; } = new float[4];
+
+    public float[] OpsinBiasesCbrt { get; set; } = new float[4];
+
+    public float[] QuantBiases { get; set; } = new float[4];
+}

src/ImageSharp/Formats/Jxl/Processing/Decoder/JxlOutputEncodingInfo.cs

@@ -0,0 +1,78 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+using System.Numerics;
+using SixLabors.ImageSharp.Formats.Jxl.Cms;
+using SixLabors.ImageSharp.Formats.Jxl.IO.Metadata;
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Decoder;
+
+// Prefer class instead of struct because it's too large for a struct
+// Note that this struct wasn't well documented, so it's not that easy to add
+// XML documentation here.
+internal sealed class JxlOutputEncodingInfo
+{
+    public JxlColorEncoding? OriginalColorEncoding { get; set; }
+
+    public float OriginalIntensityTarget { get; set; }
+
+    public JxlMatrix3x3F OriginalInverseMatrix { get; set; }
+
+    public bool DefaultTransform { get; set; }
+
+    public bool XybEncoded { get; set; }
+
+    /// <summary>
+    /// Gets or sets the requested color encoding.
+    /// </summary>
+    public JxlColorEncoding ColorEncoding { get; set; } = new();
+
+    public JxlColorEncoding LinearColorEncoding { get; set; } = new();
+
+    public bool ColorEncodingIsOriginal { get; set; }
+
+    public JxlOpsinParameters OpsinParameters { get; set; } = new();
+
+    public bool AllDefaultOpsin { get; set; }
+
+    public float InverseGamma { get; set; }
+
+    public Vector3 Luminances { get; set; }
+
+    public float DesiredIntensityTarget { get; set; }
+
+    public bool CmsSet { get; set; }
+
+    public JxlCmsInterface Cms { get; set; }
+
+    public void SetFromMetadata(JxlCodecMetadata metadata)
+    {
+        JxlImageMetadata imageMetadata = metadata.ImageMetadata ?? throw new InvalidOperationException("Missing image metadata");
+
+        this.OriginalColorEncoding = imageMetadata.ColorEncoding;
+        this.OriginalIntensityTarget = imageMetadata.IntensityTarget;
+        this.DesiredIntensityTarget = this.OriginalIntensityTarget;
+
+        JxlOpsinInverseMatrix inverseMatrix = metadata.CustomTransformData?.OpsinInverseMatrix ?? throw new InvalidOperationException("Missing Opsin inverse matrix or transform data");
+        this.OriginalInverseMatrix = inverseMatrix.InverseMatrix;
+        this.DefaultTransform = inverseMatrix.AllDefault;
+        this.XybEncoded = imageMetadata.XybEncoded;
+
+        JxlOpsinParameters parameters = this.OpsinParameters;
+
+        imageMetadata.OpsinBiases.CopyTo(parameters.OpsinBiases);
+        parameters.OpsinBiasesCbrt[0] = MathF.Cbrt(parameters.OpsinBiases[0]);
+        parameters.OpsinBiasesCbrt[1] = MathF.Cbrt(parameters.OpsinBiases[1]);
+        parameters.OpsinBiasesCbrt[2] = MathF.Cbrt(parameters.OpsinBiases[2]);
+
+        parameters.OpsinBiasesCbrt[3] = 1;
+        parameters.OpsinBiases[3] = 1;
+
+        inverseMatrix.QuantBiases.AsSpan().CopyTo(parameters.QuantBiases);
+
+        bool origOK = JxlXybDecoder.CanOutputToColorEncoding(this.OriginalColorEncoding ?? throw new InvalidCastException("Missing color encoding"));
+        bool origGrey = this.OriginalColorEncoding.IsGray;
+
+        return this.SetColorEncoding(!this.XybEncoded || origOK ? this.OriginalColorEncoding : JxlColorEncoding.LinearSrgb(origGrey));
+    }
+}

src/ImageSharp/Formats/Jxl/Processing/Decoder/JxlXybDecoder.cs

@@ -0,0 +1,180 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+using System.Numerics;
+using System.Runtime.CompilerServices;
+using System.Runtime.InteropServices;
+using SixLabors.ImageSharp.Formats.Jxl.Cms;
+using SixLabors.ImageSharp.Formats.Jxl.Memory.ImageTypes;
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Decoder;
+
+/// <summary>
+/// Decodes the XYB color format (which JPEG XL uses) into RGB.
+/// </summary>
+internal static class JxlXybDecoder
+{
+    /// <summary>
+    /// Converts XYB to RGB using SIMD, one vector at a time.
+    /// </summary>
+    /// <param name="opsinX">X channel</param>
+    /// <param name="opsinY">Y channel</param>
+    /// <param name="opsinB">B channel</param>
+    /// <param name="opsinParameters">Opsin parameters &amp; configuration</param>
+    /// <param name="linearR">Output R</param>
+    /// <param name="linearG">Output G</param>
+    /// <param name="linearB">Output B</param>
+    public static void ConvertXybToRgb(
+        Vector<float> opsinX,
+        Vector<float> opsinY,
+        Vector<float> opsinB,
+        JxlOpsinParameters opsinParameters,
+        ref Vector<float> linearR,
+        ref Vector<float> linearG,
+        ref Vector<float> linearB)
+    {
+        Vector<float> negBiasR = new(opsinParameters.OpsinBiaases[0]);
+        Vector<float> negBiasG = new(opsinParameters.OpsinBiaases[1]);
+        Vector<float> negBiasB = new(opsinParameters.OpsinBiaases[2]);
+
+        Vector<float> gammaR = opsinX + opsinY;
+        Vector<float> gammaG = opsinY - opsinX;
+        Vector<float> gammaB = opsinB;
+
+        Vector<float> gammaR2 = gammaR * gammaR;
+        Vector<float> gammaG2 = gammaG * gammaG;
+        Vector<float> gammaB2 = gammaB * gammaB;
+
+        Vector<float> mixedR = (gammaR2 * gammaR) + negBiasR;
+        Vector<float> mixedG = (gammaG2 * gammaG) + negBiasG;
+        Vector<float> mixedB = (gammaB2 * gammaB) + negBiasB;
+
+        Span<float> inverseMatrix = opsinParameters.GetInverseOpsinMatrixSpan();
+
+        linearR = LoadDuplicate128(ref inverseMatrix[0 * 4]) * mixedR;
+        linearG = LoadDuplicate128(ref inverseMatrix[3 * 4]) * mixedR;
+        linearB = LoadDuplicate128(ref inverseMatrix[6 * 4]) * mixedR;
+
+        linearR = (LoadDuplicate128(ref inverseMatrix[1 * 4]) * mixedG) + linearR;
+        linearG = (LoadDuplicate128(ref inverseMatrix[4 * 4]) * mixedG) + linearG;
+        linearB = (LoadDuplicate128(ref inverseMatrix[7 * 4]) * mixedG) + linearB;
+
+        linearR = (LoadDuplicate128(ref inverseMatrix[2 * 4]) * mixedB) + linearR;
+        linearG = (LoadDuplicate128(ref inverseMatrix[5 * 4]) * mixedB) + linearG;
+        linearB = (LoadDuplicate128(ref inverseMatrix[8 * 4]) * mixedB) + linearB;
+    }
+
+    public static bool OpsinToLinear(JxlImage3F opsin, Rectangle rect, JxlImage3F linear, JxlOpsinParameters opsinParameters)
+    {
+        if (!SameSize(rect, linear))
+        {
+            return false;
+        }
+
+        if (Vector<float>.Count < 4)
+        {
+            // TODO: support 64bit vectors or no SIMD?
+            throw new PlatformNotSupportedException("XYB to RGB conversion requires at least 128-bit SIMD");
+        }
+
+        // Reuse variables instead of creating them over
+        // and over again
+        Unsafe.SkipInit(out Vector<float> linearR);
+        Unsafe.SkipInit(out Vector<float> linearG);
+        Unsafe.SkipInit(out Vector<float> linearB);
+
+        for (int y = 0; y < rect.Height; y++)
+        {
+            ReadOnlySpan<float> rowOpsin0 = opsin.PlaneRow(rect, 0, y);
+            ReadOnlySpan<float> rowOpsin1 = opsin.PlaneRow(rect, 1, y);
+            ReadOnlySpan<float> rowOpsin2 = opsin.PlaneRow(rect, 2, y);
+
+            ref float rowOpsin0Reference = ref MemoryMarshal.GetReference(rowOpsin0);
+            ref float rowOpsin1Reference = ref MemoryMarshal.GetReference(rowOpsin1);
+            ref float rowOpsin2Reference = ref MemoryMarshal.GetReference(rowOpsin2);
+
+            Span<float> rowLinear0 = linear.PlaneRow(0, y);
+            Span<float> rowLinear1 = linear.PlaneRow(1, y);
+            Span<float> rowLinear2 = linear.PlaneRow(2, y);
+
+            ref float rowLinear0Reference = ref MemoryMarshal.GetReference(rowLinear0);
+            ref float rowLinear1Reference = ref MemoryMarshal.GetReference(rowLinear1);
+            ref float rowLinear2Reference = ref MemoryMarshal.GetReference(rowLinear2);
+
+            for (int x = 0; x < rect.Height; x += Vector<float>.Count)
+            {
+                Vector<float> inOpsinX = Vector.LoadUnsafe(ref Unsafe.Add(ref rowOpsin0Reference, x));
+                Vector<float> inOpsinY = Vector.LoadUnsafe(ref Unsafe.Add(ref rowOpsin1Reference, x));
+                Vector<float> inOpsinB = Vector.LoadUnsafe(ref Unsafe.Add(ref rowOpsin2Reference, x));
+
+                ConvertXybToRgb(inOpsinX, inOpsinY, inOpsinB, opsinParameters, ref linearR, ref linearG, ref linearB);
+
+                linearR.StoreUnsafe(ref Unsafe.Add(ref rowLinear0Reference, x));
+                linearG.StoreUnsafe(ref Unsafe.Add(ref rowLinear1Reference, x));
+                linearB.StoreUnsafe(ref Unsafe.Add(ref rowLinear2Reference, x));
+            }
+        }
+
+        return true;
+    }
+
+    /// <summary>
+    /// A SIMD utility method which reads next 128 bits
+    /// (which in this case happens to be next 4 floats),
+    /// and duplicates them to fit in the CPU vector size.
+    /// For example,
+    ///
+    ///     128 bit vectors: A B C D    (as-is)
+    ///     256 bit vectors: A B C D A B C D (duplicate once)
+    ///     512 bit vectors: A B C D A B C D A B C D A B C D (duplicate three times)
+    ///
+    /// Vector&lt;T&gt; has support for arbitrarily large
+    /// vector sizes. For example, some ARM CPUs support 2048-bit
+    /// vectors through Vector&lt;T&gt;. In that specific case, this
+    /// method can be used for future-proofing.
+    ///
+    /// Note that this method, albeit future-proof, may be considered
+    /// slow for smaller vector sizes (think CPUs with 128bit or 256bit vectors).
+    /// </summary>
+    /// <param name="reference">Reference to first element to load &amp; duplicate.</param>
+    /// <returns>Vector with first 128 bits duplicated across the vector width.</returns>
+    private static Vector<float> LoadDuplicate128(ref float reference)
+    {
+        Span<float> value = stackalloc float[Vector<float>.Count];
+        Span<float> values128 = [
+            reference,
+            Unsafe.Add(ref reference, 1),
+            Unsafe.Add(ref reference, 2),
+            Unsafe.Add(ref reference, 3)
+        ];
+
+        for (int i = 0; i < Vector<float>.Count; i += 4)
+        {
+            values128[i..].CopyTo(value[i..]);
+        }
+
+        return new(value);
+    }
+
+    public static bool CanOutputToColorEncoding(JxlColorEncoding colorEncoding)
+    {
+        if (!colorEncoding.HaveFields)
+        {
+            return false;
+        }
+
+        JxlCustomTransferFunction tf = colorEncoding.TransferFunction;
+
+        if (!tf.IsPq && !tf.IsSrgb && !tf.HaveGamma && !tf.IsLinear && !tf.IsHlg && !tf.IsDci && !tf.Is709)
+        {
+            return false;
+        }
+
+        if (colorEncoding.IsGray && colorEncoding.WhitePoint != JxlWhitePoint.D65)
+        {
+            return false;
+        }
+
+        return true;
+    }
+}

src/ImageSharp/Formats/Jxl/Processing/JxlCoefficientOrder.cs

@@ -0,0 +1,91 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+using System.Runtime.CompilerServices;
+using System.Security.Principal;
+using SixLabors.ImageSharp.Formats.Jxl.IO.Entropy;
+using SixLabors.ImageSharp.Formats.Jxl.Processing.Decoder;
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Processing;
+
+/// <summary>
+/// Helps with reordering coefficients.
+/// </summary>
+internal static class JxlCoefficientOrder
+{
+    public const int Limit = 6156;
+
+    public const int CoefficientOrderMaxSize = Limit * JxlFrameDimensions.DctBlockSize;
+
+    public const int PermutationContexts = 8;
+
+    /// <summary>
+    /// Gets the pattern which coefficients must follow to compute offsets.
+    /// </summary>
+    public static ReadOnlySpan<int> CoefficientOrderOffsets =>
+    [
+        0,    1,    2,    3,    4,    5,    6,    10,   14,   18,
+        34,   50,   66,   68,   70,   72,   76,   80,   84,   92,
+        100,  108,  172,  236,  300,  332,  364,  396,  652,  908,
+        1164, 1292, 1420, 1548, 2572, 3596, 4620, 5132, 5644, Limit
+    ];
+
+    /// <summary>
+    /// Gets the pattern which coefficients must follow to compute offsets.
+    /// </summary>
+    public static ReadOnlySpan<byte> StrategyOrder =>
+    [
+        0, 1, 1, 1, 2, 3, 4, 4, 5,  5,  6,  6,  1,  1,
+        1, 1, 1, 1, 7, 8, 8, 9, 10, 10, 11, 12, 12,
+    ];
+
+    [MethodImpl(MethodImplOptions.AggressiveInlining)]
+    private static int CoeffOrderOffset(int o, int c) => CoefficientOrderOffsets[(3 * o) + c] * JxlFrameDimensions.DctBlockSize;
+
+    [MethodImpl(MethodImplOptions.AggressiveInlining)]
+    public static uint CoeffOrderContext(uint value)
+    {
+        uint token = 0;
+        uint nbits = 0;
+        uint bits = 0;
+
+        new JxlAnsHybridUIntConfiguration(0, 0, 0).Encode(value, ref token, ref nbits, ref bits);
+
+        return Math.Min(token, PermutationContexts - 1u);
+    }
+
+    public static bool ReadPermutation(int skip, int size, Span<int> order, JxlBitReader bitReader, JxlAnsSymbolReader reader, Span<byte> contextMap)
+    {
+        Span<uint> lehmer = stackalloc uint[size];
+        lehmer.Clear();
+
+        Span<uint> temp = stackalloc uint[size * 2];
+        temp.Clear();
+
+        uint end = reader.ReadHybridUnsignedInteger(CoeffOrderContext((int)size), bitReader, contextMap) + skip;
+
+        if (end > size)
+        {
+            throw new InvalidOperationException("Invalid permutation size");
+        }
+
+        uint last = 0;
+
+        for (int i = skip; i < end; i++)
+        {
+            lehmer[i] = reader.ReadHybridUnsignedInteger(CoeffOrderContext(last), bitReader, contextMap);
+            last = lehmer[i];
+            if (lehmer[i] >= size - i)
+            {
+                throw new InvalidOperationException("Invalid lehmer code");
+            }
+        }
+
+        if (order.IsEmpty)
+        {
+            return true;
+        }
+
+        return JxlLehmerCode.DecodeLehmerCode(lehmer, temp, size, order);
+    }
+}

src/ImageSharp/Formats/Jxl/Processing/JxlConvolve.cs

@@ -0,0 +1,424 @@
+// 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.Formats.Jxl.Memory.ImageTypes;
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Processing;
+
+/// <summary>
+/// Convolution filters
+/// </summary>
+internal static class JxlConvolve
+{
+    /// <summary>
+    /// Weighted sum of 1x5 pixels around ix, iy with [wx2, wx1, wx0, wx1, wx2].
+    /// </summary>
+    public static float WeightedSumBorder(
+        JxlImageF input,
+        Func<long, long, int> wrapY,
+        long ix,
+        long iy,
+        int width,
+        int height,
+        float wx0,
+        float wx1,
+        float wx2)
+    {
+        ReadOnlySpan<float> row = input.GetRow(wrapY(iy, height));
+
+        float inM2 = row[WrapMirror(ix - 2, width)];
+        float inP2 = row[WrapMirror(ix + 2, width)];
+        float inM1 = row[WrapMirror(ix - 1, width)];
+        float inP1 = row[WrapMirror(ix + 1, width)];
+        float in00 = row[(int)ix];
+
+        float sum2 = wx2 * (inM2 + inP2);
+        float sum1 = wx1 * (inM1 + inP1);
+        float sum0 = wx0 * in00;
+
+        return sum2 + (sum1 + sum0);
+    }
+
+    public static Vector<float> WeightedSum(
+        JxlImageF input,
+        Func<long, long, int> wrapY,
+        int ix,
+        long iy,
+        int height,
+        Vector<float> wx0,
+        Vector<float> wx1,
+        Vector<float> wx2)
+    {
+        ReadOnlySpan<float> center = input.GetRow(wrapY(iy, height))[ix..];
+        ref float centerRef = ref MemoryMarshal.GetReference(center);
+
+        Vector<float> inM2 = Vector.LoadUnsafe(ref Unsafe.Subtract(ref centerRef, 2));
+        Vector<float> inP2 = Vector.LoadUnsafe(ref Unsafe.Add(ref centerRef, 2));
+        Vector<float> inM1 = Vector.LoadUnsafe(ref Unsafe.Subtract(ref centerRef, 1));
+        Vector<float> inP1 = Vector.LoadUnsafe(ref Unsafe.Add(ref centerRef, 1));
+        Vector<float> in00 = Vector.LoadUnsafe(ref centerRef);
+
+        Vector<float> sum2 = wx2 * (inM2 + inP2);
+        Vector<float> sum1 = wx1 * (inM1 + inP1);
+        Vector<float> sum0 = wx0 * in00;
+
+        return sum2 + (sum1 + sum0);
+    }
+
+    public static float Symmetric5Border(JxlImageF input, Func<long, long, int> wrapY, long ix, long iy, JxlWeightsSymmetric5 weights)
+    {
+        float w0 = weights.GetCVector()[0];
+        float w1 = weights.GetRVector()[0];
+        float w2 = weights.GetR2Vector()[0];
+        float w4 = weights.GetDVector()[0];
+        float w5 = weights.GetCVector()[0];
+        float w8 = weights.GetD2Vector()[0];
+
+        int width = input.XSize;
+        int height = input.YSize;
+
+        float sum0 = WeightedSumBorder(input, wrapY, ix, iy, width, height, w0, w1, w2)
+            + WeightedSumBorder(input, wrapY, ix, iy - 2, width, height, w2, w5, w8);
+
+        float sum1 = WeightedSumBorder(input, wrapY, ix, iy + 2, width, height, w2, w5, w8);
+
+        sum0 += WeightedSumBorder(input, wrapY, ix, iy + 1, width, height, w1, w4, w5);
+        sum1 += WeightedSumBorder(input, wrapY, ix, iy - 1, width, height, w1, w4, w5);
+
+        return sum0 + sum1;
+    }
+
+    public static void Symmetric5Interior(
+        JxlImageF image,
+        int ix,
+        Func<long, long, int> wrapY,
+        int rix,
+        long iy,
+        JxlWeightsSymmetric5 weights,
+        Span<float> rowOut)
+    {
+        Vector<float> w0 = LoadDuplicate128(weights.GetCVector());  // c
+        Vector<float> w1 = LoadDuplicate128(weights.GetRVector());  // r
+        Vector<float> w2 = LoadDuplicate128(weights.GetR2Vector()); // R
+        Vector<float> w4 = LoadDuplicate128(weights.GetDVector());  // d
+        Vector<float> w5 = LoadDuplicate128(weights.GetLVector());  // L
+        Vector<float> w8 = LoadDuplicate128(weights.GetD2Vector()); // D
+
+        int height = image.YSize;
+        Vector<float> sum0 = WeightedSum(image, wrapY, ix, iy, height, w0, w1, w2)
+            + WeightedSum(image, wrapY, ix, iy - 2, height, w2, w5, w8);
+
+        Vector<float> sum1 = WeightedSum(image, wrapY, ix, iy + 2, height, w2, w5, w8);
+
+        sum0 += WeightedSum(image, wrapY, ix, iy - 1, height, w1, w4, w5);
+        sum1 += WeightedSum(image, wrapY, ix, iy + 1, height, w1, w4, w5);
+
+        (sum0 + sum1).StoreUnsafe(ref Unsafe.Add(ref MemoryMarshal.GetReference(rowOut), rix));
+    }
+
+    public static void Symmetric5Row(
+        JxlImageF image,
+        Func<long, long, int> wrapY,
+        in Rectangle rect,
+        long iy,
+        JxlWeightsSymmetric5 weights,
+        Span<float> rowOut)
+    {
+        const int radius = 2;
+        int xEnd = rect.Right;
+
+        int rix = 0;
+        int ix = rect.X;
+
+        int n = Vector<float>.Count;
+        int alignedX = RoundUpTo(radius, n);
+
+        for (; ix < Math.Min(alignedX, xEnd); ix++, rix++)
+        {
+            rowOut[rix] = Symmetric5Border(image, wrapY, ix, iy, weights);
+        }
+
+        for (; ix + n + radius <= xEnd; ix += n, rix += n)
+        {
+            Symmetric5Interior(image, ix, wrapY, rix, iy, weights, rowOut);
+        }
+
+        for (; ix < xEnd; ix++, rix++)
+        {
+            rowOut[rix] = Symmetric5Border(image, wrapY, ix, iy, weights);
+        }
+    }
+
+    public static bool Symmetric5(
+        JxlImageF input,
+        in Rectangle rectangle,
+        JxlWeightsSymmetric5 weights,
+        JxlImageF output,
+        Rectangle outputRect)
+    {
+        if (rectangle.Width != outputRect.Width || rectangle.Height != outputRect.Height)
+        {
+            return false;
+        }
+
+        int height = rectangle.Height;
+
+        for (int riy = 0; riy < height; riy++)
+        {
+            int iy = rectangle.Y + riy;
+
+            if (iy < 2 || iy >= rectangle.Height - 2)
+            {
+                Symmetric5Row(input, WrapMirror, in rectangle, iy, weights, output.GetRow(outputRect, riy));
+            }
+            else
+            {
+                Symmetric5Row(input, in rectangle, iy, weights, output.GetRow(outputRect, riy));
+            }
+        }
+
+        return true;
+    }
+
+    public static float SlowSymmetric3Pixel(
+        JxlImageF image,
+        int x,
+        int y,
+        int width,
+        int height,
+        JxlWeightsSymmetric3 weights,
+        Func<long, long, int> wrapX,
+        Func<long, long, int> wrapY)
+    {
+        float sum = 0.0f;
+
+        float c0 = weights.GetCVector()[0];
+        float r0 = weights.GetRVector()[0];
+        float d0 = weights.GetDVector()[0];
+
+        for (int ky = -1; ky <= 1; ky++)
+        {
+            int yy = wrapY(y + ky, height);
+            ReadOnlySpan<float> row = image.GetRow(yy);
+
+            float wc = (ky == 0) ? c0 : r0;
+            float wlr = (ky == 0) ? r0 : d0;
+
+            int xm1 = wrapX(x - 1, width);
+            int xp1 = wrapX(x + 1, width);
+
+            sum += (row[x] * wc) + ((row[xm1] + row[xp1]) * wlr);
+        }
+
+        return sum;
+    }
+
+    public static void SlowSymmetric3Row(
+        JxlImageF image,
+        int y,
+        int width,
+        int height,
+        JxlWeightsSymmetric3 weights,
+        Span<float> outputRow,
+        Func<long, long, int> wrapY)
+    {
+        outputRow[0] = SlowSymmetric3Pixel(
+            image,
+            0,
+            y,
+            width,
+            height,
+            weights,
+            WrapMirror,
+            wrapY);
+
+        for (int x = 1; x < width - 1; x++)
+        {
+            outputRow[x] = SlowSymmetric3Pixel(
+                image,
+                x,
+                y,
+                width,
+                height,
+                weights,
+                WrapUnchanged,
+                wrapY);
+        }
+
+        outputRow[width - 1] = SlowSymmetric3Pixel(
+            image,
+            width - 1,
+            y,
+            width,
+            height,
+            weights,
+            WrapMirror,
+            wrapY);
+    }
+
+    public static void SlowSymmetric3(
+        JxlImageF input,
+        Rectangle rect,
+        JxlWeightsSymmetric3 weights,
+        JxlImageF output)
+    {
+        int width = rect.Width;
+        int height = rect.Height;
+
+        const int radius = 1;
+
+        for (int y = 0; y < height; y++)
+        {
+            Span<float> rowOut = output.GetRow(y);
+
+            if (y < radius || y >= height - radius)
+            {
+                SlowSymmetric3Row(
+                    input,
+                    y,
+                    width,
+                    height,
+                    weights,
+                    rowOut,
+                    WrapMirror);
+            }
+            else
+            {
+                SlowSymmetric3Row(
+                    input,
+                    y,
+                    width,
+                    height,
+                    weights,
+                    rowOut,
+                    WrapUnchanged);
+            }
+        }
+    }
+
+    public static float SlowSeparablePixel(
+        JxlImageF image,
+        Rectangle rect,
+        int x,
+        int y,
+        int radius,
+        ReadOnlySpan<float> horzWeights,
+        ReadOnlySpan<float> vertWeights)
+    {
+        int width = image.XSize;
+        int height = image.YSize;
+
+        float sum = 0;
+
+        for (int dy = -radius; dy <= radius; dy++)
+        {
+            float wy = vertWeights[Math.Abs(dy) * 4];
+            int sy = WrapMirror(rect.Y + y + dy, height);
+            ReadOnlySpan<float> row = image.GetRow(sy);
+
+            for (int dx = -radius; dx <= radius; dx++)
+            {
+                float wx = horzWeights[Math.Abs(dx) * 4];
+                int sx = WrapMirror(rect.X + x + dx, width);
+                sum += row[sx] * wx * wy;
+            }
+        }
+
+        return sum;
+    }
+
+    public static void SlowSeparable(
+        JxlImageF input,
+        Rectangle inputRect,
+        JxlWeightsSeparable5 weights,
+        JxlImageF output,
+        Rectangle outputRect,
+        int radius)
+    {
+        ReadOnlySpan<float> horz = weights.Horizontal;
+        ReadOnlySpan<float> vert = weights.Vertical;
+
+        for (int y = 0; y < inputRect.Height; y++)
+        {
+            Span<float> rowOut = output.GetRow(outputRect, y);
+
+            for (int x = 0; x < inputRect.Width; x++)
+            {
+                rowOut[x] = SlowSeparablePixel(
+                    input,
+                    inputRect,
+                    x,
+                    y,
+                    radius,
+                    horz,
+                    vert);
+            }
+        }
+    }
+
+    public static void SlowSeparable5(
+        JxlImageF input,
+        Rectangle inputRect,
+        JxlWeightsSeparable5 weights,
+        JxlImageF output,
+        Rectangle outputRect)
+        => SlowSeparable(input, inputRect, weights, output, outputRect, 2);
+
+    public static void FirstL1(ReadOnlySpan<float> c, Span<float> dst)
+    {
+        dst[0] = c[0];
+        for (int i = 1; i < dst.Length; i++)
+        {
+            dst[i] = c[i - 1];
+        }
+    }
+
+    public static void FirstL2(ReadOnlySpan<float> c, Span<float> dst)
+    {
+        dst[0] = c[1];
+        dst[1] = c[0];
+
+        for (int i = 2; i < dst.Length; i++)
+        {
+            dst[i] = c[i - 2];
+        }
+    }
+
+    /// <summary>
+    /// A SIMD utility method which takes in the 128 bit vector
+    /// and duplicates its values to fit in the CPU vector size.
+    /// For example,
+    ///
+    ///     128 bit vectors: A B C D    (as-is)
+    ///     256 bit vectors: A B C D A B C D (duplicate once)
+    ///     512 bit vectors: A B C D A B C D A B C D A B C D (duplicate three times)
+    ///
+    /// Vector&lt;T&gt; has support for arbitrarily large
+    /// vector sizes. For example, some ARM CPUs support 2048-bit
+    /// vectors through Vector&lt;T&gt;. In that specific case, this
+    /// method can be used for future-proofing.
+    ///
+    /// Note that this method, albeit future-proof, may be considered
+    /// slow for smaller vector sizes (think CPUs with 256bit vectors).
+    /// </summary>
+    /// <param name="vec">Vector to duplicate.</param>
+    /// <returns>New vector that is duplicated across the width.</returns>
+    [MethodImpl(MethodImplOptions.AggressiveInlining)]
+    private static Vector<float> LoadDuplicate128(Vector128<float> vec)
+    {
+        Span<float> value = stackalloc float[Vector<float>.Count];
+        for (int i = 0; i < Vector<float>.Count; i += 4)
+        {
+            vec.CopyTo(value[i..]);
+        }
+
+        return new(value);
+    }
+
+    [MethodImpl(MethodImplOptions.AggressiveInlining)]
+    private static int RoundUpTo(int value, int multiple) => ((value + multiple - 1) / multiple) * multiple;
+}

src/ImageSharp/Formats/Jxl/Processing/JxlDctScales.cs

@@ -5,7 +5,7 @@ namespace SixLabors.ImageSharp.Formats.Jxl.Processing;
 
 /// <summary>
 /// Read-only cosine lookups for the Discrete Cosine Transform (DCT),
-/// a mathematical function used for quantization.
+/// a mathematical function used for quantization and coefficient reordering.
 /// </summary>
 internal static class JxlDctScales
 {

src/ImageSharp/Formats/Jxl/Processing/JxlToc.cs

@@ -0,0 +1,191 @@
+// Copyright (c) Six Labors.
+// Licensed under the Six Labors Split License.
+
+using SixLabors.ImageSharp.Formats.Jxl.Fields;
+using SixLabors.ImageSharp.Formats.Jxl.Processing.Decoder;
+
+namespace SixLabors.ImageSharp.Formats.Jxl.Processing;
+
+/// <summary>
+/// Table of Contents encoding
+/// </summary>
+internal static class JxlToc
+{
+    private static readonly JxlU32Enc TocDistribution = new(
+        JxlFieldExpressions.Bits(10),
+        JxlFieldExpressions.BitsOffset(14, 2024),
+        JxlFieldExpressions.BitsOffset(22, 17408),
+        JxlFieldExpressions.BitsOffset(30, 4211712));
+
+    public static int AcGroupIndex(int pass, int group, int numGroups, int numDcGroups)
+        => 2 + numDcGroups + (pass * numGroups) + group;
+
+    public static int NumberOfTocEntries(int numGroups, int numDcGroups, int numPasses)
+    {
+        if (numGroups == 1 && numPasses == 1)
+        {
+            return 1;
+        }
+
+        return AcGroupIndex(0, 0, numGroups, numDcGroups) + (numGroups * numPasses);
+    }
+
+    private const int BitsPerByte = 8;
+    private const int MaxTocEntries = 65536;
+
+    public static bool ReadToc(
+        Configuration configuration,
+        int tocEntries,
+        JxlBitReader reader,
+        List<uint> sizes,
+        List<byte> permutation)
+    {
+        if (tocEntries > MaxTocEntries)
+        {
+            return false; // too many TOC entries
+        }
+
+        sizes.Clear();
+        sizes.Capacity = tocEntries;
+
+        for (int i = 0; i < tocEntries; i++)
+        {
+            sizes.Add(0);
+        }
+
+        if (reader.TotalBitsConsumed >= reader.TotalBytes * BitsPerByte)
+        {
+            return false; // not enough bytes
+        }
+
+        bool CheckBitBudget(int numEntries)
+        {
+            long minimalBitCost = numEntries * (2 + 10);
+            long bitBudget = reader.TotalBytes * BitsPerByte;
+            long expenses = reader.TotalBitsConsumed;
+
+            return expenses <= bitBudget &&
+                minimalBitCost <= bitBudget - expenses;
+        }
+
+        if (tocEntries <= 0)
+        {
+            return false;
+        }
+
+        if (reader.ReadBits32(1) == 1)
+        {
+            if (!CheckBitBudget(tocEntries))
+            {
+                return false;
+            }
+
+            permutation.Clear();
+
+            for (int i = 0; i < tocEntries; i++)
+            {
+                permutation.Add(default);
+            }
+
+            if (!DecodePermutation(
+                    configuration,
+                    0,
+                    tocEntries,
+                    permutation,
+                    reader))
+            {
+                return false;
+            }
+        }
+
+        if (!reader.JumpToByteBoundary())
+        {
+            return false;
+        }
+
+        if (!CheckBitBudget(tocEntries))
+        {
+            return false;
+        }
+
+        for (int i = 0; i < tocEntries; i++)
+        {
+            sizes[i] = JxlU32Coder.Read(TocDistribution, reader);
+        }
+
+        if (!reader.JumpToByteBoundary())
+        {
+            return false;
+        }
+
+        return CheckBitBudget(0);
+    }
+
+    public static bool ReadGroupOffsets(
+        Configuration configuration,
+        int tocEntries,
+        JxlBitReader reader,
+        List<ulong> offsets,
+        List<uint> sizes,
+        out ulong totalSize)
+    {
+        totalSize = 0;
+
+        List<byte> permutation = [];
+
+        if (!ReadToc(
+                configuration,
+                tocEntries,
+                reader,
+                sizes,
+                permutation))
+        {
+            return false;
+        }
+
+        offsets.Clear();
+        offsets.Capacity = tocEntries;
+
+        for (int i = 0; i < tocEntries; i++)
+        {
+            offsets.Add(0);
+        }
+
+        ulong offset = 0;
+
+        for (int i = 0; i < tocEntries; i++)
+        {
+            ulong size = sizes[i];
+
+            if (offset + size < offset)
+            {
+                return false;
+            }
+
+            offsets[i] = offset;
+            offset += size;
+        }
+
+        totalSize = offset;
+
+        if (permutation.Count != 0)
+        {
+            List<ulong> permutedOffsets = new(tocEntries);
+            List<uint> permutedSizes = new(tocEntries);
+
+            foreach (byte index in permutation)
+            {
+                permutedOffsets.Add(offsets[index]);
+                permutedSizes.Add(sizes[index]);
+            }
+
+            offsets.Clear();
+            offsets.AddRange(permutedOffsets);
+
+            sizes.Clear();
+            sizes.AddRange(permutedSizes);
+        }
+
+        return true;
+    }
+}

src/ImageSharp/Formats/Jxl/Processing/JxlXorShift.cs

@@ -7,6 +7,8 @@ namespace SixLabors.ImageSharp.Formats.Jxl.Processing;
 
 internal sealed class JxlXorShift
 {
+    public const int Generators = 8;
+
     private readonly ulong[] s0 = new ulong[8];
     private readonly ulong[] s1 = new ulong[8];