fix typo

Remove SSE tests from SConstruct.
Minimize libimagequant files.

SConstruct

@@ -479,6 +479,7 @@ opts.AddVariables(
     BoolVariable('COLOR_PRINT', 'Print build status information in color', 'True'),
     BoolVariable('BIGINT', 'Compile support for 64-bit integers in mapnik::value', 'True'),
     BoolVariable('QUIET', 'Reduce build verbosity', 'False'),
+    EnumVariable('USE_SSE', 'Enable support for SSE vector instructions','yes',['yes','no','platform_default'])
     )
 
 # variables to pickle after successful configure step
@@ -1870,6 +1871,7 @@ if not preconfigured:
     env.Prepend(CPPPATH = '#deps/mapbox/geometry/include')
     env.Prepend(CPPPATH = '#deps/mapbox/protozero/include')
     env.Prepend(CPPPATH = '#deps/mapbox/polylabel/include')
+    env.Prepend(CPPPATH = '#deps/pngquant')
     # prepend deps dir for auxillary headers
     env.Prepend(CPPPATH = '#deps')
 

deps/pngquant/COPYRIGHT

@@ -0,0 +1,36 @@
+© 1997-2002 by Greg Roelofs; based on an idea by Stefan Schneider.
+© 2009-2014 by Kornel Lesiński.
+
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or without modification,
+are permitted provided that the following conditions are met:
+
+1. Redistributions of source code must retain the above copyright notice,
+   this list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above copyright notice,
+   this list of conditions and the following disclaimer in the documentation
+   and/or other materials provided with the distribution.
+
+THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
+DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
+FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
+CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
+OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
+OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+
+- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+
+© 1989, 1991 by Jef Poskanzer.
+
+Permission to use, copy, modify, and distribute this software and its
+documentation for any purpose and without fee is hereby granted, provided
+that the above copyright notice appear in all copies and that both that
+copyright notice and this permission notice appear in supporting
+documentation.  This software is provided "as is" without express or
+implied warranty.

deps/pngquant/MANUAL.md

@@ -0,0 +1,490 @@
+# libimagequant—Image Quantization Library
+
+Small, portable C library for high-quality conversion of RGBA images to 8-bit indexed-color (palette) images.
+It's powering [pngquant2](http://pngquant.org).
+
+## License
+
+[BSD](https://raw.github.com/pornel/pngquant/master/lib/COPYRIGHT).
+It can be linked with both free and closed-source software.
+
+## Download
+
+The [library](http://pngquant.org/lib) is currently a part of the [pngquant2 project](https://github.com/pornel/pngquant/tree/master/lib).
+
+Files needed for the library are only in the `lib/` directory inside the repository (and you can ignore the rest).
+
+## Compiling and Linking
+
+The library can be linked with ANSI C and C++ programs. It has no external dependencies.
+
+To build on Unix-like systems run:
+
+    make -C lib
+
+it will create `lib/libimagequant.a` which you can link with your program.
+
+    gcc yourprogram.c /path/to/lib/libimagequant.a
+
+On BSD, use `gmake` (GNU make) rather than the native `make`.
+
+Alternatively you can compile the library with your program simply by including all `.c` files (and define `NDEBUG` to get a fast version):
+
+    gcc -std=c99 -O3 -DNDEBUG lib/*.c yourprogram.c
+
+### Compiling on Windows/Visual Studio
+
+The library can be compiled with any C compiler that has at least basic support for C99 (GCC, clang, ICC, C++ Builder, even Tiny C Compiler), but Visual Studio 2012 and older are not up to date with the 1999 C standard. There are 2 options for using `libimagequant` on Windows:
+
+ * Use Visual Studio **2013** (MSVC 18) and an [MSVC-compatible branch of the library](https://github.com/pornel/pngquant/tree/msvc/lib)
+ * Or use GCC from [MinGW](http://www.mingw.org). Use GCC to build `libimagequant.a` (using the instructions above for Unix) and add it along with `libgcc.a` (shipped with the MinGW compiler) to your VC project.
+
+## Overview
+
+The basic flow is:
+
+1. Create attributes object and configure the library.
+2. Create image object from RGBA bitmap or data source.
+3. Perform quantization (generate palette).
+4. Store remapped image and final palette.
+5. Free memory.
+
+Please note that libimagequant only handles raw uncompressed bitmaps in memory and is completely independent of any file format.
+
+<p>
+
+    #include "lib/libimagequant.h"
+
+    liq_attr *attr = liq_attr_create();
+    liq_image *image = liq_image_create_rgba(attr, bitmap, width, height, 0);
+    liq_result *res = liq_quantize_image(attr, image);
+
+    liq_write_remapped_image(res, image, bitmap, bitmap_size);
+    const liq_palette *pal = liq_get_palette(res);
+
+    // use image and palette here
+
+    liq_attr_destroy(attr);
+    liq_image_destroy(image);
+    liq_result_destroy(res);
+
+Functions returning `liq_error` return `LIQ_OK` (`0`) on success and non-zero on error.
+
+It's safe to pass `NULL` to any function accepting `liq_attr`, `liq_image`, `liq_result` (in that case the error code `LIQ_INVALID_POINTER` will be returned). These objects can be reused multiple times.
+
+There are 3 ways to create image object for quantization:
+
+  * `liq_image_create_rgba()` for simple, contiguous RGBA bitmaps (width×height×4 bytes large array).
+  * `liq_image_create_rgba_rows()` for non-contiguous RGBA bitmaps (that have padding between rows or reverse order, e.g. BMP).
+  * `liq_image_create_custom()` for RGB, ABGR, YUV and all other formats that can be converted on-the-fly to RGBA (you have to supply the conversion function).
+
+Note that "image" here means raw uncompressed pixels. If you have a compressed image file, such as PNG, you must use another library (e.g. libpng or lodepng) to decode it first.
+
+## Functions
+
+----
+
+    liq_attr* liq_attr_create(void);
+
+Returns object that will hold initial settings (attributes) for the library. The object should be freed using `liq_attr_destroy()` after it's no longer needed.
+
+Returns `NULL` in the unlikely case that the library cannot run on the current machine (e.g. the library has been compiled for SSE-capable x86 CPU and run on VIA C3 CPU).
+
+----
+
+    liq_error liq_set_max_colors(liq_attr* attr, int colors);
+
+Specifies maximum number of colors to use. The default is 256. Instead of setting a fixed limit it's better to use `liq_set_quality()`.
+
+Returns `LIQ_VALUE_OUT_OF_RANGE` if number of colors is outside the range 2-256.
+
+----
+
+    int liq_get_max_colors(liq_attr* attr);
+
+Returns the value set by `liq_set_max_colors()`.
+
+----
+
+    liq_error liq_set_quality(liq_attr* attr, int minimum, int maximum);
+
+Quality is in range `0` (worst) to `100` (best) and values are analoguous to JPEG quality (i.e. `80` is usually good enough).
+
+Quantization will attempt to use the lowest number of colors needed to achieve `maximum` quality. `maximum` value of `100` is the default and means conversion as good as possible.
+
+If it's not possible to convert the image with at least `minimum` quality (i.e. 256 colors is not enough to meet the minimum quality), then `liq_quantize_image()` will fail. The default minumum is `0` (proceeds regardless of quality).
+
+Quality measures how well the generated palette fits image given to `liq_quantize_image()`. If a different image is remapped with `liq_write_remapped_image()` then actual quality may be different.
+
+Regardless of the quality settings the number of colors won't exceed the maximum (see `liq_set_max_colors()`).
+
+Returns `LIQ_VALUE_OUT_OF_RANGE` if target is lower than minimum or any of them is outside the 0-100 range.
+Returns `LIQ_INVALID_POINTER` if `attr` appears to be invalid.
+
+    liq_attr *attr = liq_attr_create();
+    liq_set_quality(attr, 50, 80); // use quality 80 if possible. Give up if quality drops below 50.
+
+----
+
+    int liq_get_min_quality(liq_attr* attr);
+
+Returns the lower bound set by `liq_set_quality()`.
+
+----
+
+    int liq_get_max_quality(liq_attr* attr);
+
+Returns the upper bound set by `liq_set_quality()`.
+
+----
+
+    liq_image *liq_image_create_rgba(liq_attr *attr, void* bitmap, int width, int height, double gamma);
+
+Creates image object that represents a bitmap later used for quantization and remapping. The bitmap must be contiguous run of RGBA pixels (alpha is the last component, 0 = transparent, 255 = opaque).
+
+The bitmap must not be modified or freed until this object is freed with `liq_image_destroy()`. See also `liq_image_set_memory_ownership()`.
+
+`width` and `height` are dimensions in pixels. An image 10x10 pixel large will need 400-byte bitmap.
+
+`gamma` can be `0` for images with the typical 1/2.2 [gamma](http://en.wikipedia.org/wiki/Gamma_correction).
+Otherwise `gamma` must be > 0 and < 1, e.g. `0.45455` (1/2.2) or `0.55555` (1/1.8). Generated palette will use the same gamma unless `liq_set_output_gamma()` is used. If `liq_set_output_gamma` is not used, then it only affects whether brighter or darker areas of the image will get more palette colors allocated.
+
+Returns `NULL` on failure, e.g. if `bitmap` is `NULL` or `width`/`height` is <= 0.
+
+----
+
+    liq_image *liq_image_create_rgba_rows(liq_attr *attr, void* rows[], int width, int height, double gamma);
+
+Same as `liq_image_create_rgba()`, but takes array of pointers to rows in the bitmap. This allows defining bitmaps with reversed rows (like in BMP), "stride" different than width or using only fragment of a larger bitmap, etc.
+
+`rows` array must have at least `height` elements and each row must be at least `width` RGBA pixels wide.
+
+    unsigned char *bitmap = …;
+    void *rows = malloc(height * sizeof(void*));
+    int bytes_per_row = width * 4 + padding; // stride
+    for(int i=0; i < height; i++) {
+        rows[i] = bitmap + i * bytes_per_row;
+    }
+    liq_image *img = liq_image_create_rgba_rows(attr, rows, width, height, 0);
+    // …
+    liq_image_destroy(img);
+    free(rows);
+
+The row pointers and bitmap must not be modified or freed until this object is freed with `liq_image_destroy()` (you can change that with `liq_image_set_memory_ownership()`).
+
+See also `liq_image_create_rgba()` and `liq_image_create_custom()`.
+
+----
+
+    liq_result *liq_quantize_image(liq_attr *attr, liq_image *input_image);
+
+Performs quantization (palette generation) based on settings in `attr` and pixels of the image.
+
+Returns `NULL` if quantization fails, e.g. due to limit set in `liq_set_quality()`.
+
+See `liq_write_remapped_image()`.
+
+----
+
+    liq_error liq_set_dithering_level(liq_result *res, float dither_level);
+
+Enables/disables dithering in `liq_write_remapped_image()`. Dithering level must be between `0` and `1` (inclusive). Dithering level `0` enables fast non-dithered remapping. Otherwise a variation of Floyd-Steinberg error diffusion is used.
+
+Precision of the dithering algorithm depends on the speed setting, see `liq_set_speed()`.
+
+Returns `LIQ_VALUE_OUT_OF_RANGE` if the dithering level is outside the 0-1 range.
+
+----
+
+    liq_error liq_write_remapped_image(liq_result *result, liq_image *input_image, void *buffer, size_t buffer_size);
+
+Remaps the image to palette and writes its pixels to the given buffer, 1 pixel per byte. Buffer must be large enough to fit the entire image, i.e. width×height bytes large. For safety, pass size of the buffer as `buffer_size`.
+
+For best performance call `liq_get_palette()` *after* this function, as palette is improved during remapping.
+
+Returns `LIQ_BUFFER_TOO_SMALL` if given size of the buffer is not enough to fit the entire image.
+
+    int buffer_size = width*height;
+    char *buffer = malloc(buffer_size);
+    if (LIQ_OK == liq_write_remapped_image(result, input_image, buffer, buffer_size)) {
+        liq_palette *pal = liq_get_palette(result);
+        // save image
+    }
+
+See `liq_get_palette()` and `liq_write_remapped_image_rows()`.
+
+Please note that it only writes raw uncompressed pixels to memory. It does not perform any compression. If you'd like to create a PNG file then you need to pass the raw pixel data to another library, e.g. libpng or lodepng. See `rwpng.c` in `pngquant` project for an example how to do that.
+
+----
+
+    const liq_palette *liq_get_palette(liq_result *result);
+
+Returns pointer to palette optimized for image that has been quantized or remapped (final refinements are applied to the palette during remapping).
+
+It's valid to call this method before remapping, if you don't plan to remap any images or want to use same palette for multiple images.
+
+`liq_palette->count` contains number of colors (up to 256), `liq_palette->entries[n]` contains RGBA value for nth palette color.
+
+The palette is **temporary and read-only**. You must copy the palette elsewhere *before* calling `liq_result_destroy()`.
+
+Returns `NULL` on error.
+
+----
+
+    void liq_attr_destroy(liq_attr *);
+    void liq_image_destroy(liq_image *);
+    void liq_result_destroy(liq_result *);
+
+Releases memory owned by the given object. Object must not be used any more after it has been freed.
+
+Freeing `liq_result` also frees any `liq_palette` obtained from it.
+
+## Advanced Functions
+
+----
+
+    liq_error liq_set_speed(liq_attr* attr, int speed);
+
+Higher speed levels disable expensive algorithms and reduce quantization precision. The default speed is `3`. Speed `1` gives marginally better quality at significant CPU cost. Speed `10` has usually 5% lower quality, but is 8 times faster than the default.
+
+High speeds combined with `liq_set_quality()` will use more colors than necessary and will be less likely to meet minimum required quality.
+
+<table><caption>Features dependent on speed</caption>
+<tr><th>Noise-sensitive dithering</th><td>speed 1 to 5</td></tr>
+<tr><th>Forced posterization</th><td>8-10 or if image has more than million colors</td></tr>
+<tr><th>Quantization error known</th><td>1-7 or if minimum quality is set</td></tr>
+<tr><th>Additional quantization techniques</th><td>1-6</td></tr>
+</table>
+
+Returns `LIQ_VALUE_OUT_OF_RANGE` if the speed is outside the 1-10 range.
+
+----
+
+    int liq_get_speed(liq_attr* attr);
+
+Returns the value set by `liq_set_speed()`.
+
+----
+
+    liq_error liq_set_min_opacity(liq_attr* attr, int min);
+
+Alpha values higher than this will be rounded to opaque. This is a workaround for Internet Explorer 6 that truncates semitransparent values to completely transparent. The default is `255` (no change). 238 is a suggested value.
+
+Returns `LIQ_VALUE_OUT_OF_RANGE` if the value is outside the 0-255 range.
+
+----
+
+    int liq_get_min_opacity(liq_attr* attr);
+
+Returns the value set by `liq_set_min_opacity()`.
+
+----
+
+    liq_set_min_posterization(liq_attr* attr, int bits);
+
+Ignores given number of least significant bits in all channels, posterizing image to `2^bits` levels. `0` gives full quality. Use `2` for VGA or 16-bit RGB565 displays, `4` if image is going to be output on a RGB444/RGBA4444 display (e.g. low-quality textures on Android).
+
+Returns `LIQ_VALUE_OUT_OF_RANGE` if the value is outside the 0-4 range.
+
+----
+
+    int liq_get_min_posterization(liq_attr* attr);
+
+Returns the value set by `liq_set_min_posterization()`.
+
+----
+
+    liq_set_last_index_transparent(liq_attr* attr, int is_last);
+
+`0` (default) makes alpha colors sorted before opaque colors. Non-`0` mixes colors together except completely transparent color, which is moved to the end of the palette. This is a workaround for programs that blindly assume the last palette entry is transparent.
+
+----
+
+    liq_image *liq_image_create_custom(liq_attr *attr, liq_image_get_rgba_row_callback *row_callback, void *user_info, int width, int height, double gamma);
+
+<p>
+
+    void image_get_rgba_row_callback(liq_color row_out[], int row_index, int width, void *user_info) {
+        for(int column_index=0; column_index < width; column_index++) {
+            row_out[column_index] = /* generate pixel at (row_index, column_index) */;
+        }
+    }
+
+Creates image object that will use callback to read image data. This allows on-the-fly conversion of images that are not in the RGBA color space.
+
+`user_info` value will be passed to the callback. It may be useful for storing pointer to program's internal representation of the image.
+
+The callback must read/generate `row_index`-th row and write its RGBA pixels to the `row_out` array. Row `width` is given for convenience and will always equal to image width.
+
+The callback will be called multiple times for each row. Quantization and remapping require at least two full passes over image data, so caching of callback's work makes no sense — in such case it's better to convert entire image and use `liq_image_create_rgba()` instead.
+
+To use RGB image:
+
+    void rgb_to_rgba_callback(liq_color row_out[], int row_index, int width, void *user_info) {
+        unsigned char *rgb_row = ((unsigned char *)user_info) + 3*width*row_index;
+
+        for(int i=0; i < width; i++) {
+            row_out[i].r = rgb_row[i*3];
+            row_out[i].g = rgb_row[i*3+1];
+            row_out[i].b = rgb_row[i*3+2];
+            row_out[i].a = 255;
+        }
+    }
+    liq_image *img = liq_image_create_custom(attr, rgb_to_rgba_callback, rgb_bitmap, width, height, 0);
+
+The library doesn't support RGB bitmaps "natively", because supporting only single format allows compiler to inline more code, 4-byte pixel alignment is faster, and SSE instructions operate on 4 values at once, so alpha support is almost free.
+
+----
+
+    liq_error liq_image_set_memory_ownership(liq_image *image, int ownership_flags);
+
+Passes ownership of bitmap and/or rows memory to the `liq_image` object, so you don't have to free it yourself. Memory owned by the object will be freed at its discretion with `free` function specified in `liq_attr_create_with_allocator()` (by default it's stdlib's `free()`).
+
+* `LIQ_OWN_PIXELS` makes bitmap owned by the object. The bitmap will be freed automatically at any point when it's no longer needed. If you set this flag you must **not** free the bitmap yourself. If the image has been created with `liq_image_create_rgba_rows()` then the bitmap address is assumed to be the lowest address of any row.
+
+* `LIQ_OWN_ROWS` makes array of row pointers (but not bitmap pointed by these rows) owned by the object. Rows will be freed when object is deallocated. If you set this flag you must **not** free the rows array yourself. This flag is valid only if the object has been created with `liq_image_create_rgba_rows()`.
+
+These flags can be combined with binary *or*, i.e. `LIQ_OWN_PIXELS | LIQ_OWN_ROWS`.
+
+This function must not be used if the image has been created with `liq_image_create_custom()`.
+
+Returns `LIQ_VALUE_OUT_OF_RANGE` if invalid flags are specified or image is not backed by a bitmap.
+
+----
+
+    liq_error liq_write_remapped_image_rows(liq_result *result, liq_image *input_image, unsigned char **row_pointers);
+
+Similar to `liq_write_remapped_image()`. Writes remapped image, at 1 byte per pixel, to each row pointed by `row_pointers` array. The array must have at least as many elements as height of the image, and each row must have at least as many bytes as width of the image. Rows must not overlap.
+
+For best performance call `liq_get_palette()` *after* this function, as remapping may change the palette.
+
+Returns `LIQ_INVALID_POINTER` if `result` or `input_image` is `NULL`.
+
+----
+
+    double liq_get_quantization_error(liq_result *result);
+
+Returns mean square error of quantization (square of difference between pixel values in the original image and remapped image). Alpha channel and gamma correction are taken into account, so the result isn't exactly the mean square error of all channels.
+
+For most images MSE 1-5 is excellent. 7-10 is OK. 20-30 will have noticeable errors. 100 is awful.
+
+This function should be called *after* `liq_write_remapped_image()`. It may return `-1` if the value is not available (this is affected by `liq_set_speed()` and `liq_set_quality()`).
+
+----
+
+    double liq_get_quantization_quality(liq_result *result);
+
+Analoguous to `liq_get_quantization_error()`, but returns quantization error as quality value in the same 0-100 range that is used by `liq_set_quality()`.
+
+This function should be called *after* `liq_write_remapped_image()`. It may return `-1` if the value is not available (this is affected by `liq_set_speed()` and `liq_set_quality()`).
+
+This function can be used to add upper limit to quality options presented to the user, e.g.
+
+    liq_attr *attr = liq_attr_create();
+    liq_image *img = liq_image_create_rgba(…);
+    liq_result *res = liq_quantize_image(attr, img);
+    int max_attainable_quality = liq_get_quantization_quality(res);
+    printf("Please select quality between 0 and %d: ", max_attainable_quality);
+    int user_selected_quality = prompt();
+    if (user_selected_quality < max_attainable_quality) {
+        liq_set_quality(user_selected_quality, 0);
+        liq_result_destroy(res);
+        res = liq_quantize_image(attr, img);
+    }
+    liq_write_remapped_image(…);
+
+----
+
+    void liq_set_log_callback(liq_attr*, liq_log_callback_function*, void *user_info);
+
+<p>
+
+    void log_callback_function(const liq_attr*, const char *message, void *user_info) {}
+
+----
+
+    void liq_set_log_flush_callback(liq_attr*, liq_log_flush_callback_function*, void *user_info);
+<p>
+
+    void log_flush_callback_function(const liq_attr*, void *user_info) {}
+
+Sets up callback function to be called when the library reports work progress or errors. The callback must not call any library functions.
+
+`user_info` value will be passed to the callback.
+
+`NULL` callback clears the current callback.
+
+In the log callback the `message` is a zero-terminated string containing informative message to output. It is valid only until the callback returns.
+
+`liq_set_log_flush_callback()` sets up callback function that will be called after the last log callback, which can be used to flush buffers and free resources used by the log callback.
+
+----
+
+    liq_attr* liq_attr_create_with_allocator(void* (*malloc)(size_t), void (*free)(void*));
+
+Same as `liq_attr_create`, but uses given `malloc` and `free` replacements to allocate all memory used by the library.
+
+The `malloc` function must return 16-byte aligned memory on x86 (and on other architectures memory aligned for `double` and pointers). Conversely, if your stdlib's `malloc` doesn't return appropriately aligned memory, you should use this function to provide aligned replacements.
+
+----
+
+    liq_attr* liq_attr_copy(liq_attr *orig);
+
+Creates an independent copy of `liq_attr`. The copy should also be freed using `liq_attr_destroy()`.
+
+---
+
+    liq_error liq_set_output_gamma(liq_result* res, double gamma);
+
+Sets gamma correction for generated palette and remapped image. Must be > 0 and < 1, e.g. `0.45455` for gamma 1/2.2 in PNG images. By default output gamma is same as gamma of the input image.
+
+----
+
+    int liq_image_get_width(const liq_image *img);
+    int liq_image_get_height(const liq_image *img);
+    double liq_get_output_gamma(const liq_result *result);
+
+Getters for `width`, `height` and `gamma` of the input image.
+
+If the input is invalid, these all return -1.
+
+---
+
+    liq_error liq_image_add_fixed_color(liq_image* img, liq_color color);
+
+Reserves a color in the output palette created from this image. It behaves as if the given color was used in the image and was very important.
+
+RGB values of `liq_color` are assumed to have the same gamma as the image.
+
+It must be called before the image is quantized.
+
+Returns error if more than 256 colors are added. If image is quantized to fewer colors than the number of fixed colors added, then excess fixed colors will be ignored.
+
+---
+
+    int liq_version();
+
+Returns version of the library as an integer. Same as `LIQ_VERSION`. Human-readable version is defined as `LIQ_VERSION_STRING`.
+
+## Multithreading
+
+The library is stateless and doesn't use any global or thread-local storage. It doesn't use any locks.
+
+* Different threads can perform unrelated quantizations/remappings at the same time (e.g. each thread working on a different image).
+* The same `liq_attr`, `liq_result`, etc. can be accessed from different threads, but not at the same time (e.g. you can create `liq_attr` in one thread and free it in another).
+
+The library needs to sort unique colors present in the image. Although the sorting algorithm does few things to make stack usage minimal in typical cases, there is no guarantee against extremely degenerate cases, so threads should have automatically growing stack.
+
+### OpenMP
+
+The library will parallelize some operations if compiled with OpenMP.
+
+You must not increase number of maximum threads after `liq_image` has been created, as it allocates some per-thread buffers.
+
+Callback of `liq_image_create_custom()` may be called from different threads at the same time.
+
+## Acknowledgements
+
+Thanks to Irfan Skiljan for helping test the first version of the library.
+
+The library is developed by [Kornel Lesiński](mailto:%20kornel@pngquant.org).

deps/pngquant/blur.c

@@ -0,0 +1,112 @@
+
+#include "libimagequant.h"
+#include "pam.h"
+#include "blur.h"
+
+/*
+ Blurs image horizontally (width 2*size+1) and writes it transposed to dst (called twice gives 2d blur)
+ */
+static void transposing_1d_blur(unsigned char *restrict src, unsigned char *restrict dst, unsigned int width, unsigned int height, const unsigned int size)
+{
+    for(unsigned int j=0; j < height; j++) {
+        unsigned char *restrict row = src + j*width;
+
+        // accumulate sum for pixels outside line
+        unsigned int sum;
+        sum = row[0]*size;
+        for(unsigned int i=0; i < size; i++) {
+            sum += row[i];
+        }
+
+        // blur with left side outside line
+        for(unsigned int i=0; i < size; i++) {
+            sum -= row[0];
+            sum += row[i+size];
+
+            dst[i*height + j] = sum / (size*2);
+        }
+
+        for(unsigned int i=size; i < width-size; i++) {
+            sum -= row[i-size];
+            sum += row[i+size];
+
+            dst[i*height + j] = sum / (size*2);
+        }
+
+        // blur with right side outside line
+        for(unsigned int i=width-size; i < width; i++) {
+            sum -= row[i-size];
+            sum += row[width-1];
+
+            dst[i*height + j] = sum / (size*2);
+        }
+    }
+}
+
+/**
+ * Picks maximum of neighboring pixels (blur + lighten)
+ */
+LIQ_PRIVATE void liq_max3(unsigned char *src, unsigned char *dst, unsigned int width, unsigned int height)
+{
+    for(unsigned int j=0; j < height; j++) {
+        const unsigned char *row = src + j*width,
+        *prevrow = src + (j > 1 ? j-1 : 0)*width,
+        *nextrow = src + MIN(height-1,j+1)*width;
+
+        unsigned char prev,curr=row[0],next=row[0];
+
+        for(unsigned int i=0; i < width-1; i++) {
+            prev=curr;
+            curr=next;
+            next=row[i+1];
+
+            unsigned char t1 = MAX(prev,next);
+            unsigned char t2 = MAX(nextrow[i],prevrow[i]);
+            *dst++ = MAX(curr,MAX(t1,t2));
+        }
+        unsigned char t1 = MAX(curr,next);
+        unsigned char t2 = MAX(nextrow[width-1],prevrow[width-1]);
+        *dst++ = MAX(t1,t2);
+    }
+}
+
+/**
+ * Picks minimum of neighboring pixels (blur + darken)
+ */
+LIQ_PRIVATE void liq_min3(unsigned char *src, unsigned char *dst, unsigned int width, unsigned int height)
+{
+    for(unsigned int j=0; j < height; j++) {
+        const unsigned char *row = src + j*width,
+        *prevrow = src + (j > 1 ? j-1 : 0)*width,
+        *nextrow = src + MIN(height-1,j+1)*width;
+
+        unsigned char prev,curr=row[0],next=row[0];
+
+        for(unsigned int i=0; i < width-1; i++) {
+            prev=curr;
+            curr=next;
+            next=row[i+1];
+
+            unsigned char t1 = MIN(prev,next);
+            unsigned char t2 = MIN(nextrow[i],prevrow[i]);
+            *dst++ = MIN(curr,MIN(t1,t2));
+        }
+        unsigned char t1 = MIN(curr,next);
+        unsigned char t2 = MIN(nextrow[width-1],prevrow[width-1]);
+        *dst++ = MIN(t1,t2);
+    }
+}
+
+/*
+ Filters src image and saves it to dst, overwriting tmp in the process.
+ Image must be width*height pixels high. Size controls radius of box blur.
+ */
+LIQ_PRIVATE void liq_blur(unsigned char *src, unsigned char *tmp, unsigned char *dst, unsigned int width, unsigned int height, unsigned int size)
+{
+    assert(size > 0);
+    if (width < 2*size+1 || height < 2*size+1) {
+        return;
+    }
+    transposing_1d_blur(src, tmp, width, height, size);
+    transposing_1d_blur(tmp, dst, height, width, size);
+}

deps/pngquant/blur.h

@@ -0,0 +1,4 @@
+
+LIQ_PRIVATE void liq_blur(unsigned char *src, unsigned char *tmp, unsigned char *dst, unsigned int width, unsigned int height, unsigned int size);
+LIQ_PRIVATE void liq_max3(unsigned char *src, unsigned char *dst, unsigned int width, unsigned int height);
+LIQ_PRIVATE void liq_min3(unsigned char *src, unsigned char *dst, unsigned int width, unsigned int height);

deps/pngquant/libimagequant.h

@@ -0,0 +1,110 @@
+/*
+ * http://pngquant.org
+ */
+
+#ifndef LIBIMAGEQUANT_H
+#define LIBIMAGEQUANT_H
+
+#ifndef LIQ_EXPORT
+#define LIQ_EXPORT extern
+#endif
+
+#define LIQ_VERSION 20401
+#define LIQ_VERSION_STRING "2.4.1"
+
+#ifndef LIQ_PRIVATE
+#if defined(__GNUC__) || defined (__llvm__)
+#define LIQ_PRIVATE __attribute__((visibility("hidden")))
+#else
+#define LIQ_PRIVATE
+#endif
+#endif
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+#include <stddef.h>
+
+typedef struct liq_attr liq_attr;
+typedef struct liq_image liq_image;
+typedef struct liq_result liq_result;
+
+typedef struct liq_color {
+    unsigned char r, g, b, a;
+} liq_color;
+
+typedef struct liq_palette {
+    unsigned int count;
+    liq_color entries[256];
+} liq_palette;
+
+typedef enum liq_error {
+    LIQ_OK = 0,
+    LIQ_QUALITY_TOO_LOW = 99,
+    LIQ_VALUE_OUT_OF_RANGE = 100,
+    LIQ_OUT_OF_MEMORY,
+    LIQ_NOT_READY,
+    LIQ_BITMAP_NOT_AVAILABLE,
+    LIQ_BUFFER_TOO_SMALL,
+    LIQ_INVALID_POINTER,
+} liq_error;
+
+enum liq_ownership {LIQ_OWN_ROWS=4, LIQ_OWN_PIXELS=8};
+
+LIQ_EXPORT liq_attr* liq_attr_create(void);
+LIQ_EXPORT liq_attr* liq_attr_create_with_allocator(void* (*malloc)(size_t), void (*free)(void*));
+LIQ_EXPORT liq_attr* liq_attr_copy(liq_attr *orig);
+LIQ_EXPORT void liq_attr_destroy(liq_attr *attr);
+
+LIQ_EXPORT liq_error liq_set_max_colors(liq_attr* attr, int colors);
+LIQ_EXPORT int liq_get_max_colors(const liq_attr* attr);
+LIQ_EXPORT liq_error liq_set_speed(liq_attr* attr, int speed);
+LIQ_EXPORT int liq_get_speed(const liq_attr* attr);
+LIQ_EXPORT liq_error liq_set_min_opacity(liq_attr* attr, int min);
+LIQ_EXPORT int liq_get_min_opacity(const liq_attr* attr);
+LIQ_EXPORT liq_error liq_set_min_posterization(liq_attr* attr, int bits);
+LIQ_EXPORT int liq_get_min_posterization(const liq_attr* attr);
+LIQ_EXPORT liq_error liq_set_quality(liq_attr* attr, int minimum, int maximum);
+LIQ_EXPORT int liq_get_min_quality(const liq_attr* attr);
+LIQ_EXPORT int liq_get_max_quality(const liq_attr* attr);
+LIQ_EXPORT void liq_set_last_index_transparent(liq_attr* attr, int is_last);
+
+typedef void liq_log_callback_function(const liq_attr*, const char *message, void* user_info);
+typedef void liq_log_flush_callback_function(const liq_attr*, void* user_info);
+LIQ_EXPORT void liq_set_log_callback(liq_attr*, liq_log_callback_function*, void* user_info);
+LIQ_EXPORT void liq_set_log_flush_callback(liq_attr*, liq_log_flush_callback_function*, void* user_info);
+
+LIQ_EXPORT liq_image *liq_image_create_rgba_rows(liq_attr *attr, void* rows[], int width, int height, double gamma);
+LIQ_EXPORT liq_image *liq_image_create_rgba(liq_attr *attr, void* bitmap, int width, int height, double gamma);
+
+typedef void liq_image_get_rgba_row_callback(liq_color row_out[], int row, int width, void* user_info);
+LIQ_EXPORT liq_image *liq_image_create_custom(liq_attr *attr, liq_image_get_rgba_row_callback *row_callback, void* user_info, int width, int height, double gamma);
+
+LIQ_EXPORT liq_error liq_image_set_memory_ownership(liq_image *image, int ownership_flags);
+LIQ_EXPORT liq_error liq_image_add_fixed_color(liq_image *img, liq_color color);
+LIQ_EXPORT int liq_image_get_width(const liq_image *img);
+LIQ_EXPORT int liq_image_get_height(const liq_image *img);
+LIQ_EXPORT void liq_image_destroy(liq_image *img);
+
+LIQ_EXPORT liq_result *liq_quantize_image(liq_attr *options, liq_image *input_image);
+
+LIQ_EXPORT liq_error liq_set_dithering_level(liq_result *res, float dither_level);
+LIQ_EXPORT liq_error liq_set_output_gamma(liq_result* res, double gamma);
+LIQ_EXPORT double liq_get_output_gamma(const liq_result *result);
+
+LIQ_EXPORT const liq_palette *liq_get_palette(liq_result *result);
+
+LIQ_EXPORT liq_error liq_write_remapped_image(liq_result *result, liq_image *input_image, void *buffer, size_t buffer_size);
+LIQ_EXPORT liq_error liq_write_remapped_image_rows(liq_result *result, liq_image *input_image, unsigned char **row_pointers);
+
+LIQ_EXPORT double liq_get_quantization_error(liq_result *result);
+LIQ_EXPORT int liq_get_quantization_quality(liq_result *result);
+
+LIQ_EXPORT void liq_result_destroy(liq_result *);
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif

deps/pngquant/mediancut.c

@@ -0,0 +1,507 @@
+/*
+** Copyright (C) 1989, 1991 by Jef Poskanzer.
+** Copyright (C) 1997, 2000, 2002 by Greg Roelofs; based on an idea by
+**                                Stefan Schneider.
+** © 2009-2013 by Kornel Lesinski.
+**
+** Permission to use, copy, modify, and distribute this software and its
+** documentation for any purpose and without fee is hereby granted, provided
+** that the above copyright notice appear in all copies and that both that
+** copyright notice and this permission notice appear in supporting
+** documentation.  This software is provided "as is" without express or
+** implied warranty.
+*/
+
+#include <stdlib.h>
+#include <stddef.h>
+
+#include "libimagequant.h"
+#include "pam.h"
+#include "mediancut.h"
+
+#define index_of_channel(ch) (offsetof(f_pixel,ch)/sizeof(float))
+
+static f_pixel averagepixels(unsigned int clrs, const hist_item achv[], float min_opaque_val, const f_pixel center);
+
+struct box {
+    f_pixel color;
+    f_pixel variance;
+    double sum, total_error, max_error;
+    unsigned int ind;
+    unsigned int colors;
+};
+
+ALWAYS_INLINE static double variance_diff(double val, const double good_enough);
+inline static double variance_diff(double val, const double good_enough)
+{
+    val *= val;
+    if (val < good_enough*good_enough) return val*0.25;
+    return val;
+}
+
+/** Weighted per-channel variance of the box. It's used to decide which channel to split by */
+static f_pixel box_variance(const hist_item achv[], const struct box *box)
+{
+    f_pixel mean = box->color;
+    double variancea=0, variancer=0, varianceg=0, varianceb=0;
+
+    for(unsigned int i = 0; i < box->colors; ++i) {
+        f_pixel px = achv[box->ind + i].acolor;
+        double weight = achv[box->ind + i].adjusted_weight;
+        variancea += variance_diff(mean.a - px.a, 2.0/256.0)*weight;
+        variancer += variance_diff(mean.r - px.r, 1.0/256.0)*weight;
+        varianceg += variance_diff(mean.g - px.g, 1.0/256.0)*weight;
+        varianceb += variance_diff(mean.b - px.b, 1.0/256.0)*weight;
+    }
+
+    return (f_pixel){
+        .a = variancea*(4.0/16.0),
+        .r = variancer*(7.0/16.0),
+        .g = varianceg*(9.0/16.0),
+        .b = varianceb*(5.0/16.0),
+    };
+}
+
+static double box_max_error(const hist_item achv[], const struct box *box)
+{
+    f_pixel mean = box->color;
+    double max_error = 0;
+
+    for(unsigned int i = 0; i < box->colors; ++i) {
+        const double diff = colordifference(mean, achv[box->ind + i].acolor);
+        if (diff > max_error) {
+            max_error = diff;
+        }
+    }
+    return max_error;
+}
+
+ALWAYS_INLINE static double color_weight(f_pixel median, hist_item h);
+
+static inline void hist_item_swap(hist_item *l, hist_item *r)
+{
+    if (l != r) {
+        hist_item t = *l;
+        *l = *r;
+        *r = t;
+    }
+}
+
+ALWAYS_INLINE static unsigned int qsort_pivot(const hist_item *const base, const unsigned int len);
+inline static unsigned int qsort_pivot(const hist_item *const base, const unsigned int len)
+{
+    if (len < 32) {
+        return len/2;
+    }
+
+    const unsigned int aidx=8, bidx=len/2, cidx=len-1;
+    const unsigned int a=base[aidx].tmp.sort_value, b=base[bidx].tmp.sort_value, c=base[cidx].tmp.sort_value;
+    return (a < b) ? ((b < c) ? bidx : ((a < c) ? cidx : aidx ))
+                   : ((b > c) ? bidx : ((a < c) ? aidx : cidx ));
+}
+
+ALWAYS_INLINE static unsigned int qsort_partition(hist_item *const base, const unsigned int len);
+inline static unsigned int qsort_partition(hist_item *const base, const unsigned int len)
+{
+    unsigned int l = 1, r = len;
+    if (len >= 8) {
+        hist_item_swap(&base[0], &base[qsort_pivot(base,len)]);
+    }
+
+    const unsigned int pivot_value = base[0].tmp.sort_value;
+    while (l < r) {
+        if (base[l].tmp.sort_value >= pivot_value) {
+            l++;
+        } else {
+            while(l < --r && base[r].tmp.sort_value <= pivot_value) {}
+            hist_item_swap(&base[l], &base[r]);
+        }
+    }
+    l--;
+    hist_item_swap(&base[0], &base[l]);
+
+    return l;
+}
+
+/** quick select algorithm */
+static void hist_item_sort_range(hist_item *base, unsigned int len, unsigned int sort_start)
+{
+    for(;;) {
+        const unsigned int l = qsort_partition(base, len), r = l+1;
+
+        if (l > 0 && sort_start < l) {
+            len = l;
+        }
+        else if (r < len && sort_start > r) {
+            base += r; len -= r; sort_start -= r;
+        }
+        else break;
+    }
+}
+
+/** sorts array to make sum of weights lower than halfvar one side, returns edge between <halfvar and >halfvar parts of the set */
+static hist_item *hist_item_sort_halfvar(hist_item *base, unsigned int len, double *const lowervar, const double halfvar)
+{
+    do {
+        const unsigned int l = qsort_partition(base, len), r = l+1;
+
+        // check if sum of left side is smaller than half,
+        // if it is, then it doesn't need to be sorted
+        unsigned int t = 0; double tmpsum = *lowervar;
+        while (t <= l && tmpsum < halfvar) tmpsum += base[t++].color_weight;
+
+        if (tmpsum < halfvar) {
+            *lowervar = tmpsum;
+        } else {
+            if (l > 0) {
+                hist_item *res = hist_item_sort_halfvar(base, l, lowervar, halfvar);
+                if (res) return res;
+            } else {
+                // End of left recursion. This will be executed in order from the first element.
+                *lowervar += base[0].color_weight;
+                if (*lowervar > halfvar) return &base[0];
+            }
+        }
+
+        if (len > r) {
+            base += r; len -= r; // tail-recursive "call"
+        } else {
+            *lowervar += base[r].color_weight;
+            return (*lowervar > halfvar) ? &base[r] : NULL;
+        }
+    } while(1);
+}
+
+static f_pixel get_median(const struct box *b, hist_item achv[]);
+
+typedef struct {
+    unsigned int chan; float variance;
+} channelvariance;
+
+static int comparevariance(const void *ch1, const void *ch2)
+{
+    return ((const channelvariance*)ch1)->variance > ((const channelvariance*)ch2)->variance ? -1 :
+          (((const channelvariance*)ch1)->variance < ((const channelvariance*)ch2)->variance ? 1 : 0);
+}
+
+/** Finds which channels need to be sorted first and preproceses achv for fast sort */
+static double prepare_sort(struct box *b, hist_item achv[])
+{
+    /*
+     ** Sort dimensions by their variance, and then sort colors first by dimension with highest variance
+     */
+    channelvariance channels[4] = {
+        {index_of_channel(r), b->variance.r},
+        {index_of_channel(g), b->variance.g},
+        {index_of_channel(b), b->variance.b},
+        {index_of_channel(a), b->variance.a},
+    };
+
+    qsort(channels, 4, sizeof(channels[0]), comparevariance);
+
+    for(unsigned int i=0; i < b->colors; i++) {
+        const float *chans = (const float *)&achv[b->ind + i].acolor;
+        // Only the first channel really matters. When trying median cut many times
+        // with different histogram weights, I don't want sort randomness to influence outcome.
+        achv[b->ind + i].tmp.sort_value = ((unsigned int)(chans[channels[0].chan]*65535.0)<<16) |
+                                       (unsigned int)((chans[channels[2].chan] + chans[channels[1].chan]/2.0 + chans[channels[3].chan]/4.0)*65535.0);
+    }
+
+    const f_pixel median = get_median(b, achv);
+
+    // box will be split to make color_weight of each side even
+    const unsigned int ind = b->ind, end = ind+b->colors;
+    double totalvar = 0;
+    for(unsigned int j=ind; j < end; j++) totalvar += (achv[j].color_weight = color_weight(median, achv[j]));
+    return totalvar / 2.0;
+}
+
+/** finds median in unsorted set by sorting only minimum required */
+static f_pixel get_median(const struct box *b, hist_item achv[])
+{
+    const unsigned int median_start = (b->colors-1)/2;
+
+    hist_item_sort_range(&(achv[b->ind]), b->colors,
+                         median_start);
+
+    if (b->colors&1) return achv[b->ind + median_start].acolor;
+
+    // technically the second color is not guaranteed to be sorted correctly
+    // but most of the time it is good enough to be useful
+    return averagepixels(2, &achv[b->ind + median_start], 1.0, (f_pixel){0.5,0.5,0.5,0.5});
+}
+
+/*
+ ** Find the best splittable box. -1 if no boxes are splittable.
+ */
+static int best_splittable_box(struct box* bv, unsigned int boxes, const double max_mse)
+{
+    int bi=-1; double maxsum=0;
+    for(unsigned int i=0; i < boxes; i++) {
+        if (bv[i].colors < 2) {
+            continue;
+        }
+
+        // looks only at max variance, because it's only going to split by it
+        const double cv = MAX(bv[i].variance.r, MAX(bv[i].variance.g,bv[i].variance.b));
+        double thissum = bv[i].sum * MAX(bv[i].variance.a, cv);
+
+        if (bv[i].max_error > max_mse) {
+            thissum = thissum* bv[i].max_error/max_mse;
+        }
+
+        if (thissum > maxsum) {
+            maxsum = thissum;
+            bi = i;
+        }
+    }
+    return bi;
+}
+
+inline static double color_weight(f_pixel median, hist_item h)
+{
+    float diff = colordifference(median, h.acolor);
+    // if color is "good enough", don't split further
+    if (diff < 2.f/256.f/256.f) diff /= 2.f;
+    return sqrt(diff) * (sqrt(1.0+h.adjusted_weight)-1.0);
+}
+
+static void set_colormap_from_boxes(colormap *map, struct box* bv, unsigned int boxes, hist_item *achv);
+static void adjust_histogram(hist_item *achv, const colormap *map, const struct box* bv, unsigned int boxes);
+
+double box_error(const struct box *box, const hist_item achv[])
+{
+    f_pixel avg = box->color;
+
+    double total_error=0;
+    for (unsigned int i = 0; i < box->colors; ++i) {
+        total_error += colordifference(avg, achv[box->ind + i].acolor) * achv[box->ind + i].perceptual_weight;
+    }
+
+    return total_error;
+}
+
+
+static bool total_box_error_below_target(double target_mse, struct box bv[], unsigned int boxes, const histogram *hist)
+{
+    target_mse *= hist->total_perceptual_weight;
+    double total_error=0;
+
+    for(unsigned int i=0; i < boxes; i++) {
+        // error is (re)calculated lazily
+        if (bv[i].total_error >= 0) {
+            total_error += bv[i].total_error;
+        }
+        if (total_error > target_mse) return false;
+    }
+
+    for(unsigned int i=0; i < boxes; i++) {
+        if (bv[i].total_error < 0) {
+            bv[i].total_error = box_error(&bv[i], hist->achv);
+            total_error += bv[i].total_error;
+        }
+        if (total_error > target_mse) return false;
+    }
+
+    return true;
+}
+
+/*
+ ** Here is the fun part, the median-cut colormap generator.  This is based
+ ** on Paul Heckbert's paper, "Color Image Quantization for Frame Buffer
+ ** Display," SIGGRAPH 1982 Proceedings, page 297.
+ */
+LIQ_PRIVATE colormap *mediancut(histogram *hist, const float min_opaque_val, unsigned int newcolors, const double target_mse, const double max_mse, void* (*malloc)(size_t), void (*free)(void*))
+{
+    hist_item *achv = hist->achv;
+    struct box bv[newcolors];
+
+    /*
+     ** Set up the initial box.
+     */
+    bv[0].ind = 0;
+    bv[0].colors = hist->size;
+    bv[0].color = averagepixels(bv[0].colors, &achv[bv[0].ind], min_opaque_val, (f_pixel){0.5,0.5,0.5,0.5});
+    bv[0].variance = box_variance(achv, &bv[0]);
+    bv[0].max_error = box_max_error(achv, &bv[0]);
+    bv[0].sum = 0;
+    bv[0].total_error = -1;
+    for(unsigned int i=0; i < bv[0].colors; i++) bv[0].sum += achv[i].adjusted_weight;
+
+    unsigned int boxes = 1;
+
+    // remember smaller palette for fast searching
+    colormap *representative_subset = NULL;
+    unsigned int subset_size = ceilf(powf(newcolors,0.7f));
+
+    /*
+     ** Main loop: split boxes until we have enough.
+     */
+    while (boxes < newcolors) {
+
+        if (boxes == subset_size) {
+            representative_subset = pam_colormap(boxes, malloc, free);
+            set_colormap_from_boxes(representative_subset, bv, boxes, achv);
+        }
+
+        // first splits boxes that exceed quality limit (to have colors for things like odd green pixel),
+        // later raises the limit to allow large smooth areas/gradients get colors.
+        const double current_max_mse = max_mse + (boxes/(double)newcolors)*16.0*max_mse;
+        const int bi = best_splittable_box(bv, boxes, current_max_mse);
+        if (bi < 0)
+            break;        /* ran out of colors! */
+
+        unsigned int indx = bv[bi].ind;
+        unsigned int clrs = bv[bi].colors;
+
+        /*
+         Classic implementation tries to get even number of colors or pixels in each subdivision.
+
+         Here, instead of popularity I use (sqrt(popularity)*variance) metric.
+         Each subdivision balances number of pixels (popular colors) and low variance -
+         boxes can be large if they have similar colors. Later boxes with high variance
+         will be more likely to be split.
+
+         Median used as expected value gives much better results than mean.
+         */
+
+        const double halfvar = prepare_sort(&bv[bi], achv);
+        double lowervar=0;
+
+        // hist_item_sort_halfvar sorts and sums lowervar at the same time
+        // returns item to break at …minus one, which does smell like an off-by-one error.
+        hist_item *break_p = hist_item_sort_halfvar(&achv[indx], clrs, &lowervar, halfvar);
+        unsigned int break_at = MIN(clrs-1, break_p - &achv[indx] + 1);
+
+        /*
+         ** Split the box.
+         */
+        double sm = bv[bi].sum;
+        double lowersum = 0;
+        for(unsigned int i=0; i < break_at; i++) lowersum += achv[indx + i].adjusted_weight;
+
+        const f_pixel previous_center = bv[bi].color;
+        bv[bi].colors = break_at;
+        bv[bi].sum = lowersum;
+        bv[bi].color = averagepixels(bv[bi].colors, &achv[bv[bi].ind], min_opaque_val, previous_center);
+        bv[bi].total_error = -1;
+        bv[bi].variance = box_variance(achv, &bv[bi]);
+        bv[bi].max_error = box_max_error(achv, &bv[bi]);
+        bv[boxes].ind = indx + break_at;
+        bv[boxes].colors = clrs - break_at;
+        bv[boxes].sum = sm - lowersum;
+        bv[boxes].color = averagepixels(bv[boxes].colors, &achv[bv[boxes].ind], min_opaque_val, previous_center);
+        bv[boxes].total_error = -1;
+        bv[boxes].variance = box_variance(achv, &bv[boxes]);
+        bv[boxes].max_error = box_max_error(achv, &bv[boxes]);
+
+        ++boxes;
+
+        if (total_box_error_below_target(target_mse, bv, boxes, hist)) {
+            break;
+        }
+    }
+
+    colormap *map = pam_colormap(boxes, malloc, free);
+    set_colormap_from_boxes(map, bv, boxes, achv);
+
+    map->subset_palette = representative_subset;
+    adjust_histogram(achv, map, bv, boxes);
+
+    return map;
+}
+
+static void set_colormap_from_boxes(colormap *map, struct box* bv, unsigned int boxes, hist_item *achv)
+{
+    /*
+     ** Ok, we've got enough boxes.  Now choose a representative color for
+     ** each box.  There are a number of possible ways to make this choice.
+     ** One would be to choose the center of the box; this ignores any structure
+     ** within the boxes.  Another method would be to average all the colors in
+     ** the box - this is the method specified in Heckbert's paper.
+     */
+
+    for(unsigned int bi = 0; bi < boxes; ++bi) {
+        map->palette[bi].acolor = bv[bi].color;
+
+        /* store total color popularity (perceptual_weight is approximation of it) */
+        map->palette[bi].popularity = 0;
+        for(unsigned int i=bv[bi].ind; i < bv[bi].ind+bv[bi].colors; i++) {
+            map->palette[bi].popularity += achv[i].perceptual_weight;
+        }
+    }
+}
+
+/* increase histogram popularity by difference from the final color (this is used as part of feedback loop) */
+static void adjust_histogram(hist_item *achv, const colormap *map, const struct box* bv, unsigned int boxes)
+{
+    for(unsigned int bi = 0; bi < boxes; ++bi) {
+        for(unsigned int i=bv[bi].ind; i < bv[bi].ind+bv[bi].colors; i++) {
+            achv[i].adjusted_weight *= sqrt(1.0 +colordifference(map->palette[bi].acolor, achv[i].acolor)/4.0);
+            achv[i].tmp.likely_colormap_index = bi;
+        }
+    }
+}
+
+static f_pixel averagepixels(unsigned int clrs, const hist_item achv[], const float min_opaque_val, const f_pixel center)
+{
+    double r = 0, g = 0, b = 0, a = 0, new_a=0, sum = 0;
+    float maxa = 0;
+
+    // first find final opacity in order to blend colors at that opacity
+    for(unsigned int i = 0; i < clrs; ++i) {
+        const f_pixel px = achv[i].acolor;
+        new_a += px.a * achv[i].adjusted_weight;
+        sum += achv[i].adjusted_weight;
+
+        /* find if there are opaque colors, in case we're supposed to preserve opacity exactly (ie_bug) */
+        if (px.a > maxa) maxa = px.a;
+    }
+
+    if (sum) new_a /= sum;
+
+    /** if there was at least one completely opaque color, "round" final color to opaque */
+    if (new_a >= min_opaque_val && maxa >= (255.0/256.0)) new_a = 1;
+
+    sum=0;
+    // reverse iteration for cache locality with previous loop
+    for(int i = clrs-1; i >= 0; i--) {
+        double tmp, weight = 1.0f;
+        f_pixel px = achv[i].acolor;
+
+        /* give more weight to colors that are further away from average
+         this is intended to prevent desaturation of images and fading of whites
+         */
+        tmp = (center.r - px.r);
+        weight += tmp*tmp;
+        tmp = (center.g - px.g);
+        weight += tmp*tmp;
+        tmp = (center.b - px.b);
+        weight += tmp*tmp;
+
+        weight *= achv[i].adjusted_weight;
+        sum += weight;
+
+        if (px.a) {
+            px.r /= px.a;
+            px.g /= px.a;
+            px.b /= px.a;
+        }
+
+        r += px.r * new_a * weight;
+        g += px.g * new_a * weight;
+        b += px.b * new_a * weight;
+        a += new_a * weight;
+    }
+
+    if (sum) {
+    a /= sum;
+    r /= sum;
+    g /= sum;
+    b /= sum;
+    }
+
+    assert(!isnan(r) && !isnan(g) && !isnan(b) && !isnan(a));
+
+    return (f_pixel){.r=r, .g=g, .b=b, .a=a};
+}

deps/pngquant/mediancut.h

@@ -0,0 +1,2 @@
+
+LIQ_PRIVATE colormap *mediancut(histogram *hist, const float min_opaque_val, unsigned int newcolors, const double target_mse, const double max_mse, void* (*malloc)(size_t), void (*free)(void*));

deps/pngquant/mempool.c

@@ -0,0 +1,63 @@
+
+#include "libimagequant.h"
+#include "mempool.h"
+#include <stdlib.h>
+#include <stdint.h>
+#include <assert.h>
+
+#define ALIGN_MASK 15UL
+#define MEMPOOL_RESERVED ((sizeof(struct mempool)+ALIGN_MASK) & ~ALIGN_MASK)
+
+struct mempool {
+    unsigned int used, size;
+    void* (*malloc)(size_t);
+    void (*free)(void*);
+    struct mempool *next;
+};
+LIQ_PRIVATE void* mempool_create(mempool *mptr, const unsigned int size, unsigned int max_size, void* (*malloc)(size_t), void (*free)(void*))
+{
+    if (*mptr && ((*mptr)->used+size) <= (*mptr)->size) {
+        unsigned int prevused = (*mptr)->used;
+        (*mptr)->used += (size+15UL) & ~0xFUL;
+        return ((char*)(*mptr)) + prevused;
+    }
+
+    mempool old = *mptr;
+    if (!max_size) max_size = (1<<17);
+    max_size = size+ALIGN_MASK > max_size ? size+ALIGN_MASK : max_size;
+
+    *mptr = malloc(MEMPOOL_RESERVED + max_size);
+    if (!*mptr) return NULL;
+    **mptr = (struct mempool){
+        .malloc = malloc,
+        .free = free,
+        .size = MEMPOOL_RESERVED + max_size,
+        .used = sizeof(struct mempool),
+        .next = old,
+    };
+    uintptr_t mptr_used_start = (uintptr_t)(*mptr) + (*mptr)->used;
+    (*mptr)->used += (ALIGN_MASK + 1 - (mptr_used_start & ALIGN_MASK)) & ALIGN_MASK; // reserve bytes required to make subsequent allocations aligned
+    assert(!(((uintptr_t)(*mptr) + (*mptr)->used) & ALIGN_MASK));
+
+    return mempool_alloc(mptr, size, size);
+}
+
+LIQ_PRIVATE void* mempool_alloc(mempool *mptr, unsigned int size, unsigned int max_size)
+{
+    if (((*mptr)->used+size) <= (*mptr)->size) {
+        unsigned int prevused = (*mptr)->used;
+        (*mptr)->used += (size + ALIGN_MASK) & ~ALIGN_MASK;
+        return ((char*)(*mptr)) + prevused;
+    }
+
+    return mempool_create(mptr, size, max_size, (*mptr)->malloc, (*mptr)->free);
+}
+
+LIQ_PRIVATE void mempool_destroy(mempool m)
+{
+    while (m) {
+        mempool next = m->next;
+        m->free(m);
+        m = next;
+    }
+}

deps/pngquant/mempool.h

@@ -0,0 +1,13 @@
+#ifndef MEMPOOL_H
+#define MEMPOOL_H
+
+#include <stddef.h>
+
+struct mempool;
+typedef struct mempool *mempool;
+
+LIQ_PRIVATE void* mempool_create(mempool *mptr, unsigned int size, unsigned int capacity, void* (*malloc)(size_t), void (*free)(void*));
+LIQ_PRIVATE void* mempool_alloc(mempool *mptr, unsigned int size, unsigned int capacity);
+LIQ_PRIVATE void mempool_destroy(mempool m);
+
+#endif

deps/pngquant/nearest.c

@@ -0,0 +1,207 @@
+
+#include "libimagequant.h"
+#include "pam.h"
+#include "nearest.h"
+#include "mempool.h"
+#include <stdlib.h>
+
+struct sorttmp {
+    float radius;
+    unsigned int index;
+};
+
+struct head {
+    // colors less than radius away from vantage_point color will have best match in candidates
+    f_pixel vantage_point;
+    float radius;
+    unsigned int num_candidates;
+    f_pixel *candidates_color;
+    unsigned short *candidates_index;
+};
+
+struct nearest_map {
+    const colormap *map;
+    float nearest_other_color_dist[256];
+    mempool mempool;
+    struct head heads[];
+};
+
+static float distance_from_nearest_other_color(const colormap *map, const unsigned int i)
+{
+    float second_best=MAX_DIFF;
+    for(unsigned int j=0; j < map->colors; j++) {
+        if (i == j) continue;
+        float diff = colordifference(map->palette[i].acolor, map->palette[j].acolor);
+        if (diff <= second_best) {
+            second_best = diff;
+        }
+    }
+    return second_best;
+}
+
+static int compareradius(const void *ap, const void *bp)
+{
+    float a = ((const struct sorttmp*)ap)->radius;
+    float b = ((const struct sorttmp*)bp)->radius;
+    return a > b ? 1 : (a < b ? -1 : 0);
+}
+
+static struct head build_head(f_pixel px, const colormap *map, unsigned int num_candidates, mempool *m, float error_margin, bool skip_index[], unsigned int *skipped)
+{
+    struct sorttmp colors[map->colors];
+    unsigned int colorsused=0;
+
+    for(unsigned int i=0; i < map->colors; i++) {
+        if (skip_index[i]) continue; // colors in skip_index have been eliminated already in previous heads
+        colors[colorsused].index = i;
+        colors[colorsused].radius = colordifference(px, map->palette[i].acolor);
+        colorsused++;
+    }
+
+    qsort(&colors, colorsused, sizeof(colors[0]), compareradius);
+    assert(colorsused < 2 || colors[0].radius <= colors[1].radius); // closest first
+
+    num_candidates = MIN(colorsused, num_candidates);
+
+    struct head h = {
+        .candidates_color = mempool_alloc(m, num_candidates * sizeof(h.candidates_color[0]), 0),
+        .candidates_index = mempool_alloc(m, num_candidates * sizeof(h.candidates_index[0]), 0),
+        .vantage_point = px,
+        .num_candidates = num_candidates,
+    };
+    for(unsigned int i=0; i < num_candidates; i++) {
+        h.candidates_color[i] = map->palette[colors[i].index].acolor;
+        h.candidates_index[i] = colors[i].index;
+    }
+    // if all colors within this radius are included in candidates, then there cannot be any other better match
+    // farther away from the vantage point than half of the radius. Due to alpha channel must assume pessimistic radius.
+    h.radius = min_colordifference(px, h.candidates_color[num_candidates-1])/4.0f; // /4 = half of radius, but radius is squared
+
+    for(unsigned int i=0; i < num_candidates; i++) {
+        // divide again as that's matching certain subset within radius-limited subset
+        // - 1/256 is a tolerance for miscalculation (seems like colordifference isn't exact)
+        if (colors[i].radius < h.radius/4.f - error_margin) {
+            skip_index[colors[i].index]=true;
+            (*skipped)++;
+        }
+    }
+    return h;
+}
+
+static colormap *get_subset_palette(const colormap *map)
+{
+    if (map->subset_palette) {
+        return map->subset_palette;
+    }
+
+    unsigned int subset_size = (map->colors+3)/4;
+    colormap *subset_palette = pam_colormap(subset_size, map->malloc, map->free);
+
+    for(unsigned int i=0; i < subset_size; i++) {
+        subset_palette->palette[i] = map->palette[i];
+    }
+
+    return subset_palette;
+}
+
+LIQ_PRIVATE struct nearest_map *nearest_init(const colormap *map, bool fast)
+{
+    colormap *subset_palette = get_subset_palette(map);
+    const unsigned int num_vantage_points = map->colors > 16 ? MIN(map->colors/(fast ? 4 : 3), subset_palette->colors) : 0;
+    const unsigned long heads_size = sizeof(struct head) * (num_vantage_points+1); // +1 is fallback head
+
+    const unsigned long mempool_size = (sizeof(f_pixel) + sizeof(unsigned int)) * subset_palette->colors * map->colors/5 + (1<<14);
+    mempool m = NULL;
+    struct nearest_map *centroids = mempool_create(&m, sizeof(*centroids) + heads_size /* heads array is appended to it */, mempool_size, map->malloc, map->free);
+    centroids->mempool = m;
+
+    for(unsigned int i=0; i < map->colors; i++) {
+        const float dist = distance_from_nearest_other_color(map,i);
+        centroids->nearest_other_color_dist[i] = dist / 4.f; // half of squared distance
+    }
+
+    centroids->map = map;
+
+    unsigned int skipped=0;
+    assert(map->colors > 0);
+    bool skip_index[map->colors]; for(unsigned int j=0; j < map->colors; j++) skip_index[j]=false;
+
+    // floats and colordifference calculations are not perfect
+    const float error_margin = fast ? 0 : 8.f/256.f/256.f;
+    unsigned int h=0;
+    for(; h < num_vantage_points; h++) {
+        unsigned int num_candiadtes = 1+(map->colors - skipped)/((1+num_vantage_points-h)/2);
+
+        centroids->heads[h] = build_head(subset_palette->palette[h].acolor, map, num_candiadtes, &centroids->mempool, error_margin, skip_index, &skipped);
+        if (centroids->heads[h].num_candidates == 0) {
+            break;
+        }
+    }
+
+    // assumption that there is no better color within radius of vantage point color
+    // holds true only for colors within convex hull formed by palette colors.
+    // The fallback must contain all colors, since there are too many edge cases to cover.
+    if (!fast) for(unsigned int j=0; j < map->colors; j++) {
+        skip_index[j] = false;
+    }
+
+    centroids->heads[h] = build_head((f_pixel){0,0,0,0}, map, map->colors, &centroids->mempool, error_margin, skip_index, &skipped);
+    centroids->heads[h].radius = MAX_DIFF;
+
+    // get_subset_palette could have created a copy
+    if (subset_palette != map->subset_palette) {
+        pam_freecolormap(subset_palette);
+    }
+
+    return centroids;
+}
+
+LIQ_PRIVATE unsigned int nearest_search(const struct nearest_map *centroids, const f_pixel px, int likely_colormap_index, const float min_opaque_val, float *diff)
+{
+    const bool iebug = px.a > min_opaque_val;
+
+    const struct head *const heads = centroids->heads;
+
+    assert(likely_colormap_index < centroids->map->colors);
+    const float guess_diff = colordifference(centroids->map->palette[likely_colormap_index].acolor, px);
+    if (guess_diff < centroids->nearest_other_color_dist[likely_colormap_index]) {
+        if (diff) *diff = guess_diff;
+        return likely_colormap_index;
+    }
+
+    for(unsigned int i=0; /* last head will always be selected */ ; i++) {
+        float vantage_point_dist = colordifference(px, heads[i].vantage_point);
+
+        if (vantage_point_dist <= heads[i].radius) {
+            assert(heads[i].num_candidates);
+            unsigned int ind=0;
+            float dist = colordifference(px, heads[i].candidates_color[0]);
+
+            /* penalty for making holes in IE */
+            if (iebug && heads[i].candidates_color[0].a < 1) {
+                dist += 1.f/1024.f;
+            }
+
+            for(unsigned int j=1; j < heads[i].num_candidates; j++) {
+                float newdist = colordifference(px, heads[i].candidates_color[j]);
+
+                /* penalty for making holes in IE */
+                if (iebug && heads[i].candidates_color[j].a < 1) {
+                    newdist += 1.f/1024.f;
+                }
+
+                if (newdist < dist) {
+                    dist = newdist;
+                    ind = j;
+                }
+            }
+            if (diff) *diff = dist;
+            return heads[i].candidates_index[ind];
+        }
+    }
+}
+
+LIQ_PRIVATE void nearest_free(struct nearest_map *centroids)
+{
+    mempool_destroy(centroids->mempool);
+}

deps/pngquant/nearest.h

@@ -0,0 +1,8 @@
+//
+//  nearest.h
+//  pngquant
+//
+struct nearest_map;
+LIQ_PRIVATE struct nearest_map *nearest_init(const colormap *palette, const bool fast);
+LIQ_PRIVATE unsigned int nearest_search(const struct nearest_map *map, const f_pixel px, const int palette_index_guess, const float min_opaque, float *diff);
+LIQ_PRIVATE void nearest_free(struct nearest_map *map);

deps/pngquant/pam.c

@@ -0,0 +1,278 @@
+/* pam.c - pam (portable alpha map) utility library
+**
+** Copyright (C) 1989, 1991 by Jef Poskanzer.
+** Copyright (C) 1997, 2000, 2002 by Greg Roelofs; based on an idea by
+**                                Stefan Schneider.
+** © 2009-2013 by Kornel Lesinski.
+**
+** Permission to use, copy, modify, and distribute this software and its
+** documentation for any purpose and without fee is hereby granted, provided
+** that the above copyright notice appear in all copies and that both that
+** copyright notice and this permission notice appear in supporting
+** documentation.  This software is provided "as is" without express or
+** implied warranty.
+*/
+
+#include <stdlib.h>
+#include <string.h>
+
+#include "libimagequant.h"
+#include "pam.h"
+#include "mempool.h"
+
+LIQ_PRIVATE bool pam_computeacolorhash(struct acolorhash_table *acht, const rgba_pixel *const pixels[], unsigned int cols, unsigned int rows, const unsigned char *importance_map)
+{
+    const unsigned int maxacolors = acht->maxcolors, ignorebits = acht->ignorebits;
+    const unsigned int channel_mask = 255U>>ignorebits<<ignorebits;
+    const unsigned int channel_hmask = (255U>>ignorebits) ^ 0xFFU;
+    const unsigned int posterize_mask = channel_mask << 24 | channel_mask << 16 | channel_mask << 8 | channel_mask;
+    const unsigned int posterize_high_mask = channel_hmask << 24 | channel_hmask << 16 | channel_hmask << 8 | channel_hmask;
+    struct acolorhist_arr_head *const buckets = acht->buckets;
+
+    unsigned int colors = acht->colors;
+    const unsigned int hash_size = acht->hash_size;
+
+    const unsigned int stacksize = sizeof(acht->freestack)/sizeof(acht->freestack[0]);
+    struct acolorhist_arr_item **freestack = acht->freestack;
+    unsigned int freestackp=acht->freestackp;
+
+    /* Go through the entire image, building a hash table of colors. */
+    for(unsigned int row = 0; row < rows; ++row) {
+
+        float boost=1.0;
+        for(unsigned int col = 0; col < cols; ++col) {
+            if (importance_map) {
+                boost = 0.5f+ (double)*importance_map++/255.f;
+            }
+
+            // RGBA color is casted to long for easier hasing/comparisons
+            union rgba_as_int px = {pixels[row][col]};
+            unsigned int hash;
+            if (!px.rgba.a) {
+                // "dirty alpha" has different RGBA values that end up being the same fully transparent color
+                px.l=0; hash=0;
+            } else {
+                // mask posterizes all 4 channels in one go
+                px.l = (px.l & posterize_mask) | ((px.l & posterize_high_mask) >> (8-ignorebits));
+                // fancier hashing algorithms didn't improve much
+                hash = px.l % hash_size;
+            }
+
+            /* head of the hash function stores first 2 colors inline (achl->used = 1..2),
+               to reduce number of allocations of achl->other_items.
+             */
+            struct acolorhist_arr_head *achl = &buckets[hash];
+            if (achl->inline1.color.l == px.l && achl->used) {
+                achl->inline1.perceptual_weight += boost;
+                continue;
+            }
+            if (achl->used) {
+                if (achl->used > 1) {
+                    if (achl->inline2.color.l == px.l) {
+                        achl->inline2.perceptual_weight += boost;
+                        continue;
+                    }
+                    // other items are stored as an array (which gets reallocated if needed)
+                    struct acolorhist_arr_item *other_items = achl->other_items;
+                    unsigned int i = 0;
+                    for (; i < achl->used-2; i++) {
+                        if (other_items[i].color.l == px.l) {
+                            other_items[i].perceptual_weight += boost;
+                            goto continue_outer_loop;
+                        }
+                    }
+
+                    // the array was allocated with spare items
+                    if (i < achl->capacity) {
+                        other_items[i] = (struct acolorhist_arr_item){
+                            .color = px,
+                            .perceptual_weight = boost,
+                        };
+                        achl->used++;
+                        ++colors;
+                        continue;
+                    }
+
+                    if (++colors > maxacolors) {
+                        acht->colors = colors;
+                        acht->freestackp = freestackp;
+                        return false;
+                    }
+
+                    struct acolorhist_arr_item *new_items;
+                    unsigned int capacity;
+                    if (!other_items) { // there was no array previously, alloc "small" array
+                        capacity = 8;
+                        if (freestackp <= 0) {
+                            // estimate how many colors are going to be + headroom
+                            const int mempool_size = ((acht->rows + rows-row) * 2 * colors / (acht->rows + row + 1) + 1024) * sizeof(struct acolorhist_arr_item);
+                            new_items = mempool_alloc(&acht->mempool, sizeof(struct acolorhist_arr_item)*capacity, mempool_size);
+                        } else {
+                            // freestack stores previously freed (reallocated) arrays that can be reused
+                            // (all pesimistically assumed to be capacity = 8)
+                            new_items = freestack[--freestackp];
+                        }
+                    } else {
+                        // simply reallocs and copies array to larger capacity
+                        capacity = achl->capacity*2 + 16;
+                        if (freestackp < stacksize-1) {
+                            freestack[freestackp++] = other_items;
+                        }
+                        const int mempool_size = ((acht->rows + rows-row) * 2 * colors / (acht->rows + row + 1) + 32*capacity) * sizeof(struct acolorhist_arr_item);
+                        new_items = mempool_alloc(&acht->mempool, sizeof(struct acolorhist_arr_item)*capacity, mempool_size);
+                        if (!new_items) return false;
+                        memcpy(new_items, other_items, sizeof(other_items[0])*achl->capacity);
+                    }
+
+                    achl->other_items = new_items;
+                    achl->capacity = capacity;
+                    new_items[i] = (struct acolorhist_arr_item){
+                        .color = px,
+                        .perceptual_weight = boost,
+                    };
+                    achl->used++;
+                } else {
+                    // these are elses for first checks whether first and second inline-stored colors are used
+                    achl->inline2.color.l = px.l;
+                    achl->inline2.perceptual_weight = boost;
+                    achl->used = 2;
+                    ++colors;
+                }
+            } else {
+                achl->inline1.color.l = px.l;
+                achl->inline1.perceptual_weight = boost;
+                achl->used = 1;
+                ++colors;
+            }
+
+            continue_outer_loop:;
+        }
+
+    }
+    acht->colors = colors;
+    acht->cols = cols;
+    acht->rows += rows;
+    acht->freestackp = freestackp;
+    return true;
+}
+
+LIQ_PRIVATE struct acolorhash_table *pam_allocacolorhash(unsigned int maxcolors, unsigned int surface, unsigned int ignorebits, void* (*malloc)(size_t), void (*free)(void*))
+{
+    const unsigned int estimated_colors = MIN(maxcolors, surface/(ignorebits + (surface > 512*512 ? 5 : 4)));
+    const unsigned int hash_size = estimated_colors < 66000 ? 6673 : (estimated_colors < 200000 ? 12011 : 24019);
+
+    mempool m = NULL;
+    const unsigned int buckets_size = hash_size * sizeof(struct acolorhist_arr_head);
+    const unsigned int mempool_size = sizeof(struct acolorhash_table) + buckets_size + estimated_colors * sizeof(struct acolorhist_arr_item);
+    struct acolorhash_table *t = mempool_create(&m, sizeof(*t) + buckets_size, mempool_size, malloc, free);
+    if (!t) return NULL;
+    *t = (struct acolorhash_table){
+        .mempool = m,
+        .hash_size = hash_size,
+        .maxcolors = maxcolors,
+        .ignorebits = ignorebits,
+    };
+    memset(t->buckets, 0, hash_size * sizeof(struct acolorhist_arr_head));
+    return t;
+}
+
+#define PAM_ADD_TO_HIST(entry) { \
+    hist->achv[j].acolor = to_f(gamma_lut, entry.color.rgba); \
+    total_weight += hist->achv[j].adjusted_weight = hist->achv[j].perceptual_weight = MIN(entry.perceptual_weight, max_perceptual_weight); \
+    ++j; \
+}
+
+LIQ_PRIVATE histogram *pam_acolorhashtoacolorhist(const struct acolorhash_table *acht, const double gamma, void* (*malloc)(size_t), void (*free)(void*))
+{
+    histogram *hist = malloc(sizeof(hist[0]));
+    if (!hist || !acht) return NULL;
+    *hist = (histogram){
+        .achv = malloc(acht->colors * sizeof(hist->achv[0])),
+        .size = acht->colors,
+        .free = free,
+        .ignorebits = acht->ignorebits,
+    };
+    if (!hist->achv) return NULL;
+
+    float gamma_lut[256];
+    to_f_set_gamma(gamma_lut, gamma);
+
+    /* Limit perceptual weight to 1/10th of the image surface area to prevent
+       a single color from dominating all others. */
+    float max_perceptual_weight = 0.1f * acht->cols * acht->rows;
+    double total_weight = 0;
+
+    for(unsigned int j=0, i=0; i < acht->hash_size; ++i) {
+        const struct acolorhist_arr_head *const achl = &acht->buckets[i];
+        if (achl->used) {
+            PAM_ADD_TO_HIST(achl->inline1);
+
+            if (achl->used > 1) {
+                PAM_ADD_TO_HIST(achl->inline2);
+
+                for(unsigned int k=0; k < achl->used-2; k++) {
+                    PAM_ADD_TO_HIST(achl->other_items[k]);
+                }
+            }
+        }
+    }
+
+    hist->total_perceptual_weight = total_weight;
+    return hist;
+}
+
+
+LIQ_PRIVATE void pam_freeacolorhash(struct acolorhash_table *acht)
+{
+    mempool_destroy(acht->mempool);
+}
+
+LIQ_PRIVATE void pam_freeacolorhist(histogram *hist)
+{
+    hist->free(hist->achv);
+    hist->free(hist);
+}
+
+LIQ_PRIVATE colormap *pam_colormap(unsigned int colors, void* (*malloc)(size_t), void (*free)(void*))
+{
+    assert(colors > 0 && colors < 65536);
+
+    colormap *map;
+    const size_t colors_size = colors * sizeof(map->palette[0]);
+    map = malloc(sizeof(colormap) + colors_size);
+    if (!map) return NULL;
+    *map = (colormap){
+        .malloc = malloc,
+        .free = free,
+        .subset_palette = NULL,
+        .colors = colors,
+    };
+    memset(map->palette, 0, colors_size);
+    return map;
+}
+
+LIQ_PRIVATE colormap *pam_duplicate_colormap(colormap *map)
+{
+    colormap *dupe = pam_colormap(map->colors, map->malloc, map->free);
+    for(unsigned int i=0; i < map->colors; i++) {
+        dupe->palette[i] = map->palette[i];
+    }
+    if (map->subset_palette) {
+        dupe->subset_palette = pam_duplicate_colormap(map->subset_palette);
+    }
+    return dupe;
+}
+
+LIQ_PRIVATE void pam_freecolormap(colormap *c)
+{
+    if (c->subset_palette) pam_freecolormap(c->subset_palette);
+    c->free(c);
+}
+
+LIQ_PRIVATE void to_f_set_gamma(float gamma_lut[], const double gamma)
+{
+    for(int i=0; i < 256; i++) {
+        gamma_lut[i] = pow((double)i/255.0, internal_gamma/gamma);
+    }
+}
+

deps/pngquant/pam.h

@@ -0,0 +1,283 @@
+/* pam.h - pam (portable alpha map) utility library
+ **
+ ** Colormap routines.
+ **
+ ** Copyright (C) 1989, 1991 by Jef Poskanzer.
+ ** Copyright (C) 1997 by Greg Roelofs.
+ **
+ ** Permission to use, copy, modify, and distribute this software and its
+ ** documentation for any purpose and without fee is hereby granted, provided
+ ** that the above copyright notice appear in all copies and that both that
+ ** copyright notice and this permission notice appear in supporting
+ ** documentation.  This software is provided "as is" without express or
+ ** implied warranty.
+ */
+
+#ifndef PAM_H
+#define PAM_H
+
+#include <math.h>
+#include <assert.h>
+#include <stdlib.h>
+#include <stdbool.h>
+
+#ifndef MAX
+#  define MAX(a,b)  ((a) > (b)? (a) : (b))
+#  define MIN(a,b)  ((a) < (b)? (a) : (b))
+#endif
+
+#define MAX_DIFF 1e20
+
+#if __SSE__ || (_M_IX86_FP > 0)
+#  include <xmmintrin.h>
+#  ifdef _MSC_VER
+#    include <intrin.h>
+#    define SSE_ALIGN
+#  else
+#    define SSE_ALIGN __attribute__ ((aligned (16)))
+#    if defined(__i386__) && defined(__PIC__)
+#       define cpuid(func,ax,bx,cx,dx)\
+        __asm__ __volatile__ ( \
+        "push %%ebx\n" \
+        "cpuid\n" \
+        "mov %%ebx, %1\n" \
+        "pop %%ebx\n" \
+        : "=a" (ax), "=r" (bx), "=c" (cx), "=d" (dx) \
+        : "a" (func));
+#    else
+#       define cpuid(func,ax,bx,cx,dx)\
+        __asm__ __volatile__ ("cpuid":\
+        "=a" (ax), "=b" (bx), "=c" (cx), "=d" (dx) : "a" (func));
+#    endif
+#endif
+#else
+#  define SSE_ALIGN
+#endif
+
+#if defined(__GNUC__) || defined (__llvm__)
+#define ALWAYS_INLINE __attribute__((always_inline)) inline
+#define NEVER_INLINE __attribute__ ((noinline))
+#elif defined(_MSC_VER)
+#define inline __inline
+#define restrict __restrict
+#define ALWAYS_INLINE __forceinline
+#define NEVER_INLINE __declspec(noinline)
+#else
+#define ALWAYS_INLINE inline
+#define NEVER_INLINE
+#endif
+
+/* from pam.h */
+
+typedef struct {
+    unsigned char r, g, b, a;
+} rgba_pixel;
+
+typedef struct {
+    float a, r, g, b;
+} SSE_ALIGN f_pixel;
+
+static const double internal_gamma = 0.5499;
+
+LIQ_PRIVATE void to_f_set_gamma(float gamma_lut[], const double gamma);
+
+/**
+ Converts 8-bit color to internal gamma and premultiplied alpha.
+ (premultiplied color space is much better for blending of semitransparent colors)
+ */
+ALWAYS_INLINE static f_pixel to_f(const float gamma_lut[], const rgba_pixel px);
+inline static f_pixel to_f(const float gamma_lut[], const rgba_pixel px)
+{
+    float a = px.a/255.f;
+
+    return (f_pixel) {
+        .a = a,
+        .r = gamma_lut[px.r]*a,
+        .g = gamma_lut[px.g]*a,
+        .b = gamma_lut[px.b]*a,
+    };
+}
+
+inline static rgba_pixel to_rgb(const float gamma, const f_pixel px)
+{
+    if (px.a < 1.f/256.f) {
+        return (rgba_pixel){0,0,0,0};
+    }
+
+    float r = px.r / px.a,
+          g = px.g / px.a,
+          b = px.b / px.a,
+          a = px.a;
+
+    r = powf(r, gamma/internal_gamma);
+    g = powf(g, gamma/internal_gamma);
+    b = powf(b, gamma/internal_gamma);
+
+    // 256, because numbers are in range 1..255.9999… rounded down
+    r *= 256.f;
+    g *= 256.f;
+    b *= 256.f;
+    a *= 256.f;
+
+    return (rgba_pixel){
+        .r = r>=255.f ? 255 : r,
+        .g = g>=255.f ? 255 : g,
+        .b = b>=255.f ? 255 : b,
+        .a = a>=255.f ? 255 : a,
+    };
+}
+
+ALWAYS_INLINE static double colordifference_ch(const double x, const double y, const double alphas);
+inline static double colordifference_ch(const double x, const double y, const double alphas)
+{
+    // maximum of channel blended on white, and blended on black
+    // premultiplied alpha and backgrounds 0/1 shorten the formula
+    const double black = x-y, white = black+alphas;
+    return black*black + white*white;
+}
+
+ALWAYS_INLINE static float colordifference_stdc(const f_pixel px, const f_pixel py);
+inline static float colordifference_stdc(const f_pixel px, const f_pixel py)
+{
+    // px_b.rgb = px.rgb + 0*(1-px.a) // blend px on black
+    // px_b.a   = px.a   + 1*(1-px.a)
+    // px_w.rgb = px.rgb + 1*(1-px.a) // blend px on white
+    // px_w.a   = px.a   + 1*(1-px.a)
+
+    // px_b.rgb = px.rgb              // difference same as in opaque RGB
+    // px_b.a   = 1
+    // px_w.rgb = px.rgb - px.a       // difference simplifies to formula below
+    // px_w.a   = 1
+
+    // (px.rgb - px.a) - (py.rgb - py.a)
+    // (px.rgb - py.rgb) + (py.a - px.a)
+
+    const double alphas = py.a-px.a;
+    return colordifference_ch(px.r, py.r, alphas) +
+           colordifference_ch(px.g, py.g, alphas) +
+           colordifference_ch(px.b, py.b, alphas);
+}
+
+ALWAYS_INLINE static double min_colordifference_ch(const double x, const double y, const double alphas);
+inline static double min_colordifference_ch(const double x, const double y, const double alphas)
+{
+    const double black = x-y, white = black+alphas;
+    return MIN(black*black , white*white) * 2.f;
+}
+
+/* least possible difference between colors (difference varies depending on background they're blended on) */
+ALWAYS_INLINE static float min_colordifference(const f_pixel px, const f_pixel py);
+inline static float min_colordifference(const f_pixel px, const f_pixel py)
+{
+    const double alphas = py.a-px.a;
+    return min_colordifference_ch(px.r, py.r, alphas) +
+           min_colordifference_ch(px.g, py.g, alphas) +
+           min_colordifference_ch(px.b, py.b, alphas);
+}
+
+ALWAYS_INLINE static float colordifference(f_pixel px, f_pixel py);
+inline static float colordifference(f_pixel px, f_pixel py)
+{
+#if __SSE__ || (_M_IX86_FP > 0)
+    const __m128 vpx = _mm_load_ps((const float*)&px);
+    const __m128 vpy = _mm_load_ps((const float*)&py);
+
+    // y.a - x.a
+    __m128 alphas = _mm_sub_ss(vpy, vpx);
+    alphas = _mm_shuffle_ps(alphas,alphas,0); // copy first to all four
+
+    __m128 onblack = _mm_sub_ps(vpx, vpy); // x - y
+    __m128 onwhite = _mm_add_ps(onblack, alphas); // x - y + (y.a - x.a)
+
+    onblack = _mm_mul_ps(onblack, onblack);
+    onwhite = _mm_mul_ps(onwhite, onwhite);
+    const __m128 max = _mm_add_ps(onwhite, onblack);
+
+    // add rgb, not a
+    const __m128 maxhl = _mm_movehl_ps(max, max);
+    const __m128 tmp = _mm_add_ps(max, maxhl);
+    const __m128 sum = _mm_add_ss(maxhl, _mm_shuffle_ps(tmp, tmp, 1));
+
+    const float res = _mm_cvtss_f32(sum);
+    assert(fabs(res - colordifference_stdc(px,py)) < 0.001);
+    return res;
+#else
+    return colordifference_stdc(px,py);
+#endif
+}
+
+/* from pamcmap.h */
+union rgba_as_int {
+    rgba_pixel rgba;
+    unsigned int l;
+};
+
+typedef struct {
+    f_pixel acolor;
+    float adjusted_weight,   // perceptual weight changed to tweak how mediancut selects colors
+          perceptual_weight; // number of pixels weighted by importance of different areas of the picture
+
+    float color_weight;      // these two change every time histogram subset is sorted
+    union {
+        unsigned int sort_value;
+        unsigned char likely_colormap_index;
+    } tmp;
+} hist_item;
+
+typedef struct {
+    hist_item *achv;
+    void (*free)(void*);
+    double total_perceptual_weight;
+    unsigned int size;
+    unsigned int ignorebits;
+} histogram;
+
+typedef struct {
+    f_pixel acolor;
+    float popularity;
+    bool fixed; // if true it's user-supplied and must not be changed (e.g in voronoi iteration)
+} colormap_item;
+
+typedef struct colormap {
+    unsigned int colors;
+    void* (*malloc)(size_t);
+    void (*free)(void*);
+    struct colormap *subset_palette;
+    colormap_item palette[];
+} colormap;
+
+struct acolorhist_arr_item {
+    union rgba_as_int color;
+    float perceptual_weight;
+};
+
+struct acolorhist_arr_head {
+    unsigned int used, capacity;
+    struct {
+        union rgba_as_int color;
+        float perceptual_weight;
+    } inline1, inline2;
+    struct acolorhist_arr_item *other_items;
+};
+
+struct acolorhash_table {
+    struct mempool *mempool;
+    unsigned int ignorebits, maxcolors, colors, cols, rows;
+    unsigned int hash_size;
+    unsigned int freestackp;
+    struct acolorhist_arr_item *freestack[512];
+    struct acolorhist_arr_head buckets[];
+};
+
+LIQ_PRIVATE void pam_freeacolorhash(struct acolorhash_table *acht);
+LIQ_PRIVATE struct acolorhash_table *pam_allocacolorhash(unsigned int maxcolors, unsigned int surface, unsigned int ignorebits, void* (*malloc)(size_t), void (*free)(void*));
+LIQ_PRIVATE histogram *pam_acolorhashtoacolorhist(const struct acolorhash_table *acht, const double gamma, void* (*malloc)(size_t), void (*free)(void*));
+LIQ_PRIVATE bool pam_computeacolorhash(struct acolorhash_table *acht, const rgba_pixel *const pixels[], unsigned int cols, unsigned int rows, const unsigned char *importance_map);
+
+LIQ_PRIVATE void pam_freeacolorhist(histogram *h);
+
+LIQ_PRIVATE colormap *pam_colormap(unsigned int colors, void* (*malloc)(size_t), void (*free)(void*));
+LIQ_PRIVATE colormap *pam_duplicate_colormap(colormap *map);
+LIQ_PRIVATE void pam_freecolormap(colormap *c);
+
+#endif

deps/pngquant/viter.c

@@ -0,0 +1,91 @@
+
+#include "libimagequant.h"
+#include "pam.h"
+#include "viter.h"
+#include "nearest.h"
+#include <stdlib.h>
+#include <string.h>
+
+#ifdef _OPENMP
+#include <omp.h>
+#else
+#define omp_get_max_threads() 1
+#define omp_get_thread_num() 0
+#endif
+
+/*
+ * Voronoi iteration: new palette color is computed from weighted average of colors that map to that palette entry.
+ */
+LIQ_PRIVATE void viter_init(const colormap *map, const unsigned int max_threads, viter_state average_color[])
+{
+    memset(average_color, 0, sizeof(average_color[0])*(VITER_CACHE_LINE_GAP+map->colors)*max_threads);
+}
+
+LIQ_PRIVATE void viter_update_color(const f_pixel acolor, const float value, const colormap *map, unsigned int match, const unsigned int thread, viter_state average_color[])
+{
+    match += thread * (VITER_CACHE_LINE_GAP+map->colors);
+    average_color[match].a += acolor.a * value;
+    average_color[match].r += acolor.r * value;
+    average_color[match].g += acolor.g * value;
+    average_color[match].b += acolor.b * value;
+    average_color[match].total += value;
+}
+
+LIQ_PRIVATE void viter_finalize(colormap *map, const unsigned int max_threads, const viter_state average_color[])
+{
+    for (unsigned int i=0; i < map->colors; i++) {
+        double a=0, r=0, g=0, b=0, total=0;
+
+        // Aggregate results from all threads
+        for(unsigned int t=0; t < max_threads; t++) {
+            const unsigned int offset = (VITER_CACHE_LINE_GAP+map->colors) * t + i;
+
+            a += average_color[offset].a;
+            r += average_color[offset].r;
+            g += average_color[offset].g;
+            b += average_color[offset].b;
+            total += average_color[offset].total;
+        }
+
+        if (total && !map->palette[i].fixed) {
+            map->palette[i].acolor = (f_pixel){
+                .a = a / total,
+                .r = r / total,
+                .g = g / total,
+                .b = b / total,
+            };
+        } else {
+            total = i/1024.0;
+        }
+        map->palette[i].popularity = total;
+    }
+}
+
+LIQ_PRIVATE double viter_do_iteration(histogram *hist, colormap *const map, const float min_opaque_val, viter_callback callback, const bool fast_palette)
+{
+    const unsigned int max_threads = omp_get_max_threads();
+    viter_state average_color[(VITER_CACHE_LINE_GAP+map->colors) * max_threads];
+    viter_init(map, max_threads, average_color);
+    struct nearest_map *const n = nearest_init(map, fast_palette);
+    hist_item *const achv = hist->achv;
+    const int hist_size = hist->size;
+
+    double total_diff=0;
+    #pragma omp parallel for if (hist_size > 3000) \
+        schedule(static) default(none) shared(average_color,callback) reduction(+:total_diff)
+    for(int j=0; j < hist_size; j++) {
+        float diff;
+        unsigned int match = nearest_search(n, achv[j].acolor, achv[j].tmp.likely_colormap_index, min_opaque_val, &diff);
+        achv[j].tmp.likely_colormap_index = match;
+        total_diff += diff * achv[j].perceptual_weight;
+
+        viter_update_color(achv[j].acolor, achv[j].perceptual_weight, map, match, omp_get_thread_num(), average_color);
+
+        if (callback) callback(&achv[j], diff);
+    }
+
+    nearest_free(n);
+    viter_finalize(map, max_threads, average_color);
+
+    return total_diff / hist->total_perceptual_weight;
+}

deps/pngquant/viter.h

@@ -0,0 +1,19 @@
+
+#ifndef VITER_H
+#define VITER_H
+
+// Spread memory touched by different threads at least 64B apart which I assume is the cache line size. This should avoid memory write contention.
+#define VITER_CACHE_LINE_GAP ((64+sizeof(viter_state)-1)/sizeof(viter_state))
+
+typedef struct {
+    double a, r, g, b, total;
+} viter_state;
+
+typedef void (*viter_callback)(hist_item *item, float diff);
+
+LIQ_PRIVATE void viter_init(const colormap *map, const unsigned int max_threads, viter_state state[]);
+LIQ_PRIVATE void viter_update_color(const f_pixel acolor, const float value, const colormap *map, unsigned int match, const unsigned int thread, viter_state average_color[]);
+LIQ_PRIVATE void viter_finalize(colormap *map, const unsigned int max_threads, const viter_state state[]);
+LIQ_PRIVATE double viter_do_iteration(histogram *hist, colormap *const map, const float min_opaque_val, viter_callback callback, const bool fast_palette);
+
+#endif

include/mapnik/png_io.hpp

@@ -47,14 +47,20 @@ extern "C"
 namespace mapnik {
 
 struct png_options {
+
+    enum quantization_enum { HEXTREE, OCTTREE, IMGQUANT };
+
     int colors;
     int filters;
     int compression;
     int strategy;
     int trans_mode;
+    int iq_speed;
+    double iq_dither;
     double gamma;
     bool paletted;
-    bool use_hextree;
+    quantization_enum quantization;
+    bool use_miniz;
 
     png_options() :
         colors(256),
@@ -62,9 +68,12 @@ struct png_options {
         compression(Z_DEFAULT_COMPRESSION),
         strategy(Z_DEFAULT_STRATEGY),
         trans_mode(-1),
+        iq_speed(3),
+        iq_dither(-1),
         gamma(-1),
         paletted(true),
-        use_hextree(true) {}
+        quantization(HEXTREE),
+        use_miniz(false) {}
 };
 
 template <typename T>
@@ -701,6 +710,11 @@ void save_as_png8_pal(T1 & file,
     save_as_png8<T1, T2, rgba_palette>(file, image, pal, pal.palette(), pal.alpha_table(), opts);
 }
 
+template <typename T1, typename T2>
+void save_as_png8_libimagequant(T1 & file,
+                                T2 const& image,
+                                png_options const& opts);
+
 }
 
 #endif // MAPNIK_PNG_IO_HPP

src/build.py

@@ -22,6 +22,7 @@
 import os
 import sys
 import glob
+import platform
 from copy import copy
 from subprocess import Popen, PIPE
 
@@ -273,6 +274,7 @@ source = Split(
     renderer_common/pattern_alignment.cpp
     util/math.cpp
     value.cpp
+    png_io.cpp
     """
     )
 
@@ -366,6 +368,22 @@ source += Split(
     """
     )
 
+# libimagequant
+lib_env.Append(CFLAGS = "-O3 -fno-math-errno -funroll-loops -fomit-frame-pointer -std=c99")
+if env["USE_SSE"] == "yes":
+    lib_env.Append(CFLAGS="-msse -mfpmath=sse")
+    # As of GCC 4.5, 387 fp math is significantly slower in C99 mode without this.
+    # Note: CPUs without SSE2 use 387 for doubles, even when SSE fp math is set.
+    if 'gcc' in env['CC']:
+        lib_env.Append(CFLAGS='-fexcess-precision=fast')
+elif env["USE_SSE"] == "no":
+    lib_env.Append(CFLAGS="-mno-sse")
+elif env["USE_SSE"] == "platform_default":
+    # Let compiler decide
+    None
+
+source += glob.glob("../deps/pngquant/" + "*.c")
+
 if env['RUNTIME_LINK'] == "static":
     source += glob.glob('../deps/agg/src/' + '*.cpp')
 

src/image_util_png.cpp

@@ -83,8 +83,23 @@ void handle_png_options(std::string const& type,
         }
         else if (key == "m" && val)
         {
-            if (*val == "o") opts.use_hextree = false;
-            else if (*val == "h") opts.use_hextree = true;
+            if (*val == "o") opts.quantization = png_options::OCTTREE;
+            else if (*val == "h") opts.quantization = png_options::HEXTREE;
+            else if (*val == "iq") opts.quantization = png_options::IMGQUANT;
+        }
+        else if (key == "iqs")
+        {
+            if (!val || !mapnik::util::string2int(*val, opts.iq_speed) || opts.iq_speed < 1 || opts.iq_speed > 10)
+            {
+                throw image_writer_exception("invalid iq speed parameter: " + to_string(val));
+            }
+        }
+        else if (key == "iqd")
+        {
+            if (!val || !mapnik::util::string2double(*val, opts.iq_dither) || opts.iq_dither < 0 || opts.iq_dither > 1)
+            {
+                throw image_writer_exception("invalid iq dithering parameter: " + to_string(val));
+            }
         }
         else if (key == "e" && val && *val == "miniz")
         {
@@ -236,14 +251,18 @@ void process_rgba8_png_pal(T const& image,
     }
     else if (opts.paletted)
     {
-        if (opts.use_hextree)
+        if (opts.quantization == png_options::HEXTREE)
         {
             save_as_png8_hex(stream, image, opts);
         }
-        else
+        else if (opts.quantization == png_options::OCTTREE)
         {
             save_as_png8_oct(stream, image, opts);
         }
+        else
+        {
+            save_as_png8_libimagequant(stream, image, opts);
+        }
     }
     else
     {
@@ -264,14 +283,18 @@ void process_rgba8_png(T const& image,
     handle_png_options(t, opts);
     if (opts.paletted)
     {
-        if (opts.use_hextree)
+        if (opts.quantization == png_options::HEXTREE)
         {
             save_as_png8_hex(stream, image, opts);
         }
-        else
+        else if (opts.quantization == png_options::OCTTREE)
         {
             save_as_png8_oct(stream, image, opts);
         }
+        else
+        {
+            save_as_png8_libimagequant(stream, image, opts);
+        }
     }
     else
     {

src/png_io.cpp

@@ -0,0 +1,130 @@
+/*****************************************************************************
+ *
+ * This file is part of Mapnik (c++ mapping toolkit)
+ *
+ * Copyright (C) 2014 Daniel Patterson
+ *
+ * This library is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public
+ * License as published by the Free Software Foundation; either
+ * version 2.1 of the License, or (at your option) any later version.
+ *
+ * This library is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with this library; if not, write to the Free Software
+ * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA
+ *
+ *****************************************************************************/
+
+// mapnik
+#include <mapnik/png_io.hpp>
+#include <mapnik/image_view.hpp>
+
+extern "C"
+{
+#include <libimagequant.h>
+}
+
+namespace mapnik
+{
+
+  // Note: This only works with rgba8_t image types
+  template <typename T1, typename T2>
+  void save_as_png8_libimagequant(T1 & file,
+                                  T2 const& image,
+                                  png_options const& opts)
+  {
+      unsigned width = image.width();
+      unsigned height = image.height();
+
+      uint32_t const* buf[height];
+
+      // TODO: this won't work on big-endian architectures, liq expects data in RGBA byte order
+      for (size_t y = 0; y < height; ++y)
+      {
+          buf[y] = image.get_row(y);
+      }
+
+      // Set up libimagequant
+      liq_attr *attr = liq_attr_create();
+      liq_set_speed(attr, opts.iq_speed);
+      liq_image *liq_image = liq_image_create_rgba_rows(attr, (void **)buf, width, height, 0);
+
+      // Do the quantization
+      liq_result *res = liq_quantize_image(attr, liq_image);
+
+      if (opts.iq_dither != -1)
+      {
+          liq_set_dithering_level(res, opts.iq_dither);
+      }
+      if (opts.gamma > 0)
+      {
+          liq_set_output_gamma(res, opts.gamma);
+      }
+
+      // Extract the generated 8-bit imag
+      image_gray8 quantized_image(width, height);
+      liq_write_remapped_image(res, liq_image, (void *)quantized_image.data(), width*height*sizeof(gray8_t));
+
+      // Extract the generated palette and alpha vector from libimagequant
+      std::vector<mapnik::rgb> palette;
+      std::vector<unsigned> alpha;
+      const liq_palette *liq_pal = liq_get_palette(res);
+      for (unsigned i = 0; i < liq_pal->count; ++i)
+      {
+          palette.push_back({ liq_pal->entries[i].r, liq_pal->entries[i].g, liq_pal->entries[i].b });
+          alpha.push_back(liq_pal->entries[i].a);
+      }
+
+      // Encode the image, using a bit-depth depending on the number of colors in the generated
+      // palette
+      if (palette.size() > 16)  // 8-bit encoding
+      {
+          save_as_png(file, palette, quantized_image, width, height, 8, alpha, opts);
+      }
+      else if (palette.size() > 1) // up to 16 colors, 4-bit encoding
+      {
+          unsigned image_width  = ((width + 7) >> 1) & ~3U; // 4-bit image, round up to 32-bit boundary
+          unsigned image_height = height;
+          image_gray8 reduced_image(image_width, image_height);
+          for (unsigned y = 0; y < height; ++y)
+          {
+              mapnik::image_gray8::pixel_type const * row = quantized_image.get_row(y);
+              mapnik::image_gray8::pixel_type  * row_out = reduced_image.get_row(y);
+              std::uint8_t index = 0;
+              for (unsigned x = 0; x < width; ++x)
+              {
+                  index = row[x];
+                  if (x%2 == 0)
+                  {
+                      index = index<<4;
+                  }
+                  row_out[x>>1] |= index;
+              }
+          }
+          save_as_png(file, palette, reduced_image, width, height, 4, alpha, opts);
+      }
+      else // 1 color, 1-bit encoding
+      {
+          unsigned image_width  = ((width + 15) >> 3) & ~1U; // 1-bit image, round up to 16-bit boundary
+          unsigned image_height = height;
+          image_gray8 reduced_image(image_width, image_height);
+          reduced_image.set(0);
+          save_as_png(file, palette, reduced_image, width, height, 1, alpha, opts);
+      }
+
+      // Cleanup libimagequant
+      liq_attr_destroy(attr);
+      liq_image_destroy(liq_image);
+      liq_result_destroy(res);
+  }
+
+  // Both image_view and image expose a get_row(), so we can support both parameter types
+  template void save_as_png8_libimagequant(std::ostream &, mapnik::image_view<mapnik::image<mapnik::rgba8_t>> const&, png_options const&);
+  template void save_as_png8_libimagequant(std::ostream &, mapnik::image<mapnik::rgba8_t> const&, png_options const&);
+
+}