#include <stdlib.h> // malloc, EXIT_*
#include <string.h> // memset
#include <inttypes.h>
#include <math.h>
#include <png.h>

#define HILS 10
#define SMOOTHLENGTH 20

#define EXIT_PNG(F) if (!F) { \
    fprintf(stderr, "%s\n", bild.message); \
    return EXIT_FAILURE; \
}
#define u8 unsigned char

void get_hilbert(int *curve, int a)
{
    curve[0] = 0;
    curve[1] = 1;
    curve[2] = a + 1;
    curve[3] = a;
    int step = 0;
    int sidelen = 2;
    while (sidelen < a) {
        int cnt;
        switch (step) {
        case 0: // rotate right, write down
            cnt = sidelen * sidelen;
            for (int i = 0; i < cnt; ++i) {
                int x = curve[i];
                int y = x / a;
                x %= a;
                curve[cnt + i] = (x + sidelen) * a + y;
            }
            break;
        case 1: // mirror right, write right
            cnt = sidelen * sidelen * 2;
            for (int i = 0; i < cnt; ++i) {
                int x = curve[i];
                int y = x / a;
                x %= a;
                curve[2 * cnt - i - 1] = y * a + 2 * sidelen - x - 1;
            }
            sidelen *= 2;
            break;
        case 2: // rotate left, write right
            cnt = sidelen * sidelen;
            for (int i = 0; i < cnt; ++i) {
                int x = curve[i];
                int y = x / a;
                x %= a;
                curve[cnt + i] = x * a + y + sidelen;
            }
            break;
        case 3: // mirror down, write down
            cnt = sidelen * sidelen * 2;
            for (int i = 0; i < cnt; ++i) {
                int x = curve[i];
                int y = x / a;
                x %= a;
                curve[2 * cnt - i - 1] = (2 * sidelen - y - 1) * a + x;
            }
            sidelen *= 2;
            break;
        default:
            exit(EXIT_FAILURE);
        }
        step = (step + 1) % 4;
    }
}


double rv = 0.5;
double getrand()
{
    rv = 10000.0 * sqrt(rv);
    rv -= floor(rv);
    return rv;
}

int main(int argc, char **args)
{
    png_image bild;
    memset(&bild, 0, sizeof(bild));
    bild.version = PNG_IMAGE_VERSION;
    bild.version = PNG_IMAGE_VERSION;
    bild.format = PNG_FORMAT_GRAY;
    bild.flags = PNG_IMAGE_FLAG_COLORSPACE_NOT_sRGB;


    int a = 1 << HILS;

    bild.width = a;
    bild.height = a;
    u8 *pixels = malloc(a * a);


    int *hilb = malloc(a * a * sizeof(int));
    get_hilbert(hilb, a);
    //~ int *hilbinv = malloc(a * a * sizeof(int));
    //~ for (int i = 0; i < a*a; ++i) {
        //~ hilbinv[hilb[i]] = i;
    //~ }

    for (int i = 0; i < a * a; ++i) {
        pixels[i] = getrand() * 255;
    }


    float *smoothpixels = malloc(a * a * sizeof(float));
    float *data = malloc(a * a * sizeof(float));
    for (int i = 0; i < a * a; ++i) {
        data[i] = pixels[i];
    }

    for (int j = 0; j < 500; ++j) {
        for (int i = 2; i < a * a - 2; ++i) {
            if (i < SMOOTHLENGTH || i >= (a * a - SMOOTHLENGTH))
                continue;
            //~ float px = data[i];
            float v = 0.0f;
            for (int o = -SMOOTHLENGTH; o < SMOOTHLENGTH+1; ++o) {
                float otherv = pixels[i + o];
                //~ if (fabsf(otherv - px) < TRESH)
                    v += otherv;
                //~ else
                    //~ v += px;
            }
            smoothpixels[i] = v / (SMOOTHLENGTH+SMOOTHLENGTH+1);
        }

        for (int i = 0; i < SMOOTHLENGTH; ++i) {
            smoothpixels[i] = smoothpixels[SMOOTHLENGTH];
            smoothpixels[a * a - i - 1] = smoothpixels[a * a - SMOOTHLENGTH - 1];
        }

        for (int i = 0; i < a * a; ++i) {
            data[hilb[i]] = smoothpixels[i];
        }
    }

    float min = data[0];
    float max = data[0];
    for (int i = 0; i < a * a; ++i) {
        if (data[i] < min)
            min = data[i];
        else if (data[i] > max)
            max = data[i];
    }

    float f = 255.0f / (max - min);
    for (int i = 0; i < a * a; ++i) {
        pixels[i] = (data[i] - min) * f;
    }



    EXIT_PNG(png_image_write_to_stdio(&bild, stdout, 0, pixels, 0, NULL));
    return EXIT_SUCCESS;
}