#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define M_PI 3.14159265358979323846
typedef struct {
float x, y, z;
} Vec3;
Vec3 Vec3_create(float x, float y, float z) {
Vec3 v;
v.x = x;
v.y = y;
v.z = z;
return v;
}
Vec3 Vec3_add(Vec3 v1, Vec3 v2) {
Vec3 result;
result.x = v1.x + v2.x;
result.y = v1.y + v2.y;
result.z = v1.z + v2.z;
return result;
}
Vec3 Vec3_subtract(Vec3 v1, Vec3 v2) {
Vec3 result;
result.x = v1.x - v2.x;
result.y = v1.y - v2.y;
result.z = v1.z - v2.z;
return result;
}
Vec3 Vec3_multiply(Vec3 v, float scalar) {
Vec3 result;
result.x = v.x * scalar;
result.y = v.y * scalar;
result.z = v.z * scalar;
return result;
}
Vec3 Vec3_multiply_componentwise(Vec3 v1, Vec3 v2) {
Vec3 result;
result.x = v1.x * v2.x;
result.y = v1.y * v2.y;
result.z = v1.z * v2.z;
return result;
}
float Vec3_dot(Vec3 v1, Vec3 v2) {
return v1.x * v2.x + v1.y * v2.y + v1.z * v2.z;
}
float Vec3_length2(Vec3 v) {
return v.x * v.x + v.y * v.y + v.z * v.z;
}
float Vec3_length(Vec3 v) {
return sqrt(Vec3_length2(v));
}
Vec3 Vec3_normalize(Vec3 v) {
float nor2 = Vec3_length2(v);
if (nor2 > 0) {
float invNor = 1 / sqrt(nor2);
v.x *= invNor;
v.y *= invNor;
v.z *= invNor;
}
return v;
}
typedef struct {
Vec3 center;
float radius, radius2;
Vec3 surfaceColor, emissionColor;
float transparency, reflection;
} Sphere;
Sphere Sphere_create(Vec3 center, float radius, Vec3 surfaceColor, float reflection, float transparency, Vec3 emissionColor) {
Sphere s;
s.center = center;
s.radius = radius;
s.radius2 = radius * radius;
s.surfaceColor = surfaceColor;
s.emissionColor = emissionColor;
s.transparency = transparency;
s.reflection = reflection;
return s;
}
int Sphere_intersect(Sphere s, Vec3 rayorig, Vec3 raydir, float* t0, float* t1) {
Vec3 l = Vec3_subtract(s.center, rayorig);
float tca = Vec3_dot(l, raydir);
if (tca < 0) return 0;
float d2 = Vec3_dot(l, l) - tca * tca;
if (d2 > s.radius2) return 0;
float thc = sqrt(s.radius2 - d2);
*t0 = tca - thc;
*t1 = tca + thc;
return 1;
}
float mix(float a, float b, float mix) {
return b * mix + a * (1 - mix);
}
Vec3 trace(Vec3 rayorig, Vec3 raydir, Sphere* spheres, int numSpheres, int depth) {
float tnear = INFINITY;
Sphere* sphere = NULL;
for (int i = 0; i < numSpheres; ++i) {
float t0 = INFINITY, t1 = INFINITY;
if (Sphere_intersect(spheres[i], rayorig, raydir, &t0, &t1)) {
if (t0 < 0) t0 = t1;
if (t0 < tnear) {
tnear = t0;
sphere = &spheres[i];
}
}
}
if (!sphere) return Vec3_create(2, 2, 2);
Vec3 surfaceColor = Vec3_create(0, 0, 0);
Vec3 phit = Vec3_add(rayorig, Vec3_multiply(raydir, tnear));
Vec3 nhit = Vec3_subtract(phit, sphere->center);
nhit = Vec3_normalize(nhit);
float bias = 1e-4;
int inside = 0;
if (Vec3_dot(raydir, nhit) > 0) {
nhit = Vec3_multiply(nhit, -1);
inside = 1;
}
if ((sphere->transparency > 0 || sphere->reflection > 0) && depth < 5) {
float facingratio = -Vec3_dot(raydir, nhit);
float fresneleffect = mix(pow(1 - facingratio, 3), 1, 0.1);
Vec3 refldir = Vec3_subtract(raydir, Vec3_multiply(nhit, 2 * Vec3_dot(raydir, nhit)));
refldir = Vec3_normalize(refldir);
Vec3 reflection = trace(Vec3_add(phit, Vec3_multiply(nhit, bias)), refldir, spheres, numSpheres, depth + 1);
Vec3 refraction = Vec3_create(0, 0, 0);
if (sphere->transparency) {
float ior = 1.1, eta = (inside) ? ior : 1 / ior;
float cosi = -Vec3_dot(nhit, raydir);
float k = 1 - eta * eta * (1 - cosi * cosi);
Vec3 refrdir = Vec3_add(Vec3_multiply(raydir, eta), Vec3_multiply(nhit, eta * cosi - sqrt(k)));
refrdir = Vec3_normalize(refrdir);
refraction = trace(Vec3_subtract(phit, Vec3_multiply(nhit, bias)), refrdir, spheres, numSpheres, depth + 1);
}
surfaceColor = Vec3_multiply_componentwise(Vec3_add(Vec3_multiply(reflection, fresneleffect), Vec3_multiply(refraction, (1 - fresneleffect) * sphere->transparency)), sphere->surfaceColor);
} else {
for (int i = 0; i < numSpheres; ++i) {
if (spheres[i].emissionColor.x > 0) {
Vec3 transmission = Vec3_create(1, 1, 1);
Vec3 lightDirection = Vec3_subtract(spheres[i].center, phit);
lightDirection = Vec3_normalize(lightDirection);
for (int j = 0; j < numSpheres; ++j) {
if (i != j) {
float t0, t1;
if (Sphere_intersect(spheres[j], Vec3_add(phit, Vec3_multiply(nhit, bias)), lightDirection, &t0, &t1)) {
transmission = Vec3_create(0, 0, 0);
break;
}
}
}
surfaceColor = Vec3_add(surfaceColor, Vec3_multiply_componentwise(Vec3_multiply_componentwise(sphere->surfaceColor, transmission), Vec3_create(fmax(0, Vec3_dot(nhit, lightDirection)), fmax(0, Vec3_dot(nhit, lightDirection)), fmax(0, Vec3_dot(nhit, lightDirection)))));
}
}
}
return Vec3_add(surfaceColor, sphere->emissionColor);
}
void render(Sphere* spheres, int numSpheres) {
unsigned int width = 640, height = 480;
Vec3* image = (Vec3*)malloc(width * height * sizeof(Vec3));
Vec3* pixel = image;
float invWidth = 1 / (float)width;
float invHeight = 1 / (float)height;
float fov = 30;
float aspectratio = width / (float)height;
float angle = tan(M_PI * 0.5 * fov / 180.);
for (unsigned int y = 0; y < height; ++y) {
for (unsigned int x = 0; x < width; ++x, ++pixel) {
float xx = (2 * ((x + 0.5) * invWidth) - 1) * angle * aspectratio;
float yy = (1 - 2 * ((y + 0.5) * invHeight)) * angle;
Vec3 raydir = Vec3_create(xx, yy, -1);
raydir = Vec3_normalize(raydir);
*pixel = trace(Vec3_create(0, 0, 0), raydir, spheres, numSpheres, 0);
}
}
FILE* file = fopen("untitled.ppm", "wb");
fprintf(file, "P6\n%d %d\n255\n", width, height);
for (unsigned int i = 0; i < width * height; ++i) {
fprintf(file, "%c%c%c", (unsigned char)(fmin(1, image[i].x) * 255), (unsigned char)(fmin(1, image[i].y) * 255), (unsigned char)(fmin(1, image[i].z) * 255));
}
fclose(file);
free(image);
}
int main(int argc, char** argv) {
srand(13);
int numSpheres = 6;
Sphere* spheres = (Sphere*)malloc(numSpheres * sizeof(Sphere));
spheres[0] = Sphere_create(Vec3_create(0.0, -10004, -20), 10000, Vec3_create(0.20, 0.20, 0.20), 0, 0.0, Vec3_create(0, 0, 0));
spheres[1] = Sphere_create(Vec3_create(0.0, 0, -20), 4, Vec3_create(1.00, 0.32, 0.36), 1, 0.5, Vec3_create(0, 0, 0));
spheres[2] = Sphere_create(Vec3_create(5.0, -1, -15), 2, Vec3_create(0.90, 0.76, 0.46), 1, 0.0, Vec3_create(0, 0, 0));
spheres[3] = Sphere_create(Vec3_create(5.0, 0, -25), 3, Vec3_create(0.65, 0.77, 0.97), 1, 0.0, Vec3_create(0, 0, 0));
spheres[4] = Sphere_create(Vec3_create(-5.5, 0, -15), 3, Vec3_create(0.90, 0.90, 0.90), 1, 0.0, Vec3_create(0, 0, 0));
spheres[5] = Sphere_create(Vec3_create(0.0, 20, -30), 3, Vec3_create(0.00, 0.00, 0.00), 0, 0.0, Vec3_create(3, 3, 3));
render(spheres, numSpheres);
free(spheres);
return 0;
}