minor performance improvements

src/main/java/eu/jonahbauer/raytracing/Examples.java

@@ -51,14 +51,13 @@ public class Examples {
     public static @NotNull Example getSimpleScene(int height) {
         if (height <= 0) height = 675;
         return new Example(
-                new Scene(
-                        getSkyBox(),
+                new Scene(getSkyBox(), List.of(
                         new Sphere(new Vec3(0, -100.5, -1.0), 100.0, new LambertianMaterial(new Color(0.8, 0.8, 0.0))),
                         new Sphere(new Vec3(0, 0, -1.2), 0.5, new LambertianMaterial(new Color(0.1, 0.2, 0.5))),
                         new Sphere(new Vec3(-1.0, 0, -1.2), 0.5, new DielectricMaterial(1.5)),
                         new Sphere(new Vec3(-1.0, 0, -1.2), 0.4, new DielectricMaterial(1 / 1.5)),
                         new Sphere(new Vec3(1.0, 0, -1.2), 0.5, new MetallicMaterial(new Color(0.8, 0.6, 0.2), 0.0))
-                ),
+                )),
                 SimpleCamera.builder()
                         .withImage(height * 16 / 9, height)
                         .build()
@@ -119,14 +118,13 @@ public class Examples {
     public static @NotNull Example getSquares(int height) {
         if (height <= 0) height = 600;
         return new Example(
-                new Scene(
-                        getSkyBox(),
+                new Scene(getSkyBox(), List.of(
                         new Parallelogram(new Vec3(-3, -2, 5), new Vec3(0, 0, -4), new Vec3(0, 4, 0), new LambertianMaterial(new Color(1.0, 0.2, 0.2))),
                         new Parallelogram(new Vec3(-2, -2, 0), new Vec3(4, 0, 0), new Vec3(0, 4, 0), new LambertianMaterial(new Color(0.2, 1.0, 0.2))),
                         new Parallelogram(new Vec3(3, -2, 1), new Vec3(0, 0, 4), new Vec3(0, 4, 0), new LambertianMaterial(new Color(0.2, 0.2, 1.0))),
                         new Parallelogram(new Vec3(-2, 3, 1), new Vec3(4, 0, 0), new Vec3(0, 0, 4), new LambertianMaterial(new Color(1.0, 0.5, 0.0))),
                         new Parallelogram(new Vec3(-2, -3, 5), new Vec3(4, 0, 0), new Vec3(0, 0, -4), new LambertianMaterial(new Color(0.2, 0.8, 0.8)))
-                ),
+                )),
                 SimpleCamera.builder()
                         .withImage(height, height)
                         .withFieldOfView(Math.toRadians(80))
@@ -139,12 +137,12 @@ public class Examples {
     public static @NotNull Example getLight(int height) {
         if (height <= 0) height = 225;
         return new Example(
-                new Scene(
+                new Scene(List.of(
                         new Sphere(new Vec3(0, -1000, 0), 1000, new LambertianMaterial(new Color(0.2, 0.2, 0.2))),
                         new Sphere(new Vec3(0, 2, 0), 2, new LambertianMaterial(new Color(0.2, 0.2, 0.2))),
                         new Parallelogram(new Vec3(3, 1, -2), new Vec3(2, 0, 0), new Vec3(0, 2, 0), new DiffuseLight(new Color(4.0, 4.0, 4.0))),
                         new Sphere(new Vec3(0, 7, 0), 2, new DiffuseLight(new Color(4.0, 4.0, 4.0)))
-                ),
+                )),
                 SimpleCamera.builder()
                         .withImage(height * 16 / 9, height)
                         .withFieldOfView(Math.toRadians(20))
@@ -163,15 +161,16 @@ public class Examples {
         var light = new DiffuseLight(new Color(15.0, 15.0, 15.0));
 
         return new Example(
-                new Scene(
-                        new Box(
-                                new AABB(new Vec3(0, 0, 0), new Vec3(555, 555, 555)),
-                                white, white, red, green, white, null
-                        ),
+                new Scene(List.of(
+                        new Box(new Vec3(0, 0, 0), new Vec3(555, 555, 555), white, white, red, green, white, null),
                         new Parallelogram(new Vec3(343, 554, 332), new Vec3(-130, 0, 0), new Vec3(0, 0, -105), light),
-                        new Box(new Vec3(0, 0, 0), new Vec3(165, 330, 165), white).rotateY(Math.toRadians(15)).translate(new Vec3(265, 0, 295)),
-                        new Box(new Vec3(0, 0, 0), new Vec3(165, 165, 165), white).rotateY(Math.toRadians(-18)).translate(new Vec3(130, 0, 65))
-                ),
+                        new Box(new Vec3(0, 0, 0), new Vec3(165, 330, 165), white)
+                                .rotateY(Math.toRadians(15))
+                                .translate(new Vec3(265, 0, 295)),
+                        new Box(new Vec3(0, 0, 0), new Vec3(165, 165, 165), white)
+                                .rotateY(Math.toRadians(-18))
+                                .translate(new Vec3(130, 0, 65))
+                )),
                 SimpleCamera.builder()
                         .withImage(height, height)
                         .withFieldOfView(Math.toRadians(40))
@@ -189,21 +188,22 @@ public class Examples {
         var light = new DiffuseLight(new Color(7.0, 7.0, 7.0));
 
         return new Example(
-                new Scene(
-                        new Box(
-                                new AABB(new Vec3(0, 0, 0), new Vec3(555, 555, 555)),
-                                white, white, red, green, white, null
-                        ),
+                new Scene(List.of(
+                        new Box(new Vec3(0, 0, 0), new Vec3(555, 555, 555), white, white, red, green, white, null),
                         new Parallelogram(new Vec3(113, 554, 127), new Vec3(330, 0, 0), new Vec3(0, 0, 305), light),
                         new ConstantMedium(
-                                new Box(new Vec3(0, 0, 0), new Vec3(165, 330, 165), white).rotateY(Math.toRadians(15)).translate(new Vec3(265, 0, 295)),
+                                new Box(new Vec3(0, 0, 0), new Vec3(165, 330, 165), white)
+                                        .rotateY(Math.toRadians(15))
+                                        .translate(new Vec3(265, 0, 295)),
                                 0.01, new IsotropicMaterial(Color.BLACK)
                         ),
                         new ConstantMedium(
-                                new Box(new Vec3(0, 0, 0), new Vec3(165, 165, 165), white).rotateY(Math.toRadians(-18)).translate(new Vec3(130, 0, 65)),
+                                new Box(new Vec3(0, 0, 0), new Vec3(165, 165, 165), white)
+                                        .rotateY(Math.toRadians(-18))
+                                        .translate(new Vec3(130, 0, 65)),
                                 0.01, new IsotropicMaterial(Color.WHITE)
                         )
-                ),
+                )),
                 SimpleCamera.builder()
                         .withImage(height, height)
                         .withFieldOfView(Math.toRadians(40))
@@ -220,22 +220,17 @@ public class Examples {
         var white = new LambertianMaterial(new Color(.73, .73, .73));
         var green = new LambertianMaterial(new Color(.12, .45, .15));
         var light = new DiffuseLight(new Color(7.0, 7.0, 7.0));
-        var aluminum = new MetallicMaterial(new Color(0.8, 0.85, 0.88));
         var glass = new DielectricMaterial(1.5);
 
+        var room = new Box(new Vec3(0, 0, 0), new Vec3(555, 555, 555), white, white, red, green, white, null);
+        var lamp = new Parallelogram(new Vec3(343, 554, 332), new Vec3(-130, 0, 0), new Vec3(0, 0, -105), light);
+        var box = new Box(new Vec3(0, 0, 0), new Vec3(165, 330, 165), white)
+                .rotateY(Math.toRadians(15))
+                .translate(new Vec3(265, 0, 295));
+        var sphere = new Sphere(new Vec3(190, 90, 190), 90, glass);
+
         return new Example(
-                new Scene(
-                        new Box(
-                                new AABB(new Vec3(0, 0, 0), new Vec3(555, 555, 555)),
-                                white, white, red, green, white, null
-                        ),
-                        new Parallelogram(new Vec3(343, 554, 332), new Vec3(-130, 0, 0), new Vec3(0, 0, -105), light),
-                        new Box(
-                                new AABB(new Vec3(0, 0, 0), new Vec3(165, 330, 165)),
-                                white, white, white, white, white, aluminum
-                        ).rotateY(Math.toRadians(15)).translate(new Vec3(265, 0, 295)),
-                        new Sphere(new Vec3(190, 90, 190), 90, glass)
-                ),
+                new Scene(List.of(room, box), List.of(lamp, sphere)),
                 SimpleCamera.builder()
                         .withImage(height, height)
                         .withFieldOfView(Math.toRadians(40))
@@ -297,10 +292,9 @@ public class Examples {
         if (height <= 0) height = 450;
 
         return new Example(
-                new Scene(
-                        getSkyBox(),
+                new Scene(getSkyBox(), List.of(
                         new Sphere(Vec3.ZERO, 2, new LambertianMaterial(new ImageTexture("/earthmap.jpg")))
-                ),
+                )),
                 SimpleCamera.builder()
                         .withImage(height * 16 / 9, height)
                         .withFieldOfView(Math.toRadians(20))
@@ -316,11 +310,10 @@ public class Examples {
         var material = new LambertianMaterial(new PerlinTexture(4));
 
         return new Example(
-                new Scene(
-                        getSkyBox(),
+                new Scene(getSkyBox(), List.of(
                         new Sphere(new Vec3(0, -1000, 0), 1000, material),
                         new Sphere(new Vec3(0, 2, 0), 2, material)
-                ),
+                )),
                 SimpleCamera.builder()
                         .withImage(height * 16 / 9, height)
                         .withFieldOfView(Math.toRadians(20))

src/main/java/eu/jonahbauer/raytracing/Main.java

@@ -13,6 +13,8 @@ import java.nio.file.Path;
 import java.util.function.IntFunction;
 
 public class Main {
+    public static final boolean DEBUG = false;
+
     public static void main(String[] args) throws IOException {
         var config = Config.parse(args);
         var example = config.example;

src/main/java/eu/jonahbauer/raytracing/math/AABB.java

@@ -35,26 +35,49 @@ public record AABB(@NotNull Vec3 min, @NotNull Vec3 max) {
         return Optional.ofNullable(bbox);
     }
 
+    /**
+     * {@return the range of x values}
+     */
     public @NotNull Range x() {
         return new Range(min.x(), max.x());
     }
 
+    /**
+     * {@return the range of y values}
+     */
     public @NotNull Range y() {
         return new Range(min.y(), max.y());
     }
 
+    /**
+     * {@return the range of z values}
+     */
     public @NotNull Range z() {
         return new Range(min.z(), max.z());
     }
 
+    /**
+     * {@return the center of this bounding box}
+     */
     public @NotNull Vec3 center() {
         return Vec3.average(min, max, 2);
     }
 
+    /**
+     * Expands this bounding box to include the other bounding box.
+     * @param box a bounding box
+     * @return the expanded bounding box
+     */
     public @NotNull AABB expand(@NotNull AABB box) {
         return new AABB(Vec3.min(this.min, box.min), Vec3.max(this.max, box.max));
     }
 
+    /**
+     * Tests whether the {@code ray} intersects this bounding box withing the {@code range}
+     * @param ray a ray
+     * @param range a range of valid {@code t}s
+     * @return {@code true} iff the ray intersects this bounding box, {@code false} otherwise
+     */
     public boolean hit(@NotNull Ray ray, @NotNull Range range) {
         var origin = ray.origin();
         var direction = ray.direction();
@@ -85,6 +108,13 @@ public record AABB(@NotNull Vec3 min, @NotNull Vec3 max) {
         return tlmax < tumin && tumin >= range.min() && tlmax <= range.max();
     }
 
+    /**
+     * Computes the {@code t} values of the intersections of a ray with the axis-aligned planes through a point.
+     * @param corner the point
+     * @param origin the origin point of the ray
+     * @param invDirection the {@linkplain Vec3#inv() inverted} direction of the ray
+     * @return a three-element array of the {@code t} values of the intersection with the yz-, xz- and xy-plane through {@code corner}
+     */
     public static double @NotNull[] intersect(@NotNull Vec3 corner, @NotNull Vec3 origin, @NotNull Vec3 invDirection) {
         return new double[] {
                 (corner.x() - origin.x()) * invDirection.x(),

src/main/java/eu/jonahbauer/raytracing/math/Ray.java

@@ -11,16 +11,6 @@ public record Ray(@NotNull Vec3 origin, @NotNull Vec3 direction) {
     }
 
     public @NotNull Vec3 at(double t) {
-        return new Vec3(
-                Math.fma(t, direction.x(), origin.x()),
-                Math.fma(t, direction.y(), origin.y()),
-                Math.fma(t, direction.z(), origin.z())
-        );
-    }
-
-    public int vmask() {
-        return (direction().x() < 0 ? 1 : 0)
-                | (direction().y() < 0 ? 2 : 0)
-                | (direction().z() < 0 ? 4 : 0);
+        return Vec3.fma(t, direction, origin);
     }
 }

src/main/java/eu/jonahbauer/raytracing/math/Vec3.java

@@ -5,6 +5,8 @@ import org.jetbrains.annotations.NotNull;
 import java.util.Optional;
 import java.util.random.RandomGenerator;
 
+import static eu.jonahbauer.raytracing.Main.DEBUG;
+
 public record Vec3(double x, double y, double z) {
     public static final Vec3 ZERO = new Vec3(0, 0, 0);
     public static final Vec3 MAX = new Vec3(Double.MAX_VALUE, Double.MAX_VALUE, Double.MAX_VALUE);
@@ -14,43 +16,63 @@ public record Vec3(double x, double y, double z) {
     public static final Vec3 UNIT_Z = new Vec3(0, 0, 1);
 
     public Vec3 {
-        assert Double.isFinite(x) && Double.isFinite(y) && Double.isFinite(z) : "x, y and z must be finite";
-    }
-
-    /**
-     * {@return a uniformly random vector with components in the range [-1, 1)}
-     */
-    public static @NotNull Vec3 random(@NotNull RandomGenerator random) {
-        return new Vec3(
-                Math.fma(2, random.nextDouble(), -1),
-                Math.fma(2, random.nextDouble(), -1),
-                Math.fma(2, random.nextDouble(), -1)
-        );
+        if (DEBUG) {
+            if (!Double.isFinite(x) || !Double.isFinite(y) || !Double.isFinite(z)) {
+                throw new IllegalArgumentException("x, y and z must be finite");
+            }
+        }
     }
 
     /**
      * {@return a uniformly random unit vector}
      */
-    public static @NotNull Vec3 random(@NotNull RandomGenerator random, boolean unit) {
-        if (!unit) return random(random);
-        Vec3 vec;
+    public static @NotNull Vec3 random(@NotNull RandomGenerator random) {
+        double x, y, z;
         double squared;
         do {
-            vec = random(random);
-            squared = vec.squared();
+            x = Math.fma(2, random.nextDouble(), -1);
+            y = Math.fma(2, random.nextDouble(), -1);
+            z = Math.fma(2, random.nextDouble(), -1);
+            squared = x * x + y * y + z * z;
         } while (squared > 1);
-        return vec.div(Math.sqrt(squared));
+        var factor = 1 / Math.sqrt(squared);
+        return new Vec3(x * factor, y * factor, z * factor);
     }
 
+    /**
+     * {@return a uniformly random unit vector on the opposite hemisphere of the given <code>direction</code>}
+     */
+    public static @NotNull Vec3 randomOppositeHemisphere(@NotNull RandomGenerator random, @NotNull Vec3 direction) {
+        double x, y, z;
+        double squared;
+        do {
+            x = Math.fma(2, random.nextDouble(), -1);
+            y = Math.fma(2, random.nextDouble(), -1);
+            z = Math.fma(2, random.nextDouble(), -1);
+            squared = x * x + y * y + z * z;
+        } while (squared > 1 || direction.x() * x + direction.y() * y + direction.z() * z >= 0);
+        var factor = 1 / Math.sqrt(squared);
+        return new Vec3(x * factor, y * factor, z * factor);
+    }
+
+    /**
+     * Reflects a vector on the given {@code normal} vector.
+     * @param vec a vector
+     * @param normal the surface normal (must be a unit vector)
+     * @return the reflected vector
+     */
     public static @NotNull Vec3 reflect(@NotNull Vec3 vec, @NotNull Vec3 normal) {
         var factor = - 2 * normal.times(vec);
-        return new Vec3(
-                Math.fma(factor, normal.x(), vec.x()),
-                Math.fma(factor, normal.y(), vec.y()),
-                Math.fma(factor, normal.z(), vec.z())
-        );
+        return Vec3.fma(factor, normal, vec);
     }
 
+    /**
+     * Refracts a vector on the given {@code normal} vector.
+     * @param vec a vector
+     * @param normal the surface normal (must be a unit vector)
+     * @param ri the refractive index
+     * @return the refracted vector
+     */
     public static @NotNull Optional<Vec3> refract(@NotNull Vec3 vec, @NotNull Vec3 normal, double ri) {
         vec = vec.unit();
         var cosTheta = Math.min(- vec.times(normal), 1.0);
@@ -62,16 +84,35 @@ public record Vec3(double x, double y, double z) {
         return Optional.of(rOutPerp.plus(rOutParallel));
     }
 
+    /**
+     * Rotates a vector around an {@code axis}.
+     * @param vec a vector
+     * @param axis the rotation axis
+     * @param angle the angle in radians
+     * @return the rotated vector
+     */
     public static @NotNull Vec3 rotate(@NotNull Vec3 vec, @NotNull Vec3 axis, double angle) {
         Vec3 vxp = axis.cross(vec);
         Vec3 vxvxp = axis.cross(vxp);
         return vec.plus(vxp.times(Math.sin(angle))).plus(vxvxp.times(1 - Math.cos(angle)));
     }
 
+    /**
+     * {@return the euclidean distance between two vectors}
+     * @param a a vector
+     * @param b another vector
+     */
     public static double distance(@NotNull Vec3 a, @NotNull Vec3 b) {
         return a.minus(b).length();
     }
 
+    /**
+     * Computes a running average of vectors.
+     * @param current the current running average
+     * @param next the next vector
+     * @param index the one-based index of the next vector
+     * @return the new running average
+     */
     public static @NotNull Vec3 average(@NotNull Vec3 current, @NotNull Vec3 next, int index) {
         var factor = 1d / index;
         return new Vec3(
@@ -81,6 +122,11 @@ public record Vec3(double x, double y, double z) {
         );
     }
 
+    /**
+     * {@return a component-wise maximum vector}
+     * @param a a vector
+     * @param b another vector
+     */
     public static @NotNull Vec3 max(@NotNull Vec3 a, @NotNull Vec3 b) {
         return new Vec3(
                 Math.max(a.x(), b.x()),
@@ -89,6 +135,11 @@ public record Vec3(double x, double y, double z) {
         );
     }
 
+    /**
+     * {@return a component-wise minimum vector}
+     * @param a a vector
+     * @param b another vector
+     */
     public static @NotNull Vec3 min(@NotNull Vec3 a, @NotNull Vec3 b) {
         return new Vec3(
                 Math.min(a.x(), b.x()),
@@ -97,6 +148,24 @@ public record Vec3(double x, double y, double z) {
         );
     }
 
+    /**
+     * {@return <code>a * b + c</code>}
+     * @param a scalar
+     * @param b a vector
+     * @param c another vector
+     */
+    public static @NotNull Vec3 fma(double a, @NotNull Vec3 b, @NotNull Vec3 c) {
+        return new Vec3(
+                Math.fma(a, b.x(), c.x()),
+                Math.fma(a, b.y(), c.y()),
+                Math.fma(a, b.z(), c.z())
+        );
+    }
+
+    public static double tripleProduct(@NotNull Vec3 a, @NotNull Vec3 b, @NotNull Vec3 c) {
+        return a.x * b.y * c.z + a.y * b.z * c.x + a.z * b.x * c.y - c.x * b.y * a.z - c.y * b.z * a.x - c.z * b.x * a.y;
+    }
+
     public @NotNull Vec3 plus(double x, double y, double z) {
         return new Vec3(this.x + x, this.y + y, this.z + z);
     }
@@ -105,59 +174,115 @@ public record Vec3(double x, double y, double z) {
         return new Vec3(this.x - x, this.y - y, this.z - z);
     }
 
-    public @NotNull Vec3 plus(@NotNull Vec3 b) {
-        return new Vec3(this.x + b.x, this.y + b.y, this.z + b.z);
-    }
-    
-    public @NotNull Vec3 minus(@NotNull Vec3 b) {
-        return new Vec3(this.x - b.x, this.y - b.y, this.z - b.z);
-    }
-    
-    public double times(@NotNull Vec3 b) {
-        return this.x * b.x + this.y * b.y + this.z * b.z;
-    }
-    
-    public @NotNull Vec3 times(double b) {
-        return new Vec3(this.x * b, this.y * b, this.z * b);
+    /**
+     * Adds a vector to this vector
+     * @param other a vector
+     * @return the sum of this and the other vector
+     */
+    public @NotNull Vec3 plus(@NotNull Vec3 other) {
+        return new Vec3(this.x + other.x, this.y + other.y, this.z + other.z);
     }
 
+    /**
+     * Subtracts a vector from this vector
+     * @param other a vector
+     * @return the difference of this and the other vector
+     */
+    public @NotNull Vec3 minus(@NotNull Vec3 other) {
+        return new Vec3(this.x - other.x, this.y - other.y, this.z - other.z);
+    }
+
+    /**
+     * Computes the scalar product of this and another vector
+     * @param other a vector
+     * @return the scalar product
+     */
+    public double times(@NotNull Vec3 other) {
+        return this.x * other.x + this.y * other.y + this.z * other.z;
+    }
+
+    /**
+     * Multiplies this vector with a scalar
+     * @param t a scalar
+     * @return the product of this vector and the scalar
+     */
+    public @NotNull Vec3 times(double t) {
+        return new Vec3(this.x * t, this.y * t, this.z * t);
+    }
+
+    /**
+     * Negates this vector.
+     * {@return the negated vector}
+     */
     public @NotNull Vec3 neg() {
         return new Vec3(-x, -y, -z);
     }
 
+    /**
+     * Inverts each component of this vector.
+     * @return the inverted vector.
+     */
     public @NotNull Vec3 inv() {
         return new Vec3(1 / x, 1 / y, 1 / z);
     }
 
-    public @NotNull Vec3 cross(@NotNull Vec3 b) {
+    /**
+     * Computes the cross-product of this and another vector
+     * @param other a vector
+     * @return the cross-product
+     */
+    public @NotNull Vec3 cross(@NotNull Vec3 other) {
         return new Vec3(
-                this.y() * b.z() - b.y() * this.z(),
-                this.z() * b.x() - b.z() * this.x(),
-                this.x() * b.y() - b.x() * this.y()
+                Math.fma(this.y, other.z, - other.y * this.z),
+                Math.fma(this.z, other.x, - other.z * this.x),
+                Math.fma(this.x, other.y, - other.x * this.y)
         );
     }
 
-    public @NotNull Vec3 div(double b) {
-        return new Vec3(this.x / b, this.y / b, this.z / b);
+    /**
+     * Divides this vector by a scalar
+     * @param t a scalar
+     * @return this vector divided by the scalar
+     */
+    public @NotNull Vec3 div(double t) {
+        return new Vec3(this.x / t, this.y / t, this.z / t);
     }
 
+    /**
+     * {@return the squared length of this vector}
+     */
     public double squared() {
         return this.x * this.x + this.y * this.y + this.z * this.z;
     }
-    
+
+    /**
+     * {@return the length of this vector}
+     */
     public double length() {
         return Math.sqrt(squared());
     }
 
+    /**
+     * {@return whether this vector is near zero}
+     */
     public boolean isNearZero() {
         var s = 1e-8;
         return Math.abs(x) < s && Math.abs(y) < s && Math.abs(z) < s;
     }
-    
+
+    /**
+     * {@return a unit vector with the same direction as this vector}
+     */
     public @NotNull Vec3 unit() {
-        return div(length());
+        var squared = squared();
+        if (squared == 1) return this;
+        return div(Math.sqrt(squared));
     }
 
+    /**
+     * {@return the n-th component of this vector}
+     * @param axis the component index
+     */
     public double get(int axis) {
         return switch (axis) {
             case 0 -> x;
@@ -178,4 +303,9 @@ public record Vec3(double x, double y, double z) {
     public @NotNull Vec3 withZ(double z) {
         return new Vec3(x, y, z);
     }
+
+    @Override
+    public @NotNull String toString() {
+        return "(" + x + "," + y + "," + z + ")";
+    }
 }

src/main/java/eu/jonahbauer/raytracing/render/camera/SimpleCamera.java

@@ -60,7 +60,7 @@ public final class SimpleCamera implements Camera {
 
         this.pixel00 = origin.plus(direction.times(builder.focusDistance))
                 .minus(u.times(0.5 * viewportWidth)).minus(v.times(0.5 * viewportHeight))
-                .plus(pixelU.div(2)).plus(pixelV.div(2));
+                .plus(pixelU.times(.5)).plus(pixelV.times(.5));
     }
 
     /**

src/main/java/eu/jonahbauer/raytracing/render/material/MetallicMaterial.java

@@ -25,7 +25,7 @@ public record MetallicMaterial(@NotNull Texture texture, double fuzz) implements
     public @NotNull Optional<ScatterResult> scatter(@NotNull Ray ray, @NotNull HitResult hit, @NotNull RandomGenerator random) {
         var newDirection = Vec3.reflect(ray.direction(), hit.normal());
         if (fuzz > 0) {
-            newDirection = newDirection.unit().plus(Vec3.random(random, true).times(fuzz));
+            newDirection = Vec3.fma(fuzz, Vec3.random(random), newDirection.unit());
         }
         var attenuation = texture.get(hit);
         return Optional.of(new SpecularScatterResult(attenuation, new Ray(hit.position(), newDirection)));

src/main/java/eu/jonahbauer/raytracing/render/renderer/SimpleRenderer.java

@@ -16,6 +16,8 @@ import java.util.SplittableRandom;
 import java.util.random.RandomGenerator;
 import java.util.stream.IntStream;
 
+import static eu.jonahbauer.raytracing.Main.DEBUG;
+
 public final class SimpleRenderer implements Renderer {
     private final int sqrtSamplesPerPixel;
     private final int maxDepth;
@@ -62,7 +64,7 @@ public final class SimpleRenderer implements Renderer {
      * a time and updating the canvas after each sample.
      */
     private void renderIterative(@NotNull Camera camera, @NotNull Scene scene, @NotNull Canvas canvas) {
-        var random = new Random();
+        var random = new Random(0);
 
         // render one sample after the other
         int i = 0;
@@ -94,7 +96,7 @@ public final class SimpleRenderer implements Renderer {
      * per pixel and updating the canvas after each pixel.
      */
     private void renderNonIterative(@NotNull Camera camera, @NotNull Scene scene, @NotNull Canvas canvas) {
-        var splittable = new SplittableRandom();
+        var splittable = new SplittableRandom(0);
         // render one pixel after the other
         getScanlineStream(camera.getHeight(), parallel).forEach(y -> {
             var random = splittable.split();
@@ -104,6 +106,9 @@ public final class SimpleRenderer implements Renderer {
                 for (int sj = 0; sj < sqrtSamplesPerPixel; sj++) {
                     for (int si = 0; si < sqrtSamplesPerPixel; si++) {
                         var ray = camera.cast(x, y, si, sj, sqrtSamplesPerPixel, random);
+                        if (DEBUG) {
+                            System.out.println("Casting ray " + ray + " through pixel (" + x + "," + y + ") at subpixel (" + si + "," + sj + ")...");
+                        }
                         var c = getColor(scene, ray, random);
                         color = Color.average(color, c, ++i);
                     }
@@ -125,38 +130,82 @@ public final class SimpleRenderer implements Renderer {
         var attenuation = Color.WHITE;
 
         while (depth-- > 0) {
-            var optional = scene.hit(ray, new Range(0.001, Double.POSITIVE_INFINITY));
+            var optional = scene.hit(ray);
             if (optional.isEmpty()) {
-                color = Color.add(color, Color.multiply(attenuation, scene.getBackgroundColor(ray)));
+                var background = scene.getBackgroundColor(ray);
+                color = Color.add(color, Color.multiply(attenuation, background));
+                if (DEBUG) {
+                    System.out.println("  Hit background: " + background);
+                }
                 break;
             }
 
             var hit = optional.get();
+            if (DEBUG) {
+                System.out.println("  Hit " + hit.target() + " at t=" + hit.t() + " (" + hit.position() + ")");
+            }
             var material = hit.material();
             var emitted = material.emitted(hit);
+            if (DEBUG && !Color.BLACK.equals(emitted)) {
+                System.out.println("  Emitted: " + emitted);
+            }
+
             var result = material.scatter(ray, hit, random);
             color = Color.add(color, Color.multiply(attenuation, emitted));
 
-            if (result.isEmpty()) break;
+            if (result.isEmpty()) {
+                if (DEBUG) {
+                    System.out.println("  Absorbed");
+                }
+                break;
+            }
 
             switch (result.get()) {
                 case Material.SpecularScatterResult(var a, var scattered) -> {
                     attenuation = Color.multiply(attenuation, a);
                     ray = scattered;
+
+                    if (DEBUG) {
+                        System.out.println("  Specular scattering with albedo " + a);
+                    }
                 }
                 case Material.PdfScatterResult(var a, var pdf) -> {
-                    if (scene.getLights() == null) {
+                    if (scene.getTargets() == null) {
                         attenuation = Color.multiply(attenuation, a);
                         ray = new Ray(hit.position(), pdf.generate(random));
+
+                        if (DEBUG) {
+                            System.out.println("  Pdf scattering with albedo " + a);
+                        }
                     } else {
-                        var mixed = new MixtureProbabilityDensityFunction(new TargetingProbabilityDensityFunction(hit.position(), scene.getLights()), pdf);
-                        var direction = mixed.generate(random);
-                        var factor = pdf.value(direction) / mixed.value(direction);
+                        var mixed = new MixtureProbabilityDensityFunction(new TargetingProbabilityDensityFunction(hit.position(), scene.getTargets()), pdf, 0.5);
+                        var direction = mixed.generate(random).unit();
+
+                        var idealPdf = pdf.value(direction);
+                        var actualPdf = mixed.value(direction);
+                        if (actualPdf == 0) break; // when actualPdf is 0, the ray should have never been generated by mixed.generate
+
+                        var factor =  idealPdf / actualPdf;
+
                         attenuation = Color.multiply(attenuation, Color.multiply(a, factor));
                         ray = new Ray(hit.position(), direction);
+
+                        if (DEBUG) {
+                            System.out.println("  Pdf scattering with albedo " + a + " and factor " + factor);
+                        }
                     }
                 }
             }
+
+            if (DEBUG) {
+                System.out.println("  Combined color is " + color);
+                System.out.println("  Combined attenuation is " + attenuation);
+                System.out.println("  New ray is " + ray);
+            }
+        }
+
+        if (DEBUG) {
+            System.out.println("  Final color is " + color);
         }
 
         return color;

src/main/java/eu/jonahbauer/raytracing/render/renderer/pdf/CosineProbabilityDensityFunction.java

@@ -15,13 +15,13 @@ public record CosineProbabilityDensityFunction(@NotNull Vec3 normal) implements
 
     @Override
     public double value(@NotNull Vec3 direction) {
-        var cos = normal.times(direction.unit());
+        var cos = normal.times(direction);
         return Math.max(0, cos / Math.PI);
     }
 
     @Override
     public @NotNull Vec3 generate(@NotNull RandomGenerator random) {
-        var out = normal().plus(Vec3.random(random, true));
+        var out = normal().plus(Vec3.random(random));
         return out.isNearZero() ? normal() : out;
     }
 }

src/main/java/eu/jonahbauer/raytracing/render/renderer/pdf/MixtureProbabilityDensityFunction.java

@@ -6,6 +6,13 @@ import org.jetbrains.annotations.NotNull;
 import java.util.Objects;
 import java.util.random.RandomGenerator;
 
+/**
+ * Mixes between two probability density functions (pdf) using a weight. When the weight is closer to zero, the
+ * influence of the second pdf is stronger. When the weight is closer to one, the influence of the first pdf is stronger.
+ * @param a the first probability density function
+ * @param b the second probability density function
+ * @param weight a weight in the range [0, 1]
+ */
 public record MixtureProbabilityDensityFunction(
         @NotNull ProbabilityDensityFunction a,
         @NotNull ProbabilityDensityFunction b,
@@ -23,7 +30,9 @@ public record MixtureProbabilityDensityFunction(
 
     @Override
     public double value(@NotNull Vec3 direction) {
-        return weight * a.value(direction) + (1 - weight) * b.value(direction);
+        var v = a.value(direction);
+        var w = b.value(direction);
+        return Math.fma(weight, v, Math.fma(-weight, w, w));
     }
 
     @Override

src/main/java/eu/jonahbauer/raytracing/render/renderer/pdf/ProbabilityDensityFunction.java

@@ -5,7 +5,21 @@ import org.jetbrains.annotations.NotNull;
 
 import java.util.random.RandomGenerator;
 
+/**
+ * A probability density function used for sampling random directions when scattering a ray.
+ */
 public interface ProbabilityDensityFunction {
+
+    /**
+     * {@return the value of this probability density function at the given point}
+     * @param direction the direction
+     */
     double value(@NotNull Vec3 direction);
+
+    /**
+     * Generates a random direction based on this probability density function.
+     * @param random a random number generator
+     * @return the random direction
+     */
     @NotNull Vec3 generate(@NotNull RandomGenerator random);
 }

src/main/java/eu/jonahbauer/raytracing/render/renderer/pdf/SphereProbabilityDensityFunction.java

@@ -5,6 +5,9 @@ import org.jetbrains.annotations.NotNull;
 
 import java.util.random.RandomGenerator;
 
+/**
+ * A probability density function sampling the sphere uniformly.
+ */
 public record SphereProbabilityDensityFunction() implements ProbabilityDensityFunction {
 
     @Override
@@ -14,6 +17,6 @@ public record SphereProbabilityDensityFunction() implements ProbabilityDensityFu
 
     @Override
     public @NotNull Vec3 generate(@NotNull RandomGenerator random) {
-        return Vec3.random(random, true);
+        return Vec3.random(random);
     }
 }

src/main/java/eu/jonahbauer/raytracing/render/renderer/pdf/TargetingProbabilityDensityFunction.java

@@ -4,25 +4,38 @@ import eu.jonahbauer.raytracing.math.Vec3;
 import eu.jonahbauer.raytracing.scene.Target;
 import org.jetbrains.annotations.NotNull;
 
+import java.util.ArrayList;
+import java.util.List;
 import java.util.Objects;
 import java.util.random.RandomGenerator;
 
 /**
  * A probability density function targeting a target.
+ * @see Target
  */
-public record TargetingProbabilityDensityFunction(@NotNull Vec3 origin, @NotNull Target target) implements ProbabilityDensityFunction {
-    public TargetingProbabilityDensityFunction {
-        Objects.requireNonNull(origin, "origin");
-        Objects.requireNonNull(target, "target");
+public final class TargetingProbabilityDensityFunction implements ProbabilityDensityFunction {
+    private final @NotNull Vec3 origin;
+    private final @NotNull List<@NotNull Target> targets;
+
+    public TargetingProbabilityDensityFunction(@NotNull Vec3 origin, @NotNull List<@NotNull Target> targets) {
+        this.origin = Objects.requireNonNull(origin, "origin");
+        this.targets = new ArrayList<>(targets);
     }
 
     @Override
     public double value(@NotNull Vec3 direction) {
-        return target.getProbabilityDensity(origin, direction);
+        var weight = 1d / targets.size();
+        var sum = 0.0;
+
+        for (var target : targets) {
+            sum = Math.fma(weight, target.getProbabilityDensity(origin, direction), sum);
+        }
+
+        return sum;
     }
 
     @Override
     public @NotNull Vec3 generate(@NotNull RandomGenerator random) {
-        return target.getTargetingDirection(origin, random);
+        return targets.get(random.nextInt(targets.size())).getTargetingDirection(origin, random);
     }
 }

src/main/java/eu/jonahbauer/raytracing/render/texture/CheckerTexture.java

@@ -5,7 +5,6 @@ import org.jetbrains.annotations.NotNull;
 
 public record CheckerTexture(double scale, @NotNull Texture even, @NotNull Texture odd) implements Texture {
 
-
     @Override
     public @NotNull Color get(double u, double v, @NotNull Vec3 p) {
         var x = (int) Math.floor(p.x() / scale);

src/main/java/eu/jonahbauer/raytracing/render/texture/Color.java

@@ -1,11 +1,15 @@
 package eu.jonahbauer.raytracing.render.texture;
 
+import eu.jonahbauer.raytracing.math.Ray;
 import eu.jonahbauer.raytracing.math.Vec3;
+import eu.jonahbauer.raytracing.scene.SkyBox;
 import org.jetbrains.annotations.NotNull;
 
 import java.util.Random;
 
-public record Color(double r, double g, double b) implements Texture {
+import static eu.jonahbauer.raytracing.Main.DEBUG;
+
+public record Color(double r, double g, double b) implements Texture, SkyBox {
     public static final @NotNull Color BLACK = new Color(0.0, 0.0, 0.0);
     public static final @NotNull Color WHITE = new Color(1.0, 1.0, 1.0);
     public static final @NotNull Color RED = new Color(1.0, 0.0, 0.0);
@@ -82,6 +86,17 @@ public record Color(double r, double g, double b) implements Texture {
         this(red / 255f, green / 255f, blue / 255f);
     }
 
+    public Color {
+        if (DEBUG) {
+            if (!Double.isFinite(r) || !Double.isFinite(g) || !Double.isFinite(b)) {
+                throw new IllegalArgumentException("r, g and b must be finite");
+            }
+            if (r < 0 || g < 0 || b < 0) {
+                throw new IllegalArgumentException("r, g and b must be non-negative");
+            }
+        }
+    }
+
     public int red() {
         return toInt(r);
     }
@@ -99,6 +114,11 @@ public record Color(double r, double g, double b) implements Texture {
         return this;
     }
 
+    @Override
+    public @NotNull Color getColor(@NotNull Ray ray) {
+        return this;
+    }
+
     @Override
     public boolean isUVRequired() {
         return false;

src/main/java/eu/jonahbauer/raytracing/render/texture/PerlinTexture.java

@@ -53,7 +53,7 @@ public final class PerlinTexture implements Texture {
         this.mask = count - 1;
         this.randvec = new Vec3[count];
         for (int i = 0; i < count; i++) {
-            this.randvec[i] = Vec3.random(random, true);
+            this.randvec[i] = Vec3.random(random);
         }
         this.permX = generatePerm(count, random);
         this.permY = generatePerm(count, random);

src/main/java/eu/jonahbauer/raytracing/render/texture/Texture.java

@@ -5,12 +5,27 @@ import eu.jonahbauer.raytracing.scene.HitResult;
 import org.jetbrains.annotations.NotNull;
 
 public interface Texture {
+    /**
+     * {@return the color of <code>this</code> texture for a hit}
+     */
     default @NotNull Color get(@NotNull HitResult hit) {
         return get(hit.u(), hit.v(), hit.position());
     }
 
+    /**
+     * {@return the color of <code>this</code> texture at the specified position}
+     * @param u the texture u coordinate
+     * @param v the texture v coordinate
+     * @param p the position
+     */
     @NotNull Color get(double u, double v, @NotNull Vec3 p);
 
+    /**
+     * Returns whether {@link #get(double, double, Vec3)} uses the {@code u} and/or {@code v} parameters.
+     * When a texture indicates that the {@code u} and {@code v} coordinates are not required, the calculation may be
+     * skipped and {@link Double#NaN} will be passed.
+     * @return whether {@link #get(double, double, Vec3)} uses the {@code u} and/or {@code v} parameters
+     */
     default boolean isUVRequired() {
         return true;
     }

src/main/java/eu/jonahbauer/raytracing/scene/HitResult.java

@@ -1,13 +1,27 @@
 package eu.jonahbauer.raytracing.scene;
 
+import eu.jonahbauer.raytracing.math.Range;
+import eu.jonahbauer.raytracing.math.Ray;
 import eu.jonahbauer.raytracing.math.Vec3;
 import eu.jonahbauer.raytracing.render.material.Material;
+import eu.jonahbauer.raytracing.render.texture.Texture;
 import org.jetbrains.annotations.NotNull;
 
 import java.util.Objects;
 
+/**
+ * The result of a {@linkplain Hittable#hit(Ray, Range) hit}.
+ * @param t the {@code t} value at which the hit occurs
+ * @param position the position of the hit
+ * @param normal the surface normal at the hit position
+ * @param target the hit target (for debug purposes only)
+ * @param material the material of the surface
+ * @param u the texture u coordinate (or {@code Double.NaN} if the {@linkplain Material#texture() material's texture} does {@linkplain Texture#isUVRequired() not depend} on the uv-coordinates)
+ * @param v the texture v coordinate (or {@code Double.NaN} if the {@linkplain Material#texture() material's texture} does {@linkplain Texture#isUVRequired() not depend} on the uv-coordinates)
+ * @param isFrontFace whether the front or the back of the surface was it
+ */
 public record HitResult(
-        double t, @NotNull Vec3 position, @NotNull Vec3 normal,
+        double t, @NotNull Vec3 position, @NotNull Vec3 normal, @NotNull Hittable target,
         @NotNull Material material, double u, double v, boolean isFrontFace
 ) implements Comparable<HitResult> {
     public HitResult {
@@ -16,7 +30,7 @@ public record HitResult(
     }
 
     public @NotNull HitResult withPositionAndNormal(@NotNull Vec3 position, @NotNull Vec3 normal) {
-        return new HitResult(t, position, normal, material, u, v, isFrontFace);
+        return new HitResult(t, position, normal, target, material, u, v, isFrontFace);
     }
 
     @Override

src/main/java/eu/jonahbauer/raytracing/scene/Hittable.java

@@ -11,11 +11,25 @@ import org.jetbrains.annotations.NotNull;
 import java.util.Optional;
 
 public interface Hittable {
+    @NotNull Range FORWARD = new Range(0.001, Double.POSITIVE_INFINITY);
 
     /**
-     * {@return the value <code>t</code> such that <code>ray.at(t)</code> is the intersection of this shaped closest to
-     * the ray origin, or <code>Double.NaN</code> if the ray does not intersect this shape}
+     * @see #hit(Ray, Range)
+     */
+    default @NotNull Optional<HitResult> hit(@NotNull Ray ray) {
+        return hit(ray, FORWARD);
+    }
+
+    /**
+     * Tests whether the {@code ray} intersects {@code this} hittable.
+     * <p>
+     * The second parameter {@code range} allows the implementation to skip unnecessary calculations if it can
+     * determine that a hit (if any) will fall outside the valid range of {@code t}s. The returned hit may still be
+     * outside the valid range and has to be checked by the caller.
      * @param ray a ray
+     * @param range the range of valid {@code t}s
+     * @return the result of the hit test, containing (among others) the value {@code t} such that {@code ray.at(t)} is
+     *         a point on {@code this} hittable
      */
     @NotNull Optional<HitResult> hit(@NotNull Ray ray, @NotNull Range range);
 
@@ -25,10 +39,10 @@ public interface Hittable {
     @NotNull AABB getBoundingBox();
 
     default @NotNull Hittable translate(@NotNull Vec3 offset) {
-        return new Translate(this, offset);
+        return Translate.create(this, offset);
     }
 
     default @NotNull Hittable rotateY(double angle) {
-        return new RotateY(this, angle);
+        return RotateY.create(this, angle);
     }
 }

src/main/java/eu/jonahbauer/raytracing/scene/Scene.java

@@ -8,6 +8,7 @@ import eu.jonahbauer.raytracing.scene.util.HittableCollection;
 import org.jetbrains.annotations.NotNull;
 import org.jetbrains.annotations.Nullable;
 
+import java.util.ArrayList;
 import java.util.List;
 import java.util.Objects;
 
@@ -15,33 +16,28 @@ public final class Scene extends HittableCollection {
     private final @NotNull HittableCollection objects;
     private final @NotNull SkyBox background;
 
-    private final @Nullable Target light;
+    private final @Nullable List<@NotNull Target> targets;
 
     public Scene(@NotNull List<? extends @NotNull Hittable> objects) {
-        this(Color.BLACK, objects);
+        this(objects, null);
     }
 
-    public Scene(@NotNull Color background, @NotNull List<? extends @NotNull Hittable> objects) {
-        this(SkyBox.solid(background), objects);
+    public Scene(@NotNull List<? extends @NotNull Hittable> objects, @Nullable List<? extends @NotNull Target> targets) {
+        this(Color.BLACK, objects, targets);
     }
 
     public Scene(@NotNull SkyBox background, @NotNull List<? extends @NotNull Hittable> objects) {
-        this.objects = new HittableBinaryTree(objects);
+        this(background, objects, null);
+    }
+
+    public Scene(@NotNull SkyBox background, @NotNull List<? extends @NotNull Hittable> objects, @Nullable List<? extends @NotNull Target> targets) {
+        var list = new ArrayList<Hittable>(objects.size() + (targets != null ? targets.size() : 0));
+        list.addAll(objects);
+        if (targets != null) list.addAll(targets);
+        this.objects = new HittableBinaryTree(list);
         this.background = Objects.requireNonNull(background);
 
-        this.light = (Target) objects.get(1);
-    }
-
-    public Scene(@NotNull Hittable @NotNull... objects) {
-        this(List.of(objects));
-    }
-
-    public Scene(@NotNull Color background, @NotNull Hittable @NotNull... objects) {
-        this(background, List.of(objects));
-    }
-
-    public Scene(@NotNull SkyBox background, @NotNull Hittable @NotNull... objects) {
-        this(background, List.of(objects));
+        this.targets = targets != null ? List.copyOf(targets) : null;
     }
 
     @Override
@@ -54,8 +50,8 @@ public final class Scene extends HittableCollection {
         return objects.getBoundingBox();
     }
 
-    public @Nullable Target getLights() {
-        return light;
+    public @Nullable List<@NotNull Target> getTargets() {
+        return targets;
     }
 
     public @NotNull Color getBackgroundColor(@NotNull Ray ray) {

src/main/java/eu/jonahbauer/raytracing/scene/SkyBox.java

@@ -18,8 +18,4 @@ public interface SkyBox {
             return Color.lerp(bottom, top, alt / Math.PI + 0.5);
         };
     }
-
-    static @NotNull SkyBox solid(@NotNull Color color) {
-        return _ -> color;
-    }
 }

src/main/java/eu/jonahbauer/raytracing/scene/hittable2d/Hittable2D.java

@@ -52,10 +52,40 @@ public abstract class Hittable2D implements Hittable {
 
         var frontFace = denominator < 0;
         return Optional.of(new HitResult(
-                t, position, frontFace ? normal : normal.neg(),
+                t, position, frontFace ? normal : normal.neg(), this,
                 material, alpha, beta, frontFace
         ));
     }
 
+    protected double hit0(@NotNull Ray ray, @NotNull Range range) {
+        var denominator = ray.direction().times(normal);
+        if (Math.abs(denominator) < 1e-8) return Double.NaN; // parallel
+
+        var t = (d - ray.origin().times(normal)) / denominator;
+        if (!range.surrounds(t)) return Double.NaN;
+
+        var position = ray.at(t);
+        var p = position.minus(origin);
+
+        var alpha = Vec3.tripleProduct(w, p, v);
+        var beta = Vec3.tripleProduct(w, u, p);
+        if (!isInterior(alpha, beta)) return Double.NaN;
+
+        return t;
+    }
+
+    protected @NotNull Vec3 get(double alpha, double beta) {
+        return new Vec3(
+                Math.fma(beta, v.x(), Math.fma(alpha, u.x(), origin.x())),
+                Math.fma(beta, v.y(), Math.fma(alpha, u.y(), origin.y())),
+                Math.fma(beta, v.z(), Math.fma(alpha, u.z(), origin.z()))
+        );
+    }
+
     protected abstract boolean isInterior(double alpha, double beta);
+
+    @Override
+    public @NotNull String toString() {
+        return this.getClass().getSimpleName() + "(origin=" + origin + ", u=" + u + ", v=" + v + ")";
+    }
 }

src/main/java/eu/jonahbauer/raytracing/scene/hittable2d/Parallelogram.java

@@ -1,7 +1,6 @@
 package eu.jonahbauer.raytracing.scene.hittable2d;
 
 import eu.jonahbauer.raytracing.math.AABB;
-import eu.jonahbauer.raytracing.math.Range;
 import eu.jonahbauer.raytracing.math.Ray;
 import eu.jonahbauer.raytracing.math.Vec3;
 import eu.jonahbauer.raytracing.render.material.Material;
@@ -31,14 +30,14 @@ public final class Parallelogram extends Hittable2D implements Target {
 
     @Override
     public double getProbabilityDensity(@NotNull Vec3 origin, @NotNull Vec3 direction) {
-        var result = hit(new Ray(origin, direction), new Range(0.001, Double.POSITIVE_INFINITY));
-        if (result.isEmpty()) return 0;
+        if (Double.isNaN(hit0(new Ray(origin, direction), FORWARD))) return 0;
 
-        var a = this.origin;
-        var b = this.origin.plus(u);
-        var c = this.origin.plus(v);
-        var d = b.plus(v);
-        var angle = PdfUtil.getSolidAngle(origin, a, b, d) + PdfUtil.getSolidAngle(origin, c, b, d);
+        var o = this.origin.minus(origin);
+        var a = o.unit();
+        var b = o.plus(u).unit();
+        var c = o.plus(v).unit();
+        var d = o.plus(u).plus(v).unit();
+        var angle = PdfUtil.getSolidAngle(a, b, d) + PdfUtil.getSolidAngle(c, b, d);
         return 1 / angle;
     }
 
@@ -46,6 +45,6 @@ public final class Parallelogram extends Hittable2D implements Target {
     public @NotNull Vec3 getTargetingDirection(@NotNull Vec3 origin, @NotNull RandomGenerator random) {
         var alpha = random.nextDouble();
         var beta = random.nextDouble();
-        return this.origin.plus(u.times(alpha)).plus(v.times(beta)).minus(origin);
+        return get(alpha, beta).minus(origin);
     }
 }

src/main/java/eu/jonahbauer/raytracing/scene/hittable3d/Box.java

@@ -28,6 +28,15 @@ public final class Box implements Hittable, Target {
         this(box, Objects.requireNonNull(material, "material"), material, material, material, material, material);
     }
 
+    public Box(
+            @NotNull Vec3 a, @NotNull Vec3 b,
+            @Nullable Material top, @Nullable Material bottom,
+            @Nullable Material left, @Nullable Material right,
+            @Nullable Material front, @Nullable Material back
+    ) {
+        this(new AABB(a, b), top, bottom, left, right, front, back);
+    }
+
     public Box(
             @NotNull AABB box,
             @Nullable Material top, @Nullable Material bottom,
@@ -110,7 +119,7 @@ public final class Box implements Hittable, Target {
         var uv = material.texture().isUVRequired();
         var u = uv ? side.getTextureU(box, position) : Double.NaN;
         var v = uv ? side.getTextureV(box, position) : Double.NaN;
-        return Optional.of(new HitResult(t, position, normal, material, u, v, frontFace));
+        return Optional.of(new HitResult(t, position, normal, this, material, u, v, frontFace));
     }
 
     @Override
@@ -121,7 +130,7 @@ public final class Box implements Hittable, Target {
     @Override
     public double getProbabilityDensity(@NotNull Vec3 origin, @NotNull Vec3 direction) {
         if (contains(origin)) return 1 / (4 * Math.PI);
-        if (hit(new Ray(origin, direction), new Range(0.001, Double.POSITIVE_INFINITY)).isEmpty()) return 0;
+        if (hit(new Ray(origin, direction)).isEmpty()) return 0;
 
         var solidAngle = 0d;
         for (var s : Side.values()) {
@@ -134,7 +143,7 @@ public final class Box implements Hittable, Target {
 
     @Override
     public @NotNull Vec3 getTargetingDirection(@NotNull Vec3 origin, @NotNull RandomGenerator random) {
-        if (contains(origin)) return Vec3.random(random, true);
+        if (contains(origin)) return Vec3.random(random);
 
         // determine sides facing the origin and their solid angles
         int visible = 0;
@@ -173,6 +182,11 @@ public final class Box implements Hittable, Target {
                && box.min().z() < point.z() && point.z() < box.max().z();
     }
 
+    @Override
+    public @NotNull String toString() {
+        return "Box(min=" + box.min() + ", max=" + box.max() + ")";
+    }
+
     private enum Side {
         NEG_X(Vec3.UNIT_X.neg()),
         NEG_Y(Vec3.UNIT_Y.neg()),
@@ -279,12 +293,12 @@ public final class Box implements Hittable, Target {
          * {@return the solid angle covered by <code>this</code> side of the <code>box</code> when viewed from <code>pos</code>}
          */
         public double getSolidAngle(@NotNull AABB box, @NotNull Vec3 pos) {
-            var a = get(box, 0, 0);
-            var b = get(box, 0, 1);
-            var c = get(box, 1, 1);
-            var d = get(box, 1, 0);
+            var a = get(box, 0, 0).minus(pos).unit();
+            var b = get(box, 0, 1).minus(pos).unit();
+            var c = get(box, 1, 1).minus(pos).unit();
+            var d = get(box, 1, 0).minus(pos).unit();
 
-            return PdfUtil.getSolidAngle(pos, a, b, d) + PdfUtil.getSolidAngle(pos, c, b, d);
+            return PdfUtil.getSolidAngle(a, b, d) + PdfUtil.getSolidAngle(c, b, d);
         }
     }
 }

src/main/java/eu/jonahbauer/raytracing/scene/hittable3d/ConstantMedium.java

@@ -32,7 +32,7 @@ public record ConstantMedium(@NotNull Hittable boundary, double density, @NotNul
         if (hitDistance > distance) return Optional.empty();
 
         var t = tmin + hitDistance / length;
-        return Optional.of(new HitResult(t, ray.at(t), Vec3.UNIT_X, material, 0, 0, true)); // arbitrary normal, u, v and isFrontFace
+        return Optional.of(new HitResult(t, ray.at(t), Vec3.UNIT_X, this, material, 0, 0, true)); // arbitrary normal, u, v and isFrontFace
     }
 
     @Override

src/main/java/eu/jonahbauer/raytracing/scene/hittable3d/Sphere.java

@@ -21,12 +21,18 @@ public final class Sphere implements Hittable, Target {
 
     private final @NotNull AABB bbox;
 
+    private final @NotNull Vec3 normalizedCenter;
+    private final double invRadius;
+
     public Sphere(@NotNull Vec3 center, double radius, @NotNull Material material) {
         this.center = Objects.requireNonNull(center, "center");
         this.material = Objects.requireNonNull(material, "material");
         if (radius <= 0 || !Double.isFinite(radius)) throw new IllegalArgumentException("radius must be positive");
         this.radius = radius;
 
+        this.invRadius = 1 / radius;
+        this.normalizedCenter = this.center.times(-this.invRadius);
+
         this.bbox = new AABB(
                 center.minus(radius, radius, radius),
                 center.plus(radius, radius, radius)
@@ -35,23 +41,11 @@ public final class Sphere implements Hittable, Target {
 
     @Override
     public @NotNull Optional<HitResult> hit(@NotNull Ray ray, @NotNull Range range) {
-        var oc = ray.origin().minus(center);
-
-        var a = ray.direction().squared();
-        var h = ray.direction().times(oc);
-        var c = oc.squared() - radius * radius;
-
-        var discriminant = h * h - a * c;
-        if (discriminant < 0) return Optional.empty();
-
-        var sd = Math.sqrt(discriminant);
-
-        double t = (- h - sd) / a;
-        if (!range.surrounds(t)) t = (- h + sd) / a;
-        if (!range.surrounds(t)) return Optional.empty();
+        var t = hit0(ray, range);
+        if (Double.isNaN(t)) return Optional.empty();
 
         var position = ray.at(t);
-        var normal = position.minus(center).div(radius);
+        var normal = Vec3.fma(invRadius, position, normalizedCenter);
         var frontFace = normal.times(ray.direction()) < 0;
 
         double u;
@@ -67,11 +61,30 @@ public final class Sphere implements Hittable, Target {
         }
 
         return Optional.of(new HitResult(
-                t, position, frontFace ? normal : normal.neg(),
+                t, position, frontFace ? normal : normal.neg(), this,
                 material, u, v, frontFace
         ));
     }
 
+    private double hit0(@NotNull Ray ray, @NotNull Range range) {
+        var oc = ray.origin().minus(center);
+
+        var a = ray.direction().squared();
+        var h = ray.direction().times(oc);
+        var c = oc.squared() - radius * radius;
+
+        var discriminant = h * h - a * c;
+        if (discriminant < 0) return Double.NaN;
+
+        var sd = Math.sqrt(discriminant);
+
+        double t = (- h - sd) / a;
+        if (!range.surrounds(t)) t = (- h + sd) / a;
+        if (!range.surrounds(t)) return Double.NaN;
+
+        return t;
+    }
+
     @Override
     public @NotNull AABB getBoundingBox() {
         return bbox;
@@ -79,7 +92,7 @@ public final class Sphere implements Hittable, Target {
 
     @Override
     public double getProbabilityDensity(@NotNull Vec3 origin, @NotNull Vec3 direction) {
-        if (hit(new Ray(origin, direction), new Range(0.001, Double.POSITIVE_INFINITY)).isEmpty()) return 0;
+        if (Double.isNaN(hit0(new Ray(origin, direction), FORWARD))) return 0;
 
         var cos = Math.sqrt(1 - radius * radius / (center.minus(origin).squared()));
         var solidAngle = 2 * Math.PI * (1 - cos);
@@ -89,12 +102,16 @@ public final class Sphere implements Hittable, Target {
     @Override
     public @NotNull Vec3 getTargetingDirection(@NotNull Vec3 origin, @NotNull RandomGenerator random) {
         var direction = center.minus(origin);
+        var out = Vec3.randomOppositeHemisphere(random, direction);
+        return new Vec3(
+                Math.fma(radius, out.x(), center.x() - origin.x()),
+                Math.fma(radius, out.y(), center.y() - origin.y()),
+                Math.fma(radius, out.z(), center.z() - origin.z())
+        );
+    }
 
-        Vec3 target;
-        do {
-            target = Vec3.random(random, true);
-        } while (target.times(direction) >= 0);
-
-        return target.times(radius).plus(center).minus(origin);
+    @Override
+    public @NotNull String toString() {
+        return "Sphere(center=" + center + ", radius=" + radius + ")";
     }
 }

src/main/java/eu/jonahbauer/raytracing/scene/transform/RotateY.java

@@ -10,13 +10,21 @@ import org.jetbrains.annotations.NotNull;
 
 import java.util.random.RandomGenerator;
 
-public final class RotateY extends Transform {
+public sealed class RotateY extends Transform {
     private final double cos;
     private final double sin;
 
     private final @NotNull AABB bbox;
 
-    public RotateY(@NotNull Hittable object, double angle) {
+    public static @NotNull RotateY create(@NotNull Hittable object, double angle) {
+        if (object instanceof Target) {
+            return new RotateYTarget(object, angle);
+        } else {
+            return new RotateY(object, angle);
+        }
+    }
+
+    private RotateY(@NotNull Hittable object, double angle) {
         super(object);
         this.cos = Math.cos(angle);
         this.sin = Math.sin(angle);
@@ -48,7 +56,7 @@ public final class RotateY extends Transform {
     }
 
     @Override
-    protected @NotNull Ray transform(@NotNull Ray ray) {
+    protected final @NotNull Ray transform(@NotNull Ray ray) {
         var origin = ray.origin();
         var direction = ray.direction();
 
@@ -58,7 +66,7 @@ public final class RotateY extends Transform {
     }
 
     @Override
-    protected @NotNull HitResult transform(@NotNull HitResult result) {
+    protected final @NotNull HitResult transform(@NotNull HitResult result) {
         var position = result.position();
         var newPosition = untransform(position);
 
@@ -68,11 +76,11 @@ public final class RotateY extends Transform {
         return result.withPositionAndNormal(newPosition, newNormal);
     }
 
-    private @NotNull Vec3 transform(@NotNull Vec3 vec) {
+    protected final @NotNull Vec3 transform(@NotNull Vec3 vec) {
         return new Vec3(cos * vec.x() - sin * vec.z(), vec.y(), sin * vec.x() + cos * vec.z());
     }
 
-    private @NotNull Vec3 untransform(@NotNull Vec3 vec) {
+    protected final @NotNull Vec3 untransform(@NotNull Vec3 vec) {
         return new Vec3(cos * vec.x() + sin * vec.z(), vec.y(), - sin * vec.x() + cos * vec.z());
     }
 
@@ -82,14 +90,24 @@ public final class RotateY extends Transform {
     }
 
     @Override
-    public double getProbabilityDensity(@NotNull Vec3 origin, @NotNull Vec3 direction) {
-        if (!(object instanceof Target target)) throw new UnsupportedOperationException();
-        return target.getProbabilityDensity(transform(origin), transform(direction));
+    public @NotNull String toString() {
+        return object + " rotated by " + Math.toDegrees(Math.atan2(sin, cos)) + "° around the y axis";
     }
 
-    @Override
-    public @NotNull Vec3 getTargetingDirection(@NotNull Vec3 origin, @NotNull RandomGenerator random) {
-        if (!(object instanceof Target target)) throw new UnsupportedOperationException();
-        return untransform(target.getTargetingDirection(transform(origin), random));
+    private static final class RotateYTarget extends RotateY implements Target {
+
+        private RotateYTarget(@NotNull Hittable object, double angle) {
+            super(object, angle);
+        }
+
+        @Override
+        public double getProbabilityDensity(@NotNull Vec3 origin, @NotNull Vec3 direction) {
+            return ((Target) object).getProbabilityDensity(transform(origin), transform(direction));
+        }
+
+        @Override
+        public @NotNull Vec3 getTargetingDirection(@NotNull Vec3 origin, @NotNull RandomGenerator random) {
+            return untransform(((Target) object).getTargetingDirection(transform(origin), random));
+        }
     }
 }

src/main/java/eu/jonahbauer/raytracing/scene/transform/Transform.java

@@ -12,7 +12,7 @@ import java.util.Objects;
 import java.util.Optional;
 import java.util.random.RandomGenerator;
 
-public abstract class Transform implements Hittable, Target {
+public abstract class Transform implements Hittable {
     protected final @NotNull Hittable object;
 
     protected Transform(@NotNull Hittable object) {

src/main/java/eu/jonahbauer/raytracing/scene/transform/Translate.java

@@ -10,11 +10,19 @@ import org.jetbrains.annotations.NotNull;
 
 import java.util.random.RandomGenerator;
 
-public final class Translate extends Transform {
-    private final @NotNull Vec3 offset;
+public sealed class Translate extends Transform {
+    protected final @NotNull Vec3 offset;
     private final @NotNull AABB bbox;
 
-    public Translate(@NotNull Hittable object, @NotNull Vec3 offset) {
+    public static @NotNull Translate create(@NotNull Hittable object, @NotNull Vec3 offset) {
+        if (object instanceof Target) {
+            return new TranslateTarget(object, offset);
+        } else {
+            return new Translate(object, offset);
+        }
+    }
+
+    private Translate(@NotNull Hittable object, @NotNull Vec3 offset) {
         super(object);
         this.offset = offset;
 
@@ -26,29 +34,40 @@ public final class Translate extends Transform {
     }
 
     @Override
-    protected @NotNull Ray transform(@NotNull Ray ray) {
+    protected final @NotNull Ray transform(@NotNull Ray ray) {
         return new Ray(ray.origin().minus(offset), ray.direction());
     }
 
     @Override
-    protected @NotNull HitResult transform(@NotNull HitResult result) {
+    protected final @NotNull HitResult transform(@NotNull HitResult result) {
         return result.withPositionAndNormal(result.position().plus(offset), result.normal());
     }
 
     @Override
-    public @NotNull AABB getBoundingBox() {
+    public final @NotNull AABB getBoundingBox() {
         return bbox;
     }
 
     @Override
-    public double getProbabilityDensity(@NotNull Vec3 origin, @NotNull Vec3 direction) {
-        if (!(object instanceof Target target)) throw new UnsupportedOperationException();
-        return target.getProbabilityDensity(origin.minus(offset), direction);
+    public @NotNull String toString() {
+        return object + " translated by " + offset;
     }
 
-    @Override
-    public @NotNull Vec3 getTargetingDirection(@NotNull Vec3 origin, @NotNull RandomGenerator random) {
-        if (!(object instanceof Target target)) throw new UnsupportedOperationException();
-        return target.getTargetingDirection(origin.minus(offset), random);
+    private static final class TranslateTarget extends Translate implements Target {
+        private TranslateTarget(@NotNull Hittable object, @NotNull Vec3 offset) {
+            super(object, offset);
+        }
+
+        @Override
+        public double getProbabilityDensity(@NotNull Vec3 origin, @NotNull Vec3 direction) {
+            if (!(object instanceof Target target)) throw new UnsupportedOperationException();
+            return target.getProbabilityDensity(origin.minus(offset), direction);
+        }
+
+        @Override
+        public @NotNull Vec3 getTargetingDirection(@NotNull Vec3 origin, @NotNull RandomGenerator random) {
+            if (!(object instanceof Target target)) throw new UnsupportedOperationException();
+            return target.getTargetingDirection(origin.minus(offset), random);
+        }
     }
 }

src/main/java/eu/jonahbauer/raytracing/scene/util/PdfUtil.java

@@ -8,11 +8,12 @@ public final class PdfUtil {
         throw new UnsupportedOperationException();
     }
 
-    public static double getSolidAngle(@NotNull Vec3 o, @NotNull Vec3 a, @NotNull Vec3 b, @NotNull Vec3 c) {
-        var i = a.minus(o).unit();
-        var j = b.minus(o).unit();
-        var k = c.minus(o).unit();
-
-        return 2 * Math.atan(Math.abs(i.times(j.cross(k))) / (1 + i.times(j) + j.times(k) + k.times(i)));
+    /**
+     * {@return the solid angle of the triangle abc as seen from the origin} The vectors {@code a}, {@code b} and {@code c}
+     * must be unit vectors.
+     */
+    public static double getSolidAngle(@NotNull Vec3 a, @NotNull Vec3 b, @NotNull Vec3 c) {
+        var angle = 2 * Math.atan(Math.abs(Vec3.tripleProduct(a, b, c)) / (1 + a.times(b) + b.times(c) + c.times(a)));
+        return angle < 0 ? 2 * Math.PI + angle : angle;
     }
 }