< Summary

Information
Class: JeremyAnsel.ColorQuant.WuAlphaColorQuantizer
Assembly: JeremyAnsel.ColorQuant
File(s): File 1: https://raw.githubusercontent.com/JeremyAnsel/JeremyAnsel.ColorQuant/6d79217e72af9e3af1a8a29c606732adac1e8d87/JeremyAnsel.ColorQuant/JeremyAnsel.ColorQuant/WuAlphaColorQuantizer.cs
File 2: https://raw.githubusercontent.com/JeremyAnsel/JeremyAnsel.ColorQuant/6d79217e72af9e3af1a8a29c606732adac1e8d87/JeremyAnsel.ColorQuant/JeremyAnsel.ColorQuant/WuAlphaColorQuantizer2.cs
Line coverage
97%
Covered lines: 695
Uncovered lines: 16
Coverable lines: 711
Total lines: 1599
Line coverage: 97.7%
Branch coverage
93%
Covered branches: 176
Total branches: 188
Branch coverage: 93.6%
Method coverage

Feature is only available for sponsors

Upgrade to PRO version

Metrics

MethodBranch coverage Crap Score Cyclomatic complexity Line coverage
File 1: .ctor()100%11100%
File 1: Quantize(...)100%11100%
File 1: Quantize(...)100%66100%
File 1: GetIndex(...)100%11100%
File 1: Volume(...)100%11100%
File 1: Bottom(...)80%5596.96%
File 1: Top(...)80%5596.96%
File 1: Clear()100%11100%
File 1: Build3DHistogram(...)100%22100%
File 1: Get3DMoments()100%88100%
File 1: Variance(...)100%11100%
File 1: Maximize(...)100%88100%
File 1: Cut(...)100%2525100%
File 1: Mark(...)100%88100%
File 1: BuildCube(...)87.5%161696%
File 1: GenerateResult(...)83.33%6680.95%
File 2: .ctor()100%11100%
File 2: Quantize(...)100%11100%
File 2: Quantize(...)100%66100%
File 2: GetIndex(...)100%11100%
File 2: Volume(...)100%11100%
File 2: Volume(...)100%11100%
File 2: Bottom(...)80%5596.96%
File 2: Bottom(...)80%5596.96%
File 2: Top(...)80%5596.96%
File 2: Top(...)80%5596.96%
File 2: Clear()100%11100%
File 2: Build3DHistogram(...)100%22100%
File 2: Get3DMoments()100%88100%
File 2: Variance(...)100%11100%
File 2: Maximize(...)100%88100%
File 2: Cut(...)100%2525100%
File 2: Mark(...)100%88100%
File 2: BuildCube(...)87.5%161696%
File 2: GenerateResult(...)83.33%6681.81%

File(s)

https://raw.githubusercontent.com/JeremyAnsel/JeremyAnsel.ColorQuant/6d79217e72af9e3af1a8a29c606732adac1e8d87/JeremyAnsel.ColorQuant/JeremyAnsel.ColorQuant/WuAlphaColorQuantizer.cs

#LineLine coverage
 1// <copyright file="WuAlphaColorQuantizer.cs" company="Jérémy Ansel">
 2// Copyright (c) 2014-2019 Jérémy Ansel
 3// </copyright>
 4// <license>
 5// Licensed under the MIT license. See LICENSE.txt
 6// </license>
 7
 8#if !NET8_0_OR_GREATER
 9namespace JeremyAnsel.ColorQuant
 10{
 11    using System;
 12    using System.Diagnostics.CodeAnalysis;
 13
 14    /// <summary>
 15    /// A Wu's color quantizer with alpha channel.
 16    /// </summary>
 17    /// <remarks>
 18    /// <para>
 19    /// Based on C Implementation of Xiaolin Wu's Color Quantizer (v. 2)
 20    /// (see Graphics Gems volume II, pages 126-133)
 21    /// (<see href="http://www.ece.mcmaster.ca/~xwu/cq.c"/>).
 22    /// </para>
 23    /// <para>
 24    /// Algorithm: Greedy orthogonal bipartition of RGB space for variance
 25    /// minimization aided by inclusion-exclusion tricks.
 26    /// For speed no nearest neighbor search is done. Slightly
 27    /// better performance can be expected by more sophisticated
 28    /// but more expensive versions.
 29    /// </para>
 30    /// </remarks>
 31    [SuppressMessage("Microsoft.Naming", "CA1709:IdentifiersShouldBeCasedCorrectly", MessageId = "Wu", Justification = "
 32    public sealed class WuAlphaColorQuantizer : IColorQuantizer
 33    {
 34        /// <summary>
 35        /// The maximum color count for the quantize.
 36        /// </summary>
 37        public const int MaxColors = 256;
 38
 39        /// <summary>
 40        /// The index bits.
 41        /// </summary>
 42        private const int IndexBits = 6;
 43
 44        /// <summary>
 45        /// The index alpha bits.
 46        /// </summary>
 47        private const int IndexAlphaBits = 4;
 48
 49        /// <summary>
 50        /// The index count.
 51        /// </summary>
 52        private const int IndexCount = (1 << IndexBits) + 1;
 53
 54        /// <summary>
 55        /// The index alpha count.
 56        /// </summary>
 57        private const int IndexAlphaCount = (1 << IndexAlphaBits) + 1;
 58
 59        /// <summary>
 60        /// The table length.
 61        /// </summary>
 62        private const int TableLength = IndexCount * IndexCount * IndexCount * IndexAlphaCount;
 63
 64        /// <summary>
 65        /// Moment of <c>P(c)</c>.
 66        /// </summary>
 167        private readonly long[] vwt = new long[TableLength];
 68
 69        /// <summary>
 70        /// Moment of <c>r*P(c)</c>.
 71        /// </summary>
 172        private readonly long[] vmr = new long[TableLength];
 73
 74        /// <summary>
 75        /// Moment of <c>g*P(c)</c>.
 76        /// </summary>
 177        private readonly long[] vmg = new long[TableLength];
 78
 79        /// <summary>
 80        /// Moment of <c>b*P(c)</c>.
 81        /// </summary>
 182        private readonly long[] vmb = new long[TableLength];
 83
 84        /// <summary>
 85        /// Moment of <c>a*P(c)</c>.
 86        /// </summary>
 187        private readonly long[] vma = new long[TableLength];
 88
 89        /// <summary>
 90        /// Moment of <c>c^2*P(c)</c>.
 91        /// </summary>
 192        private readonly double[] m2 = new double[TableLength];
 93
 94        /// <summary>
 95        /// Color space tag.
 96        /// </summary>
 197        private readonly byte[] tag = new byte[TableLength];
 98
 99        /// <summary>
 100        /// Quantizes an image.
 101        /// </summary>
 102        /// <param name="image">The image (ARGB).</param>
 103        /// <returns>The result.</returns>
 104        public ColorQuantizerResult Quantize(byte[] image)
 105        {
 14106            return this.Quantize(image, 256);
 107        }
 108
 109        /// <summary>
 110        /// Quantizes an image.
 111        /// </summary>
 112        /// <param name="image">The image (ARGB).</param>
 113        /// <param name="colorCount">The color count.</param>
 114        /// <returns>The result.</returns>
 115        public ColorQuantizerResult Quantize(byte[] image, int colorCount)
 116        {
 17117            if (image == null)
 118            {
 2119                throw new ArgumentNullException(nameof(image));
 120            }
 121
 15122            if (colorCount < 1 || colorCount > 256)
 123            {
 2124                throw new ArgumentOutOfRangeException(nameof(colorCount));
 125            }
 126
 13127            this.Clear();
 128
 13129            this.Build3DHistogram(image);
 13130            this.Get3DMoments();
 131
 132            AlphaBox[] cube;
 13133            this.BuildCube(out cube, ref colorCount);
 134
 13135            return this.GenerateResult(image, colorCount, cube);
 136        }
 137
 138        /// <summary>
 139        /// Gets an index.
 140        /// </summary>
 141        /// <param name="r">The red value.</param>
 142        /// <param name="g">The green value.</param>
 143        /// <param name="b">The blue value.</param>
 144        /// <param name="a">The alpha value.</param>
 145        /// <returns>The index.</returns>
 146        private static int GetIndex(int r, int g, int b, int a)
 147        {
 333090916148            return (r << ((IndexBits * 2) + IndexAlphaBits))
 333090916149                + (r << (IndexBits + IndexAlphaBits + 1))
 333090916150                + (g << (IndexBits + IndexAlphaBits))
 333090916151                + (r << (IndexBits * 2))
 333090916152                + (r << (IndexBits + 1))
 333090916153                + (g << IndexBits)
 333090916154                + ((r + g + b) << IndexAlphaBits)
 333090916155                + r + g + b + a;
 156        }
 157
 158        /// <summary>
 159        /// Computes sum over a AlphaBox of any given statistic.
 160        /// </summary>
 161        /// <param name="cube">The cube.</param>
 162        /// <param name="moment">The moment.</param>
 163        /// <returns>The result.</returns>
 164        private static double Volume(AlphaBox cube, long[] moment)
 165        {
 18570166            return moment[cube.Index1111]
 18570167                - moment[cube.Index1110]
 18570168                - moment[cube.Index1101]
 18570169                + moment[cube.Index1100]
 18570170                - moment[cube.Index1011]
 18570171                + moment[cube.Index1010]
 18570172                + moment[cube.Index1001]
 18570173                - moment[cube.Index1000]
 18570174                - moment[cube.Index0111]
 18570175                + moment[cube.Index0110]
 18570176                + moment[cube.Index0101]
 18570177                - moment[cube.Index0100]
 18570178                + moment[cube.Index0011]
 18570179                - moment[cube.Index0010]
 18570180                - moment[cube.Index0001]
 18570181                + moment[cube.Index0000];
 182        }
 183
 184        /// <summary>
 185        /// Computes part of Volume(cube, moment) that doesn't depend on r1, g1, or b1 (depending on direction).
 186        /// </summary>
 187        /// <param name="cube">The cube.</param>
 188        /// <param name="direction">The direction.</param>
 189        /// <param name="moment">The moment.</param>
 190        /// <returns>The result.</returns>
 191        private static long Bottom(AlphaBox cube, int direction, long[] moment)
 192        {
 193            switch (direction)
 194            {
 195                // Red
 196                case 3:
 6675197                    return -moment[cube.Index0111]
 6675198                        + moment[cube.Index0110]
 6675199                        + moment[cube.Index0101]
 6675200                        - moment[cube.Index0100]
 6675201                        + moment[cube.Index0011]
 6675202                        - moment[cube.Index0010]
 6675203                        - moment[cube.Index0001]
 6675204                        + moment[cube.Index0000];
 205
 206                // Green
 207                case 2:
 6675208                    return -moment[cube.Index1011]
 6675209                        + moment[cube.Index1010]
 6675210                        + moment[cube.Index1001]
 6675211                        - moment[cube.Index1000]
 6675212                        + moment[cube.Index0011]
 6675213                        - moment[cube.Index0010]
 6675214                        - moment[cube.Index0001]
 6675215                        + moment[cube.Index0000];
 216
 217                // Blue
 218                case 1:
 6675219                    return -moment[cube.Index1101]
 6675220                        + moment[cube.Index1100]
 6675221                        + moment[cube.Index1001]
 6675222                        - moment[cube.Index1000]
 6675223                        + moment[cube.Index0101]
 6675224                        - moment[cube.Index0100]
 6675225                        - moment[cube.Index0001]
 6675226                        + moment[cube.Index0000];
 227
 228                // Alpha
 229                case 0:
 6675230                    return -moment[cube.Index1110]
 6675231                        + moment[cube.Index1100]
 6675232                        + moment[cube.Index1010]
 6675233                        - moment[cube.Index1000]
 6675234                        + moment[cube.Index0110]
 6675235                        - moment[cube.Index0100]
 6675236                        - moment[cube.Index0010]
 6675237                        + moment[cube.Index0000];
 238
 239                default:
 0240                    throw new ArgumentOutOfRangeException(nameof(direction));
 241            }
 242        }
 243
 244        /// <summary>
 245        /// Computes remainder of Volume(cube, moment), substituting position for r1, g1, or b1 (depending on direction)
 246        /// </summary>
 247        /// <param name="cube">The cube.</param>
 248        /// <param name="direction">The direction.</param>
 249        /// <param name="position">The position.</param>
 250        /// <param name="moment">The moment.</param>
 251        /// <returns>The result.</returns>
 252        private static long Top(AlphaBox cube, int direction, int position, long[] moment)
 253        {
 254            switch (direction)
 255            {
 256                // Red
 257                case 3:
 49581258                    return moment[GetIndex(position, cube.G1, cube.B1, cube.A1)]
 49581259                        - moment[GetIndex(position, cube.G1, cube.B1, cube.A0)]
 49581260                        - moment[GetIndex(position, cube.G1, cube.B0, cube.A1)]
 49581261                        + moment[GetIndex(position, cube.G1, cube.B0, cube.A0)]
 49581262                        - moment[GetIndex(position, cube.G0, cube.B1, cube.A1)]
 49581263                        + moment[GetIndex(position, cube.G0, cube.B1, cube.A0)]
 49581264                        + moment[GetIndex(position, cube.G0, cube.B0, cube.A1)]
 49581265                        - moment[GetIndex(position, cube.G0, cube.B0, cube.A0)];
 266
 267                // Green
 268                case 2:
 68541269                    return moment[GetIndex(cube.R1, position, cube.B1, cube.A1)]
 68541270                        - moment[GetIndex(cube.R1, position, cube.B1, cube.A0)]
 68541271                        - moment[GetIndex(cube.R1, position, cube.B0, cube.A1)]
 68541272                        + moment[GetIndex(cube.R1, position, cube.B0, cube.A0)]
 68541273                        - moment[GetIndex(cube.R0, position, cube.B1, cube.A1)]
 68541274                        + moment[GetIndex(cube.R0, position, cube.B1, cube.A0)]
 68541275                        + moment[GetIndex(cube.R0, position, cube.B0, cube.A1)]
 68541276                        - moment[GetIndex(cube.R0, position, cube.B0, cube.A0)];
 277
 278                // Blue
 279                case 1:
 75741280                    return moment[GetIndex(cube.R1, cube.G1, position, cube.A1)]
 75741281                        - moment[GetIndex(cube.R1, cube.G1, position, cube.A0)]
 75741282                        - moment[GetIndex(cube.R1, cube.G0, position, cube.A1)]
 75741283                        + moment[GetIndex(cube.R1, cube.G0, position, cube.A0)]
 75741284                        - moment[GetIndex(cube.R0, cube.G1, position, cube.A1)]
 75741285                        + moment[GetIndex(cube.R0, cube.G1, position, cube.A0)]
 75741286                        + moment[GetIndex(cube.R0, cube.G0, position, cube.A1)]
 75741287                        - moment[GetIndex(cube.R0, cube.G0, position, cube.A0)];
 288
 289                // Alpha
 290                case 0:
 23365291                    return moment[GetIndex(cube.R1, cube.G1, cube.B1, position)]
 23365292                        - moment[GetIndex(cube.R1, cube.G1, cube.B0, position)]
 23365293                        - moment[GetIndex(cube.R1, cube.G0, cube.B1, position)]
 23365294                        + moment[GetIndex(cube.R1, cube.G0, cube.B0, position)]
 23365295                        - moment[GetIndex(cube.R0, cube.G1, cube.B1, position)]
 23365296                        + moment[GetIndex(cube.R0, cube.G1, cube.B0, position)]
 23365297                        + moment[GetIndex(cube.R0, cube.G0, cube.B1, position)]
 23365298                        - moment[GetIndex(cube.R0, cube.G0, cube.B0, position)];
 299
 300                default:
 0301                    throw new ArgumentOutOfRangeException(nameof(direction));
 302            }
 303        }
 304
 305        /// <summary>
 306        /// Clears the tables.
 307        /// </summary>
 308        private void Clear()
 309        {
 13310            this.vwt.AsSpan().Clear();
 13311            this.vmr.AsSpan().Clear();
 13312            this.vmg.AsSpan().Clear();
 13313            this.vmb.AsSpan().Clear();
 13314            this.vma.AsSpan().Clear();
 13315            this.m2.AsSpan().Clear();
 13316            this.tag.AsSpan().Clear();
 13317        }
 318
 319        /// <summary>
 320        /// Builds a 3-D color histogram of <c>counts, r/g/b, c^2</c>.
 321        /// </summary>
 322        /// <param name="image">The image.</param>
 323        private void Build3DHistogram(byte[] image)
 324        {
 167775262325            for (int i = 0; i < image.Length; i += 4)
 326            {
 83887618327                int a = image[i + 3];
 83887618328                int r = image[i + 2];
 83887618329                int g = image[i + 1];
 83887618330                int b = image[i];
 331
 83887618332                int inr = r >> (8 - IndexBits);
 83887618333                int ing = g >> (8 - IndexBits);
 83887618334                int inb = b >> (8 - IndexBits);
 83887618335                int ina = a >> (8 - IndexAlphaBits);
 336
 83887618337                int ind = GetIndex(inr + 1, ing + 1, inb + 1, ina + 1);
 338
 83887618339                this.vwt[ind]++;
 83887618340                this.vmr[ind] += r;
 83887618341                this.vmg[ind] += g;
 83887618342                this.vmb[ind] += b;
 83887618343                this.vma[ind] += a;
 83887618344                this.m2[ind] += (r * r) + (g * g) + (b * b) + (a * a);
 345            }
 13346        }
 347
 348        /// <summary>
 349        /// Converts the histogram into moments so that we can rapidly calculate
 350        /// the sums of the above quantities over any desired AlphaBox.
 351        /// </summary>
 352        private void Get3DMoments()
 353        {
 13354            Span<long> volume = stackalloc long[IndexCount * IndexAlphaCount];
 13355            Span<long> volumeR = stackalloc long[IndexCount * IndexAlphaCount];
 13356            Span<long> volumeG = stackalloc long[IndexCount * IndexAlphaCount];
 13357            Span<long> volumeB = stackalloc long[IndexCount * IndexAlphaCount];
 13358            Span<long> volumeA = stackalloc long[IndexCount * IndexAlphaCount];
 13359            Span<double> volume2 = stackalloc double[IndexCount * IndexAlphaCount];
 360
 13361            Span<long> area = stackalloc long[IndexAlphaCount];
 13362            Span<long> areaR = stackalloc long[IndexAlphaCount];
 13363            Span<long> areaG = stackalloc long[IndexAlphaCount];
 13364            Span<long> areaB = stackalloc long[IndexAlphaCount];
 13365            Span<long> areaA = stackalloc long[IndexAlphaCount];
 13366            Span<double> area2 = stackalloc double[IndexAlphaCount];
 367
 1690368            for (int r = 1; r < IndexCount; r++)
 369            {
 832370                volume.Clear();
 832371                volumeR.Clear();
 832372                volumeG.Clear();
 832373                volumeB.Clear();
 832374                volumeA.Clear();
 832375                volume2.Clear();
 376
 108160377                for (int g = 1; g < IndexCount; g++)
 378                {
 53248379                    area.Clear();
 53248380                    areaR.Clear();
 53248381                    areaG.Clear();
 53248382                    areaB.Clear();
 53248383                    areaA.Clear();
 53248384                    area2.Clear();
 385
 6922240386                    for (int b = 1; b < IndexCount; b++)
 387                    {
 3407872388                        long line = 0;
 3407872389                        long lineR = 0;
 3407872390                        long lineG = 0;
 3407872391                        long lineB = 0;
 3407872392                        long lineA = 0;
 3407872393                        double line2 = 0;
 394
 115867648395                        for (int a = 1; a < IndexAlphaCount; a++)
 396                        {
 54525952397                            int ind1 = GetIndex(r, g, b, a);
 398
 54525952399                            line += this.vwt[ind1];
 54525952400                            lineR += this.vmr[ind1];
 54525952401                            lineG += this.vmg[ind1];
 54525952402                            lineB += this.vmb[ind1];
 54525952403                            lineA += this.vma[ind1];
 54525952404                            line2 += this.m2[ind1];
 405
 54525952406                            area[a] += line;
 54525952407                            areaR[a] += lineR;
 54525952408                            areaG[a] += lineG;
 54525952409                            areaB[a] += lineB;
 54525952410                            areaA[a] += lineA;
 54525952411                            area2[a] += line2;
 412
 54525952413                            int inv = (b * IndexAlphaCount) + a;
 414
 54525952415                            volume[inv] += area[a];
 54525952416                            volumeR[inv] += areaR[a];
 54525952417                            volumeG[inv] += areaG[a];
 54525952418                            volumeB[inv] += areaB[a];
 54525952419                            volumeA[inv] += areaA[a];
 54525952420                            volume2[inv] += area2[a];
 421
 54525952422                            int ind2 = ind1 - GetIndex(1, 0, 0, 0);
 423
 54525952424                            this.vwt[ind1] = this.vwt[ind2] + volume[inv];
 54525952425                            this.vmr[ind1] = this.vmr[ind2] + volumeR[inv];
 54525952426                            this.vmg[ind1] = this.vmg[ind2] + volumeG[inv];
 54525952427                            this.vmb[ind1] = this.vmb[ind2] + volumeB[inv];
 54525952428                            this.vma[ind1] = this.vma[ind2] + volumeA[inv];
 54525952429                            this.m2[ind1] = this.m2[ind2] + volume2[inv];
 430                        }
 431                    }
 432                }
 433            }
 13434        }
 435
 436        /// <summary>
 437        /// Computes the weighted variance of a AlphaBox.
 438        /// </summary>
 439        /// <param name="cube">The cube.</param>
 440        /// <returns>The result.</returns>
 441        private double Variance(AlphaBox cube)
 442        {
 1577443            double dr = Volume(cube, this.vmr);
 1577444            double dg = Volume(cube, this.vmg);
 1577445            double db = Volume(cube, this.vmb);
 1577446            double da = Volume(cube, this.vma);
 447
 1577448            double xx =
 1577449                this.m2[cube.Index1111]
 1577450                - this.m2[cube.Index1110]
 1577451                - this.m2[cube.Index1101]
 1577452                + this.m2[cube.Index1100]
 1577453                - this.m2[cube.Index1011]
 1577454                + this.m2[cube.Index1010]
 1577455                + this.m2[cube.Index1001]
 1577456                - this.m2[cube.Index1000]
 1577457                - this.m2[cube.Index0111]
 1577458                + this.m2[cube.Index0110]
 1577459                + this.m2[cube.Index0101]
 1577460                - this.m2[cube.Index0100]
 1577461                + this.m2[cube.Index0011]
 1577462                - this.m2[cube.Index0010]
 1577463                - this.m2[cube.Index0001]
 1577464                + this.m2[cube.Index0000];
 465
 1577466            return xx - (((dr * dr) + (dg * dg) + (db * db) + (da * da)) / Volume(cube, this.vwt));
 467        }
 468
 469        /// <summary>
 470        /// We want to minimize the sum of the variances of two sub-AlphaBoxes.
 471        /// The sum(c^2) terms can be ignored since their sum over both sub-AlphaBoxes
 472        /// is the same (the sum for the whole AlphaBox) no matter where we split.
 473        /// The remaining terms have a minus sign in the variance formula,
 474        /// so we drop the minus sign and maximize the sum of the two terms.
 475        /// </summary>
 476        /// <param name="cube">The cube.</param>
 477        /// <param name="direction">The direction.</param>
 478        /// <param name="first">The first position.</param>
 479        /// <param name="last">The last position.</param>
 480        /// <param name="cut">The cutting point.</param>
 481        /// <param name="wholeR">The whole red.</param>
 482        /// <param name="wholeG">The whole green.</param>
 483        /// <param name="wholeB">The whole blue.</param>
 484        /// <param name="wholeA">The whole alpha.</param>
 485        /// <param name="wholeW">The whole weight.</param>
 486        /// <returns>The result.</returns>
 487        private double Maximize(AlphaBox cube, int direction, int first, int last, out int cut, double wholeR, double wh
 488        {
 5340489            long baseR = Bottom(cube, direction, this.vmr);
 5340490            long baseG = Bottom(cube, direction, this.vmg);
 5340491            long baseB = Bottom(cube, direction, this.vmb);
 5340492            long baseA = Bottom(cube, direction, this.vma);
 5340493            long baseW = Bottom(cube, direction, this.vwt);
 494
 5340495            double max = 0.0;
 5340496            cut = -1;
 497
 361312498            for (int i = first; i < last; i++)
 499            {
 175316500                double halfW1 = baseW + Top(cube, direction, i, this.vwt);
 501
 175316502                if (halfW1 == 0)
 503                {
 504                    continue;
 505                }
 506
 141073507                double halfW2 = wholeW - halfW1;
 508
 141073509                if (halfW2 == 0)
 510                {
 511                    continue;
 512                }
 513
 10478514                double halfR1 = baseR + Top(cube, direction, i, this.vmr);
 10478515                double halfG1 = baseG + Top(cube, direction, i, this.vmg);
 10478516                double halfB1 = baseB + Top(cube, direction, i, this.vmb);
 10478517                double halfA1 = baseA + Top(cube, direction, i, this.vma);
 518
 10478519                double temp = ((halfR1 * halfR1) + (halfG1 * halfG1) + (halfB1 * halfB1) + (halfA1 * halfA1)) / halfW1;
 520
 10478521                double halfR2 = wholeR - halfR1;
 10478522                double halfG2 = wholeG - halfG1;
 10478523                double halfB2 = wholeB - halfB1;
 10478524                double halfA2 = wholeA - halfA1;
 525
 10478526                temp += ((halfR2 * halfR2) + (halfG2 * halfG2) + (halfB2 * halfB2) + (halfA2 * halfA2)) / halfW2;
 527
 10478528                if (temp > max)
 529                {
 2870530                    max = temp;
 2870531                    cut = i;
 532                }
 533            }
 534
 5340535            return max;
 536        }
 537
 538        /// <summary>
 539        /// Cuts a AlphaBox.
 540        /// </summary>
 541        /// <param name="set1">The first set.</param>
 542        /// <param name="set2">The second set.</param>
 543        /// <returns>Returns a value indicating whether the AlphaBox has been split.</returns>
 544        private bool Cut(AlphaBox set1, AlphaBox set2)
 545        {
 1335546            double wholeR = Volume(set1, this.vmr);
 1335547            double wholeG = Volume(set1, this.vmg);
 1335548            double wholeB = Volume(set1, this.vmb);
 1335549            double wholeA = Volume(set1, this.vma);
 1335550            double wholeW = Volume(set1, this.vwt);
 551
 552            int cutr;
 553            int cutg;
 554            int cutb;
 555            int cuta;
 556
 1335557            double maxr = this.Maximize(set1, 3, set1.R0 + 1, set1.R1, out cutr, wholeR, wholeG, wholeB, wholeA, wholeW)
 1335558            double maxg = this.Maximize(set1, 2, set1.G0 + 1, set1.G1, out cutg, wholeR, wholeG, wholeB, wholeA, wholeW)
 1335559            double maxb = this.Maximize(set1, 1, set1.B0 + 1, set1.B1, out cutb, wholeR, wholeG, wholeB, wholeA, wholeW)
 1335560            double maxa = this.Maximize(set1, 0, set1.A0 + 1, set1.A1, out cuta, wholeR, wholeG, wholeB, wholeA, wholeW)
 561
 562            int dir;
 563
 1335564            if ((maxr >= maxg) && (maxr >= maxb) && (maxr >= maxa))
 565            {
 815566                dir = 3;
 567
 815568                if (cutr < 0)
 569                {
 546570                    return false;
 571                }
 572            }
 520573            else if ((maxg >= maxr) && (maxg >= maxb) && (maxg >= maxa))
 574            {
 160575                dir = 2;
 576            }
 360577            else if ((maxb >= maxr) && (maxb >= maxg) && (maxb >= maxa))
 578            {
 194579                dir = 1;
 580            }
 581            else
 582            {
 166583                dir = 0;
 584            }
 585
 789586            set2.R1 = set1.R1;
 789587            set2.G1 = set1.G1;
 789588            set2.B1 = set1.B1;
 789589            set2.A1 = set1.A1;
 590
 591            switch (dir)
 592            {
 593                // Red
 594                case 3:
 269595                    set2.R0 = set1.R1 = cutr;
 269596                    set2.G0 = set1.G0;
 269597                    set2.B0 = set1.B0;
 269598                    set2.A0 = set1.A0;
 269599                    break;
 600
 601                // Green
 602                case 2:
 160603                    set2.G0 = set1.G1 = cutg;
 160604                    set2.R0 = set1.R0;
 160605                    set2.B0 = set1.B0;
 160606                    set2.A0 = set1.A0;
 160607                    break;
 608
 609                // Blue
 610                case 1:
 194611                    set2.B0 = set1.B1 = cutb;
 194612                    set2.R0 = set1.R0;
 194613                    set2.G0 = set1.G0;
 194614                    set2.A0 = set1.A0;
 194615                    break;
 616
 617                // Alpha
 618                case 0:
 166619                    set2.A0 = set1.A1 = cuta;
 166620                    set2.R0 = set1.R0;
 166621                    set2.G0 = set1.G0;
 166622                    set2.B0 = set1.B0;
 623                    break;
 624            }
 625
 789626            set1.Update();
 789627            set2.Update();
 628
 789629            return true;
 630        }
 631
 632        /// <summary>
 633        /// Marks a color space tag.
 634        /// </summary>
 635        /// <param name="cube">The cube.</param>
 636        /// <param name="label">A label.</param>
 637        private void Mark(AlphaBox cube, byte label)
 638        {
 47428639            for (int r = cube.R0 + 1; r <= cube.R1; r++)
 640            {
 1553152641                for (int g = cube.G0 + 1; g <= cube.G1; g++)
 642                {
 25624576643                    for (int b = cube.B0 + 1; b <= cube.B1; b++)
 644                    {
 133169152645                        for (int a = cube.A0 + 1; a <= cube.A1; a++)
 646                        {
 54525952647                            this.tag[GetIndex(r, g, b, a)] = label;
 648                        }
 649                    }
 650                }
 651            }
 802652        }
 653
 654        /// <summary>
 655        /// Builds the cube.
 656        /// </summary>
 657        /// <param name="cube">The cube.</param>
 658        /// <param name="colorCount">The color count.</param>
 659        private void BuildCube(out AlphaBox[] cube, ref int colorCount)
 660        {
 13661            cube = new AlphaBox[colorCount];
 13662            double[] vv = new double[colorCount];
 663
 6682664            for (int i = 0; i < colorCount; i++)
 665            {
 3328666                cube[i] = new AlphaBox();
 667            }
 668
 13669            cube[0].R0 = cube[0].G0 = cube[0].B0 = cube[0].A0 = 0;
 13670            cube[0].R1 = cube[0].G1 = cube[0].B1 = IndexCount - 1;
 13671            cube[0].A1 = IndexAlphaCount - 1;
 13672            cube[0].Update();
 673
 13674            int next = 0;
 675
 2670676            for (int i = 1; i < colorCount; i++)
 677            {
 1335678                if (this.Cut(cube[next], cube[i]))
 679                {
 789680                    vv[next] = cube[next].Volume > 1 ? this.Variance(cube[next]) : 0.0;
 789681                    vv[i] = cube[i].Volume > 1 ? this.Variance(cube[i]) : 0.0;
 682                }
 683                else
 684                {
 546685                    vv[next] = 0.0;
 546686                    i--;
 687                }
 688
 1335689                next = 0;
 690
 1335691                double temp = vv[0];
 166158692                for (int k = 1; k <= i; k++)
 693                {
 81744694                    if (vv[k] > temp)
 695                    {
 1259696                        temp = vv[k];
 1259697                        next = k;
 698                    }
 699                }
 700
 1335701                if (temp <= 0.0)
 702                {
 13703                    colorCount = i + 1;
 13704                    break;
 705                }
 706            }
 0707        }
 708
 709        /// <summary>
 710        /// Generates the quantized result.
 711        /// </summary>
 712        /// <param name="image">The image.</param>
 713        /// <param name="colorCount">The color count.</param>
 714        /// <param name="cube">The cube.</param>
 715        /// <returns>The result.</returns>
 716        private ColorQuantizerResult GenerateResult(byte[] image, int colorCount, AlphaBox[] cube)
 717        {
 13718            var quantizedImage = new ColorQuantizerResult(image.Length / 4, colorCount);
 719
 1630720            for (int k = 0; k < colorCount; k++)
 721            {
 802722                this.Mark(cube[k], (byte)k);
 723
 802724                double weight = Volume(cube[k], this.vwt);
 725
 802726                if (weight != 0)
 727                {
 802728                    quantizedImage.Palette[(k * 4) + 3] = (byte)(Volume(cube[k], this.vma) / weight);
 802729                    quantizedImage.Palette[(k * 4) + 2] = (byte)(Volume(cube[k], this.vmr) / weight);
 802730                    quantizedImage.Palette[(k * 4) + 1] = (byte)(Volume(cube[k], this.vmg) / weight);
 802731                    quantizedImage.Palette[k * 4] = (byte)(Volume(cube[k], this.vmb) / weight);
 732                }
 733                else
 734                {
 0735                    quantizedImage.Palette[(k * 4) + 3] = 0xff;
 0736                    quantizedImage.Palette[(k * 4) + 2] = 0;
 0737                    quantizedImage.Palette[(k * 4) + 1] = 0;
 0738                    quantizedImage.Palette[k * 4] = 0;
 739                }
 740            }
 741
 167775262742            for (int i = 0; i < image.Length / 4; i++)
 743            {
 83887618744                int a = image[(i * 4) + 3] >> (8 - IndexAlphaBits);
 83887618745                int r = image[(i * 4) + 2] >> (8 - IndexBits);
 83887618746                int g = image[(i * 4) + 1] >> (8 - IndexBits);
 83887618747                int b = image[i * 4] >> (8 - IndexBits);
 748
 83887618749                int ind = GetIndex(r + 1, g + 1, b + 1, a + 1);
 750
 83887618751                quantizedImage.Bytes[i] = this.tag[ind];
 752            }
 753
 13754            return quantizedImage;
 755        }
 756    }
 757}
 758#endif
 759

https://raw.githubusercontent.com/JeremyAnsel/JeremyAnsel.ColorQuant/6d79217e72af9e3af1a8a29c606732adac1e8d87/JeremyAnsel.ColorQuant/JeremyAnsel.ColorQuant/WuAlphaColorQuantizer2.cs

#LineLine coverage
 1// <copyright file="WuAlphaColorQuantizer2.cs" company="Jérémy Ansel">
 2// Copyright (c) 2014-2019 Jérémy Ansel
 3// </copyright>
 4// <license>
 5// Licensed under the MIT license. See LICENSE.txt
 6// </license>
 7
 8#if NET8_0_OR_GREATER
 9namespace JeremyAnsel.ColorQuant
 10{
 11    using System;
 12    using System.Diagnostics.CodeAnalysis;
 13    using System.Runtime.Intrinsics;
 14
 15    /// <summary>
 16    /// A Wu's color quantizer with alpha channel.
 17    /// </summary>
 18    /// <remarks>
 19    /// <para>
 20    /// Based on C Implementation of Xiaolin Wu's Color Quantizer (v. 2)
 21    /// (see Graphics Gems volume II, pages 126-133)
 22    /// (<see href="http://www.ece.mcmaster.ca/~xwu/cq.c"/>).
 23    /// </para>
 24    /// <para>
 25    /// Algorithm: Greedy orthogonal bipartition of RGB space for variance
 26    /// minimization aided by inclusion-exclusion tricks.
 27    /// For speed no nearest neighbor search is done. Slightly
 28    /// better performance can be expected by more sophisticated
 29    /// but more expensive versions.
 30    /// </para>
 31    /// </remarks>
 32    [SuppressMessage("Microsoft.Naming", "CA1709:IdentifiersShouldBeCasedCorrectly", MessageId = "Wu", Justification = "
 33    public sealed class WuAlphaColorQuantizer : IColorQuantizer
 34    {
 35        /// <summary>
 36        /// The maximum color count for the quantize.
 37        /// </summary>
 38        public const int MaxColors = 256;
 39
 40        /// <summary>
 41        /// The index bits.
 42        /// </summary>
 43        private const int IndexBits = 6;
 44
 45        /// <summary>
 46        /// The index alpha bits.
 47        /// </summary>
 48        private const int IndexAlphaBits = 4;
 49
 50        /// <summary>
 51        /// The index count.
 52        /// </summary>
 53        private const int IndexCount = (1 << IndexBits) + 1;
 54
 55        /// <summary>
 56        /// The index alpha count.
 57        /// </summary>
 58        private const int IndexAlphaCount = (1 << IndexAlphaBits) + 1;
 59
 60        /// <summary>
 61        /// The table length.
 62        /// </summary>
 63        private const int TableLength = IndexCount * IndexCount * IndexCount * IndexAlphaCount;
 64
 65        /// <summary>
 66        /// Moment of <c>P(c)</c>.
 67        /// </summary>
 168        private readonly long[] vwt = new long[TableLength];
 69
 170        private readonly Vector256<long>[] vmc = new Vector256<long>[TableLength];
 71
 72        /// <summary>
 73        /// Moment of <c>c^2*P(c)</c>.
 74        /// </summary>
 175        private readonly double[] m2 = new double[TableLength];
 76
 77        /// <summary>
 78        /// Color space tag.
 79        /// </summary>
 180        private readonly byte[] tag = new byte[TableLength];
 81
 82        /// <summary>
 83        /// Quantizes an image.
 84        /// </summary>
 85        /// <param name="image">The image (ARGB).</param>
 86        /// <returns>The result.</returns>
 87        public ColorQuantizerResult Quantize(byte[] image)
 88        {
 1489            return this.Quantize(image, 256);
 90        }
 91
 92        /// <summary>
 93        /// Quantizes an image.
 94        /// </summary>
 95        /// <param name="image">The image (ARGB).</param>
 96        /// <param name="colorCount">The color count.</param>
 97        /// <returns>The result.</returns>
 98        public ColorQuantizerResult Quantize(byte[] image, int colorCount)
 99        {
 17100            if (image == null)
 101            {
 2102                throw new ArgumentNullException(nameof(image));
 103            }
 104
 15105            if (colorCount < 1 || colorCount > 256)
 106            {
 2107                throw new ArgumentOutOfRangeException(nameof(colorCount));
 108            }
 109
 13110            this.Clear();
 111
 13112            this.Build3DHistogram(image);
 13113            this.Get3DMoments();
 114
 115            AlphaBox[] cube;
 13116            this.BuildCube(out cube, ref colorCount);
 117
 13118            return this.GenerateResult(image, colorCount, cube);
 119        }
 120
 121        /// <summary>
 122        /// Gets an index.
 123        /// </summary>
 124        /// <param name="r">The red value.</param>
 125        /// <param name="g">The green value.</param>
 126        /// <param name="b">The blue value.</param>
 127        /// <param name="a">The alpha value.</param>
 128        /// <returns>The index.</returns>
 129        private static int GetIndex(int r, int g, int b, int a)
 130        {
 332839444131            return (r << ((IndexBits * 2) + IndexAlphaBits))
 332839444132                + (r << (IndexBits + IndexAlphaBits + 1))
 332839444133                + (g << (IndexBits + IndexAlphaBits))
 332839444134                + (r << (IndexBits * 2))
 332839444135                + (r << (IndexBits + 1))
 332839444136                + (g << IndexBits)
 332839444137                + ((r + g + b) << IndexAlphaBits)
 332839444138                + r + g + b + a;
 139        }
 140
 141        /// <summary>
 142        /// Computes sum over a box of any given statistic.
 143        /// </summary>
 144        /// <param name="cube">The cube.</param>
 145        /// <param name="moment">The moment.</param>
 146        /// <returns>The result.</returns>
 147        private static long Volume(AlphaBox cube, long[] moment)
 148        {
 3714149            return moment[cube.Index1111]
 3714150                - moment[cube.Index1110]
 3714151                - moment[cube.Index1101]
 3714152                + moment[cube.Index1100]
 3714153                - moment[cube.Index1011]
 3714154                + moment[cube.Index1010]
 3714155                + moment[cube.Index1001]
 3714156                - moment[cube.Index1000]
 3714157                - moment[cube.Index0111]
 3714158                + moment[cube.Index0110]
 3714159                + moment[cube.Index0101]
 3714160                - moment[cube.Index0100]
 3714161                + moment[cube.Index0011]
 3714162                - moment[cube.Index0010]
 3714163                - moment[cube.Index0001]
 3714164                + moment[cube.Index0000];
 165        }
 166
 167        /// <summary>
 168        /// Computes sum over a box of any given statistic.
 169        /// </summary>
 170        /// <param name="cube">The cube.</param>
 171        /// <param name="moment">The moment.</param>
 172        /// <returns>The result.</returns>
 173        private static Vector256<long> Volume(AlphaBox cube, Vector256<long>[] moment)
 174        {
 3714175            return moment[cube.Index1111]
 3714176                - moment[cube.Index1110]
 3714177                - moment[cube.Index1101]
 3714178                + moment[cube.Index1100]
 3714179                - moment[cube.Index1011]
 3714180                + moment[cube.Index1010]
 3714181                + moment[cube.Index1001]
 3714182                - moment[cube.Index1000]
 3714183                - moment[cube.Index0111]
 3714184                + moment[cube.Index0110]
 3714185                + moment[cube.Index0101]
 3714186                - moment[cube.Index0100]
 3714187                + moment[cube.Index0011]
 3714188                - moment[cube.Index0010]
 3714189                - moment[cube.Index0001]
 3714190                + moment[cube.Index0000];
 191        }
 192
 193        /// <summary>
 194        /// Computes part of Volume(cube, moment) that doesn't depend on r1, g1, or b1 (depending on direction).
 195        /// </summary>
 196        /// <param name="cube">The cube.</param>
 197        /// <param name="direction">The direction.</param>
 198        /// <param name="moment">The moment.</param>
 199        /// <returns>The result.</returns>
 200        private static long Bottom(AlphaBox cube, int direction, long[] moment)
 201        {
 202            switch (direction)
 203            {
 204                // Red
 205                case 3:
 1335206                    return -moment[cube.Index0111]
 1335207                        + moment[cube.Index0110]
 1335208                        + moment[cube.Index0101]
 1335209                        - moment[cube.Index0100]
 1335210                        + moment[cube.Index0011]
 1335211                        - moment[cube.Index0010]
 1335212                        - moment[cube.Index0001]
 1335213                        + moment[cube.Index0000];
 214
 215                // Green
 216                case 2:
 1335217                    return -moment[cube.Index1011]
 1335218                        + moment[cube.Index1010]
 1335219                        + moment[cube.Index1001]
 1335220                        - moment[cube.Index1000]
 1335221                        + moment[cube.Index0011]
 1335222                        - moment[cube.Index0010]
 1335223                        - moment[cube.Index0001]
 1335224                        + moment[cube.Index0000];
 225
 226                // Blue
 227                case 1:
 1335228                    return -moment[cube.Index1101]
 1335229                        + moment[cube.Index1100]
 1335230                        + moment[cube.Index1001]
 1335231                        - moment[cube.Index1000]
 1335232                        + moment[cube.Index0101]
 1335233                        - moment[cube.Index0100]
 1335234                        - moment[cube.Index0001]
 1335235                        + moment[cube.Index0000];
 236
 237                // Alpha
 238                case 0:
 1335239                    return -moment[cube.Index1110]
 1335240                        + moment[cube.Index1100]
 1335241                        + moment[cube.Index1010]
 1335242                        - moment[cube.Index1000]
 1335243                        + moment[cube.Index0110]
 1335244                        - moment[cube.Index0100]
 1335245                        - moment[cube.Index0010]
 1335246                        + moment[cube.Index0000];
 247
 248                default:
 0249                    throw new ArgumentOutOfRangeException(nameof(direction));
 250            }
 251        }
 252
 253        /// <summary>
 254        /// Computes part of Volume(cube, moment) that doesn't depend on r1, g1, or b1 (depending on direction).
 255        /// </summary>
 256        /// <param name="cube">The cube.</param>
 257        /// <param name="direction">The direction.</param>
 258        /// <param name="moment">The moment.</param>
 259        /// <returns>The result.</returns>
 260        private static Vector256<long> Bottom(AlphaBox cube, int direction, Vector256<long>[] moment)
 261        {
 262            switch (direction)
 263            {
 264                // Red
 265                case 3:
 1335266                    return -moment[cube.Index0111]
 1335267                        + moment[cube.Index0110]
 1335268                        + moment[cube.Index0101]
 1335269                        - moment[cube.Index0100]
 1335270                        + moment[cube.Index0011]
 1335271                        - moment[cube.Index0010]
 1335272                        - moment[cube.Index0001]
 1335273                        + moment[cube.Index0000];
 274
 275                // Green
 276                case 2:
 1335277                    return -moment[cube.Index1011]
 1335278                        + moment[cube.Index1010]
 1335279                        + moment[cube.Index1001]
 1335280                        - moment[cube.Index1000]
 1335281                        + moment[cube.Index0011]
 1335282                        - moment[cube.Index0010]
 1335283                        - moment[cube.Index0001]
 1335284                        + moment[cube.Index0000];
 285
 286                // Blue
 287                case 1:
 1335288                    return -moment[cube.Index1101]
 1335289                        + moment[cube.Index1100]
 1335290                        + moment[cube.Index1001]
 1335291                        - moment[cube.Index1000]
 1335292                        + moment[cube.Index0101]
 1335293                        - moment[cube.Index0100]
 1335294                        - moment[cube.Index0001]
 1335295                        + moment[cube.Index0000];
 296
 297                // Alpha
 298                case 0:
 1335299                    return -moment[cube.Index1110]
 1335300                        + moment[cube.Index1100]
 1335301                        + moment[cube.Index1010]
 1335302                        - moment[cube.Index1000]
 1335303                        + moment[cube.Index0110]
 1335304                        - moment[cube.Index0100]
 1335305                        - moment[cube.Index0010]
 1335306                        + moment[cube.Index0000];
 307
 308                default:
 0309                    throw new ArgumentOutOfRangeException(nameof(direction));
 310            }
 311        }
 312
 313        /// <summary>
 314        /// Computes remainder of Volume(cube, moment), substituting position for r1, g1, or b1 (depending on direction)
 315        /// </summary>
 316        /// <param name="cube">The cube.</param>
 317        /// <param name="direction">The direction.</param>
 318        /// <param name="position">The position.</param>
 319        /// <param name="moment">The moment.</param>
 320        /// <returns>The result.</returns>
 321        private static long Top(AlphaBox cube, int direction, int position, long[] moment)
 322        {
 323            switch (direction)
 324            {
 325                // Red
 326                case 3:
 41673327                    return moment[GetIndex(position, cube.G1, cube.B1, cube.A1)]
 41673328                        - moment[GetIndex(position, cube.G1, cube.B1, cube.A0)]
 41673329                        - moment[GetIndex(position, cube.G1, cube.B0, cube.A1)]
 41673330                        + moment[GetIndex(position, cube.G1, cube.B0, cube.A0)]
 41673331                        - moment[GetIndex(position, cube.G0, cube.B1, cube.A1)]
 41673332                        + moment[GetIndex(position, cube.G0, cube.B1, cube.A0)]
 41673333                        + moment[GetIndex(position, cube.G0, cube.B0, cube.A1)]
 41673334                        - moment[GetIndex(position, cube.G0, cube.B0, cube.A0)];
 335
 336                // Green
 337                case 2:
 57609338                    return moment[GetIndex(cube.R1, position, cube.B1, cube.A1)]
 57609339                        - moment[GetIndex(cube.R1, position, cube.B1, cube.A0)]
 57609340                        - moment[GetIndex(cube.R1, position, cube.B0, cube.A1)]
 57609341                        + moment[GetIndex(cube.R1, position, cube.B0, cube.A0)]
 57609342                        - moment[GetIndex(cube.R0, position, cube.B1, cube.A1)]
 57609343                        + moment[GetIndex(cube.R0, position, cube.B1, cube.A0)]
 57609344                        + moment[GetIndex(cube.R0, position, cube.B0, cube.A1)]
 57609345                        - moment[GetIndex(cube.R0, position, cube.B0, cube.A0)];
 346
 347                // Blue
 348                case 1:
 58761349                    return moment[GetIndex(cube.R1, cube.G1, position, cube.A1)]
 58761350                        - moment[GetIndex(cube.R1, cube.G1, position, cube.A0)]
 58761351                        - moment[GetIndex(cube.R1, cube.G0, position, cube.A1)]
 58761352                        + moment[GetIndex(cube.R1, cube.G0, position, cube.A0)]
 58761353                        - moment[GetIndex(cube.R0, cube.G1, position, cube.A1)]
 58761354                        + moment[GetIndex(cube.R0, cube.G1, position, cube.A0)]
 58761355                        + moment[GetIndex(cube.R0, cube.G0, position, cube.A1)]
 58761356                        - moment[GetIndex(cube.R0, cube.G0, position, cube.A0)];
 357
 358                // Alpha
 359                case 0:
 17273360                    return moment[GetIndex(cube.R1, cube.G1, cube.B1, position)]
 17273361                        - moment[GetIndex(cube.R1, cube.G1, cube.B0, position)]
 17273362                        - moment[GetIndex(cube.R1, cube.G0, cube.B1, position)]
 17273363                        + moment[GetIndex(cube.R1, cube.G0, cube.B0, position)]
 17273364                        - moment[GetIndex(cube.R0, cube.G1, cube.B1, position)]
 17273365                        + moment[GetIndex(cube.R0, cube.G1, cube.B0, position)]
 17273366                        + moment[GetIndex(cube.R0, cube.G0, cube.B1, position)]
 17273367                        - moment[GetIndex(cube.R0, cube.G0, cube.B0, position)];
 368
 369                default:
 0370                    throw new ArgumentOutOfRangeException(nameof(direction));
 371            }
 372        }
 373
 374        /// <summary>
 375        /// Computes remainder of Volume(cube, moment), substituting position for r1, g1, or b1 (depending on direction)
 376        /// </summary>
 377        /// <param name="cube">The cube.</param>
 378        /// <param name="direction">The direction.</param>
 379        /// <param name="position">The position.</param>
 380        /// <param name="moment">The moment.</param>
 381        /// <returns>The result.</returns>
 382        private static Vector256<long> Top(AlphaBox cube, int direction, int position, Vector256<long>[] moment)
 383        {
 384            switch (direction)
 385            {
 386                // Red
 387                case 3:
 1977388                    return moment[GetIndex(position, cube.G1, cube.B1, cube.A1)]
 1977389                        - moment[GetIndex(position, cube.G1, cube.B1, cube.A0)]
 1977390                        - moment[GetIndex(position, cube.G1, cube.B0, cube.A1)]
 1977391                        + moment[GetIndex(position, cube.G1, cube.B0, cube.A0)]
 1977392                        - moment[GetIndex(position, cube.G0, cube.B1, cube.A1)]
 1977393                        + moment[GetIndex(position, cube.G0, cube.B1, cube.A0)]
 1977394                        + moment[GetIndex(position, cube.G0, cube.B0, cube.A1)]
 1977395                        - moment[GetIndex(position, cube.G0, cube.B0, cube.A0)];
 396
 397                // Green
 398                case 2:
 2733399                    return moment[GetIndex(cube.R1, position, cube.B1, cube.A1)]
 2733400                        - moment[GetIndex(cube.R1, position, cube.B1, cube.A0)]
 2733401                        - moment[GetIndex(cube.R1, position, cube.B0, cube.A1)]
 2733402                        + moment[GetIndex(cube.R1, position, cube.B0, cube.A0)]
 2733403                        - moment[GetIndex(cube.R0, position, cube.B1, cube.A1)]
 2733404                        + moment[GetIndex(cube.R0, position, cube.B1, cube.A0)]
 2733405                        + moment[GetIndex(cube.R0, position, cube.B0, cube.A1)]
 2733406                        - moment[GetIndex(cube.R0, position, cube.B0, cube.A0)];
 407
 408                // Blue
 409                case 1:
 4245410                    return moment[GetIndex(cube.R1, cube.G1, position, cube.A1)]
 4245411                        - moment[GetIndex(cube.R1, cube.G1, position, cube.A0)]
 4245412                        - moment[GetIndex(cube.R1, cube.G0, position, cube.A1)]
 4245413                        + moment[GetIndex(cube.R1, cube.G0, position, cube.A0)]
 4245414                        - moment[GetIndex(cube.R0, cube.G1, position, cube.A1)]
 4245415                        + moment[GetIndex(cube.R0, cube.G1, position, cube.A0)]
 4245416                        + moment[GetIndex(cube.R0, cube.G0, position, cube.A1)]
 4245417                        - moment[GetIndex(cube.R0, cube.G0, position, cube.A0)];
 418
 419                // Alpha
 420                case 0:
 1523421                    return moment[GetIndex(cube.R1, cube.G1, cube.B1, position)]
 1523422                        - moment[GetIndex(cube.R1, cube.G1, cube.B0, position)]
 1523423                        - moment[GetIndex(cube.R1, cube.G0, cube.B1, position)]
 1523424                        + moment[GetIndex(cube.R1, cube.G0, cube.B0, position)]
 1523425                        - moment[GetIndex(cube.R0, cube.G1, cube.B1, position)]
 1523426                        + moment[GetIndex(cube.R0, cube.G1, cube.B0, position)]
 1523427                        + moment[GetIndex(cube.R0, cube.G0, cube.B1, position)]
 1523428                        - moment[GetIndex(cube.R0, cube.G0, cube.B0, position)];
 429
 430                default:
 0431                    throw new ArgumentOutOfRangeException(nameof(direction));
 432            }
 433        }
 434
 435        /// <summary>
 436        /// Clears the tables.
 437        /// </summary>
 438        private void Clear()
 439        {
 13440            this.vwt.AsSpan().Clear();
 13441            this.vmc.AsSpan().Clear();
 13442            this.m2.AsSpan().Clear();
 13443            this.tag.AsSpan().Clear();
 13444        }
 445
 446        /// <summary>
 447        /// Builds a 3-D color histogram of <c>counts, r/g/b, c^2</c>.
 448        /// </summary>
 449        /// <param name="image">The image.</param>
 450        private void Build3DHistogram(byte[] image)
 451        {
 167775262452            for (int i = 0; i < image.Length; i += 4)
 453            {
 83887618454                int a = image[i + 3];
 83887618455                int r = image[i + 2];
 83887618456                int g = image[i + 1];
 83887618457                int b = image[i];
 458
 83887618459                int inr = r >> (8 - IndexBits);
 83887618460                int ing = g >> (8 - IndexBits);
 83887618461                int inb = b >> (8 - IndexBits);
 83887618462                int ina = a >> (8 - IndexAlphaBits);
 463
 83887618464                int ind = GetIndex(inr + 1, ing + 1, inb + 1, ina + 1);
 465
 83887618466                this.vwt[ind]++;
 83887618467                this.vmc[ind] += Vector256.Create(r, g, b, a);
 83887618468                this.m2[ind] += (r * r) + (g * g) + (b * b) + (a * a);
 469            }
 13470        }
 471
 472        /// <summary>
 473        /// Converts the histogram into moments so that we can rapidly calculate
 474        /// the sums of the above quantities over any desired box.
 475        /// </summary>
 476        private void Get3DMoments()
 477        {
 13478            Span<long> volume = stackalloc long[IndexCount * IndexAlphaCount];
 13479            Span<Vector256<long>> volumeC = stackalloc Vector256<long>[IndexCount * IndexAlphaCount];
 13480            Span<double> volume2 = stackalloc double[IndexCount * IndexAlphaCount];
 481
 13482            Span<long> area = stackalloc long[IndexAlphaCount];
 13483            Span<Vector256<long>> areaC = stackalloc Vector256<long>[IndexAlphaCount];
 13484            Span<double> area2 = stackalloc double[IndexAlphaCount];
 485
 1690486            for (int r = 1; r < IndexCount; r++)
 487            {
 832488                volume.Clear();
 832489                volumeC.Clear();
 832490                volume2.Clear();
 491
 108160492                for (int g = 1; g < IndexCount; g++)
 493                {
 53248494                    area.Clear();
 53248495                    areaC.Clear();
 53248496                    area2.Clear();
 497
 6922240498                    for (int b = 1; b < IndexCount; b++)
 499                    {
 3407872500                        long line = 0;
 3407872501                        Vector256<long> lineC = Vector256<long>.Zero;
 3407872502                        double line2 = 0;
 503
 115867648504                        for (int a = 1; a < IndexAlphaCount; a++)
 505                        {
 54525952506                            int ind1 = GetIndex(r, g, b, a);
 507
 54525952508                            line += this.vwt[ind1];
 54525952509                            lineC += this.vmc[ind1];
 54525952510                            line2 += this.m2[ind1];
 511
 54525952512                            area[a] += line;
 54525952513                            areaC[a] += lineC;
 54525952514                            area2[a] += line2;
 515
 54525952516                            int inv = (b * IndexAlphaCount) + a;
 517
 54525952518                            volume[inv] += area[a];
 54525952519                            volumeC[inv] += areaC[a];
 54525952520                            volume2[inv] += area2[a];
 521
 54525952522                            int ind2 = ind1 - GetIndex(1, 0, 0, 0);
 523
 54525952524                            this.vwt[ind1] = this.vwt[ind2] + volume[inv];
 54525952525                            this.vmc[ind1] = this.vmc[ind2] + volumeC[inv];
 54525952526                            this.m2[ind1] = this.m2[ind2] + volume2[inv];
 527                        }
 528                    }
 529                }
 530            }
 13531        }
 532
 533        /// <summary>
 534        /// Computes the weighted variance of a box.
 535        /// </summary>
 536        /// <param name="cube">The cube.</param>
 537        /// <returns>The result.</returns>
 538        private double Variance(AlphaBox cube)
 539        {
 1577540            Vector256<long> d = Volume(cube, this.vmc);
 541
 1577542            double xx =
 1577543                this.m2[cube.Index1111]
 1577544                - this.m2[cube.Index1110]
 1577545                - this.m2[cube.Index1101]
 1577546                + this.m2[cube.Index1100]
 1577547                - this.m2[cube.Index1011]
 1577548                + this.m2[cube.Index1010]
 1577549                + this.m2[cube.Index1001]
 1577550                - this.m2[cube.Index1000]
 1577551                - this.m2[cube.Index0111]
 1577552                + this.m2[cube.Index0110]
 1577553                + this.m2[cube.Index0101]
 1577554                - this.m2[cube.Index0100]
 1577555                + this.m2[cube.Index0011]
 1577556                - this.m2[cube.Index0010]
 1577557                - this.m2[cube.Index0001]
 1577558                + this.m2[cube.Index0000];
 559
 1577560            Vector256<double> dv = Vector256.Create((double)d[0], (double)d[1], (double)d[2], (double)d[3]);
 1577561            return xx - (Vector256.Dot(dv, dv) / Volume(cube, this.vwt));
 562        }
 563
 564        /// <summary>
 565        /// We want to minimize the sum of the variances of two sub-boxes.
 566        /// The sum(c^2) terms can be ignored since their sum over both sub-boxes
 567        /// is the same (the sum for the whole box) no matter where we split.
 568        /// The remaining terms have a minus sign in the variance formula,
 569        /// so we drop the minus sign and maximize the sum of the two terms.
 570        /// </summary>
 571        /// <param name="cube">The cube.</param>
 572        /// <param name="direction">The direction.</param>
 573        /// <param name="first">The first position.</param>
 574        /// <param name="last">The last position.</param>
 575        /// <param name="cut">The cutting point.</param>
 576        /// <param name="wholeC">The whole color.</param>
 577        /// <param name="wholeW">The whole weight.</param>
 578        /// <returns>The result.</returns>
 579        private double Maximize(AlphaBox cube, int direction, int first, int last, out int cut, in Vector256<long> whole
 580        {
 5340581            Vector256<long> baseC = Bottom(cube, direction, this.vmc);
 5340582            long baseW = Bottom(cube, direction, this.vwt);
 583
 5340584            double max = 0.0;
 5340585            cut = -1;
 586
 361312587            for (int i = first; i < last; i++)
 588            {
 175316589                long halfW1 = baseW + Top(cube, direction, i, this.vwt);
 590
 175316591                if (halfW1 == 0)
 592                {
 593                    continue;
 594                }
 595
 141073596                long halfW2 = wholeW - halfW1;
 597
 141073598                if (halfW2 == 0)
 599                {
 600                    continue;
 601                }
 602
 10478603                Vector256<long> halfC1 = baseC + Top(cube, direction, i, this.vmc);
 10478604                Vector256<double> dv1 = Vector256.Create((double)halfC1[0], (double)halfC1[1], (double)halfC1[2], (doubl
 10478605                double temp = Vector256.Dot(dv1, dv1) / halfW1;
 606
 10478607                Vector256<long> halfC2 = wholeC - halfC1;
 10478608                Vector256<double> dv2 = Vector256.Create((double)halfC2[0], (double)halfC2[1], (double)halfC2[2], (doubl
 10478609                temp += Vector256.Dot(dv2, dv2) / halfW2;
 610
 10478611                if (temp > max)
 612                {
 2870613                    max = temp;
 2870614                    cut = i;
 615                }
 616            }
 617
 5340618            return max;
 619        }
 620
 621        /// <summary>
 622        /// Cuts a box.
 623        /// </summary>
 624        /// <param name="set1">The first set.</param>
 625        /// <param name="set2">The second set.</param>
 626        /// <returns>Returns a value indicating whether the box has been split.</returns>
 627        private bool Cut(AlphaBox set1, AlphaBox set2)
 628        {
 1335629            Vector256<long> wholeC = Volume(set1, this.vmc);
 1335630            long wholeW = Volume(set1, this.vwt);
 631
 632            int cutr;
 633            int cutg;
 634            int cutb;
 635            int cuta;
 636
 1335637            double maxr = this.Maximize(set1, 3, set1.R0 + 1, set1.R1, out cutr, wholeC, wholeW);
 1335638            double maxg = this.Maximize(set1, 2, set1.G0 + 1, set1.G1, out cutg, wholeC, wholeW);
 1335639            double maxb = this.Maximize(set1, 1, set1.B0 + 1, set1.B1, out cutb, wholeC, wholeW);
 1335640            double maxa = this.Maximize(set1, 0, set1.A0 + 1, set1.A1, out cuta, wholeC, wholeW);
 641
 642            int dir;
 643
 1335644            if ((maxr >= maxg) && (maxr >= maxb) && (maxr >= maxa))
 645            {
 815646                dir = 3;
 647
 815648                if (cutr < 0)
 649                {
 546650                    return false;
 651                }
 652            }
 520653            else if ((maxg >= maxr) && (maxg >= maxb) && (maxg >= maxa))
 654            {
 160655                dir = 2;
 656            }
 360657            else if ((maxb >= maxr) && (maxb >= maxg) && (maxb >= maxa))
 658            {
 194659                dir = 1;
 660            }
 661            else
 662            {
 166663                dir = 0;
 664            }
 665
 789666            set2.R1 = set1.R1;
 789667            set2.G1 = set1.G1;
 789668            set2.B1 = set1.B1;
 789669            set2.A1 = set1.A1;
 670
 671            switch (dir)
 672            {
 673                // Red
 674                case 3:
 269675                    set2.R0 = set1.R1 = cutr;
 269676                    set2.G0 = set1.G0;
 269677                    set2.B0 = set1.B0;
 269678                    set2.A0 = set1.A0;
 269679                    break;
 680
 681                // Green
 682                case 2:
 160683                    set2.G0 = set1.G1 = cutg;
 160684                    set2.R0 = set1.R0;
 160685                    set2.B0 = set1.B0;
 160686                    set2.A0 = set1.A0;
 160687                    break;
 688
 689                // Blue
 690                case 1:
 194691                    set2.B0 = set1.B1 = cutb;
 194692                    set2.R0 = set1.R0;
 194693                    set2.G0 = set1.G0;
 194694                    set2.A0 = set1.A0;
 194695                    break;
 696
 697                // Alpha
 698                case 0:
 166699                    set2.A0 = set1.A1 = cuta;
 166700                    set2.R0 = set1.R0;
 166701                    set2.G0 = set1.G0;
 166702                    set2.B0 = set1.B0;
 703                    break;
 704            }
 705
 789706            set1.Update();
 789707            set2.Update();
 708
 789709            return true;
 710        }
 711
 712        /// <summary>
 713        /// Marks a color space tag.
 714        /// </summary>
 715        /// <param name="cube">The cube.</param>
 716        /// <param name="label">A label.</param>
 717        private void Mark(AlphaBox cube, byte label)
 718        {
 47428719            for (int r = cube.R0 + 1; r <= cube.R1; r++)
 720            {
 1553152721                for (int g = cube.G0 + 1; g <= cube.G1; g++)
 722                {
 25624576723                    for (int b = cube.B0 + 1; b <= cube.B1; b++)
 724                    {
 133169152725                        for (int a = cube.A0 + 1; a <= cube.A1; a++)
 726                        {
 54525952727                            this.tag[GetIndex(r, g, b, a)] = label;
 728                        }
 729                    }
 730                }
 731            }
 802732        }
 733
 734        /// <summary>
 735        /// Builds the cube.
 736        /// </summary>
 737        /// <param name="cube">The cube.</param>
 738        /// <param name="colorCount">The color count.</param>
 739        private void BuildCube(out AlphaBox[] cube, ref int colorCount)
 740        {
 13741            cube = new AlphaBox[colorCount];
 13742            Span<double> vv = stackalloc double[colorCount];
 743
 6682744            for (int i = 0; i < colorCount; i++)
 745            {
 3328746                cube[i] = new AlphaBox();
 747            }
 748
 13749            cube[0].R0 = cube[0].G0 = cube[0].B0 = cube[0].A0 = 0;
 13750            cube[0].R1 = cube[0].G1 = cube[0].B1 = IndexCount - 1;
 13751            cube[0].A1 = IndexAlphaCount - 1;
 13752            cube[0].Update();
 753
 13754            int next = 0;
 755
 2670756            for (int i = 1; i < colorCount; i++)
 757            {
 1335758                if (this.Cut(cube[next], cube[i]))
 759                {
 789760                    vv[next] = cube[next].Volume > 1 ? this.Variance(cube[next]) : 0.0;
 789761                    vv[i] = cube[i].Volume > 1 ? this.Variance(cube[i]) : 0.0;
 762                }
 763                else
 764                {
 546765                    vv[next] = 0.0;
 546766                    i--;
 767                }
 768
 1335769                next = 0;
 770
 1335771                double temp = vv[0];
 166158772                for (int k = 1; k <= i; k++)
 773                {
 81744774                    if (vv[k] > temp)
 775                    {
 1259776                        temp = vv[k];
 1259777                        next = k;
 778                    }
 779                }
 780
 1335781                if (temp <= 0.0)
 782                {
 13783                    colorCount = i + 1;
 13784                    break;
 785                }
 786            }
 0787        }
 788
 789        /// <summary>
 790        /// Generates the quantized result.
 791        /// </summary>
 792        /// <param name="image">The image.</param>
 793        /// <param name="colorCount">The color count.</param>
 794        /// <param name="cube">The cube.</param>
 795        /// <returns>The result.</returns>
 796        private ColorQuantizerResult GenerateResult(byte[] image, int colorCount, AlphaBox[] cube)
 797        {
 13798            var quantizedImage = new ColorQuantizerResult(image.Length / 4, colorCount);
 799
 1630800            for (int k = 0; k < colorCount; k++)
 801            {
 802802                this.Mark(cube[k], (byte)k);
 803
 802804                double weight = Volume(cube[k], this.vwt);
 805
 802806                if (weight != 0)
 807                {
 802808                    Vector256<long> color = Volume(cube[k], this.vmc);
 802809                    quantizedImage.Palette[(k * 4) + 3] = (byte)(color[3] / weight);
 802810                    quantizedImage.Palette[(k * 4) + 2] = (byte)(color[0] / weight);
 802811                    quantizedImage.Palette[(k * 4) + 1] = (byte)(color[1] / weight);
 802812                    quantizedImage.Palette[k * 4] = (byte)(color[2] / weight);
 813                }
 814                else
 815                {
 0816                    quantizedImage.Palette[(k * 4) + 3] = 0xff;
 0817                    quantizedImage.Palette[(k * 4) + 2] = 0;
 0818                    quantizedImage.Palette[(k * 4) + 1] = 0;
 0819                    quantizedImage.Palette[k * 4] = 0;
 820                }
 821            }
 822
 167775262823            for (int i = 0; i < image.Length / 4; i++)
 824            {
 83887618825                int a = image[(i * 4) + 3] >> (8 - IndexAlphaBits);
 83887618826                int r = image[(i * 4) + 2] >> (8 - IndexBits);
 83887618827                int g = image[(i * 4) + 1] >> (8 - IndexBits);
 83887618828                int b = image[i * 4] >> (8 - IndexBits);
 829
 83887618830                int ind = GetIndex(r + 1, g + 1, b + 1, a + 1);
 831
 83887618832                quantizedImage.Bytes[i] = this.tag[ind];
 833            }
 834
 13835            return quantizedImage;
 836        }
 837    }
 838}
 839#endif
 840
Generated by: ReportGenerator 5.4.12.0
9/4/2025 - 6:55:30 PM
GitHub | reportgenerator.io

Methods/Properties

.ctor()
Quantize(System.Byte[])
Quantize(System.Byte[],System.Int32)
GetIndex(System.Int32,System.Int32,System.Int32,System.Int32)
Volume(JeremyAnsel.ColorQuant.AlphaBox,System.Int64[])
Bottom(JeremyAnsel.ColorQuant.AlphaBox,System.Int32,System.Int64[])
Top(JeremyAnsel.ColorQuant.AlphaBox,System.Int32,System.Int32,System.Int64[])
Clear()
Build3DHistogram(System.Byte[])
Get3DMoments()
Variance(JeremyAnsel.ColorQuant.AlphaBox)
Maximize(JeremyAnsel.ColorQuant.AlphaBox,System.Int32,System.Int32,System.Int32,System.Int32&,System.Double,System.Double,System.Double,System.Double,System.Double)
Cut(JeremyAnsel.ColorQuant.AlphaBox,JeremyAnsel.ColorQuant.AlphaBox)
Mark(JeremyAnsel.ColorQuant.AlphaBox,System.Byte)
BuildCube(JeremyAnsel.ColorQuant.AlphaBox[]&,System.Int32&)
GenerateResult(System.Byte[],System.Int32,JeremyAnsel.ColorQuant.AlphaBox[])
.ctor()
Quantize(System.Byte[])
Quantize(System.Byte[],System.Int32)
GetIndex(System.Int32,System.Int32,System.Int32,System.Int32)
Volume(JeremyAnsel.ColorQuant.AlphaBox,System.Int64[])
Volume(JeremyAnsel.ColorQuant.AlphaBox,System.Runtime.Intrinsics.Vector256`1<System.Int64>[])
Bottom(JeremyAnsel.ColorQuant.AlphaBox,System.Int32,System.Int64[])
Bottom(JeremyAnsel.ColorQuant.AlphaBox,System.Int32,System.Runtime.Intrinsics.Vector256`1<System.Int64>[])
Top(JeremyAnsel.ColorQuant.AlphaBox,System.Int32,System.Int32,System.Int64[])
Top(JeremyAnsel.ColorQuant.AlphaBox,System.Int32,System.Int32,System.Runtime.Intrinsics.Vector256`1<System.Int64>[])
Clear()
Build3DHistogram(System.Byte[])
Get3DMoments()
Variance(JeremyAnsel.ColorQuant.AlphaBox)
Maximize(JeremyAnsel.ColorQuant.AlphaBox,System.Int32,System.Int32,System.Int32,System.Int32&,System.Runtime.Intrinsics.Vector256`1<System.Int64>&,System.Int64)
Cut(JeremyAnsel.ColorQuant.AlphaBox,JeremyAnsel.ColorQuant.AlphaBox)
Mark(JeremyAnsel.ColorQuant.AlphaBox,System.Byte)
BuildCube(JeremyAnsel.ColorQuant.AlphaBox[]&,System.Int32&)
GenerateResult(System.Byte[],System.Int32,JeremyAnsel.ColorQuant.AlphaBox[])