Add bokeh, aperture ad focus_dist.

3 files changed, 138 insertions, 23 deletions

diff --git a/aa.h b/aa.h
index a1bb51e..2be5f2c 100644
--- a/aa.h
+++ b/aa.h
@@ -26,6 +26,8 @@ class aa_viewport {
     hitlist &world;
     unsigned samples;
     int threads;
+    double aperture;
+    double focus_dist;
 
 public:
 
@@ -35,6 +37,8 @@ public:
                 uint16_t image_height,
                 V screen_hw,
                 V screen_hh,
+                double aperture,
+                double focus_dist,
                 hitlist &world,
                 unsigned samples,
                 int threads = -1) :
@@ -46,7 +50,9 @@ public:
             screen_hh(screen_hh),
             world(world),
             samples(samples),
-            threads((threads > 0) ? threads : (int) std::thread::hardware_concurrency()) {
+            threads((threads > 0) ? threads : (int) std::thread::hardware_concurrency()),
+            aperture{aperture},
+            focus_dist{focus_dist} {
         assert(samples >= 1);
         assert(std::abs(1.0 * image_width / image_height - 1.0 * screen_hw / screen_hh) < 1e-8);
     }
@@ -56,6 +62,8 @@ public:
                 uint16_t image_width,
                 uint16_t image_height,
                 double fov_h,
+                double aperture,
+                double focus_dist,
                 hitlist &world,
                 unsigned samples,
                 int threads = -1) :
@@ -67,7 +75,9 @@ public:
             screen_hh{screen_hw * ((double) image_height / image_width)},
             world(world),
             samples(samples),
-            threads((threads > 0) ? threads : (int) std::thread::hardware_concurrency()) {
+            threads((threads > 0) ? threads : (int) std::thread::hardware_concurrency()),
+            aperture{aperture},
+            focus_dist{focus_dist} {
         assert(samples >= 1);
     }
 
@@ -94,10 +104,11 @@ public:
                         ctx.cxyz, ctx.screen_center,
                         ctx.image_width, ctx.image_height,
                         ctx.screen_hw, ctx.screen_hh,
-                        ctx.world
+                        ctx.aperture, ctx.focus_dist, ctx.world
                 };
                 bias_ctx bc{task.bias_seed};
-                images[tid] = vp.render(task.diffuse_seed, bc);
+                bokeh_ctx bokeh{task.diffuse_seed + 6543210987ULL};
+                images[tid] = vp.render(task.diffuse_seed, bc, bokeh);
             }, s_render_task{
                     .bias_seed=seedgen(), .diffuse_seed=seedgen()
             });
diff --git a/main_simple_scanner.cpp b/main_simple_scanner.cpp
index c5be80c..16a39a0 100644
--- a/main_simple_scanner.cpp
+++ b/main_simple_scanner.cpp
@@ -21,9 +21,10 @@
 // Select the scene to render
 //#define SCENE_DIFFUSE
 //#define SCENE_REFLECT
-#define SCENE_DIALECT
+//#define SCENE_DIALECT
 //#define SCENE_FOV
 //#define SCENE_FREECAM
+#define SCENE_DOF
 
 #ifdef SCENE_FOV
 #define NO_DEFAULT_CAM
@@ -33,13 +34,19 @@
 #define NO_DEFAULT_CAM
 #endif
 
+#ifdef SCENE_DOF
+#define NO_DEFAULT_CAM
+#endif
+
 static constexpr uint64_t default_diffuse_seed = 123456789012345678ULL;
 
 // T: color depth, V: pos
 template<typename T, typename V>
 void generate_image(uint16_t image_width, uint16_t image_height, double viewport_width, double focal_length,
-                    double sphere_z, double sphere_r, unsigned samples, const std::string &caption = "",
-                    unsigned caption_scale = 1) {
+                    double sphere_z, double sphere_r, unsigned samples,
+                    double aperture, double focus_dist,
+                    const std::string &caption = "", unsigned caption_scale = 1) {
+    std::cerr << "Aperture: " << aperture << std::endl;
     if (samples == 1) {
         std::cerr << "Antialiasing is disabled." << std::endl;
     } else {
@@ -56,8 +63,12 @@ void generate_image(uint16_t image_width, uint16_t image_height, double viewport
     basic_viewport<T, V> vp_noaa{
             vec3<V>::zero(), // camera position as the coordinate origin
             vec3d{0, 0, -focal_length},
-            image_width, image_height,
-            viewport_width / 2.0, ((double) image_height / image_width) * viewport_width / 2.0,
+            image_width,
+            image_height,
+            viewport_width / 2.0,
+            ((double) image_height / image_width) * viewport_width / 2.0,
+            aperture,
+            focus_dist,
             world
     };
     ////////////////
@@ -69,6 +80,7 @@ void generate_image(uint16_t image_width, uint16_t image_height, double viewport
             vec3d{0, 0, -focal_length},
             image_width, image_height,
             viewport_width / 2.0, ((double) image_height / image_width) * viewport_width / 2.0,
+            aperture, focus_dist,
             world, samples
     };
     ////////////////
@@ -84,6 +96,8 @@ void generate_image(uint16_t image_width, uint16_t image_height, double viewport
             {0, 0, 0},
             image_width, image_height,
             fov_h,
+            aperture,
+            focus_dist,
             world
     };
     ////////////////
@@ -95,6 +109,8 @@ void generate_image(uint16_t image_width, uint16_t image_height, double viewport
             {0, 0, 0},
             image_width, image_height,
             fov_h,
+            aperture,
+            focus_dist,
             world, samples
     };
     ////////////////
@@ -113,6 +129,8 @@ void generate_image(uint16_t image_width, uint16_t image_height, double viewport
             {0,0,-1},
             image_width, image_height,
             fov_h,
+            aperture,
+            focus_dist,
             world
     };
     ////////////////
@@ -124,6 +142,51 @@ void generate_image(uint16_t image_width, uint16_t image_height, double viewport
             {0,0,-1},
             image_width, image_height,
             fov_h,
+            aperture,
+            focus_dist,
+            world, samples
+    };
+    ////////////////
+
+    material_diffuse_lambertian m_ground{{0.8, 0.8, 0.0}};
+    material_diffuse_lambertian m_ball_center{{0.7, 0.3, 0.3}};
+    material_dielectric m_ball_left{1.5};
+    material_reflective m_ball_right{{0.8, 0.6, 0.2}};
+    // the earth
+    world.add_object(std::make_shared<sphere>(vec3d{0.0, -100.5, -1.0}, 100.0, m_ground));
+    // three balls
+    world.add_object(std::make_shared<sphere>(vec3d{-1.0, 0.0, -1.0}, 0.5, m_ball_left));
+    world.add_object(std::make_shared<sphere>(vec3d{-1.0, 0.0, -1.0}, -0.45, m_ball_left));
+    world.add_object(std::make_shared<sphere>(vec3d{0.0, 0.0, -1.0}, 0.5, m_ball_center));
+    world.add_object(std::make_shared<sphere>(vec3d{1.0, 0.0, -1.0}, 0.5, m_ball_right));
+#endif
+#ifdef SCENE_DOF
+    focus_dist = 1.7320508075688772;
+    // This scene needs custom viewport setting
+    const auto fov_h = 121.28449291441746 * (M_PI / 180.0);
+    ////////////////
+    // noaa rendering
+    bias_ctx no_bias{};
+    basic_viewport<T, V> vp_noaa{
+            {-2, 2, 1}, // camera position as the coordinate origin
+            {-1, 1, 0},
+            image_width, image_height,
+            fov_h,
+            aperture,
+            focus_dist,
+            world
+    };
+    ////////////////
+
+    ////////////////
+    // aa rendering
+    aa_viewport<T, V> vp_aa{
+            {-2, 2, 1}, // camera position as the coordinate origin
+            {-1, 1, 0},
+            image_width, image_height,
+            fov_h,
+            aperture,
+            focus_dist,
             world, samples
     };
     ////////////////
@@ -173,8 +236,9 @@ void generate_image(uint16_t image_width, uint16_t image_height, double viewport
 #endif
     timer tm;
     std::cerr << "Rendering..." << std::endl;
+    bokeh_ctx bokeh{12345678ULL};
     tm.start_measure();
-    auto image = ((samples == 1) ? vp_noaa.render(default_diffuse_seed, no_bias) : vp_aa.render());
+    auto image = ((samples == 1) ? vp_noaa.render(default_diffuse_seed, no_bias, bokeh) : vp_aa.render());
     tm.stop_measure();
     std::cerr << "Applying gamma2..." << std::endl;
     tm.start_measure();
@@ -194,8 +258,8 @@ void generate_image(uint16_t image_width, uint16_t image_height, double viewport
 }
 
 int main(int argc, char **argv) {
-    if (argc != 8 && argc != 9) {
-        printf("Usage: %s <image_width> <image_height> <viewport_width> <focal_length> <sphere_z> <sphere_r> <samples> [<caption>]\n",
+    if (argc != 10 && argc != 11) {
+        printf("Usage: %s <image_width> <image_height> <viewport_width> <focal_length> <sphere_z> <sphere_r> <aperture> <focus_dist> <samples> [<caption>]\n",
                argv[0]);
         return 0;
     }
@@ -205,15 +269,16 @@ int main(int argc, char **argv) {
     std::cerr << "Notice: assertion is disabled." << std::endl;
 #endif
     std::string iw{argv[1]}, ih{argv[2]}, vw{argv[3]}, fl{argv[4]},
-            sz{argv[5]}, sr{argv[6]}, sp{argv[7]}, cap{};
-    if (argc == 9) {
+            sz{argv[5]}, sr{argv[6]}, aper{argv[7]}, fd{argv[8]}, sp{argv[9]}, cap{};
+    if (argc == 11) {
         // with caption
-        cap = std::string{argv[8]};
+        cap = std::string{argv[10]};
     }
     const auto image_width = std::stoul(iw);
     generate_image<uint16_t, double>(image_width, std::stoul(ih),
                                      std::stod(vw), std::stod(fl),
                                      std::stod(sz), std::stod(sr),
-                                     std::stoul(sp), cap,
-                                     std::max((int) (1.0 * image_width * 0.010 / 8), 1));
+                                     std::stoul(sp), std::stod(aper),
+                                     std::stod(fd),
+                                     cap, std::max((int) (1.0 * image_width * 0.010 / 8), 1));
 }
\ No newline at end of file
diff --git a/viewport.h b/viewport.h
index ea6c94b..cd7e2eb 100644
--- a/viewport.h
+++ b/viewport.h
@@ -42,6 +42,33 @@ public:
     }
 };
 
+// Generate a circle of confusion of circle shape. (can be extended to any shape in the future)
+class bokeh_ctx {
+    bool enabled;
+    std::mt19937_64 mt;
+    std::uniform_real_distribution<double> uni{-1.0, 1.0};
+
+public:
+    bokeh_ctx() : enabled{false} {}
+
+    bokeh_ctx(uint64_t seed) : enabled{true}, mt{seed} {}
+
+    void operator()(double &x, double &y) {
+        if (!enabled) {
+            x = 0;
+            y = 0;
+            return;
+        }
+        double x_, y_;
+        do {
+            x_ = uni(mt);
+            y_ = uni(mt);
+        } while (x_ * x_ + y_ * y_ >= 1);
+        x = x_;
+        y = y_;
+    }
+};
+
 // TODO rename to camera
 // Single sampled viewport which supports bias sampling
 // U: color depth, V: pos
@@ -59,6 +86,8 @@ class basic_viewport {
 //    double focus_length; // distance between the focus point and the image screen
     hitlist &world;
     vec3<V> vup{0, 1, 0}; // vector determine the camera rotating
+    double aperture; // radius ratio of the aperture
+    double focus_dist;
 
     inline void check_vup() const {
         // vup must not be parallel with screen_center-cxyz
@@ -69,22 +98,25 @@ public:
 
     basic_viewport(const vec3<V> &cxyz, const vec3<V> &screen_center,
                    uint32_t image_width, uint32_t image_height,
-                   double fov_h, hitlist &world) :
+                   double fov_h, double aperture, double focus_dist, hitlist &world) :
             cxyz{cxyz}, screen_center{screen_center}, image_width{image_width}, image_height{image_height},
             screen_hw{(cxyz - screen_center).norm() * tan((double) fov_h / 2.0)},
             screen_hh{screen_hw * ((double) image_height / image_width)},
-            world{world} {
+            world{world},
+            aperture{aperture}, focus_dist{focus_dist} {
         check_vup();
     }
 
     basic_viewport(const vec3<V> &cxyz, const vec3<V> &screen_center,
                    uint32_t image_width, uint32_t image_height,
                    double screen_hw, double screen_hh,
+                   double aperture, double focus_dist,
                    hitlist &world) :
             cxyz{cxyz}, screen_center{screen_center}, image_width{image_width}, image_height{image_height},
             screen_hw{screen_hw},
             screen_hh{screen_hh},
-            world{world} {
+            world{world},
+            aperture{aperture}, focus_dist{focus_dist} {
         assert(std::abs(1.0 * image_width / image_height - 1.0 * screen_hw / screen_hh) < 1e-8);
         check_vup();
     }
@@ -95,7 +127,7 @@ public:
      * @param by bias on y axis (0.0 <= by < 1.0)
      * @return
      */
-    bitmap<U> render(uint64_t diffuse_seed, bias_ctx &bias
+    bitmap<U> render(uint64_t diffuse_seed, bias_ctx &bias, bokeh_ctx &bokeh
             /* by putting thread-specific parameters in call argument list, make users convenient*/) const {
         // The implementation keep all mutable state in local stack,
         // keeping the class immutable and thread-safe.
@@ -106,10 +138,12 @@ public:
         const int img_hw = image_width / 2, img_hh = image_height / 2;
         // screen plane is determined by coord system x`Vy`, where V is screen_center
         // for variable name we let u := x`, v := y`
-        const auto u = cross(r, vup).unit_vec() * screen_hw, v = cross(u, r).unit_vec() * screen_hh;
+        const auto u0 = cross(r, vup).unit_vec(), v0 = cross(u0, r).unit_vec();
+        const auto u = u0 * screen_hw, v = v0 * screen_hh;
         assert(dot(r, u) < 1e-8);
         assert(dot(r, v) < 1e-8);
         assert(dot(u, v) < 1e-8);
+        const V pof_scale = (V) 1.0 + (V) focus_dist / r.norm();
         // iterate over every pixel on the image
         for (int j = -img_hh + 1; j <= img_hh; ++j) { // axis y, transformation is needed
             for (int i = -img_hw; i < img_hw; ++i) { // axis x
@@ -122,7 +156,12 @@ public:
                 const auto off_v = (1.0 * j + by) / img_hh;
                 const auto off = off_u * u + off_v * v; // offset on screen plane
                 const auto dir = r + off; // direction vector from camera to current pixel on screen
-                ray3d ray{cxyz, dir}; // from camera to pixel (on the viewport)
+                const auto dir_pof = dir * pof_scale; // difference from camera to point on focus plane
+                const auto pof = cxyz + dir_pof; // point on focus plane, the destination
+                double bokeh_u, bokeh_v;
+                bokeh(bokeh_u, bokeh_v);
+                const auto source = cxyz + (aperture * bokeh_u) * u0 + (aperture * bokeh_v) * v0;
+                ray3d ray{source, pof - source}; // from camera to pixel (on the viewport)
                 const auto pixel = world.color<U>(ray, ruvg);
                 const auto x_ = i + img_hw, y_ = -j + img_hh;
                 image.set(x_, y_, pixel);