diff --git a/kd_tree.cpp b/kd_tree.cpp index 6054e4e..5fb2c15 100644 --- a/kd_tree.cpp +++ b/kd_tree.cpp @@ -136,11 +136,9 @@ kdTree::~kdTree(){} void kdTree::addPhoton(Photon p) { - Photons.push_back(p); } - void kdTree::createNodeKdTree(treeNode** node, std::vector & originalData , int* xyz, int* yzx, int* zxy, superKey key, int begin, int end, int* xyz_2, int* yzx_2, int* zxy_2) { if(end - begin < 2) @@ -400,9 +398,9 @@ size_t kdTree::getNumPhotons() { return Photons.size(); } -void kdTree::save_photon_list() const { - cout << "Writing photons to \x1b[1;33mphotons.txt\x1b[m" << endl; - ofstream ofs("photons.txt", ios::out); +void kdTree::save_photon_list(const char * file_name) const { + cout << "Writing photons to \x1b[1;33m" << file_name << "\x1b[m" << endl; + ofstream ofs(file_name, ios::out); float r, g, b; for (std::vector::const_iterator it = Photons.begin(); it != Photons.end(); it++) { rgbe2float(r, g, b, (*it).radiance); diff --git a/kd_tree.hpp b/kd_tree.hpp index 14f88dc..06605a7 100644 --- a/kd_tree.hpp +++ b/kd_tree.hpp @@ -165,7 +165,7 @@ public: void printTree(); std::vector findInRange (Vec3 min, Vec3 max) const; void find_by_distance(std::vector & found, const glm::vec3 & point, const glm::vec3 & normal, const float distance, const unsigned int max) const; - void save_photon_list() const; + void save_photon_list(const char * file_name) const; size_t getNumPhotons(); private: diff --git a/main.cpp b/main.cpp index b160b35..aa997f9 100644 --- a/main.cpp +++ b/main.cpp @@ -49,6 +49,7 @@ typedef enum TRACERS { NONE, WHITTED, MONTE_CARLO, JENSEN } tracer_t; static char * g_input_file = NULL; static char * g_photons_file = NULL; +static char * g_caustics_file = NULL; static char * g_out_file_name = NULL; static int g_samples = 25; static float g_fov = 45.0f; @@ -113,17 +114,22 @@ int main(int argc, char ** argv) { } else if(g_tracer == JENSEN) { cout << "Using " << ANSI_BOLD_YELLOW << "Jensen's photon mapping" << ANSI_RESET_STYLE << " with ray tracing." << endl; p_tracer = new PhotonTracer(g_max_depth, g_p_sample_radius); - if (g_photons_file == NULL) { + if (g_photons_file == NULL && g_caustics_file == NULL) { cout << "Building global photon map with " << ANSI_BOLD_YELLOW << g_photons / 2 << ANSI_RESET_STYLE << " primary photons per light source." << endl; - //p_tracer->photon_tracing(scn, g_photons / 2); + p_tracer->photon_tracing(scn, g_photons / 2); cout << "Building caustics photon map with " << ANSI_BOLD_YELLOW << g_photons / 2 << ANSI_RESET_STYLE << " primary photons per light source." << endl; p_tracer->photon_tracing(scn, g_photons / 2, true); p_tracer->build_photon_map(); } else { + if (g_photons_file != g_caustics_file) { + cerr << "Must specify both a photon map file and a caustics file." << endl; + return EXIT_FAILURE; + } p_tracer->build_photon_map(g_photons_file); + p_tracer->build_photon_map(g_caustics_file, true); } tracer = static_cast(p_tracer); - + } else { cerr << "Must specify a ray tracer with \"-t\"." << endl; print_usage(argv); @@ -257,7 +263,7 @@ void parse_args(int argc, char ** const argv) { exit(EXIT_FAILURE); } - while((opt = getopt(argc, argv, "-:t:s:w:f:o:r:g:e:p:h:k:")) != -1) { + while((opt = getopt(argc, argv, "-:t:s:w:f:o:r:g:e:p:h:k:c:")) != -1) { switch (opt) { case 1: g_input_file = (char *)malloc((strlen(optarg) + 1) * sizeof(char)); @@ -374,6 +380,11 @@ void parse_args(int argc, char ** const argv) { strcpy(g_photons_file, optarg); break; + case 'c': + g_caustics_file = (char *)malloc((strlen(optarg) + 1) * sizeof(char)); + strcpy(g_caustics_file, optarg); + break; + case ':': cerr << "Option \"-" << static_cast(optopt) << "\" requires an argument." << endl; print_usage(argv); diff --git a/photon_tracer.cpp b/photon_tracer.cpp index 9f79e32..808c476 100644 --- a/photon_tracer.cpp +++ b/photon_tracer.cpp @@ -217,6 +217,7 @@ void PhotonTracer::photon_tracing(Scene * s, const size_t n_photons_per_ligth, c if (al != NULL) { l_sample = al->sample_at_surface(); s_normal = al->normal_at_last_sample(); + l_sample = l_sample + (BIAS * s_normal); if (!specular || (specular && spec_figures.size() == 0)) { // Generate photon from light source in random direction. @@ -225,39 +226,30 @@ void PhotonTracer::photon_tracing(Scene * s, const size_t n_photons_per_ligth, c h_sample = normalize(sample_hemisphere(r1, r2)); rotate_sample(h_sample, s_normal); } else { - // Generate photon from light source in direction of specular reflective objects. - h_sample = spec_figures[p % spec_figures.size()]->sample_at_surface(); - h_sample = normalize(h_sample - l_sample); + // Generate photon from light source in the direction of specular reflective objects. + h_sample = normalize(spec_figures[p % spec_figures.size()]->sample_at_surface() - l_sample); } // Create the primary photon. - ls = Vec3(l_sample.x, l_sample.y, l_sample.z); - dir = Vec3(h_sample.x, h_sample.y, h_sample.z); power = (al->m_figure->m_mat->m_emission / static_cast(n_photons_per_ligth)); - ph = Photon(ls, dir, power.r, power.g, power.b, 1.0f); - -#pragma omp critical - { - m_photon_map.addPhoton(ph); - } } else if (pl != NULL) { l_sample = glm::vec3(pl->m_position.x, pl->m_position.y, pl->m_position.z); if (!specular || (specular && spec_figures.size() == 0)) { - + h_sample = normalize(sample_sphere(l_sample, 1.0f) - l_sample); } else { - // Generate photon from light source in direction of specular reflective objects. - h_sample = spec_figures[p % spec_figures.size()]->sample_at_surface(); - h_sample = normalize(h_sample - l_sample); + // Generate photon from light source in the direction of specular reflective objects. + h_sample = normalize(spec_figures[p % spec_figures.size()]->sample_at_surface() - l_sample); } - ls = Vec3(l_sample.x, l_sample.y, l_sample.z); - dir = Vec3(h_sample.x, h_sample.y, h_sample.z); power = (pl->m_diffuse / static_cast(n_photons_per_ligth)); - ph = Photon(ls, dir, power.r, power.g, power.b, 1.0f); } + ls = Vec3(l_sample.x, l_sample.y, l_sample.z); + dir = Vec3(h_sample.x, h_sample.y, h_sample.z); + ph = Photon(ls, dir, power.r, power.g, power.b, 1.0f); + trace_photon(ph, s, 0); #pragma omp atomic @@ -267,9 +259,12 @@ void PhotonTracer::photon_tracing(Scene * s, const size_t n_photons_per_ligth, c cout << "\r" << setw(3) << static_cast((static_cast(current) / static_cast(total)) * 100.0) << "% done."; } cout << endl; + + cout << "Generated " << ANSI_BOLD_YELLOW << m_photon_map.getNumPhotons() << ANSI_RESET_STYLE << " total photons." << endl; + m_photon_map.save_photon_list(specular ? "caustics.txt" : "photons.txt"); } -void PhotonTracer::build_photon_map(const char * photons_file) { +void PhotonTracer::build_photon_map(const char * photons_file, const bool caustics) { Photon ph; float x, y, z, dx, dy, dz, r, g, b, rc; ifstream ifs(photons_file, ios::in); @@ -289,14 +284,15 @@ void PhotonTracer::build_photon_map(const char * photons_file) { ifs.close(); - build_photon_map(); + build_photon_map(caustics); } -void PhotonTracer::build_photon_map() { - cout << "Generated " << ANSI_BOLD_YELLOW << m_photon_map.getNumPhotons() << ANSI_RESET_STYLE << " total photons." << endl; - m_photon_map.save_photon_list(); +void PhotonTracer::build_photon_map(const bool caustics) { cout << "Building photon map Kd-tree." << endl; - m_photon_map.buildKdTree(); + if (!caustics) + m_photon_map.buildKdTree(); + else + m_caustics_map.buildKdTree(); } void PhotonTracer::trace_photon(Photon & ph, Scene * s, const unsigned int rec_level) { @@ -328,6 +324,11 @@ void PhotonTracer::trace_photon(Photon & ph, Scene * s, const unsigned int rec_l // Store the diffuse photon and trace. if (!_f->m_mat->m_refract){ +#pragma omp critical + { + m_photon_map.addPhoton(ph); + } + r1 = random01(); r2 = random01(); sample = sample_hemisphere(r1, r2); @@ -338,24 +339,19 @@ void PhotonTracer::trace_photon(Photon & ph, Scene * s, const unsigned int rec_l p_pos = Vec3(i_pos.x, i_pos.y, i_pos.z); p_dir = Vec3(sample.x, sample.y, sample.z); photon = Photon(p_pos, p_dir, color.r, color.g, color.b, ph.ref_index); -#pragma omp critical - { - m_photon_map.addPhoton(photon); - } if (rec_level < m_max_depth) trace_photon(photon, s, rec_level + 1); - } - // Trace the reflected photon. - if (!_f->m_mat->m_refract && _f->m_mat->m_rho > 0.0f && rec_level < m_max_depth) { - color = (_f->m_mat->m_rho) * vec3(red, green, blue); - i_pos += n * BIAS; - p_pos = Vec3(i_pos.x, i_pos.y, i_pos.z); - ph_dir = normalize(reflect(vec3(ph.direction.x, ph.direction.y, ph.direction.z), n)); - p_dir = Vec3(ph_dir.x, ph_dir.y, ph_dir.z); - photon = Photon(p_pos, p_dir, color.r, color.g, color.b, ph.ref_index); - trace_photon(photon, s, rec_level + 1); + if (_f->m_mat->m_rho > 0.0f && rec_level < m_max_depth) { + color = (_f->m_mat->m_rho) * vec3(red, green, blue); + i_pos += n * BIAS; + p_pos = Vec3(i_pos.x, i_pos.y, i_pos.z); + ph_dir = normalize(reflect(vec3(ph.direction.x, ph.direction.y, ph.direction.z), n)); + p_dir = Vec3(ph_dir.x, ph_dir.y, ph_dir.z); + photon = Photon(p_pos, p_dir, color.r, color.g, color.b, ph.ref_index); + trace_photon(photon, s, rec_level + 1); + } } else if (_f->m_mat->m_refract && rec_level < m_max_depth) { // If the material has transmission enabled, calculate the Fresnel term. diff --git a/photon_tracer.hpp b/photon_tracer.hpp index 3dd6364..96afc8f 100644 --- a/photon_tracer.hpp +++ b/photon_tracer.hpp @@ -14,11 +14,12 @@ public: virtual vec3 trace_ray(Ray & r, Scene * s, unsigned int rec_level) const; void photon_tracing(Scene * s, const size_t n_photons_per_ligth = 10000, const bool specular = false); - void build_photon_map(const char * photons_file); - void build_photon_map(); + void build_photon_map(const char * photons_file, const bool caustics = false); + void build_photon_map(const bool caustics = false); private: float m_h_radius; kdTree m_photon_map; + kdTree m_caustics_map; void trace_photon(Photon & ph, Scene * s, const unsigned int rec_level); }; diff --git a/sampling.cpp b/sampling.cpp index 2be6a39..20f049b 100644 --- a/sampling.cpp +++ b/sampling.cpp @@ -70,3 +70,19 @@ void rotate_sample(vec3 & sample, const vec3 & n) { sample.x * nb.y + sample.y * n.y + sample.z * nt.y, sample.x * nb.z + sample.y * n.z + sample.z * nt.z); } + +vec3 sample_sphere(const vec3 center, const float radius) { + float theta; + float u, sqrt1muu, x, y, z; + + // Sampling formula from Wolfram Mathworld: + // http://mathworld.wolfram.com/SpherePointPicking.html + theta = random01()* (2.0f * pi()); + u = (random01() * 2.0f) - 1.0f; + sqrt1muu = glm::sqrt(1.0f - (u * u)); + x = radius * sqrt1muu * cos(theta); + y = radius * sqrt1muu * sin(theta); + z = radius * u; + + return vec3(vec3(x, y, z) + center); +} diff --git a/sampling.hpp b/sampling.hpp index 3860112..aea0eaa 100644 --- a/sampling.hpp +++ b/sampling.hpp @@ -15,5 +15,6 @@ extern vec2 sample_pixel(int i, int j, float w, float h, float a_ratio, float fo extern void create_coords_system(const vec3 &n, vec3 &nt, vec3 &nb); extern vec3 sample_hemisphere(const float r1, float r2); extern void rotate_sample(vec3 & sample, const vec3 & n); +extern vec3 sample_sphere(const vec3 center, const float radius); #endif diff --git a/scenes/scene9.json b/scenes/scene9.json new file mode 100644 index 0000000..a98d993 --- /dev/null +++ b/scenes/scene9.json @@ -0,0 +1,93 @@ +{ + "camera": { + "eye": [0.0, 0.0, 1.0], + "look": [0.0, 0.0, -1.0], + "left": [-1.0, 0.0, 0.0] + }, + + "point_light": { + "position": [0.0, 0.9, -1.0] + }, + + "sphere": { + "position": [-0.4, -0.75, -0.65], + "radius": 0.25, + "material": { + "diffuse": [1.0, 1.0, 1.0], + "rho": 0.4 + } + }, + + "sphere": { + "position": [-0.75, -0.5, -1.5], + "radius": 0.5, + "material": { + "diffuse": [0.0, 0.0, 0.0], + "rho": 1.0 + } + }, + + "sphere": { + "position": [1.0, -0.5, -1.1], + "radius": 0.5, + "material": { + "diffuse": [1.0, 1.0, 0.0], + "transmissive": true, + "ref_index": 1.33 + } + }, + + "plane": { + "position": [0.0, -1.0, 0.0], + "normal": [0.0, 1.0, 0.0], + "material": { + "diffuse": [1.0, 1.0, 1.0], + "specular": [0.0, 0.0, 0.0] + } + }, + + "plane": { + "position": [-2.0, 0.0, 0.0], + "normal": [1.0, 0.0, 0.0], + "material": { + "diffuse": [1.0, 0.0, 0.0], + "specular": [0.0, 0.0, 0.0] + } + }, + + "plane": { + "position": [2.0, 0.0, 0.0], + "normal": [-1.0, 0.0, 0.0], + "material": { + "diffuse": [0.0, 0.0, 1.0], + "specular": [0.0, 0.0, 0.0] + } + }, + + "plane": { + "position": [0.0, 1.0, 0.0], + "normal": [0.0, -1.0, 0.0], + "material": { + "diffuse": [0.0, 1.0, 1.0], + "specular": [0.0, 0.0, 0.0] + } + }, + + "plane": { + "position": [0.0, 0.0, -2.0], + "normal": [0.0, 0.0, 1.0], + "material": { + "diffuse": [1.0, 0.0, 1.0], + "specular": [0.0, 0.0, 0.0] + } + }, + + "plane": { + "position": [0.0, 0.0, 1.1], + "normal": [0.0, 0.0, -1.0], + "material": { + "diffuse": [1.0, 1.0, 0.0], + "specular": [0.0, 0.0, 0.0] + } + } +} diff --git a/sphere.cpp b/sphere.cpp index 951d829..84de230 100644 --- a/sphere.cpp +++ b/sphere.cpp @@ -44,19 +44,7 @@ vec3 Sphere::normal_at_int(Ray & r, float & t) const { } vec3 Sphere::sample_at_surface() const { - float theta; - float u, sqrt1muu, x, y, z; - - // Sampling formula from Wolfram Mathworld: - // http://mathworld.wolfram.com/SpherePointPicking.html - theta = random01()* (2.0f * pi()); - u = (random01() * 2.0f) - 1.0f; - sqrt1muu = glm::sqrt(1.0f - (u * u)); - x = m_radius * sqrt1muu * cos(theta); - y = m_radius * sqrt1muu * sin(theta); - z = m_radius * u; - - return vec3(x, y, z) + m_center; + return sample_sphere(m_center, m_radius); } void Sphere::calculate_inv_area() {