Buddhabrot Renderer

buddha.cpp

#include <omp.h>
#include <vector>
#include <complex>
#include <iostream>
#include <fstream>
#include <sstream>
#include <exception>
#include <string>
#include <cstdint>

struct Image
{
	int32_t w, h;
	double* pixels;

	Image(int width, int height)
		: w(width)
		, h(height)
	{
		pixels = new double[w*h];
		for (int i = 0; i < w*h; ++i)
			pixels[i] = 0;
	}

	~Image()
	{
		delete[] pixels;
	}

	inline void inc(int x, int y, double v)
	{
		if (x >= 0 && x < w && y >= 0 && y < h)
			pixels[x*h + y] += v;
	}
};

int buddha(std::vector<std::complex<double>>& visited, std::complex<double> c, int iters)
{
	auto z = c;
	int i = 0;
	for (visited.reserve(iters); i < iters && z.real()*z.real() + z.imag()*z.imag() <= 4.; ++i, z = z*z + c)
		visited.push_back(z);
	visited.push_back(z);
	return i;
}

void brot(Image& img, double x0, double x1, double y0, double y1, int max_iter)
{
	double dx = (x1 - x0) / double(img.w), dy = (y1 - y0) / double(img.h);

	#pragma omp parallel for schedule(dynamic) num_threads(8)
	for (int x = 0; x < img.w; ++x)
	for (int y = 0; y < img.h; ++y)
	{
		std::vector<std::complex<double>> points;
		int iters = buddha(points, std::complex<double>(x * dx + x0, y * dy + y0), max_iter);

		if (iters < max_iter)
		{
			for (const auto &p : points)
			{
				double u  = (p.real() - x0) / dx, v  = (p.imag() - y0) / dy;
				int    s  = int(u),               t  = int(v);
				double du = u - double(s),        dv = v - double(t);

				#pragma omp critical
				{
					img.inc(s,     t,     (1. - du) * (1. - dv));
					img.inc(s,     t + 1, (1. - du) * dv);
					img.inc(s + 1, t,            du * (1. - dv));
					img.inc(s + 1, t + 1,        du * dv);
				}
			}
		}
	}
}

template<typename A, typename B>
A convertTo(const B& v)
{
	std::stringstream ss;
	ss << v;
	A r;
	ss >> r;
	return r;
}

int usage(const std::string& progname)
{
	std::cout << "Usage: " << progname << " [-w width] [-h height] [-i iters] [-r real_low real_high imag_low imag_high] outfile.raw" << std::endl;
	return -1;
}

int main(int argc, char* argv[])
{
	// standard values
	int32_t width=1920, height=1080, iters=25000;
	double y0 = -1.8, y1 = 1., x0 = -2.4, x1 = 2.4;
	std::string outfile;

	// parse command line arguments
	std::vector<std::string> args(argv, argv + argc);
	try
	{
		for (size_t i = 1; i < args.size(); ++i)
		{
			if (args[i] == "--help")
				return usage(args[0]);
			else if (args[i] == "-w" || args[i] == "--width")
				width = convertTo<int>(args.at(++i));
			else if (args[i] == "-h" || args[i] == "--height")
				height = convertTo<int>(args.at(++i));
			else if (args[i] == "-i" || args[i] == "--iters")
				iters = convertTo<int>(args.at(++i));
			else if (args[i] == "-r" || args[i] == "--rect")
			{
				x0 = convertTo<double>(args.at(++i));
				x1 = convertTo<double>(args.at(++i));
				y0 = convertTo<double>(args.at(++i));
				y1 = convertTo<double>(args.at(++i));
			}
			else
			{
				if (!outfile.empty()) return usage(args[0]);
				outfile = args[i];
			}
		}
	}
	catch (const std::out_of_range&)
	{
		return usage(args[0]);
	}
	if (outfile.empty())
		return usage(args[0]);

	// actual work
	std::cout << "Baking at " << width << "x" << height << " range (" << x0 << "," << y0 << ")x(" << x1 << "," << y1 << ") for " << iters << "..." << std::endl;
	Image img(height, width);
	brot(img, y0,y1, x0,x1, iters);
	
	// save file
	std::cout << "Serving to " << outfile << "..." << std::endl;
	std::ofstream of(outfile.c_str(), std::ios::binary);
	if (!of.good())
	{
		std::cout << "Error: Cannot open " << outfile << " for writing" << std::endl;
		return -2;
	}

	of.write((char*)(&width), sizeof(int32_t));
	of.write((char*)(&height), sizeof(int32_t));
	of.write((char*)(img.pixels), sizeof(double)*img.w*img.h);

	return 0;
}