Added directional (distant light) emitter plugin

docs/conf.py

@@ -59,6 +59,7 @@ def replacement(self, node):
 .. |string| replace:: :paramtype:`string`
 .. |bsdf| replace:: :paramtype:`bsdf`
 .. |point| replace:: :paramtype:`point`
+.. |vector| replace:: :paramtype:`vector`
 .. |transform| replace:: :paramtype:`transform`
 
 .. |enoki| replace:: :monosp:`enoki`

src/emitters/CMakeLists.txt

@@ -1,9 +1,10 @@
 set(MTS_PLUGIN_PREFIX "emitters")
 
-add_plugin(area      area.cpp)
-add_plugin(point     point.cpp)
-add_plugin(constant  constant.cpp)
-add_plugin(envmap    envmap.cpp)
+add_plugin(area        area.cpp)
+add_plugin(point       point.cpp)
+add_plugin(constant    constant.cpp)
+add_plugin(envmap      envmap.cpp)
+add_plugin(directional directional.cpp)
 
 # Register the test directory
 add_tests(${CMAKE_CURRENT_SOURCE_DIR}/tests)

src/emitters/directional.cpp

@@ -0,0 +1,167 @@
+#include <mitsuba/core/bsphere.h>
+#include <mitsuba/core/properties.h>
+#include <mitsuba/core/warp.h>
+#include <mitsuba/render/emitter.h>
+#include <mitsuba/render/scene.h>
+#include <mitsuba/render/texture.h>
+
+NAMESPACE_BEGIN(mitsuba)
+
+/**!
+
+.. _emitter-distant:
+
+Distant directional emitter (:monosp:`directional`)
+---------------------------------------------------
+
+.. pluginparameters::
+
+    * - irradiance
+      - |spectrum|
+      - Spectral irradiance, which corresponds to the amount of spectral power
+        per unit area received by a hypothetical surface normal to the specified
+        direction.
+    * - to_world
+      - |transform|
+      - Emitter-to-world transformation matrix.
+    * - direction
+      - |vector|
+      - Alternative (and exclusive) to `to_world`. Direction towards which the
+        emitter is radiating in world coordinates.
+
+This emitter plugin implements a distant directional source which radiates a
+specified power per unit area along a fixed direction. By default, the emitter
+radiates in the direction of the positive Z axis, i.e. :math:`(0, 0, 1)`.
+
+*/
+
+template <typename Float, typename Spectrum>
+class DirectionalEmitter final : public Emitter<Float, Spectrum> {
+public:
+    MTS_IMPORT_BASE(Emitter, m_flags, m_world_transform)
+    MTS_IMPORT_TYPES(Scene, Texture)
+
+    DirectionalEmitter(const Properties &props) : Base(props) {
+        /* Until `set_scene` is called, we have no information
+           about the scene and default to the unit bounding sphere. */
+        m_bsphere = ScalarBoundingSphere3f(ScalarPoint3f(0.f), 1.f);
+
+        if (props.has_property("direction")) {
+            if (props.has_property("to_world"))
+                Throw("Only one of the parameters 'direction' and 'to_world' "
+                      "can be specified at the same time!'");
+
+            ScalarVector3f direction(normalize(props.vector3f("direction")));
+            auto [up, unused] = coordinate_system(direction);
+
+            m_world_transform =
+                new AnimatedTransform(ScalarTransform4f::look_at(
+                    ScalarPoint3f(0.0f), ScalarPoint3f(direction), up));
+        }
+
+        m_flags      = EmitterFlags::Infinite | EmitterFlags::DeltaDirection;
+        m_irradiance = props.texture<Texture>("irradiance", Texture::D65(1.f));
+    }
+
+    void set_scene(const Scene *scene) override {
+        m_bsphere = scene->bbox().bounding_sphere();
+        m_bsphere.radius =
+            max(math::RayEpsilon<Float>,
+                m_bsphere.radius * (1.f + math::RayEpsilon<Float>) );
+    }
+
+    Spectrum eval(const SurfaceInteraction3f & /*si*/,
+                  Mask /*active*/) const override {
+        return 0.f;
+    }
+
+    std::pair<Ray3f, Spectrum> sample_ray(Float time, Float wavelength_sample,
+                                          const Point2f &spatial_sample,
+                                          const Point2f & /*direction_sample*/,
+                                          Mask active) const override {
+        MTS_MASKED_FUNCTION(ProfilerPhase::EndpointSampleRay, active);
+
+        // 1. Sample spectrum
+        auto [wavelengths, weight] =
+            sample_wavelength<Float, Spectrum>(wavelength_sample);
+
+        // 2. Sample spatial component
+        const Transform4f &trafo = m_world_transform->eval(time, active);
+        Point2f p = warp::square_to_uniform_disk_concentric(spatial_sample);
+        Vector3f perp_offset = trafo.transform_affine(
+            Vector3f{ p.x(), p.y(), 0.f } * m_bsphere.radius);
+
+        // 3. Sample directional component
+        Vector3f d = trafo.transform_affine(Vector3f{ 0.f, 0.f, 1.f });
+
+        return std::make_pair(
+            Ray3f(m_bsphere.center - d * m_bsphere.radius + perp_offset, d,
+                  time, wavelengths),
+            unpolarized<Spectrum>(weight) *
+                (4.f * sqr(math::Pi<Float> * m_bsphere.radius)));
+    }
+
+    std::pair<DirectionSample3f, Spectrum>
+    sample_direction(const Interaction3f &it, const Point2f & /*sample*/,
+                     Mask active) const override {
+        MTS_MASKED_FUNCTION(ProfilerPhase::EndpointSampleDirection, active);
+
+        Vector3f d = m_world_transform->eval(it.time, active)
+                         .transform_affine(Vector3f{ 0.f, 0.f, 1.f });
+        Float dist = 2.f * m_bsphere.radius;
+
+        DirectionSample3f ds;
+        ds.p      = it.p - d * dist;
+        ds.n      = d;
+        ds.uv     = Point2f(0.f);
+        ds.time   = it.time;
+        ds.pdf    = 1.f;
+        ds.delta  = true;
+        ds.object = this;
+        ds.d      = -d;
+        ds.dist   = dist;
+
+        SurfaceInteraction3f si = zero<SurfaceInteraction3f>();
+        si.wavelengths          = it.wavelengths;
+
+        // No need to divide by the PDF here (always equal to 1.f)
+        return std::make_pair(
+            ds, unpolarized<Spectrum>(m_irradiance->eval(si, active)));
+    }
+
+    Float pdf_direction(const Interaction3f & /*it*/,
+                        const DirectionSample3f & /*ds*/,
+                        Mask /*active*/) const override {
+        return 0.f;
+    }
+
+    ScalarBoundingBox3f bbox() const override {
+        /* This emitter does not occupy any particular region
+        of space, return an invalid bounding box */
+        return ScalarBoundingBox3f();
+    }
+
+    void traverse(TraversalCallback *callback) override {
+        callback->put_object("irradiance", m_irradiance.get());
+    }
+
+    std::string to_string() const override {
+        std::ostringstream oss;
+        oss << "DirectionalEmitter[" << std::endl
+            << "  irradiance = " << string::indent(m_irradiance) << ","
+            << std::endl
+            << "  bsphere = " << m_bsphere << "," << std::endl
+            << "]";
+        return oss.str();
+    }
+
+    MTS_DECLARE_CLASS()
+
+protected:
+    ref<Texture> m_irradiance;
+    ScalarBoundingSphere3f m_bsphere;
+};
+
+MTS_IMPLEMENT_CLASS_VARIANT(DirectionalEmitter, Emitter)
+MTS_EXPORT_PLUGIN(DirectionalEmitter, "Distant emitter")
+NAMESPACE_END(mitsuba)

src/emitters/tests/test_constant.py

@@ -95,6 +95,6 @@ def test03_sample_direction(variant_packet_spectral):
     assert ek.allclose(emitter.pdf_direction(it, ds), InvFourPi)
     assert ek.allclose(ds.time, it.time)
 
-    # Evalutate the spectrum (divide by the pdf)
+    # Evaluate the spectrum (divide by the pdf)
     spec = spectrum.eval(it) / warp.square_to_uniform_sphere_pdf(ds.d)
     assert ek.allclose(res, spec)

src/emitters/tests/test_directional.py

@@ -0,0 +1,154 @@
+import mitsuba
+import pytest
+import enoki as ek
+
+
+xml_spectrum = {
+    "d65": """
+        <spectrum version="2.0.0" name="irradiance" type="d65"/>
+        """,
+    "regular": """
+        <spectrum version="2.0.0" name="irradiance" type="regular">
+           <float name="lambda_min" value="500"/>
+           <float name="lambda_max" value="600"/>
+           <string name="values" value="1, 2"/>
+        </spectrum>
+        """,
+}
+
+xml_spectrum_keys = list(set(xml_spectrum.keys()) - {"null"})
+
+
+def make_spectrum(spectrum_key="d65"):
+    from mitsuba.core.xml import load_string
+
+    spectrum = load_string(xml_spectrum[spectrum_key])
+    expanded = spectrum.expand()
+
+    if len(expanded) == 1:
+        spectrum = expanded[0]
+
+    return spectrum
+
+
+def make_emitter(direction=None, spectrum_key="d65"):
+    from mitsuba.core.xml import load_string
+
+    if direction is None:
+        xml_direction = ""
+    else:
+        if type(direction) is not str:
+            direction = ",".join([str(x) for x in direction])
+        xml_direction = \
+            """<vector name="direction" value="{}"/>""".format(direction)
+
+    return load_string("""
+        <emitter version="2.0.0" type="directional">
+            {d}
+            {s}
+        </emitter>
+    """.format(d=xml_direction, s=xml_spectrum[spectrum_key]))
+
+
+def test_construct(variant_scalar_rgb):
+    # Test if the emitter can be constructed as intended
+    emitter = make_emitter()
+    assert not emitter.bbox().valid()  # Degenerate bounding box
+    assert ek.allclose(
+        emitter.world_transform().eval(0.).matrix,
+        [[1, 0, 0, 0],
+         [0, 1, 0, 0],
+         [0, 0, 1, 0],
+         [0, 0, 0, 1]]
+    )  # Identity transform matrix by default
+
+    # Check transform setup correctness
+    emitter = make_emitter(direction=[0, 0, -1])
+    assert ek.allclose(
+        emitter.world_transform().eval(0.).matrix,
+        [[0, 1, 0, 0],
+         [1, 0, 0, 0],
+         [0, 0, -1, 0],
+         [0, 0, 0, 1]]
+    )
+
+
+@pytest.mark.parametrize("spectrum_key", xml_spectrum_keys)
+def test_eval(variant_scalar_spectral, spectrum_key):
+    # Check correctness of the eval() method
+    from mitsuba.core import Vector3f
+    from mitsuba.render import SurfaceInteraction3f
+
+    direction = Vector3f([0, 0, -1])
+    emitter = make_emitter(direction, spectrum_key)
+    spectrum = make_spectrum(spectrum_key)
+
+    # Incident direction in the illuminated direction
+    wi = [0, 0, 1]
+    it = SurfaceInteraction3f()
+    it.p = [0, 0, 0]
+    it.wi = wi
+    assert ek.allclose(emitter.eval(it), 0.)
+
+    # Incident direction off the illuminated direction
+    wi = [0, 0, 1.1]
+    it = SurfaceInteraction3f()
+    it.p = [0, 0, 0]
+    it.wi = wi
+    assert ek.allclose(emitter.eval(it), 0.)
+
+
+@pytest.mark.parametrize("spectrum_key", xml_spectrum_keys)
+@pytest.mark.parametrize("direction", [[0, 0, -1], [1, 1, 1], [0, 0, 1]])
+def test_sample_direction(variant_scalar_spectral, spectrum_key, direction):
+    # Check correctness of sample_direction() and pdf_direction() methods
+
+    from mitsuba.render import SurfaceInteraction3f
+    from mitsuba.core import Vector3f
+
+    direction = Vector3f(direction)
+    emitter = make_emitter(direction, spectrum_key)
+    spectrum = make_spectrum(spectrum_key)
+
+    it = SurfaceInteraction3f.zero()
+    # Some position inside the unit sphere (i.e. within the emitter's default bounding sphere)
+    it.p = [-0.5, 0.3, -0.1]
+    it.time = 1.0
+
+    # Sample direction
+    samples = [0.85, 0.13]
+    ds, res = emitter.sample_direction(it, samples)
+
+    # Direction should point *towards* the illuminated direction
+    assert ek.allclose(ds.d, -direction / ek.norm(direction))
+    assert ek.allclose(ds.pdf, 1.)
+    assert ek.allclose(emitter.pdf_direction(it, ds), 0.)
+    assert ek.allclose(ds.time, it.time)
+
+    # Check spectrum (no attenuation vs distance)
+    spec = spectrum.eval(it)
+    assert ek.allclose(res, spec)
+
+
+@pytest.mark.parametrize("spatial_sample", [[0.85, 0.13], [0.16, 0.50], [0.00, 1.00]])
+@pytest.mark.parametrize("direction", [[0, 0, -1], [1, 1, 1], [0, 0, 1]])
+def test_sample_ray(variant_scalar_rgb, spatial_sample, direction):
+    import enoki as ek
+    from mitsuba.core import Vector2f, Vector3f
+
+    emitter = make_emitter(direction=direction)
+    direction = Vector3f(direction)
+
+    time = 1.0
+    wavelength_sample = 0.3
+    directional_sample = [0.3, 0.2]
+
+    ray, wavelength = emitter.sample_ray(
+        time, wavelength_sample, spatial_sample, directional_sample)
+
+    # Check that ray direction is what is expected
+    assert ek.allclose(ray.d, direction / ek.norm(direction))
+
+    # Check that ray origin is outside of bounding sphere
+    # Bounding sphere is centered at world origin and has radius 1 without scene
+    assert ek.norm(ray.o) >= 1.

src/librender/tests/test_renders.py

@@ -202,6 +202,7 @@ def test_render(variants_all, scene_fname):
             img, var_img = bitmap_extract(bmp)
 
             # Write rendered image to a file
+            os.makedirs(dirname(ref_fname), exist_ok=True)
             xyz_to_rgb_bmp(img).write(ref_fname)
             print('Saved rendered image to: ' + ref_fname)