2019-01-22 10:16:05 +01:00
|
|
|
#include <limits>
|
2019-01-20 12:13:53 +01:00
|
|
|
#include <fstream>
|
|
|
|
#include <vector>
|
2019-01-21 21:50:33 +03:00
|
|
|
#include <algorithm>
|
2022-02-21 21:49:18 +01:00
|
|
|
#include <cassert>
|
|
|
|
#include <cmath>
|
|
|
|
|
|
|
|
struct vec3 {
|
|
|
|
float x=0, y=0, z=0;
|
|
|
|
float& operator[](const size_t i) { assert(i<3); return i==0 ? x : (1==i ? y : z); }
|
|
|
|
const float& operator[](const size_t i) const { assert(i<3); return i==0 ? x : (1==i ? y : z); }
|
|
|
|
vec3 operator*(const float v) const { return {x*v, y*v, z*v}; }
|
|
|
|
float operator*(const vec3& v) const { return x*v.x + y*v.y + z*v.z; }
|
|
|
|
vec3 operator+(const vec3& v) const { return {x+v.x, y+v.y, z+v.z}; }
|
|
|
|
vec3 operator-(const vec3& v) const { return {x-v.x, y-v.y, z-v.z}; }
|
|
|
|
vec3 operator-() const { return {-x, -y, -z}; }
|
|
|
|
float norm() const { return std::sqrt(x*x+y*y+z*z); }
|
|
|
|
vec3 normalized() const { return (*this)*(1.f/norm()); }
|
|
|
|
};
|
|
|
|
|
|
|
|
vec3 cross(const vec3 v1, const vec3 v2) {
|
|
|
|
return { v1.y*v2.z - v1.z*v2.y, v1.z*v2.x - v1.x*v2.z, v1.x*v2.y - v1.y*v2.x };
|
|
|
|
}
|
2019-01-20 12:13:53 +01:00
|
|
|
|
2019-01-20 12:16:45 +01:00
|
|
|
struct Light {
|
2021-05-26 15:58:11 +02:00
|
|
|
vec3 position;
|
2019-01-20 12:16:45 +01:00
|
|
|
float intensity;
|
|
|
|
};
|
|
|
|
|
2019-01-20 12:16:23 +01:00
|
|
|
struct Material {
|
2021-05-26 15:58:11 +02:00
|
|
|
float refractive_index = 1;
|
2022-02-21 21:49:18 +01:00
|
|
|
float albedo[4] = {1,0,0,0};
|
2021-05-26 15:58:11 +02:00
|
|
|
vec3 diffuse_color = {0,0,0};
|
|
|
|
float specular_exponent = 0;
|
2019-01-20 12:16:23 +01:00
|
|
|
};
|
|
|
|
|
2019-01-20 12:15:18 +01:00
|
|
|
struct Sphere {
|
2021-05-26 15:58:11 +02:00
|
|
|
vec3 center;
|
2019-01-20 12:15:18 +01:00
|
|
|
float radius;
|
2019-01-20 12:16:23 +01:00
|
|
|
Material material;
|
2019-01-20 12:15:18 +01:00
|
|
|
};
|
|
|
|
|
2022-02-21 21:21:26 +01:00
|
|
|
static const Material ivory = {1.0, {0.6, 0.3, 0.1, 0.0}, {0.4, 0.4, 0.3}, 50.};
|
|
|
|
static const Material glass = {1.5, {0.0, 0.5, 0.1, 0.8}, {0.6, 0.7, 0.8}, 125.};
|
|
|
|
static const Material red_rubber = {1.0, {0.9, 0.1, 0.0, 0.0}, {0.3, 0.1, 0.1}, 10.};
|
|
|
|
static const Material mirror = {1.0, {0.0, 10.0, 0.8, 0.0}, {1.0, 1.0, 1.0}, 1425.};
|
|
|
|
|
|
|
|
static const std::vector<Sphere> spheres = {
|
|
|
|
Sphere{vec3{-3, 0, -16}, 2, ivory},
|
|
|
|
Sphere{vec3{-1.0, -1.5, -12}, 2, glass},
|
|
|
|
Sphere{vec3{ 1.5, -0.5, -18}, 3, red_rubber},
|
|
|
|
Sphere{vec3{ 7, 5, -18}, 4, mirror}
|
|
|
|
};
|
|
|
|
|
|
|
|
static const std::vector<Light> lights = {
|
|
|
|
{{-20, 20, 20}, 1.5},
|
|
|
|
{{ 30, 50, -25}, 1.8},
|
|
|
|
{{ 30, 20, 30}, 1.7}
|
|
|
|
};
|
|
|
|
|
2021-05-26 15:58:11 +02:00
|
|
|
bool ray_sphere_intersect(const vec3 &orig, const vec3 &dir, const Sphere &s, float &t0) {
|
|
|
|
vec3 L = s.center - orig;
|
|
|
|
float tca = L*dir;
|
|
|
|
float d2 = L*L - tca*tca;
|
|
|
|
if (d2 > s.radius*s.radius) return false;
|
2022-02-21 21:21:26 +01:00
|
|
|
float thc = std::sqrt(s.radius*s.radius - d2);
|
2021-05-26 15:58:11 +02:00
|
|
|
t0 = tca - thc;
|
|
|
|
float t1 = tca + thc;
|
|
|
|
if (t0 < 1e-3) t0 = t1; // offset the original point to avoid occlusion by the object itself
|
|
|
|
if (t0 < 1e-3) return false;
|
|
|
|
return true;
|
|
|
|
}
|
|
|
|
|
|
|
|
vec3 reflect(const vec3 &I, const vec3 &N) {
|
2019-01-20 12:17:06 +01:00
|
|
|
return I - N*2.f*(I*N);
|
|
|
|
}
|
|
|
|
|
2021-05-26 15:58:11 +02:00
|
|
|
vec3 refract(const vec3 &I, const vec3 &N, const float eta_t, const float eta_i=1.f) { // Snell's law
|
2019-01-20 12:18:08 +01:00
|
|
|
float cosi = - std::max(-1.f, std::min(1.f, I*N));
|
2019-01-24 10:06:56 +01:00
|
|
|
if (cosi<0) return refract(I, -N, eta_i, eta_t); // if the ray comes from the inside the object, swap the air and the media
|
|
|
|
float eta = eta_i / eta_t;
|
2019-01-20 12:18:08 +01:00
|
|
|
float k = 1 - eta*eta*(1 - cosi*cosi);
|
2021-05-26 15:58:11 +02:00
|
|
|
return k<0 ? vec3{1,0,0} : I*eta + N*(eta*cosi - std::sqrt(k)); // k<0 = total reflection, no ray to refract. I refract it anyways, this has no physical meaning
|
2019-01-20 12:18:08 +01:00
|
|
|
}
|
|
|
|
|
2022-02-21 21:21:26 +01:00
|
|
|
bool scene_intersect(const vec3 &orig, const vec3 &dir, vec3 &hit, vec3 &N, Material &material) {
|
2019-01-20 12:16:23 +01:00
|
|
|
float spheres_dist = std::numeric_limits<float>::max();
|
2021-05-26 15:58:11 +02:00
|
|
|
for (const Sphere &s : spheres) {
|
2019-01-20 12:16:23 +01:00
|
|
|
float dist_i;
|
2021-05-26 15:58:11 +02:00
|
|
|
if (ray_sphere_intersect(orig, dir, s, dist_i) && dist_i < spheres_dist) {
|
2019-01-20 12:16:23 +01:00
|
|
|
spheres_dist = dist_i;
|
|
|
|
hit = orig + dir*dist_i;
|
2022-02-21 21:49:18 +01:00
|
|
|
N = (hit - s.center).normalized();
|
2021-05-26 15:58:11 +02:00
|
|
|
material = s.material;
|
2019-01-20 12:16:23 +01:00
|
|
|
}
|
|
|
|
}
|
2019-01-20 12:27:10 +01:00
|
|
|
|
|
|
|
float checkerboard_dist = std::numeric_limits<float>::max();
|
2021-05-26 15:58:11 +02:00
|
|
|
if (std::abs(dir.y)>1e-3) { // avoid division by zero
|
2019-01-20 12:27:10 +01:00
|
|
|
float d = -(orig.y+4)/dir.y; // the checkerboard plane has equation y = -4
|
2021-05-26 15:58:11 +02:00
|
|
|
vec3 pt = orig + dir*d;
|
|
|
|
if (d>1e-3 && fabs(pt.x)<10 && pt.z<-10 && pt.z>-30 && d<spheres_dist) {
|
2019-01-20 12:27:10 +01:00
|
|
|
checkerboard_dist = d;
|
|
|
|
hit = pt;
|
2021-05-26 15:58:11 +02:00
|
|
|
N = vec3{0,1,0};
|
|
|
|
material.diffuse_color = (int(.5*hit.x+1000) + int(.5*hit.z)) & 1 ? vec3{.3, .3, .3} : vec3{.3, .2, .1};
|
2019-01-20 12:27:10 +01:00
|
|
|
}
|
|
|
|
}
|
|
|
|
return std::min(spheres_dist, checkerboard_dist)<1000;
|
2019-01-20 12:16:23 +01:00
|
|
|
}
|
|
|
|
|
2022-02-21 21:21:26 +01:00
|
|
|
vec3 cast_ray(const vec3 &orig, const vec3 &dir, size_t depth=0) {
|
2021-05-26 15:58:11 +02:00
|
|
|
vec3 point, N;
|
2019-01-20 12:16:23 +01:00
|
|
|
Material material;
|
|
|
|
|
2022-02-21 21:21:26 +01:00
|
|
|
if (depth>4 || !scene_intersect(orig, dir, point, N, material))
|
2021-05-26 15:58:11 +02:00
|
|
|
return vec3{0.2, 0.7, 0.8}; // background color
|
2019-01-20 12:15:18 +01:00
|
|
|
|
2022-02-21 21:49:18 +01:00
|
|
|
vec3 reflect_dir = reflect(dir, N).normalized();
|
|
|
|
vec3 refract_dir = refract(dir, N, material.refractive_index).normalized();
|
2022-02-21 21:21:26 +01:00
|
|
|
vec3 reflect_color = cast_ray(point, reflect_dir, depth + 1);
|
|
|
|
vec3 refract_color = cast_ray(point, refract_dir, depth + 1);
|
2019-01-20 12:17:45 +01:00
|
|
|
|
2019-01-20 12:17:06 +01:00
|
|
|
float diffuse_light_intensity = 0, specular_light_intensity = 0;
|
2021-05-26 15:58:11 +02:00
|
|
|
for (const Light light : lights) {
|
2022-02-21 21:49:18 +01:00
|
|
|
vec3 light_dir = (light.position - point).normalized();
|
2021-05-26 15:58:11 +02:00
|
|
|
|
|
|
|
vec3 shadow_pt, trashnrm;
|
|
|
|
Material trashmat;
|
2022-02-21 21:21:26 +01:00
|
|
|
if (scene_intersect(point, light_dir, shadow_pt, trashnrm, trashmat) &&
|
2021-05-26 15:58:11 +02:00
|
|
|
(shadow_pt-point).norm() < (light.position-point).norm()) // checking if the point lies in the shadow of the light
|
2019-01-20 12:17:25 +01:00
|
|
|
continue;
|
2019-01-20 12:17:06 +01:00
|
|
|
|
2021-05-26 15:58:11 +02:00
|
|
|
diffuse_light_intensity += light.intensity * std::max(0.f, light_dir*N);
|
|
|
|
specular_light_intensity += std::pow(std::max(0.f, -reflect(-light_dir, N)*dir), material.specular_exponent)*light.intensity;
|
2019-01-20 12:16:45 +01:00
|
|
|
}
|
2021-05-26 15:58:11 +02:00
|
|
|
return material.diffuse_color * diffuse_light_intensity * material.albedo[0] + vec3{1., 1., 1.}*specular_light_intensity * material.albedo[1] + reflect_color*material.albedo[2] + refract_color*material.albedo[3];
|
2019-01-20 12:15:18 +01:00
|
|
|
}
|
|
|
|
|
2022-02-21 21:21:26 +01:00
|
|
|
int main() {
|
|
|
|
const int width = 1024;
|
|
|
|
const int height = 768;
|
|
|
|
const float fov = M_PI/3.;
|
2021-05-26 15:58:11 +02:00
|
|
|
std::vector<vec3> framebuffer(width*height);
|
2022-02-21 21:21:26 +01:00
|
|
|
#pragma omp parallel for
|
2019-01-24 10:06:56 +01:00
|
|
|
for (size_t j = 0; j<height; j++) { // actual rendering loop
|
2019-01-20 12:13:53 +01:00
|
|
|
for (size_t i = 0; i<width; i++) {
|
2019-01-24 10:06:56 +01:00
|
|
|
float dir_x = (i + 0.5) - width/2.;
|
|
|
|
float dir_y = -(j + 0.5) + height/2.; // this flips the image at the same time
|
|
|
|
float dir_z = -height/(2.*tan(fov/2.));
|
2022-02-21 21:49:18 +01:00
|
|
|
framebuffer[i+j*width] = cast_ray(vec3{0,0,0}, vec3{dir_x, dir_y, dir_z}.normalized());
|
2019-01-20 12:13:53 +01:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
std::ofstream ofs; // save the framebuffer to file
|
2021-05-26 15:58:11 +02:00
|
|
|
ofs.open("./out.ppm", std::ios::binary);
|
2019-01-20 12:13:53 +01:00
|
|
|
ofs << "P6\n" << width << " " << height << "\n255\n";
|
2021-05-26 15:58:11 +02:00
|
|
|
for (vec3 &c : framebuffer) {
|
2019-01-20 12:17:06 +01:00
|
|
|
float max = std::max(c[0], std::max(c[1], c[2]));
|
|
|
|
if (max>1) c = c*(1./max);
|
2021-05-26 15:58:11 +02:00
|
|
|
ofs << (char)(255 * c[0]) << (char)(255 * c[1]) << (char)(255 * c[2]);
|
2019-01-20 12:13:53 +01:00
|
|
|
}
|
|
|
|
ofs.close();
|
|
|
|
return 0;
|
|
|
|
}
|
2019-01-23 22:41:38 +01:00
|
|
|
|