+#include <stdio.h>
+#include <math.h>
+#include <png.h>
+#include <SDL/SDL_image.h>
+#include <vector>
+#include <algorithm>
+
+#define NUM_CLUSTERS 16
+
+using namespace std;
+
+struct Pixel {
+ float r, g, b, a;
+
+ Pixel(float r, float g, float b, float a) : r(r), g(g), b(b), a(a) {}
+ Pixel() {}
+};
+
+struct CompareByRGBA {
+ bool operator() (const Pixel &a, const Pixel &b) const {
+ if (a.a != b.a)
+ return (a.a < b.a);
+ if (a.r != b.r)
+ return (a.r < b.r);
+ if (a.g != b.g)
+ return (a.g < b.g);
+ return (a.b < b.b);
+ }
+};
+
+void write_png(const char *filename, unsigned char *screenbuf, int width, int height)
+{
+ FILE *fp = fopen(filename, "wb");
+ png_structp png_ptr = png_create_write_struct(PNG_LIBPNG_VER_STRING, NULL, NULL, NULL);
+ png_infop info_ptr = png_create_info_struct(png_ptr);
+
+ if (setjmp(png_jmpbuf(png_ptr))) {
+ fclose(fp);
+ fprintf(stderr, "Write to %s failed; exiting.\n", filename);
+ exit(1);
+ }
+
+ png_set_IHDR(png_ptr, info_ptr, width, height, 8, PNG_COLOR_TYPE_RGB_ALPHA, PNG_INTERLACE_NONE, PNG_COMPRESSION_TYPE_DEFAULT, PNG_FILTER_TYPE_DEFAULT);
+
+ png_bytep *row_pointers = new png_bytep[height];
+ for (int y = 0; y < height; ++y) {
+ row_pointers[y] = screenbuf + (y * width) * 4;
+ }
+
+ png_init_io(png_ptr, fp);
+ png_set_rows(png_ptr, info_ptr, row_pointers);
+ png_write_png(png_ptr, info_ptr, PNG_TRANSFORM_BGR, NULL);
+ png_destroy_write_struct(&png_ptr, &info_ptr);
+ fclose(fp);
+
+ delete[] row_pointers;
+}
+
+// [0.0, 1.0)
+float uniform_rand()
+{
+ return rand() / float(RAND_MAX + 1.0);
+}
+
+Pixel random_pixel()
+{
+ float r = uniform_rand();
+ float g = uniform_rand();
+ float b = uniform_rand();
+ float a = uniform_rand();
+
+ r *= a;
+ g *= a;
+ b *= a;
+
+ return Pixel(r, g, b, a);
+}
+
+int assign_to_cluster(const Pixel &pixel, const Pixel *center)
+{
+ int best_cluster = -1;
+ float best_dist = HUGE_VAL;
+
+ for (int k = 0; k < NUM_CLUSTERS; ++k) {
+ // TODO: Better perceptual distance than just nonlinear RGB space.
+ float dist =
+ (pixel.r - center[k].r) * (pixel.r - center[k].r) +
+ (pixel.g - center[k].g) * (pixel.g - center[k].g) +
+ (pixel.b - center[k].b) * (pixel.b - center[k].b) +
+ 3.0 * (pixel.a - center[k].a) * (pixel.a - center[k].a);
+ if (dist < best_dist) {
+ best_dist = dist;
+ best_cluster = k;
+ }
+ }
+
+ return best_cluster;
+}
+
+int main(void)
+{
+ SDL_Surface *foo = IMG_Load("/home/sesse/clusterfuck/overlay-shadow.png");
+
+ vector<Pixel> pixels;
+ for (int y = 0; y < foo->h; ++y) {
+ uint8_t *ptr = ((uint8_t *)foo->pixels + foo->pitch * y);;
+ for (int x = 0; x < foo->w; ++x) {
+ float r = *ptr++ * (1.0/255.0);
+ float g = *ptr++ * (1.0/255.0);
+ float b = *ptr++ * (1.0/255.0);
+ float a = *ptr++ * (1.0/255.0);
+
+ // Premultiply alpha.
+ r *= a;
+ g *= a;
+ b *= a;
+
+ pixels.push_back(Pixel(r, g, b, a));
+ }
+ }
+
+ // K-means: Initialization
+ Pixel center[NUM_CLUSTERS];
+ for (int i = 0; i < NUM_CLUSTERS; ++i) {
+ center[i] = random_pixel();
+ }
+
+ // K-means: Iteration
+ for (int i = 0; i < 1000; ++i) {
+ Pixel new_center[NUM_CLUSTERS];
+ int num_center[NUM_CLUSTERS];
+
+ sort(center, center + NUM_CLUSTERS, CompareByRGBA());
+
+ for (int j = 0; j < NUM_CLUSTERS; ++j) {
+ printf("%2d: %.0f, %.0f, %.0f, %.0f\n", j,
+ 255.0 * center[j].r,
+ 255.0 * center[j].g,
+ 255.0 * center[j].b,
+ 255.0 * center[j].a);
+ }
+ printf("\n");
+
+ for (int j = 0; j < NUM_CLUSTERS; ++j) {
+ new_center[j] = Pixel(0, 0, 0, 0);
+ num_center[j] = 0;
+ }
+
+ for (unsigned j = 0; j < pixels.size(); ++j) {
+ int best_cluster = assign_to_cluster(pixels[j], center);
+
+ ++num_center[best_cluster];
+ new_center[best_cluster].r += pixels[j].r;
+ new_center[best_cluster].g += pixels[j].g;
+ new_center[best_cluster].b += pixels[j].b;
+ new_center[best_cluster].a += pixels[j].a;
+ }
+
+ for (int j = 0; j < NUM_CLUSTERS; ++j) {
+ if (num_center[j] == 0) {
+ center[j] = random_pixel();
+ } else {
+ center[j].r = new_center[j].r / num_center[j];
+ center[j].g = new_center[j].g / num_center[j];
+ center[j].b = new_center[j].b / num_center[j];
+ center[j].a = new_center[j].a / num_center[j];
+ }
+ }
+ }
+
+ // Write PNG.
+ unsigned char *out = new unsigned char[foo->w * foo->h * 4];
+ unsigned char *ptr = out;
+ for (unsigned i = 0; i < pixels.size(); ++i) {
+ Pixel q = pixels[i];
+
+#if 0
+ // TPDF dither. (TODO: How do you do dither with non-uniform quantization?)
+ q.r += (uniform_rand() + uniform_rand()) * (1.0 / 8.0);
+ q.g += (uniform_rand() + uniform_rand()) * (1.0 / 8.0);
+ q.b += (uniform_rand() + uniform_rand()) * (1.0 / 8.0);
+
+ printf("POST: %f %f %f\n", q.r, q.g, q.b);
+#endif
+
+ const Pixel &p = center[assign_to_cluster(q, center)];
+ if (p.a < 1e-6) {
+ *ptr++ = 0;
+ *ptr++ = 0;
+ *ptr++ = 0;
+ *ptr++ = 0;
+ } else {
+ *ptr++ = lrintf((p.r/p.a) * 255.0);
+ *ptr++ = lrintf((p.g/p.a) * 255.0);
+ *ptr++ = lrintf((p.b/p.a) * 255.0);
+ *ptr++ = lrintf(p.a * 255.0);
+ }
+ }
+ write_png("indexed.png", out, foo->w, foo->h);
+
+ // Write indexed stuff.
+ FILE *fp = fopen("/home/sesse/logo.indexed", "wb");
+ for (unsigned i = 0; i < NUM_CLUSTERS; ++i) {
+ const Pixel& p = center[i];
+ fputc(lrintf((p.r/p.a) * 255.0), fp);
+ fputc(lrintf((p.g/p.a) * 255.0), fp);
+ fputc(lrintf((p.b/p.a) * 255.0), fp);
+ fputc(lrintf(p.a * 255.0), fp);
+ }
+ for (unsigned i = 0; i < pixels.size(); ++i) {
+ fputc(assign_to_cluster(pixels[i], center), fp);
+ }
+ fclose(fp);
+}