#include <u.h>
#include <libc.h>
#include <draw.h>
#include <geometry.h>
Image *px;
double res;
int xo, yo;
Point3 sphere = {0, 1.9, 6, 0.9}; /* w = radius */
Point3 camera = {0, 2.3, 0, 0};
Point3 light = {5, 15, 0, 0};
double
distance(Point3 p){
	double d;
	d = vec3len(subpt3(p, sphere)) - sphere.w;
	return p.y < d ? p.y : d;
}
Point3
march(Point3 p, Point3 v){
	double d, tot;
	int i;
	for(i = tot = 0; i < 100; i++){
		tot += d = distance(p);
		p = addpt3(p, mulpt3(v, d));
		if(tot >= 100.0 || d < 0.01)
			break;
	}
	p.w = tot;
	return p;
}
int
pixel(int x, int y){
	Point3 p, norm, lux;
	double d;
	p = march(camera, normvec3(Pt3((x - xo) / res, (yo - y) / res, 1, 0)));
	if(p.w > 100.0)
		return 0;
	lux = normvec3(subpt3(light, p));
	d = distance(p);
	norm.x = (distance(Pt3(p.x + 0.01, p.y, p.z, 0)) - d) / 0.01;
	norm.y = (distance(Pt3(p.x, p.y + 0.01, p.z, 0)) - d) / 0.01;
	norm.z = (distance(Pt3(p.x, p.y, p.z + 0.01, 0)) - d) / 0.01;
	d = dotvec3(lux, norm) * 255;
	norm = march(addpt3(p, mulpt3(norm, 0.05)), lux);
	if(norm.w + 0.05 < vec3len(subpt3(light, p)))
		d *= 0.17;
	return d < 0 ? 0 : d;
}
void
main(int, char**){
	int x, y;
	uchar hue;
	if(initdraw(nil, nil, nil) < 0)
		sysfatal("initdraw: %r");
	px = allocimage(display, Rect(0, 0, 1, 1), GREY8, 1, DNofill);
	if(px == nil)
		sysfatal("allocimage: %r");
	res = Dx(screen->r) < Dy(screen->r) ? Dx(screen->r) : Dy(screen->r);
	xo = screen->r.min.x + screen->r.max.x >> 1;
	yo = screen->r.min.y + screen->r.max.y >> 1;
	for(x = screen->r.min.x; x < screen->r.max.x; x++){
		for(y = screen->r.min.y; y < screen->r.max.y; y++){
			hue = pixel(x, y);
			loadimage(px, px->r, &hue, 1);
			draw(screen, Rect(x, y, x + 1, y + 1), px, nil, ZP);
		}
		flushimage(display, 1);
	}
	sleep(9999999);
}