Diamond-square algorithm. NEVER USE THIS CODE IN PRODUCTION, PLEASE. IT IS AWFUL

diamond-square.cpp

// Diamond-square algorithm

// MIT License
// 
// Copyright (c) 2019 Alexander Korzun
// 
// Permission is hereby granted, free of charge, to any person obtaining a copy
// of this software and associated documentation files (the "Software"), to deal
// in the Software without restriction, including without limitation the rights
// to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the Software is
// furnished to do so, subject to the following conditions:
// 
// The above copyright notice and this permission notice shall be included in all
// copies or substantial portions of the Software.
// 
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
// OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
// SOFTWARE.

// Compile me with `c++ $(pkg-config --cflags --libs libpng) -std=c++17 diamond-square.cpp -o diamond-square`
// Requirements: libpng, libpng++

#include <iostream>
#include <vector>
#include <random>
#include <cmath>
#include <cassert>
#include <queue>
#include <png++/png.hpp>


using namespace std;
using namespace png;


using Matrix = vector<vector<double>>;

constexpr double v2 = sqrt(2.0);


const auto SIZE = 513;


double random_float()
{
    static auto gen = mt19937(random_device()());
    static auto dist = uniform_real_distribution<double>(-1.0, 1.0);
    return dist(gen);
}


enum Type
{
    DIAMOND,
    SQUARE,
};


struct Point
{
    int x;
    int y;

    Point operator+(const Point& other) const noexcept
    {
        return {x + other.x, y + other.y};
    }

    Point operator-(const Point& other) const noexcept
    {
        return {x - other.x, y - other.y};
    }

    Point operator*(int k) const noexcept
    {
        return {x * k, y * k};
    }

    Point operator/(int k) const noexcept
    {
        return {x / k, y / k};
    }
};


struct Query
{
    Type type;
    Point p1;
    Point p2;
    Point p3;
    Point p4;
};


bool in(const Point& p)
{
    return p.x >= 0 && p.x < SIZE && p.y >= 0 && p.y < SIZE;
}


double average4(Matrix& v, const Point& a, const Point& b, const Point& c, const Point& d)
{
    double sum = 0;
    int count = 0;
    if (in(a)) {
        ++count;
        sum += v[a.y][a.x];
    }
    if (in(b)) {
        ++count;
        sum += v[b.y][b.x];
    }
    if (in(c)) {
        ++count;
        sum += v[c.y][c.x];
    }
    if (in(d)) {
        ++count;
        sum += v[d.y][d.x];
    }
    assert(count > 0);
    return sum / count;
}


void diamond(Matrix& v, queue<Query>& q, const Query& qq, vector<vector<char>>& used)
{
    auto lt = qq.p1;
    auto rt = qq.p2;
    auto lb = qq.p3;
    auto rb = qq.p4;
    auto center = (lt + rb) / 2;
    if (!in(center) || used[center.y][center.x]) {
        return;
    }
    auto len = (rt - lt).x;
    auto left = center - Point{len, 0};
    auto right = center + Point{len, 0};
    auto top = center - Point{0, len};
    auto bottom = center + Point{0, len};
    auto base = average4(v, lt, rt, lb, rb);
    v[center.y][center.x] = base + random_float() * len / SIZE;
    used[center.y][center.x] = 1;
    q.push({
        SQUARE,
        lt,
        left,
        center,
        lb,
    });
    q.push({
        SQUARE,
        top,
        lt,
        rt,
        center,
    });
    q.push({
        SQUARE,
        rt,
        center,
        right,
        rb,
    });
    q.push({
        SQUARE,
        center,
        lb,
        rb,
        bottom,
    });
}


void square(Matrix& v, queue<Query>& q, const Query& qq, vector<vector<char>>& used)
{
    auto top = qq.p1;
    auto left = qq.p2;
    auto right = qq.p3;
    auto bottom = qq.p4;
    auto center = (left + right) / 2;
    if (!in(center) || used[center.y][center.x]) {
        return;
    }
    auto len = (center - left).x;
    auto lt = center + Point{-len, -len};
    auto rt = center + Point{len, -len};
    auto lb = center + Point{-len, len};
    auto rb = center + Point{len, len};
    auto base = average4(v, left, right, top, bottom);
    v[center.y][center.x] = base + random_float() * len / SIZE / v2;
    used[center.y][center.x] = 1;
    if (len > 1) {
        q.push({
            DIAMOND,
            lt,
            top,
            left,
            center,
        });
        q.push({
            DIAMOND,
            top,
            rt,
            center,
            right,
        });
        q.push({
            DIAMOND,
            left,
            center,
            lb,
            bottom,
        });
        q.push({
            DIAMOND,
            center,
            right,
            bottom,
            rb,
        });
    }
}


void dsq(Matrix& v)
{
    vector<vector<char>> used(SIZE, vector<char>(SIZE, 0));
    queue<Query> q;
    q.push({
        DIAMOND,
        {0, 0},
        {SIZE-1, 0},
        {0, SIZE-1},
        {SIZE-1, SIZE-1},
    });

    while (!q.empty()) {
        auto qq = q.front();
        q.pop();
        if (qq.type == DIAMOND) {
            diamond(v, q, qq, used);
        } else {
            square(v, q, qq, used);
        }
    }
}


void norm(Matrix& v)
{
    double minv = 1e100;
    double maxv = -1e100;

    for (auto& row : v) {
        for (double x : row) {
            minv = min(minv, x);
            maxv = max(maxv, x);
        }
    }

    for (auto& row : v) {
        for (double& x : row) {
            x -= minv;
            x /= (maxv - minv);
        }
    }
}


int main()
{
    Matrix v(SIZE, vector<double>(SIZE, 0));
    v[0][0] = random_float();
    v[0][SIZE-1] = random_float();
    v[SIZE-1][0] = random_float();
    v[SIZE-1][SIZE-1] = random_float();

    dsq(v);
    norm(v);

    image<rgb_pixel> img(SIZE, SIZE);
    for (int i = 0; i < SIZE; ++i) {
        for (int j = 0; j < SIZE; ++j) {
            auto color = int(v[i][j] * 255);
            img[i][j] = rgb_pixel(color, color, color);
        }
    }

    img.write("output.png");
}