diff --git a/src/ImageSharp/Formats/Jxl/Memory/JxlImage3{T}.cs b/src/ImageSharp/Formats/Jxl/Memory/JxlImage3{T}.cs index 91aed2e97a..4442eb6495 100644 --- a/src/ImageSharp/Formats/Jxl/Memory/JxlImage3{T}.cs +++ b/src/ImageSharp/Formats/Jxl/Memory/JxlImage3{T}.cs @@ -94,6 +94,12 @@ internal class JxlImage3 : IDisposable DebugGuard.MustBeLessThan(y, this.YSize, nameof(y)); } + /// + /// Returns the rectangle for this image bounds. + /// + /// A rectangle with x,y=0,0 width,height=XSize,YSize. + public Rectangle GetRectangle() => new(0, 0, this.XSize, this.YSize); + public void Dispose() { foreach (JxlPlane plane in this.planes) diff --git a/src/ImageSharp/Formats/Jxl/Processing/Encoder/Huffman/JxlHuffmanEncoder.cs b/src/ImageSharp/Formats/Jxl/Processing/Encoder/Huffman/JxlHuffmanEncoder.cs new file mode 100644 index 0000000000..7caf854655 --- /dev/null +++ b/src/ImageSharp/Formats/Jxl/Processing/Encoder/Huffman/JxlHuffmanEncoder.cs @@ -0,0 +1,226 @@ +// Copyright (c) Six Labors. +// Licensed under the Six Labors Split License. + +using SixLabors.ImageSharp.Common.Helpers; + +namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Encoder.Huffman; + +/// +/// Derives & writes Huffman codes. +/// +internal static class JxlHuffmanEncoder +{ + private const int CodeLengthCodes = 18; + + private static ReadOnlySpan StorageOrder => [1, 2, 3, 4, 0, 5, 17, 6, 16, 7, 8, 9, 10, 11, 12, 13, 14, 15]; + + private static ReadOnlySpan HuffmanBitLengthHuffmanCodeSymbols => [0, 7, 3, 2, 1, 15]; + + private static ReadOnlySpan HuffmanBitLengthHuffmanCodeBitLengths => [2, 4, 3, 2, 2, 4]; + + public static void StoreHuffmanTreeOfHuffmanTreeToBitMask(int numCodes, Span codeLengthBitDepth, JxlBitWriter writer) + { + int codesToStore = CodeLengthCodes; + if (numCodes > 1) + { + for (; codesToStore > 0; codesToStore--) + { + if (codeLengthBitDepth[StorageOrder[codesToStore - 1]] != 0) + { + break; + } + } + } + + int skipSome = 0; + if (codeLengthBitDepth[StorageOrder[0]] == 0 && codeLengthBitDepth[StorageOrder[1]] == 0) + { + skipSome = 2; // skips two + if (codeLengthBitDepth[StorageOrder[2]] == 0) + { + skipSome = 3; // skips three + } + } + + writer.Write(2, skipSome); + + for (int i = skipSome; i < codesToStore; ++i) + { + int l = codeLengthBitDepth[StorageOrder[i]]; + writer.Write(HuffmanBitLengthHuffmanCodeBitLengths[l], HuffmanBitLengthHuffmanCodeSymbols[l]); + } + } + + public static void StoreHuffmanTreeToBitMask(int huffmanTreeSize, Span huffmanTree, Span huffmanTreeExtraBits, Span codeLengthBitDepth, Span codeLengthBitDepthSymbols, JxlBitWriter writer) + { + for (int i = 0; i < huffmanTreeSize; ++i) + { + int ix = huffmanTree[i]; + writer.Write(codeLengthBitDepth[ix], codeLengthBitDepthSymbols[ix]); + DebugGuard.MustBeLessThan(ix, 17, nameof(ix)); + + // Extra bits + // + // Micro optimization: + // Original: + // switch (ix) + // { + // case 16: + // writer->Write(2, huffman_tree_extra_bits[i]); + // break; + // case 17: + // writer->Write(3, huffman_tree_extra_bits[i]); + // break; + // default: + // // no-op + // break; + // } + if ((ix & 16) != 0) + { + writer.Write(2 + (ix & 1), huffmanTreeExtraBits[i]); + } + } + } + + public static void StoreSimpleHuffmanTree(Span depths, InlineArray4 symbols, int numSymbols, int maxBits, JxlBitWriter writer) + { + writer.Write(2, 1); + writer.Write(2, numSymbols - 1); + + for (int i = 0; i < numSymbols; i++) + { + for (int j = i + 1; j < numSymbols; j++) + { + if (depths[symbols[j]] < depths[symbols[i]]) + { + RuntimeUtility.Swap(ref symbols[j], ref symbols[i]); + } + } + } + + if (numSymbols == 2) + { + writer.Write(maxBits, symbols[0]); + writer.Write(maxBits, symbols[1]); + } + else if (numSymbols == 3) + { + writer.Write(maxBits, symbols[0]); + writer.Write(maxBits, symbols[1]); + writer.Write(maxBits, symbols[2]); + } + else + { + writer.Write(maxBits, symbols[0]); + writer.Write(maxBits, symbols[1]); + writer.Write(maxBits, symbols[2]); + writer.Write(maxBits, symbols[3]); + writer.Write(1, depths[symbols[0]] == 1 ? 1 : 0); + } + } + + public static void StoreHuffmanTree(Span depths, int num, JxlBitWriter writer) + { + Span arena = stackalloc byte[2 * num]; + Span huffmanTree = arena; + Span huffmanTreeExtraBits = arena[num..]; + int huffmanTreeSize = 0; + JxlHuffmanTree.WriteHuffmanTree(depths, num, ref huffmanTreeSize, huffmanTree, huffmanTreeExtraBits); + + Span huffmanTreeHistogram = stackalloc int[CodeLengthCodes]; + huffmanTreeHistogram.Clear(); + + for (int i = 0; i < huffmanTreeSize; ++i) + { + huffmanTreeHistogram[huffmanTree[i]]++; + } + + int numCodes = 0; + int code = 0; + + for (int i = 0; i < CodeLengthCodes; ++i) + { + if (huffmanTreeHistogram[i] != 0) + { + if (numCodes == 0) + { + code = i; + numCodes = 1; + } + else if (numCodes == 1) + { + numCodes = 2; + break; + } + } + } + + Span codeLengthBitDepth = stackalloc byte[CodeLengthCodes]; + Span codeLengthBitDepthSymbols = stackalloc short[CodeLengthCodes]; + codeLengthBitDepth.Clear(); + codeLengthBitDepthSymbols.Clear(); + + JxlHuffmanTree.CreateHuffmanTree(huffmanTreeHistogram, CodeLengthCodes, 5, codeLengthBitDepth); + JxlHuffmanTree.ConvertBitDepthsToSymbols(codeLengthBitDepth, CodeLengthCodes, codeLengthBitDepthSymbols); + + StoreHuffmanTreeOfHuffmanTreeToBitMask(numCodes, codeLengthBitDepth, writer); + + if (numCodes == 1) + { + codeLengthBitDepth[code] = 0; + } + + StoreHuffmanTreeToBitMask(huffmanTreeSize, huffmanTree, huffmanTreeExtraBits, codeLengthBitDepth, codeLengthBitDepthSymbols, writer); + } + + public static void BuildAndStoreHuffmanTree(Span histogram, int length, Span depth, Span bits, JxlBitWriter writer) + { + int count = 0; + InlineArray4 s4 = default; + + for (int i = 0; i < length; i++) + { + if (histogram[i] != 0) + { + if (count < 4) + { + s4[count] = i; + } + else if (count > 4) + { + break; + } + + count++; + } + } + + int maxBitsCounter = length - 1; + int maxBits = 0; + + while (maxBitsCounter != 0) + { + maxBitsCounter >>= 1; + ++maxBits; + } + + if (count <= 1) + { + writer.Write(4, 1); + writer.Write(maxBits, s4[0]); + return; + } + + JxlHuffmanTree.CreateHuffmanTree(histogram, length, 15, depth); + JxlHuffmanTree.ConvertBitDepthsToSymbols(depth, length, bits); + + if (count <= 4) + { + StoreSimpleHuffmanTree(depth, s4, count, maxBits, writer); + } + else + { + StoreHuffmanTree(depth, length, writer); + } + } +} diff --git a/src/ImageSharp/Formats/Jxl/Processing/Encoder/Huffman/JxlHuffmanTree.cs b/src/ImageSharp/Formats/Jxl/Processing/Encoder/Huffman/JxlHuffmanTree.cs new file mode 100644 index 0000000000..0e3faa65aa --- /dev/null +++ b/src/ImageSharp/Formats/Jxl/Processing/Encoder/Huffman/JxlHuffmanTree.cs @@ -0,0 +1,411 @@ +// Copyright (c) Six Labors. +// Licensed under the Six Labors Split License. + +using System.Buffers; +using System.Numerics.Tensors; +using System.Runtime.CompilerServices; +using SixLabors.ImageSharp.Common.Helpers; + +namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Encoder.Huffman; + +/// +/// Node of a Huffman tree. +/// +internal struct JxlHuffmanTree(int count, short left, short right) +{ + public int TotalCount = count; + + /// + /// Index of the left node of the tree. + /// + public short IndexLeft = left; + + /// + /// Index of the right node of the tree. If it's missing + /// then this is the value of the node. + /// + public short IndexRightOrValue = right; + + /// + /// Gets a lookup table with pre-reversed 4-bit values. + /// This lookup is used by . + /// + private static ReadOnlySpan ReverseLookup => + [ + 0x0, 0x8, 0x4, 0xc, 0x2, 0xa, 0x6, 0xe, + 0x1, 0x9, 0x5, 0xd, 0x3, 0xb, 0x7, 0xf + ]; + + public static void SetDepth(ref JxlHuffmanTree p, Span pool, Span depth, byte level) + { + if (p.IndexLeft >= 0) + { + level++; + SetDepth(ref pool[p.IndexLeft], pool, depth, level); + SetDepth(ref pool[p.IndexRightOrValue], pool, depth, level); + } + else + { + depth[p.IndexRightOrValue] = level; + } + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static int Compare(JxlHuffmanTree v0, JxlHuffmanTree v1) + { + if (v0.TotalCount != v1.TotalCount) + { + return v0.TotalCount.CompareTo(v1.TotalCount); + } + + return v0.IndexRightOrValue.CompareTo(v1.IndexRightOrValue); + } + + public static void CreateHuffmanTree(Span data, int length, int treeLimit, Span depth) + { + JxlHuffmanTree[]? pool = null; + + // This is basically the equivalent of List but is fixed-size. + // We don't need an entire collection on the heap. + int desiredTreeItems = (2 * length) + 1; + Span tree = + desiredTreeItems <= 256 + ? stackalloc JxlHuffmanTree[256].Slice(0, desiredTreeItems) + : pool = ArrayPool.Shared.Rent(desiredTreeItems); + + // Number of items in our "fixed List". So in other to + // add to the tree we do 'tree[treeRef++] = ...'. + int treeRef = 0; + + for (int countLimit = 1; ; countLimit *= 2) + { + tree.Clear(); // Always clear on every iteration + + int i = length; + for (; i != 0;) + { + --i; + if (data[i] != 0) + { + int count = Math.Max(data[i], countLimit - 1); + tree[treeRef++] = new JxlHuffmanTree(count, -1, (short)i); + } + } + + if (treeRef == 1) + { + // Fake value; will be fixed on upper level. + depth[tree[0].IndexRightOrValue] = 1; + break; + } + + tree.Sort(Compare); + + JxlHuffmanTree sentinel = new(int.MaxValue, -1, -1); + tree[treeRef++] = sentinel; + tree[treeRef++] = sentinel; // We do this twice, yes + + i = 0; + int j = treeRef + 1; + + for (int k = treeRef - 1; k != 0; --k) + { + int left; + int right; + + if (tree[i].TotalCount <= tree[j].TotalCount) + { + left = i; + i++; + } + else + { + left = j; + j++; + } + + if (tree[i].TotalCount <= tree[j].TotalCount) + { + right = i; + i++; + } + else + { + right = j; + j++; + } + + int j_end = treeRef - 1; + ref JxlHuffmanTree currTree = ref tree[j_end]; + currTree.TotalCount = tree[left].TotalCount + tree[right].TotalCount; + currTree.IndexLeft = (short)left; + currTree.IndexRightOrValue = (short)right; + + tree[treeRef++] = sentinel; + } + + SetDepth(ref tree[(2 * treeRef) - 1], tree, depth, 0); + + if (TensorPrimitives.Max((ReadOnlySpan)depth[..length]) <= treeLimit) + { + break; + } + } + + // Don't forget to return the pooled array + if (pool is not null) + { + ArrayPool.Shared.Return(pool); + } + } + + public static void Reverse(Span v, int start, int end) + { + end--; + while (start < end) + { + RuntimeUtility.Swap(ref v[end], ref v[start]); + start++; + end++; + } + } + + public static void WriteHuffmanTreeRepetitions(byte previousValue, byte value, int repetitions, ref int treeSize, Span tree, Span extraBitsData) + { + DebugGuard.MustBeGreaterThan(repetitions, 0, nameof(repetitions)); + + if (previousValue != value) + { + tree[treeSize] = value; + extraBitsData[treeSize] = 0; + treeSize++; + repetitions--; + } + + if (repetitions == 7) + { + tree[treeSize] = value; + extraBitsData[treeSize] = 0; + treeSize++; + } + + if (repetitions < 3) + { + for (int i = 0; i < repetitions; ++i) + { + tree[treeSize] = value; + extraBitsData[treeSize] = 0; + treeSize++; + } + } + else + { + repetitions -= 3; + int start = treeSize; + while (true) + { + tree[treeSize] = 16; + extraBitsData[treeSize] = (byte)(repetitions & 0x3); + treeSize++; + repetitions >>= 2; + + if (repetitions == 0) + { + break; + } + + repetitions--; + } + + Reverse(tree, start, treeSize); + Reverse(extraBitsData, start, treeSize); + } + } + + public static void WriteHuffmanTreeRepetitionsZeros(int repetitions, ref int treeSize, Span tree, Span extraBitsData) + { + if (repetitions == 11) + { + tree[treeSize] = 0; + extraBitsData[treeSize] = 0; + treeSize++; + repetitions--; + } + + if (repetitions < 3) + { + for (int i = 0; i < repetitions; ++i) + { + tree[treeSize] = 0; + extraBitsData[treeSize] = 0; + treeSize++; + } + } + else + { + repetitions -= 3; + int start = treeSize; + + while (true) + { + tree[treeSize] = 17; + extraBitsData[treeSize] = (byte)(repetitions & 0x7); + treeSize++; + repetitions >>= 3; + + if (repetitions == 0) + { + break; + } + + repetitions--; + } + + Reverse(tree, start, treeSize); + Reverse(extraBitsData, start, treeSize); + } + } + + // Decides whether or not to use Run Length Encoding (RLE). + // Basically that's where, for example, when we have a + // string of repetitive letters "aaaaaa", instead of encoding + // them all separately, it encodes "a times 6". + public static void DecideOverRleUse(Span depth, int length, ref bool useRleForNonZero, ref bool useRleForZero) + { + int totalRepsZero = 0; + int totalRepsNonZero = 0; + int countRepsZero = 1; + int countRepsNonZero = 1; + + for (int i = 0; i < length;) + { + byte value = depth[i]; + int reps = 1; + + for (int k = i + 1; k < length && depth[k] == value; k++) + { + reps++; + } + + if (reps >= 3 && value == 0) + { + totalRepsZero += reps; + countRepsZero++; + } + + if (reps >= 4 && value != 0) + { + totalRepsNonZero += reps; + countRepsNonZero++; + } + + i += reps; + } + + useRleForNonZero = totalRepsNonZero > countRepsNonZero * 2; + useRleForZero = totalRepsZero > countRepsZero * 2; + } + + public static void WriteHuffmanTree(Span depth, int length, ref int treeSize, Span tree, Span extraBitsData) + { + byte previousValue = 8; + int newLength = length; + + for (int i = 0; i < length; i++) + { + if (depth[length - i - 1] == 0) + { + newLength--; + } + else + { + break; + } + } + + bool useRleForNonZeroes = false; + bool useRleForZero = false; + + if (length > 50) + { + DecideOverRleUse(depth, newLength, ref useRleForNonZeroes, ref useRleForZero); + } + + for (int i = 0; i < newLength;) + { + byte value = depth[i]; + int reps = 1; + + if ((value != 0 && useRleForNonZeroes) || (value == 0 && useRleForZero)) + { + for (int k = i + 1; k < newLength && depth[k] == value; k++) + { + reps++; + } + } + + if (value == 0) + { + WriteHuffmanTreeRepetitionsZeros(reps, ref treeSize, tree, extraBitsData); + } + else + { + WriteHuffmanTreeRepetitions(previousValue, value, reps, ref treeSize, tree, extraBitsData); + previousValue = value; + } + + i += reps; + } + } + + public static short ReverseBits(int numBits, short bits) + { + int result = ReverseLookup[bits & 0xf]; + + for (int i = 4; i < numBits; i += 4) + { + result <<= 4; + bits = (short)(bits >> 4); + result |= ReverseLookup[bits & 0xf]; + } + + result >>= -numBits & 0x3; + + return (short)result; + } + + public static void ConvertBitDepthsToSymbols(Span depth, int len, Span bits) + { + // In Brotli, all bit depths are [1..15] + // 0 bit depth means that the symbol does not exist. + const int maxBits = 16; // 0..15 are values for bits + + Span blCount = stackalloc short[maxBits]; + blCount.Clear(); // explicitly cleared from reference + + for (int i = 0; i < len; i++) + { + blCount[depth[i]]++; + } + + blCount[0] = 0; + + Span nextCode = stackalloc short[maxBits]; // not cleared in reference + nextCode[0] = 0; + + int code = 0; + for (int i = 1; i < maxBits; ++i) + { + code = (code + blCount[i - 1]) << 1; + nextCode[i] = (short)code; + } + + for (int i = 0; i < len; ++i) + { + if (depth[i] != 0) + { + bits[i] = ReverseBits(depth[i], nextCode[depth[i]]++); + } + } + } +} diff --git a/src/ImageSharp/Formats/Jxl/Processing/Encoder/JxlGaborish.cs b/src/ImageSharp/Formats/Jxl/Processing/Encoder/JxlGaborish.cs new file mode 100644 index 0000000000..191e2c5a4a --- /dev/null +++ b/src/ImageSharp/Formats/Jxl/Processing/Encoder/JxlGaborish.cs @@ -0,0 +1,69 @@ +// Copyright (c) Six Labors. +// Licensed under the Six Labors Split License. + +using System.Runtime.CompilerServices; +using SixLabors.ImageSharp.Formats.Jxl.Memory.ImageTypes; +using SixLabors.ImageSharp.Formats.Jxl.Processing.Primitives; + +namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Encoder; + +/// +/// Gaborish transform +/// +internal static class JxlGaborish +{ + private static ReadOnlySpan GaborishLookup => [ + -0.09495815671340026f, -0.041031725066768575f, 0.013710004822696948f, + 0.006510206083837737f, -0.0014789063378272242f]; + + public static void InverseGaborish(Configuration configuration, JxlImage3F inOut, Rectangle rect, InlineArray3 mul) + { + InlineArray3 weights = default; + + for (int i = 0; i < 3; ++i) + { + double sum = 1.0 + (mul[i] * 4 * ((GaborishLookup[0] + GaborishLookup[1]) + (GaborishLookup[2] + GaborishLookup[4]) + (2 * GaborishLookup[3]))); + sum = Math.Max(sum, 1e-5); // if (sum < 1e-5) sum = 1e-5 + + float normalize = (float)(1.0f / sum); + float normalizeMul = mul[i] * normalize; + + weights[i] = new JxlWeightsSymmetric5() + { + C = JxlWeightsSymmetric5.CreateVector4(normalize), + R = JxlWeightsSymmetric5.CreateVector4(normalizeMul * GaborishLookup[0]), + R2 = JxlWeightsSymmetric5.CreateVector4(normalizeMul * GaborishLookup[2]), + D = JxlWeightsSymmetric5.CreateVector4(normalizeMul * GaborishLookup[1]), + D2 = JxlWeightsSymmetric5.CreateVector4(normalizeMul * GaborishLookup[4]), + L = JxlWeightsSymmetric5.CreateVector4(normalizeMul * GaborishLookup[3]) + }; + } + + using JxlImageF temp = new(configuration, inOut.Plane(2).XSize, inOut.Plane(2).YSize); + + if (!JxlImageOperations.CopyImage(inOut.Plane(2), temp)) + { + throw new InvalidOperationException("Image copying failed"); + } + + Rectangle xRect = RectangleUtils.Extend(rect, 3, inOut.GetRectangle()); + + if (!JxlConvolve.Symmetric5(inOut.Plane(0), xRect, ref weights[0], inOut.Plane(2), xRect)) + { + throw new InvalidOperationException("Symmetric5 convolution failed"); + } + + if (!JxlConvolve.Symmetric5(inOut.Plane(1), xRect, ref weights[1], inOut.Plane(0), xRect)) + { + throw new InvalidOperationException("Symmetric5 convolution failed"); + } + + if (!JxlConvolve.Symmetric5(temp, xRect, ref weights[2], inOut.Plane(1), xRect)) + { + throw new InvalidOperationException("Symmetric5 convolution failed"); + } + + inOut.Plane(0).Swap(inOut.Plane(1)); + inOut.Plane(0).Swap(inOut.Plane(2)); + } +} diff --git a/src/ImageSharp/Formats/Jxl/Processing/Encoder/JxlGammaCorrect.cs b/src/ImageSharp/Formats/Jxl/Processing/Encoder/JxlGammaCorrect.cs new file mode 100644 index 0000000000..1b60ed38c8 --- /dev/null +++ b/src/ImageSharp/Formats/Jxl/Processing/Encoder/JxlGammaCorrect.cs @@ -0,0 +1,51 @@ +// Copyright (c) Six Labors. +// Licensed under the Six Labors Split License. + +using System.Runtime.CompilerServices; + +namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Encoder; + +internal static class JxlGammaCorrect +{ + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static double SRgb8ToLinearDirect(double srgb) + { + if (srgb <= 0.0) + { + return 0.0; + } + + if (srgb <= 0.04045) + { + return srgb / 12.92; + } + + if (srgb >= 1.0) + { + return 1.0; + } + + return Math.Pow((srgb + 0.055) / 1.055, 2.4); + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + public static double LinearToSRgb8Direct(double linear) + { + if (linear <= 0.0) + { + return 0.0; + } + + if (linear >= 1.0) + { + return 1.0; + } + + if (linear <= 0.0031308) + { + return linear * 12.92; + } + + return (Math.Pow(linear, 1.0 / 2.4) * 1.055) - 0.055; + } +} diff --git a/src/ImageSharp/Formats/Jxl/Processing/Encoder/Noise/JxlLossFunction.cs b/src/ImageSharp/Formats/Jxl/Processing/Encoder/Noise/JxlLossFunction.cs new file mode 100644 index 0000000000..e338862bb0 --- /dev/null +++ b/src/ImageSharp/Formats/Jxl/Processing/Encoder/Noise/JxlLossFunction.cs @@ -0,0 +1,63 @@ +// Copyright (c) Six Labors. +// Licensed under the Six Labors Split License. + +using SixLabors.ImageSharp.Formats.Jxl.Processing.Noise; + +namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Encoder.Noise; + +internal sealed class JxlLossFunction(ReadOnlyMemory noiseLevels) +{ + public double Compute(Span w, Span df, bool skipRegularization = false) + { + const double reg = 0.005; + const double asym = 1.1; + + double lossFunction = 0; + + w.Clear(); + + ReadOnlySpan levels = noiseLevels.Span; + + for (int i = 0; i < levels.Length; i++) + { + JxlNoiseLevel nl = levels[i]; + + JxlNoiseIndexAndFraction pos = JxlNoiseHelper.IndexAndFraction(nl.Intensity); + + double low = w[pos.Index]; + double hi = w[pos.Index + 1]; + double val = (low * (1.0f - pos.Fraction)) + (hi * pos.Fraction); + double dist = val - nl.NoiseLevel; + + if (dist > 0) + { + lossFunction += asym * dist * dist; + df[pos.Index] -= asym * (1.0f - pos.Fraction) * dist; + df[pos.Index + 1] -= asym * pos.Fraction * dist; + } + else + { + lossFunction += dist * dist; + df[pos.Index] -= (1.0f - pos.Fraction) * dist; + df[pos.Index + 1] -= pos.Fraction * dist; + } + } + + if (skipRegularization) + { + return lossFunction; + } + + int levelsSize = levels.Length; + + for (int i = 0; i + 1 < w.Length; i++) + { + double diff = w[i] - w[i + 1]; + lossFunction += reg * levelsSize * diff * diff; + df[i] -= reg * diff * levelsSize; + df[i + 1] += reg * diff * levelsSize; + } + + return lossFunction; + } +} diff --git a/src/ImageSharp/Formats/Jxl/Processing/Encoder/Noise/JxlNoiseEncoder.cs b/src/ImageSharp/Formats/Jxl/Processing/Encoder/Noise/JxlNoiseEncoder.cs new file mode 100644 index 0000000000..026ee73f2f --- /dev/null +++ b/src/ImageSharp/Formats/Jxl/Processing/Encoder/Noise/JxlNoiseEncoder.cs @@ -0,0 +1,66 @@ +// Copyright (c) Six Labors. +// Licensed under the Six Labors Split License. + +using System.Buffers; +using System.Numerics.Tensors; +using SixLabors.ImageSharp.Formats.Jxl.Memory.ImageTypes; + +namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Encoder.Noise; + +/// +/// Noise functions for encoder. +/// +internal static class JxlNoiseEncoder +{ + public static float GetScoreSumsOfAbsoluteDifferences(JxlImage3F opsin, int x, int y, int blockSize) + { + const int smallBlockSizeX = 3; + const int smallBlockSizeY = 4; + + int numSAD = (blockSize - smallBlockSizeX) * (blockSize - smallBlockSizeY); + int counter = 0; + const int offset = 2; + + float[]? pooled = null; + + Span sad = numSAD <= 128 + ? stackalloc float[128].Slice(0, numSAD) + : pooled = ArrayPool.Shared.Rent(numSAD); + + for (int yBl = 0; yBl + smallBlockSizeY < blockSize; ++yBl) + { + for (int xBl = 0; xBl + smallBlockSizeX < blockSize; ++xBl) + { + float sadSum = 0; + + for (int cy = 0; cy < smallBlockSizeY; ++cy) + { + for (int cx = 0; cx < smallBlockSizeX; ++cx) + { + float wnd = 0.5f * (opsin.PlaneRow(1, y + yBl + cy)[x + xBl + cx] + opsin.PlaneRow(0, y + yBl + cy)[x + xBl + cx]); + float center = 0.5f * (opsin.PlaneRow(1, y + offset + cy)[x + offset + cx] + opsin.PlaneRow(0, y + offset + cy)[x + offset + cx]); + sadSum += MathF.Abs(center - wnd); + } + } + + sad[counter++] = sadSum; + } + } + + int samples = numSAD / 2; + + // As with ROAD (rank order absolute distance), we keep the smallest half of + // the values in SAD (we use here the more robust patch SAD instead of + // absolute single-pixel differences). + sad.Sort(); + + float totalSadSum = TensorPrimitives.Sum(sad); + + if (pooled is not null) + { + ArrayPool.Shared.Return(pooled); + } + + return totalSadSum / samples; + } +} diff --git a/src/ImageSharp/Formats/Jxl/Processing/Encoder/Noise/JxlNoiseHistogram.cs b/src/ImageSharp/Formats/Jxl/Processing/Encoder/Noise/JxlNoiseHistogram.cs new file mode 100644 index 0000000000..5cacaa5d02 --- /dev/null +++ b/src/ImageSharp/Formats/Jxl/Processing/Encoder/Noise/JxlNoiseHistogram.cs @@ -0,0 +1,87 @@ +// Copyright (c) Six Labors. +// Licensed under the Six Labors Split License. + +using System.Numerics; +using System.Numerics.Tensors; +using System.Runtime.CompilerServices; + +namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Encoder.Noise; + +internal sealed class JxlNoiseHistogram +{ + private const int Bins = 256; + + private readonly uint[] bins = new uint[Bins]; + + public int Mode + { + get + { + int maxIdx = 0; + + for (int i = 0; i < Bins; i++) + { + if (this.bins[i] > this.bins[maxIdx]) + { + maxIdx = i; + } + } + + return maxIdx; + } + } + + /// + /// Gets the Inter-quartile range. + /// + public double Iqr => this.Quantile(0.75) - this.Quantile(0.25); + + public void Increment(float x) => this.bins[Index(x)]++; + + public uint Get(float x) => this.bins[Index(x)]; + + public uint Bin(int bin) => this.bins[bin]; + + public double Quantile(double q01) + { + long total = 1 + TensorPrimitives.Sum((ReadOnlySpan)this.bins.AsSpan()); + long target = (long)q01 * total; + long sum = 0; + int i = 0; + + for (; i < Bins; i++) + { + sum += this.bins[i]; + + if (sum == target) + { + return i + 0.5; + } + + if (sum > target) + { + break; + } + } + + int next = i + 1; + + while (next < Bins && this.bins[next] == 0) + { + next++; + } + + double excess = target - sum; + double weightNext = this.bins[Index(next)] / excess; + + return ClampX((next * weightNext) + (i * (1.0 - weightNext))); + } + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private static T ClampX(T x) + where T : unmanaged, INumber + => T.Clamp(x, T.Zero, T.CreateSaturating(Bins - 1)); + + [MethodImpl(MethodImplOptions.AggressiveInlining)] + private static int Index(float x) => ClampX((int)x); +} diff --git a/src/ImageSharp/Formats/Jxl/Processing/JxlConvolve.cs b/src/ImageSharp/Formats/Jxl/Processing/JxlConvolve.cs index ebe7f5195d..2baad82d6d 100644 --- a/src/ImageSharp/Formats/Jxl/Processing/JxlConvolve.cs +++ b/src/ImageSharp/Formats/Jxl/Processing/JxlConvolve.cs @@ -5,6 +5,7 @@ using System.Numerics; using System.Runtime.CompilerServices; using System.Runtime.InteropServices; using System.Runtime.Intrinsics; +using SixLabors.ImageSharp.Formats.Jxl.Memory; using SixLabors.ImageSharp.Formats.Jxl.Memory.ImageTypes; using SixLabors.ImageSharp.Formats.Jxl.Processing.Primitives; @@ -70,7 +71,7 @@ internal static class JxlConvolve return sum2 + (sum1 + sum0); } - public static float Symmetric5Border(JxlImageF input, Func wrapY, long ix, long iy, JxlWeightsSymmetric5 weights) + public static float Symmetric5Border(JxlImageF input, Func wrapY, long ix, long iy, ref JxlWeightsSymmetric5 weights) { float w0 = weights.GetCVector()[0]; float w1 = weights.GetRVector()[0]; @@ -99,7 +100,7 @@ internal static class JxlConvolve Func wrapY, int rix, long iy, - JxlWeightsSymmetric5 weights, + ref JxlWeightsSymmetric5 weights, Span rowOut) { Vector w0 = LoadDuplicate128(weights.GetCVector()); // c @@ -126,7 +127,7 @@ internal static class JxlConvolve Func wrapY, in Rectangle rect, long iy, - JxlWeightsSymmetric5 weights, + ref JxlWeightsSymmetric5 weights, Span rowOut) { const int radius = 2; @@ -140,25 +141,25 @@ internal static class JxlConvolve for (; ix < Math.Min(alignedX, xEnd); ix++, rix++) { - rowOut[rix] = Symmetric5Border(image, wrapY, ix, iy, weights); + rowOut[rix] = Symmetric5Border(image, wrapY, ix, iy, ref weights); } for (; ix + n + radius <= xEnd; ix += n, rix += n) { - Symmetric5Interior(image, ix, wrapY, rix, iy, weights, rowOut); + Symmetric5Interior(image, ix, wrapY, rix, iy, ref weights, rowOut); } for (; ix < xEnd; ix++, rix++) { - rowOut[rix] = Symmetric5Border(image, wrapY, ix, iy, weights); + rowOut[rix] = Symmetric5Border(image, wrapY, ix, iy, ref weights); } } public static bool Symmetric5( - JxlImageF input, + JxlPlane input, in Rectangle rectangle, - JxlWeightsSymmetric5 weights, - JxlImageF output, + ref JxlWeightsSymmetric5 weights, + JxlPlane output, Rectangle outputRect) { if (rectangle.Width != outputRect.Width || rectangle.Height != outputRect.Height) diff --git a/src/ImageSharp/Formats/Jxl/Processing/Primitives/JxlWeightsSymmetric5.cs b/src/ImageSharp/Formats/Jxl/Processing/Primitives/JxlWeightsSymmetric5.cs index df5fa66b7b..7fb8977a50 100644 --- a/src/ImageSharp/Formats/Jxl/Processing/Primitives/JxlWeightsSymmetric5.cs +++ b/src/ImageSharp/Formats/Jxl/Processing/Primitives/JxlWeightsSymmetric5.cs @@ -6,89 +6,98 @@ using System.Runtime.Intrinsics; namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Primitives; -internal sealed class JxlWeightsSymmetric5 +internal struct JxlWeightsSymmetric5 { - private InlineArray4 c; + public InlineArray4 C; - private InlineArray4 r; + public InlineArray4 R; - private InlineArray4 r2; + public InlineArray4 R2; - private InlineArray4 d; + public InlineArray4 D; - private InlineArray4 d2; + public InlineArray4 D2; - private InlineArray4 l; + public InlineArray4 L; - public Vector128 GetCVector() + public static InlineArray4 CreateVector4(float x) { - ref float first = ref Unsafe.AsRef(in this.c[0]); + InlineArray4 array = default; + + array[0] = array[1] = array[2] = array[3] = x; + + return array; + } + + public readonly Vector128 GetCVector() + { + ref float first = ref Unsafe.AsRef(in this.C[0]); return Vector128.LoadUnsafe(ref first); } - public Vector128 GetRVector() + public readonly Vector128 GetRVector() { - ref float first = ref Unsafe.AsRef(in this.r[0]); + ref float first = ref Unsafe.AsRef(in this.R[0]); return Vector128.LoadUnsafe(ref first); } - public Vector128 GetR2Vector() + public readonly Vector128 GetR2Vector() { - ref float first = ref Unsafe.AsRef(in this.r2[0]); + ref float first = ref Unsafe.AsRef(in this.R2[0]); return Vector128.LoadUnsafe(ref first); } - public Vector128 GetDVector() + public readonly Vector128 GetDVector() { - ref float first = ref Unsafe.AsRef(in this.d[0]); + ref float first = ref Unsafe.AsRef(in this.D[0]); return Vector128.LoadUnsafe(ref first); } - public Vector128 GetD2Vector() + public readonly Vector128 GetD2Vector() { - ref float first = ref Unsafe.AsRef(in this.d2[0]); + ref float first = ref Unsafe.AsRef(in this.D2[0]); return Vector128.LoadUnsafe(ref first); } - public Vector128 GetLVector() + public readonly Vector128 GetLVector() { - ref float first = ref Unsafe.AsRef(in this.l[0]); + ref float first = ref Unsafe.AsRef(in this.L[0]); return Vector128.LoadUnsafe(ref first); } - public void SetC(Vector128 vec) + public readonly void SetC(Vector128 vec) { - ref float first = ref Unsafe.AsRef(in this.c[0]); + ref float first = ref Unsafe.AsRef(in this.C[0]); vec.StoreUnsafe(ref first); } - public void SetD(Vector128 vec) + public readonly void SetD(Vector128 vec) { - ref float first = ref Unsafe.AsRef(in this.d[0]); + ref float first = ref Unsafe.AsRef(in this.D[0]); vec.StoreUnsafe(ref first); } - public void SetD2(Vector128 vec) + public readonly void SetD2(Vector128 vec) { - ref float first = ref Unsafe.AsRef(in this.d2[0]); + ref float first = ref Unsafe.AsRef(in this.D2[0]); vec.StoreUnsafe(ref first); } - public void SetR(Vector128 vec) + public readonly void SetR(Vector128 vec) { - ref float first = ref Unsafe.AsRef(in this.r[0]); + ref float first = ref Unsafe.AsRef(in this.R[0]); vec.StoreUnsafe(ref first); } - public void SetR2(Vector128 vec) + public readonly void SetR2(Vector128 vec) { - ref float first = ref Unsafe.AsRef(in this.r2[0]); + ref float first = ref Unsafe.AsRef(in this.R2[0]); vec.StoreUnsafe(ref first); } - public void SetL(Vector128 vec) + public readonly void SetL(Vector128 vec) { - ref float first = ref Unsafe.AsRef(in this.l[0]); + ref float first = ref Unsafe.AsRef(in this.L[0]); vec.StoreUnsafe(ref first); } } diff --git a/src/ImageSharp/Formats/Jxl/Processing/Primitives/RectangleUtils.cs b/src/ImageSharp/Formats/Jxl/Processing/Primitives/RectangleUtils.cs new file mode 100644 index 0000000000..072b51983f --- /dev/null +++ b/src/ImageSharp/Formats/Jxl/Processing/Primitives/RectangleUtils.cs @@ -0,0 +1,25 @@ +// Copyright (c) Six Labors. +// Licensed under the Six Labors Split License. + +namespace SixLabors.ImageSharp.Formats.Jxl.Processing.Primitives; + +internal static class RectangleUtils +{ + public static int X0(in Rectangle rect) => rect.X; + + public static int Y0(in Rectangle rect) => rect.Y; + + public static int X1(in Rectangle rect) => rect.X + rect.Width; + + public static int Y1(in Rectangle rect) => rect.Y + rect.Height; + + public static Rectangle Extend(Rectangle curr, int border, Rectangle parent) + { + int newX0 = X0(in curr) > X0(in parent) + border ? X0(in curr) - border : X0(in parent); + int newY0 = Y0(in curr) > Y0(in parent) + border ? Y0(in curr) - border : Y0(in parent); + int newX1 = X1(in curr) + border > X1(in parent) ? X1(in parent) : X1(in curr) + border; + int newY1 = Y1(in curr) + border > Y1(in parent) ? Y1(in parent) : Y1(in curr) + border; + + return new(newX0, newY0, newX1 - newX0, newY1 - newY0); + } +} diff --git a/src/ImageSharp/Formats/Jxl/Processing/RenderPipeline/Epf0Stage.cs b/src/ImageSharp/Formats/Jxl/Processing/RenderPipeline/Epf0Stage.cs index 5df7525ea9..28f3dcb34d 100644 --- a/src/ImageSharp/Formats/Jxl/Processing/RenderPipeline/Epf0Stage.cs +++ b/src/ImageSharp/Formats/Jxl/Processing/RenderPipeline/Epf0Stage.cs @@ -19,6 +19,11 @@ internal sealed class Epf0Stage : RenderPipelineStageBase [0, 1], [0, 2], [1, -1], [1, 0], [1, 1], [2, 0] ]; + private static readonly int[][] PlusOffsets = + [ + [0, 0], [-1, 0], [0, -1], [1, 0], [0, 1] + ]; + private readonly JxlLoopFilter loopFilter; private readonly JxlImageF sigma; @@ -29,12 +34,13 @@ internal sealed class Epf0Stage : RenderPipelineStageBase this.Settings = RenderPipelineStageConfiguration.CreateSymmetricBorderOnly(3); } + /// public override string Name => "EPF0"; [MethodImpl(MethodImplOptions.AggressiveInlining)] public static void AddPixel( int row, - InlineArray7>> rows, + InlineArray3>> rows, int x, Vector256 sad, Vector256 inverseSigma, @@ -54,6 +60,7 @@ internal sealed class Epf0Stage : RenderPipelineStageBase bOut += (weight * cb) + bOut; } + /// public override void ProcessRow(Buffer2D> inputRows, Buffer2D> outputRows, int xExtraLeft, int xExtraRight, int width, int xPos, int yPos) { Span> sads = stackalloc Vector256[16].Slice(0, 12); @@ -103,6 +110,60 @@ internal sealed class Epf0Stage : RenderPipelineStageBase Vector256 vsm = Vector256.Create(sadMul[ix..]); Vector256 inverseSigma = Vector256.Create(rowSigma[bx]) * vsm; + + sads.Clear(); + + for (int c = 0; c < 3; c++) + { + Vector256 scale = Vector256.Create(this.loopFilter.EpfChannelScale[c]); + + for (int i = 0; i < 12; i++) + { + Vector256 sad = Vector256.Zero; + + foreach (Span offset in PlusOffsets) + { + Vector256 r11 = Vector256.Create((ReadOnlySpan)rows[c][3 + offset[0]][(x + offset[1])..].Span); + Vector256 c11 = Vector256.Create((ReadOnlySpan)rows[c][3 + SadOffsets[i][0] + offset[0]][(x + SadOffsets[i][1] + offset[1])..].Span); + sad += Vector256.Abs(r11 - c11); + } + + sads[i] = (sad * scale) + sads[i]; + } + } + + Vector256 xCC = Vector256.Create((ReadOnlySpan)rows[0][3 + 0][x..].Span); + Vector256 yCC = Vector256.Create((ReadOnlySpan)rows[1][3 + 0][x..].Span); + Vector256 bCC = Vector256.Create((ReadOnlySpan)rows[2][3 + 0][x..].Span); + + Vector256 w = Vector256.One; + Vector256 X = xCC; + Vector256 Y = yCC; + Vector256 B = bCC; + + for (int i = 0; i < 12; i++) + { + AddPixel(SadOffsets[i][0], rows, x + SadOffsets[i][1], sads[i], inverseSigma, ref X, ref Y, ref B, ref w); + } + + Vector256 inverseW = Vector256.One / w; + + (X * inverseW).CopyTo(GetOutputRow(outputRows, 0, 0)[x..]); + (Y * inverseW).CopyTo(GetOutputRow(outputRows, 1, 0)[x..]); + (B * inverseW).CopyTo(GetOutputRow(outputRows, 2, 0)[x..]); + } + } + + /// + public override RenderPipelineChannelMode GetChannelMode(int channel) + { + if (channel < 3) + { + return RenderPipelineChannelMode.InOut; + } + else + { + return RenderPipelineChannelMode.Ignored; } } }