patrick/ray-tracing-fsharp
1
0
Fork 0
mirror of https://github.com/Smaug123/ray-tracing-fsharp synced 2025-10-13 07:48:40 +00:00
Code Issues Packages Projects Releases Wiki Activity

Perf (#2)

Browse Source
This commit is contained in:
Patrick Stevens
2021-04-07 18:53:51 +01:00
committed by GitHub
parent eefaa92d0e
commit a001bab92b
24 changed files with 458 additions and 753 deletions
Download Patch File Download Diff File Expand all files Collapse all files

1
RayTracing.Test/RayTracing.Test.fsproj
View File

@@ -8,7 +8,6 @@
<Compile Include="TestUtils.fs" />
<Compile Include="TestPpmOutput.fs" />
<Compile Include="TestRay.fs" />
<Compile Include="TestRational.fs" />
<Compile Include="TestSphereIntersection.fs" />
<EmbeddedResource Include="PpmOutputExample.txt" />
<Compile Include="TestPixel.fs" />

4
RayTracing.Test/TestPixel.fs
View File

@@ -13,7 +13,7 @@ module TestPixel =
let ``Average of one pixel`` () =
let property (p1 : byte) (p2 : byte) (p3 : byte) : bool =
let pixel = { Red = p1 ; Green = p2 ; Blue = p3 }
Pixel.average [ pixel ]
Pixel.average [| pixel |]
|> (=) pixel
Check.QuickThrowOnFailure property
@@ -78,7 +78,7 @@ module TestPixel =
values
|> List.map (fun (a, b, c) -> { Pixel.Red = a ; Green = b ; Blue = c })
let avg = Pixel.average pixels
let avg = Pixel.average (Array.ofList pixels)
avg.Green = (pixels |> List.map (fun i -> float i.Green) |> List.average |> Math.Round |> byte)
&& avg.Red = (pixels |> List.map (fun i -> float i.Red) |> List.average |> Math.Round |> byte)

10
RayTracing.Test/TestPlane.fs
View File

@@ -9,12 +9,12 @@ module TestPlane =
[<Test>]
let ``Orthogonalise does make orthogonal vectors`` () =
let property (p : Plane<float>) : bool =
let orth = Plane.orthonormalise Num.float p |> Option.get
let property (p : Plane) : bool =
let orth = Plane.orthonormalise p |> Option.get
let v1, v2 = Plane.basis orth
Num.float.Equal (Vector.dot Num.float v1.Vector v2.Vector) Num.float.Zero
&& Num.float.Equal (Vector.dot Num.float v1.Vector v1.Vector) Num.float.One
&& Num.float.Equal (Vector.dot Num.float v2.Vector v2.Vector) Num.float.One
Float.equal (UnitVector.dot (Ray.vector v1) (Ray.vector v2)) 0.0
&& Float.equal (UnitVector.dot (Ray.vector v1) (Ray.vector v1)) 1.0
&& Float.equal (UnitVector.dot (Ray.vector v2) (Ray.vector v2)) 1.0
property
|> Prop.forAll (Arb.fromGen TestUtils.planeGen)

34
RayTracing.Test/TestRational.fs
View File

@@ -1,34 +0,0 @@
namespace RayTracing.Test
open NUnit.Framework
open FsCheck
open RayTracing
[<TestFixture>]
module TestRational =
[<Test>]
let ``ofInt compares correctly`` () =
let property (i : int) (j : int) : bool =
let i1 = Rational.ofInt i
let j1 = Rational.ofInt j
if i1 < j1 then i < j
elif i1 > j1 then i > j
else i = j
Check.QuickThrowOnFailure property
[<Test>]
let ``Addition preserves comparison`` () =
let property (i : Rational, j : Rational, k : Rational) : bool =
if i < j then
Rational.add i k < Rational.add j k
elif i = j then
Rational.add i k = Rational.add j k
else
Rational.add i k > Rational.add j k
property
|> Prop.forAll (Arb.fromGen (Gen.three TestUtils.rationalGen))
|> Check.QuickThrowOnFailure

43
RayTracing.Test/TestRay.fs
View File

@@ -10,12 +10,11 @@ module TestRay =
[<Test>]
let ``Walking along two parallel rays maintains the same vector difference`` () =
let property
(num : Num<'a>)
(
((originX : 'a, originY : 'a, originZ : 'a),
(origin2X : 'a, origin2Y : 'a, origin2Z : 'a),
(rayX : 'a, rayY : 'a, rayZ : 'a)),
magnitude : 'a
(((originX : float, originY : float, originZ : float),
(origin2X : float, origin2Y : float, origin2Z : float)),
vector : UnitVector),
magnitude : float
)
: bool
=
@@ -24,40 +23,36 @@ module TestRay =
let origin2 =
[| origin2X; origin2Y; origin2Z |] |> Point
let vector = Vector [| rayX; rayY; rayZ |]
let ray = { Origin = origin1; Vector = vector }
let ray2 = { Origin = origin2; Vector = vector }
let output = Ray.walkAlong num ray magnitude
let output = Ray.walkAlong (Ray.make origin1 vector) magnitude
let output2 =
Ray.walkAlong num ray2 magnitude
let output2 = Ray.walkAlong (Ray.make origin2 vector) magnitude
let actual =
Point.difference num output output2
Point.difference output output2
let expected =
Point.difference num origin1 origin2
Point.difference origin1 origin2
Vector.equal num actual expected
Vector.equal actual expected
let gen : Gen<float> =
Arb.generate<NormalFloat>
|> Gen.map NormalFloat.op_Explicit
let gen =
Gen.zip (Gen.three (Gen.three gen)) gen
Gen.zip (Gen.zip (Gen.two (Gen.three gen)) TestUtils.unitVectorGen) gen
property Num.float
property
|> Prop.forAll (Arb.fromGen gen)
|> Check.QuickThrowOnFailure
[<Test>]
let ``walkAlong walks the right distance`` () =
let property (ray : Ray<float>, distance : float) =
let walked = Ray.walkAlong Num.float ray distance
Point.difference Num.float walked ray.Origin
|> Vector.normSquared Num.float
|> Num.float.Equal (distance * distance)
let property (ray : Ray, distance : float) =
let walked = Ray.walkAlong ray distance
Point.difference walked (Ray.origin ray)
|> Vector.normSquared
|> Float.equal (distance * distance)
property
|> Prop.forAll (Arb.fromGen (Gen.zip TestUtils.rayGen (Arb.generate<NormalFloat> |> Gen.map NormalFloat.op_Explicit)))
@@ -65,9 +60,9 @@ module TestRay =
[<Test>]
let ``walkAlong stays on the ray`` () =
let property (ray : Ray<float>, distance : float) =
let walked = Ray.walkAlong Num.float ray distance
Ray.liesOn Num.float walked ray
let property (ray : Ray, distance : float) =
let walked = Ray.walkAlong ray distance
Ray.liesOn walked ray
property
|> Prop.forAll (Arb.fromGen (Gen.zip TestUtils.rayGen (Arb.generate<NormalFloat> |> Gen.map NormalFloat.op_Explicit)))

18
RayTracing.Test/TestSphere.fs
View File

@@ -9,33 +9,31 @@ module TestSphere =
[<Test>]
let ``Point at distance r from centre lies on sphere`` () =
let property (centre : Point<float>, radius : float, point : Point<float>) : bool =
let property (centre : Point, radius : float, point : Point) : bool =
let radius = abs radius
let sphere = Sphere.make Num.float (SphereStyle.PureReflection (1.0, Colour.White)) centre radius
Sphere.liesOn Num.float point sphere
let sphere = Sphere.make (SphereStyle.PureReflection (1.0, Colour.White)) centre radius
Sphere.liesOn point sphere
let gen : Gen<Point<float> * float * Point<float>> =
let gen : Gen<Point * float * Point> =
gen {
let! centre = TestUtils.pointGen
let! radius = Arb.generate<NormalFloat> |> Gen.map NormalFloat.op_Explicit
let! theta =
Arb.generate<NormalFloat>
|> Gen.map NormalFloat.op_Explicit
|> Gen.map Radian
let! phi =
Arb.generate<NormalFloat>
|> Gen.map NormalFloat.op_Explicit
|> Gen.map Radian
let surfacePoint =
[|
radius * Num.float.Cos phi * Num.float.Sin theta
radius * Num.float.Sin phi * Num.float.Sin theta
radius * Num.float.Cos theta
radius * cos phi * sin theta
radius * sin phi * sin theta
radius * cos theta
|]
|> Point
|> Point.difference Num.float centre
|> Point.difference centre
|> fun (Vector v) -> Point v
return centre, radius, surfacePoint
}

34
RayTracing.Test/TestSphereIntersection.fs
View File

@@ -9,28 +9,28 @@ open RayTracing
module TestSphereIntersection =
let sphere : Gen<Sphere<float>> =
let sphere : Gen<Sphere> =
gen {
let! origin = TestUtils.pointGen
let! radius = Arb.generate<NormalFloat>
return Sphere.make Num.float (SphereStyle.LightSource Colour.White) origin radius.Get
return Sphere.make (SphereStyle.LightSource Colour.White) origin radius.Get
}
[<Test>]
let ``Intersection of sphere and ray does lie on both`` () =
let property (ray : Ray<float>, sphere : Sphere<float>) : bool =
let intersections = Sphere.intersections Num.float sphere ray Colour.White
let property (ray : Ray, sphere : Sphere) : bool =
let intersections = Sphere.intersections sphere ray Colour.White
intersections
|> Seq.forall (fun (p, _, _) ->
let rayOk = Ray.liesOn Num.float p ray
let sphereOk = Sphere.liesOn Num.float p sphere
let rayOk = Ray.liesOn p ray
let sphereOk = Sphere.liesOn p sphere
rayOk && sphereOk
)
&&
intersections
|> Array.map (fun (intersection, _, _) -> Vector.normSquared Num.float (Point.difference Num.float ray.Origin intersection))
|> Array.map (fun (intersection, _, _) -> Vector.normSquared (Point.difference (Ray.origin ray) intersection))
|> Seq.pairwise
|> Seq.forall (fun (i, j) -> Num.float.Compare i j = Less)
|> Seq.forall (fun (i, j) -> Float.compare i j = Less)
property
|> Prop.forAll (Arb.fromGen (Gen.zip TestUtils.rayGen sphere))
@@ -39,24 +39,22 @@ module TestSphereIntersection =
[<Test>]
let ``Intersection of sphere and ray does lie on both, case 1`` () =
let ray =
{
Origin = Point [|1.462205539; -4.888279676; 7.123293244|]
Vector = Vector [|-9.549697616; 4.400018428; 10.41024923|]
}
let sphere = Sphere.make Num.float (SphereStyle.PureReflection (1.0, Colour.White)) (Point [|-5.688391601; -5.360125644; 9.074300761|]) 8.199747973
Ray.make' (Point [|1.462205539; -4.888279676; 7.123293244|]) (Vector [|-9.549697616; 4.400018428; 10.41024923|])
|> Option.get
let sphere = Sphere.make (SphereStyle.PureReflection (1.0, Colour.White)) (Point [|-5.688391601; -5.360125644; 9.074300761|]) 8.199747973
let intersections = Sphere.intersections Num.float sphere ray Colour.White
let intersections = Sphere.intersections sphere ray Colour.White
intersections
|> Array.map (fun (intersection, _, _) -> Vector.normSquared Num.float (Point.difference Num.float ray.Origin intersection))
|> Array.map (fun (intersection, _, _) -> Vector.normSquared (Point.difference (Ray.origin ray) intersection))
|> Seq.pairwise
|> Seq.forall (fun (i, j) -> Num.float.Compare i j = Less)
|> Seq.forall (fun (i, j) -> Float.compare i j = Less)
|> shouldEqual true
intersections
|> Seq.forall (fun (p, _, _) -> Ray.liesOn Num.float p ray)
|> Seq.forall (fun (p, _, _) -> Ray.liesOn p ray)
|> shouldEqual true
intersections
|> Seq.forall (fun (p, _, _) -> Sphere.liesOn Num.float p sphere)
|> Seq.forall (fun (p, _, _) -> Sphere.liesOn p sphere)
|> shouldEqual true

40
RayTracing.Test/TestUtils.fs
View File

@@ -1,6 +1,5 @@
namespace RayTracing.Test
open System.Numerics
open RayTracing
open System.IO
open System.Reflection
@@ -29,27 +28,6 @@ module TestUtils =
use reader = new StreamReader (stream)
reader.ReadToEnd().Replace ("\r\n", "\n")
let rationalGen : Gen<Rational> =
gen {
let! i = Gen.choose (-100, 100)
let! sign = Gen.choose (0, 1)
let! j = Gen.choose (1, 100)
return Rational.Make (BigInteger i) (if sign = 0 then BigInteger j else BigInteger(-j))
}
let rec algebraicGen () : Gen<Algebraic> =
[
rationalGen |> Gen.map Algebraic.ofRational
Gen.two (algebraicGen ()) |> Gen.map Algebraic.Sum
Gen.two (algebraicGen ()) |> Gen.map Algebraic.Times
// TODO make this nonnegative
algebraicGen () |> Gen.map Algebraic.Sqrt
// TODO make this nonzero
algebraicGen () |> Gen.map Algebraic.Reciprocal
// TODO more of these
]
|> Gen.oneof
let floatGen = Arb.generate<NormalFloat> |> Gen.map NormalFloat.op_Explicit
let pointGen =
@@ -60,17 +38,23 @@ module TestUtils =
Gen.three Arb.generate<NormalFloat>
|> Gen.map (fun (i, j, k) -> Vector [| i.Get ; j.Get ; k.Get |])
let rayGen : Gen<Ray<float>> =
let unitVectorGen =
vectorGen
|> Gen.filter (fun i -> Vector.normSquared i > 0.0)
|> Gen.map Vector.unitise
|> Gen.map Option.get
let rayGen : Gen<Ray> =
gen {
let! origin = pointGen
let! direction = vectorGen
return { Origin = origin ; Vector = direction }
let! direction = unitVectorGen
return Ray.make origin direction
}
let planeGen =
gen {
let! origin = pointGen
let! v1 = vectorGen
let! v2 = vectorGen
return Plane.makeSpannedBy { Origin = origin ; Vector = v1 } { Origin = origin ; Vector = v2 }
let! v1 = unitVectorGen
let! v2 = unitVectorGen
return Plane.makeSpannedBy (Ray.make origin v1) (Ray.make origin v2)
}

32
RayTracing/Algebraic.fs
View File

@@ -1,32 +0,0 @@
namespace RayTracing
type Algebraic =
| Sqrt of Algebraic
| Rational of Rational
| Sum of Algebraic * Algebraic
| Times of Algebraic * Algebraic
| Subtract of Algebraic * Algebraic
| Negate of Algebraic
| Reciprocal of Algebraic
module Algebraic =
let ofInt (i : int) = Rational (Rational.ofInt i)
let ofRational (r : Rational) = Rational r
let add (a1 : Algebraic) (a2 : Algebraic) : Algebraic =
Sum (a1, a2)
let times (a1 : Algebraic) (a2 : Algebraic) : Algebraic =
Times (a1, a2)
let sqrt (a1 : Algebraic) : Algebraic =
Sqrt a1
let negate (a1 : Algebraic) : Algebraic =
Negate a1
let reciprocal (a1 : Algebraic) : Algebraic =
Reciprocal a1
let equal (a1 : Algebraic) (a2 : Algebraic) : bool =
failwith "TODO"

54
RayTracing/Camera.fs
View File

@@ -1,66 +1,58 @@
namespace RayTracing
type Camera<'a> =
type Camera =
{
Num : Num<'a>
/// How tall is our viewport?
ViewportHeight : 'a
ViewportHeight : float
/// How wide is our viewport?
ViewportWidth : 'a
ViewportWidth : float
/// In which direction is the camera pointing?
View : Ray<'a>
View : Ray
/// What is the orientation of the imaginary plane
/// onto which we're collecting the pixels of the result?
/// This is normal to View and to ViewportYAxis, and its
/// origin is at distance FocalLength from View.Origin.
ViewportXAxis : Ray<'a>
ViewportXAxis : Ray
/// What is the orientation of the imaginary plane
/// onto which we're collecting the pixels of the result?
/// This is normal to View and to ViewportXAxis, and its
/// origin is at distance FocalLength from View.Origin.
ViewportYAxis : Ray<'a>
ViewportYAxis : Ray
/// How far away from the camera is the imaginary plane
/// onto which we're collecting the pixels of the result?
FocalLength : 'a
FocalLength : float
/// How many samples will we take per pixel, for anti-aliasing?
SamplesPerPixel : int
}
[<RequireQualifiedAccess>]
module Camera =
let makeBasic<'a>
(n : Num<'a>)
(focalLength : 'a)
(aspectRatio : 'a)
(origin : Point<'a>)
: Camera<'a>
let makeBasic
(focalLength : float)
(aspectRatio : float)
(origin : Point)
: Camera
=
let height = n.Double n.One
let height = 2.0
let basis = UnitVector.basis 3
let view =
{
Origin = origin
Vector = Vector [| n.Zero ; n.Zero ; n.One |]
}
basis.[2]
|> Ray.make origin
let xAxis =
{
Origin = origin
Vector = Vector [| n.One ; n.Zero ; n.Zero |]
}
basis.[0]
|> Ray.make origin
let yAxis =
{
Origin = origin
Vector = Vector [| n.Zero ; n.One ; n.Zero |]
}
basis.[1]
|> Ray.make origin
{
Num = n
FocalLength = focalLength
ViewportHeight = height
ViewportWidth = n.Times aspectRatio height
ViewportWidth = aspectRatio * height
View = view
ViewportXAxis =
Ray.parallelTo (Ray.walkAlong n view focalLength) xAxis
Ray.parallelTo (Ray.walkAlong view focalLength) xAxis
ViewportYAxis =
Ray.parallelTo (Ray.walkAlong n view focalLength) yAxis
Ray.parallelTo (Ray.walkAlong view focalLength) yAxis
SamplesPerPixel = 10
}

50
RayTracing/Domain.fs
View File

@@ -14,53 +14,3 @@ module Image =
let colCount (Image i) : int = i.[0].Length
[<RequireQualifiedAccess>]
module Vector =
let dot<'a> (num : Num<'a>) (p1 : Vector<'a>) (p2 : Vector<'a>) : 'a =
match p1, p2 with
| Vector p1, Vector p2 ->
let mutable answer = num.Zero
for i in 0..p1.Length - 1 do
answer <- num.Add answer (num.Times p1.[i] p2.[i])
answer
let scale<'a> (num : Num<'a>) (scale : 'a) (vec : Vector<'a>) : Vector<'a> =
match vec with
| Vector vec ->
vec
|> Array.map (fun i -> num.Times scale i)
|> Vector
let difference<'a> (num : Num<'a>) (v1 : Vector<'a>) (v2 : Vector<'a>) : Vector<'a> =
match v1, v2 with
| Vector v1, Vector v2 ->
Array.zip v1 v2
|> Array.map (fun (a, b) -> num.Subtract a b)
|> Vector
let unitise<'a> (num : Num<'a>) (vec : Vector<'a>) : Vector<'a> option =
let dot = dot num vec vec
match num.Compare dot num.Zero with
| Equal -> None
| _ ->
let factor = dot |> num.Reciprocal |> num.Sqrt
scale num factor vec
|> Some
let normSquared<'a> (num : Num<'a>) (vec : Vector<'a>) : 'a =
dot num vec vec
let equal<'a> (num : Num<'a>) (v1 : Vector<'a>) (v2 : Vector<'a>) : bool =
match v1, v2 with
| Vector p1, Vector p2 ->
Array.zip p1 p2
|> Array.forall (fun (a, b) -> num.Equal a b)
let rec randomUnit<'a> (num : Num<'a>) (rand : Random) (dimension : int) : Vector<'a> =
let vector =
Array.init dimension (fun _ -> num.Subtract (num.TimesInteger 2 (num.RandomBetween01 rand)) num.One)
|> Vector
|> unitise num
match vector with
| None -> randomUnit num rand dimension
| Some result -> result

37
RayTracing/Float.fs Normal file
View File

@@ -0,0 +1,37 @@
namespace RayTracing
open System
type Comparison =
| Greater
| Equal
| Less
[<RequireQualifiedAccess>]
module Float =
let Pi = acos -1.0
let tolerance = 0.00000001
let inline equal (a : float) (b : float) : bool =
abs (a - b) < tolerance
let inline random (rand : Random) : float =
float (abs (rand.Next ())) / float Int32.MaxValue
let inline compare (a : float) (b : float) : Comparison =
if abs (a - b) < tolerance then Comparison.Equal
elif a < b then Comparison.Less
else Comparison.Greater
let sortInPlaceBy<'b> (proj : 'b -> float) (a : 'b array) : 'b array =
for i in 0..a.Length - 2 do
for j in i+1..a.Length - 1 do
match compare (proj a.[i]) (proj a.[j]) with
| Greater ->
let tmp = a.[j]
a.[j] <- a.[i]
a.[i] <- tmp
| _ -> ()
a

97
RayTracing/InfinitePlane.fs
View File

@@ -2,25 +2,25 @@ namespace RayTracing
open System
type InfinitePlaneStyle<'a> =
type InfinitePlaneStyle =
/// An emitter of light.
| LightSource of Pixel
/// Perfect reflection, as you would see from a smooth flat metal surface.
/// Albedo must be between 0 and 1.
| PureReflection of albedo : 'a * colour : Pixel
| PureReflection of albedo : float * colour : Pixel
/// An ideal matte (diffusely-reflecting) surface: apparent brightness of the
/// surface is the same regardless of the angle of view.
/// Albedo must be between 0 and 1.
| LambertReflection of albedo : 'a * colour : Pixel * Random
| LambertReflection of albedo : float * colour : Pixel * Random
type InfinitePlane<'a> =
type InfinitePlane =
{
Normal : Vector<'a>
Point : Point<'a>
Normal : UnitVector
Point : Point
/// If an incoming ray has the given colour, and hits the
/// given point (which is guaranteed to be on the surface),
/// what colour ray does it output and in what direction?
Reflection : Ray<'a> -> Pixel -> Point<'a> -> Ray<'a> option * Pixel
Reflection : Ray -> Pixel -> Point -> Ray option * Pixel
}
[<RequireQualifiedAccess>]
@@ -31,33 +31,32 @@ module InfinitePlane =
/// Does not return any intersections which are behind us.
/// If the plane is made of a material which does not re-emit light, you'll
/// get a None for the outgoing ray.
let intersections<'a>
(num : Num<'a>)
(plane : InfinitePlane<'a>)
(ray : Ray<'a>)
let intersections
(plane : InfinitePlane)
(ray : Ray)
(incomingColour : Pixel)
: (Point<'a> * Ray<'a> option * Pixel) array
: (Point * Ray option * Pixel) array
=
// ((ray.Origin - plane.Point) + t ray.Vector) . plane.Normal = 0
let rayVec = ray.Vector |> Vector.unitise num |> Option.get
let denominator = Vector.dot num plane.Normal rayVec
if num.Equal denominator num.Zero then [||]
let rayVec = Ray.vector ray
let denominator = UnitVector.dot plane.Normal rayVec
if Float.equal denominator 0.0 then [||]
else
let t = num.Divide (Vector.dot num (Point.difference num plane.Point ray.Origin) plane.Normal) denominator
match num.Compare t num.Zero with
// TODO I flipped the args in this dot
let t = (UnitVector.dot' plane.Normal (Point.difference plane.Point (Ray.origin ray))) / denominator
match Float.compare t 0.0 with
| Greater ->
let strikePoint = Ray.walkAlong num { Origin = ray.Origin ; Vector = rayVec } t
let strikePoint = Ray.walkAlong ray t
let outgoing, newColour = plane.Reflection ray incomingColour strikePoint
[| strikePoint, outgoing, newColour |]
| _ -> [||]
let reflection<'a>
(num : Num<'a>)
(style : InfinitePlaneStyle<'a>)
(pointOnPlane : Point<'a>)
(normal : Vector<'a>)
: Ray<'a> -> Pixel -> Point<'a> -> Ray<'a> option * Pixel
let reflection
(style : InfinitePlaneStyle)
(pointOnPlane : Point)
(normal : UnitVector)
: Ray -> Pixel -> Point -> Ray option * Pixel
=
fun incomingRay incomingColour strikePoint ->
match style with
@@ -66,52 +65,46 @@ module InfinitePlane =
| InfinitePlaneStyle.LambertReflection (albedo, colour, rand) ->
let outgoing =
{
Origin = strikePoint
Vector =
let (Point pointOnPlane) = pointOnPlane
let sphereCentre = Ray.walkAlong num { Origin = strikePoint ; Vector = normal } num.One
let offset = Vector.randomUnit num rand pointOnPlane.Length
let target = Ray.walkAlong num { Origin = sphereCentre ; Vector = offset } num.One
Point.difference num target strikePoint
}
let (Point pointOnPlane) = pointOnPlane
let sphereCentre = Ray.walkAlong (Ray.make strikePoint normal) 1.0
let offset = UnitVector.random rand pointOnPlane.Length
let target = Ray.walkAlong (Ray.make sphereCentre offset) 1.0
Point.difference target strikePoint
|> Ray.make' strikePoint
let newColour = Pixel.combine incomingColour colour
Some outgoing, Pixel.darken num newColour albedo
outgoing, Pixel.darken newColour albedo
| InfinitePlaneStyle.PureReflection (albedo, colour) ->
let plane =
Plane.makeSpannedBy { Origin = strikePoint ; Vector = normal } incomingRay
|> Plane.orthonormalise num
Plane.makeSpannedBy (Ray.make strikePoint normal) incomingRay
|> Plane.orthonormalise
let outgoing =
match plane with
| None ->
// Incoming ray is directly along the normal
{
Origin = strikePoint
Vector = incomingRay.Vector |> Vector.scale num (num.Negate num.One)
}
Ray.flip incomingRay
|> Ray.parallelTo strikePoint
|> Some
| Some plane ->
// Incoming ray is (plane1.ray) plane1 + (plane2.ray) plane2
// We want the reflection in the normal, so need (plane1.ray) plane1 - (plane2.ray) plane2
let normalComponent = (Vector.dot num plane.V1 incomingRay.Vector)
let tangentComponent = num.Negate (Vector.dot num plane.V2 incomingRay.Vector)
{
Origin = strikePoint
Vector =
Ray.walkAlong num { Origin = Ray.walkAlong num { Origin = plane.Point ; Vector = plane.V1 } normalComponent ; Vector = plane.V2 } tangentComponent
|> Point.difference num strikePoint
}
let normalComponent = (UnitVector.dot plane.V1 (Ray.vector incomingRay))
let tangentComponent = - (UnitVector.dot plane.V2 (Ray.vector incomingRay))
tangentComponent
|> Ray.walkAlong (Ray.make (Ray.walkAlong (Ray.make plane.Point plane.V1) normalComponent) plane.V2)
|> Point.difference strikePoint
|> Ray.make' strikePoint
let newColour = Pixel.combine incomingColour colour
let darkened = Pixel.darken num newColour albedo
Some outgoing, darkened
let darkened = Pixel.darken newColour albedo
outgoing, darkened
let make<'a> (num : Num<'a>) (style : InfinitePlaneStyle<'a>) (pointOnPlane : Point<'a>) (normal : Vector<'a>) : InfinitePlane<'a> =
let make (style : InfinitePlaneStyle) (pointOnPlane : Point) (normal : UnitVector) : InfinitePlane =
{
Point = pointOnPlane
Normal = normal
Reflection = reflection num style pointOnPlane normal
Reflection = reflection style pointOnPlane normal
}

109
RayTracing/Num.fs
View File

@@ -1,109 +0,0 @@
namespace RayTracing
open System
type 'a Radian =
| Radian of 'a
type Comparison =
| Greater
| Equal
| Less
type Num<'a> =
{
Add : 'a -> 'a -> 'a
Times : 'a -> 'a -> 'a
Negate : 'a -> 'a
Reciprocal : 'a -> 'a
Zero : 'a
Compare : 'a -> 'a -> Comparison
Sqrt : 'a -> 'a
Equal : 'a -> 'a -> bool
TimesInteger : int -> 'a -> 'a
DivideInteger : 'a -> int -> 'a
One : 'a
RandomBetween01 : Random -> 'a
ArcCos : 'a -> 'a Radian
// arctan(second / first)
ArcTan2 : 'a -> 'a -> 'a Radian
Cos : 'a Radian -> 'a
Sin : 'a Radian -> 'a
Round : 'a -> int
}
member this.Double (x : 'a) : 'a = this.Add x x
member this.Subtract (x : 'a) (y : 'a) : 'a = this.Add x (this.Negate y)
member this.Divide (x : 'a) (y : 'a) : 'a = this.Times x (this.Reciprocal y)
member this.Pi =
let (Radian t) = this.ArcCos (this.Negate this.One)
t
[<RequireQualifiedAccess>]
module Num =
let float : Num<float> =
let tolerance = 0.00000001
{
Add = (+)
Times = (*)
Negate = fun x -> -x
Zero = 0.0
Reciprocal = fun i -> 1.0 / i
Compare =
fun a b ->
if abs (a - b) < tolerance then Comparison.Equal
elif a < b then Comparison.Less
else Comparison.Greater
Sqrt = sqrt
Equal = fun a b -> abs (a - b) < tolerance
TimesInteger = fun a b -> float a * b
DivideInteger = fun a b -> a / float b
One = 1.0
RandomBetween01 = fun rand -> float (abs (rand.Next ())) / float Int32.MaxValue
ArcCos = acos >> Radian
ArcTan2 = fun x -> atan2 x >> Radian
Sin = fun (Radian r) -> sin r
Cos = fun (Radian r) -> cos r
Round = fun i -> Math.Round i |> int
}
let algebraic : Num<Algebraic> =
{
Add = Algebraic.add
Times = Algebraic.times
Negate = Algebraic.negate
Zero = Algebraic.ofInt 0
Reciprocal = Algebraic.reciprocal
Compare =
fun a b ->
if a < b then Comparison.Less
elif a = b then Comparison.Equal
else Comparison.Greater
Sqrt = Algebraic.sqrt
Equal = Algebraic.equal
TimesInteger = fun _ _ -> failwith ""
DivideInteger = fun _ _ -> failwith ""
One = Algebraic.ofInt 1
RandomBetween01 = fun _ -> failwith ""
ArcCos = fun _ -> failwith ""
ArcTan2 = fun _ -> failwith ""
Cos = fun _ -> failwith ""
Sin = fun _ -> failwith ""
Round = fun _ -> failwith ""
}
let sortInPlaceBy<'a, 'b> (num : 'a Num) (proj : 'b -> 'a) (a : 'b array) : 'b array =
for i in 0..a.Length - 2 do
for j in i+1..a.Length - 1 do
match num.Compare (proj a.[i]) (proj a.[j]) with
| Greater ->
let tmp = a.[j]
a.[j] <- a.[i]
a.[i] <- tmp
| _ -> ()
a
[<RequireQualifiedAccess>]
module Radian =
let add<'a> (n : Num<'a>) (Radian r1) (Radian r2) = n.Add r1 r2 |> Radian

31
RayTracing/Pixel.fs
View File

@@ -45,19 +45,15 @@ module Colour =
[<RequireQualifiedAccess>]
module Pixel =
let average (s : Pixel seq) : Pixel =
use e = s.GetEnumerator ()
if not (e.MoveNext ()) then failwith "Input sequence was empty when averaging pixels"
let mutable count = 1
let mutable r = e.Current.Red |> float
let mutable g = e.Current.Green |> float
let mutable b = e.Current.Blue |> float
while e.MoveNext () do
count <- count + 1
r <- r + float e.Current.Red
g <- g + float e.Current.Green
b <- b + float e.Current.Blue
let count = float count
let average (s : Pixel []) : Pixel =
let mutable r = s.[0].Red |> float
let mutable g = s.[0].Green |> float
let mutable b = s.[0].Blue |> float
for i in 1..s.Length - 1 do
r <- r + float s.[i].Red
g <- g + float s.[i].Green
b <- b + float s.[i].Blue
let count = s.Length |> float
{
Red = byte (Math.Round (r / count))
Green = byte (Math.Round (g / count))
@@ -71,9 +67,10 @@ module Pixel =
Blue = (int p1.Blue * int p2.Blue) / 255 |> byte
}
let darken<'a> (num : Num<'a>) (p : Pixel) (albedo : 'a) : Pixel =
/// albedo should be between 0 and 1.
let darken (p : Pixel) (albedo : float) : Pixel =
{
Red = num.TimesInteger (int p.Red) albedo |> num.Round |> byte
Green = num.TimesInteger (int p.Green) albedo |> num.Round |> byte
Blue = num.TimesInteger (int p.Blue) albedo |> num.Round |> byte
Red = (float p.Red) * albedo |> Math.Round |> byte
Green = (float p.Green) * albedo |> Math.Round |> byte
Blue = (float p.Blue) * albedo |> Math.Round |> byte
}

41
RayTracing/Plane.fs
View File

@@ -1,47 +1,46 @@
namespace RayTracing
/// A plane spanned by two rays from a common origin.
type 'a Plane =
type Plane =
private
{
V1 : 'a Vector
V2 : 'a Vector
Point : 'a Point
V1 : UnitVector
V2 : UnitVector
Point : Point
}
type 'a OrthonormalPlane =
type OrthonormalPlane =
{
V1 : 'a Vector
V2 : 'a Vector
Point : 'a Point
V1 : UnitVector
V2 : UnitVector
Point : Point
}
[<RequireQualifiedAccess>]
module Plane =
let orthonormalise<'a> (num : 'a Num) (plane : 'a Plane) : 'a OrthonormalPlane option =
let v1 = Vector.unitise num plane.V1 |> Option.get
let coeff = Vector.dot num v1 plane.V2
let orthonormalise (plane : Plane) : OrthonormalPlane option =
let coeff = UnitVector.dot plane.V1 plane.V2
let vec2 =
Vector.difference num plane.V2 (Vector.scale num coeff v1)
|> Vector.unitise num
UnitVector.difference' plane.V2 (UnitVector.scale coeff plane.V1)
|> Vector.unitise
match vec2 with
| None -> None
| Some v2 ->
{
V1 = v1
V1 = plane.V1
V2 = v2
Point = plane.Point
}
|> Some
let makeSpannedBy<'a> (r1 : 'a Ray) (r2 : 'a Ray) : 'a Plane =
let makeSpannedBy (r1 : Ray) (r2 : Ray) : Plane =
{
V1 = r1.Vector
V2 = r2.Vector
Point = r1.Origin
V1 = Ray.vector r1
V2 = Ray.vector r2
Point = Ray.origin r1
}
let basis<'a> (plane : 'a OrthonormalPlane) : 'a Ray * 'a Ray =
{ Origin = plane.Point ; Vector = plane.V1 }, { Origin = plane.Point ; Vector = plane.V2 }
let basis (plane : OrthonormalPlane) : Ray * Ray =
Ray.make plane.Point plane.V1,
Ray.make plane.Point plane.V2

115
RayTracing/Point.fs
View File

@@ -1,44 +1,103 @@
namespace RayTracing
open System
/// An n-dimensional point.
/// We don't let you compare these for equality, because floats are hard.
[<NoEquality ; NoComparison>]
type Point<'a> = Point of 'a array
type Point = Point of float array
[<NoEquality ; NoComparison>]
type Vector<'a> = Vector of 'a array
type Vector = Vector of float array
type UnitVector = UnitVector of Vector
[<RequireQualifiedAccess>]
module Vector =
let dot (p1 : Vector) (p2 : Vector) : float =
match p1, p2 with
| Vector p1, Vector p2 ->
let mutable answer = 0.0
for i in 0..p1.Length - 1 do
answer <- answer + (p1.[i] * p2.[i])
answer
let scale (scale : float) (vec : Vector) : Vector =
match vec with
| Vector vec ->
vec
|> Array.map (fun i -> scale * i)
|> Vector
let difference (v1 : Vector) (v2 : Vector) : Vector =
match v1, v2 with
| Vector v1, Vector v2 ->
let answer = Array.zeroCreate v1.Length
for i in 0..answer.Length - 1 do
answer.[i] <- v1.[i] - v2.[i]
answer
|> Vector
let unitise (vec : Vector) : UnitVector option =
let dot = dot vec vec
if Float.equal dot 0.0 then None else
let factor = 1.0 / sqrt dot
scale factor vec
|> UnitVector
|> Some
let normSquared (vec : Vector) : float =
dot vec vec
let equal (v1 : Vector) (v2 : Vector) : bool =
match v1, v2 with
| Vector p1, Vector p2 ->
let rec go (i : int) =
if i >= p1.Length then true else
if Float.equal p1.[i] p2.[i] then go (i + 1) else false
go 0
[<RequireQualifiedAccess>]
module UnitVector =
let rec random (rand : Random) (dimension : int) : UnitVector =
let vector =
Array.init dimension (fun _ -> (2.0 * Float.random rand) - 1.0)
|> Vector
|> Vector.unitise
match vector with
| None -> random rand dimension
| Some result -> result
let dot (UnitVector a) (UnitVector b) = Vector.dot a b
let dot' (UnitVector a) (b : Vector) = Vector.dot a b
let difference (UnitVector v1) (UnitVector v2) = Vector.difference v1 v2
let difference' (UnitVector v1) (v2 : Vector) = Vector.difference v1 v2
let scale (scale : float) (UnitVector vec) = Vector.scale scale vec
let basis (dimension : int) : UnitVector [] =
Array.init dimension (fun i ->
Array.init dimension (fun j ->
if i = j then 1.0 else 0.0
)
|> Vector
|> UnitVector
)
[<RequireQualifiedAccess>]
module Point =
let difference<'a> (num : Num<'a>) (p1 : Point<'a>) (p2 : Point<'a>) : Vector<'a> =
let difference (p1 : Point) (p2 : Point) : Vector =
match p1, p2 with
| Point p1, Point p2 ->
Array.zip p1 p2
|> Array.map (fun (a, b) -> num.Subtract a b)
Array.init p1.Length (fun i ->
p1.[i] - p2.[i]
)
|> Vector
let sum<'a> (num : Num<'a>) (p1 : Point<'a>) (p2 : Point<'a>) : Point<'a> =
let equal (p1 : Point) (p2 : Point) : bool =
match p1, p2 with
| Point p1, Point p2 ->
Array.zip p1 p2
|> Array.map (fun (a, b) -> num.Add a b)
|> Point
let normSquared<'a> (num : Num<'a>) (p : Point<'a>) : 'a =
match p with
| Point p ->
p
|> Array.fold (fun s p -> num.Add (num.Times p p) s) num.Zero
let equal<'a> (num : Num<'a>) (p1 : Point<'a>) (p2 : Point<'a>) : bool =
match p1, p2 with
| Point p1, Point p2 ->
Array.zip p1 p2
|> Array.forall (fun (a, b) -> num.Equal a b)
let add<'a> (num : Num<'a>) (v1 : Point<'a>) (v2 : Point<'a>) : Point<'a> =
match v1, v2 with
| Point v1, Point v2 ->
Array.zip v1 v2
|> Array.map (fun (a, b) -> num.Add a b)
|> Point
let rec go (i : int) : bool =
if i >= p1.Length then true else
if Float.equal p1.[i] p2.[i] then go (i + 1) else false
go 0

124
RayTracing/Rational.fs
View File

@@ -1,124 +0,0 @@
namespace RayTracing
open System.Numerics
[<Struct>]
[<CustomEquality; CustomComparison>]
type Rational =
private
{
Num : BigInteger
Denom : BigInteger
IsNormal : bool
}
static member Numerator (r : Rational) : BigInteger = r.Num
static member Denominator (r : Rational) : BigInteger = r.Denom + BigInteger.One
static member Normalise (r : Rational) : Rational =
if r.IsNormal then
r
else
let rec gcd (a : BigInteger) (b : BigInteger) : BigInteger =
if a.Sign = -1 then -gcd (-a) b
elif b.Sign = -1 then -gcd a (-b)
elif a.IsZero then b
elif b.IsZero then a
else if a > b then gcd b (a % b)
elif a = b then a
else gcd b a
let gcd =
gcd (Rational.Numerator r) (Rational.Denominator r)
{ Rational.Make (Rational.Numerator r / gcd) (Rational.Denominator r / gcd) with
IsNormal = true
}
member this.Normalise () = Rational.Normalise this
static member Make (num : BigInteger) (denom : BigInteger) : Rational =
if denom.IsZero then
failwith "Invalid zero denominator"
elif denom.Sign = -1 then
{
Num = -num
Denom = (-denom) - BigInteger.One
IsNormal = false
}
else
{
Num = num
Denom = denom - BigInteger.One
IsNormal = false
}
override this.Equals (other : obj) : bool =
match other with
| :? Rational as other ->
printfn "%+A" other
match this.Normalise () with
| { Num = num; Denom = denom } ->
match other.Normalise () with
| { Num = numOther; Denom = denomOther } ->
printfn "%+A %+A %+A %+A" numOther num denom denomOther
numOther = num && denom = denomOther
| _ -> failwith "how did you do this"
override this.GetHashCode () : int =
let n = this.Normalise ()
hash (n.Num, n.Denom)
interface System.IComparable<Rational> with
member this.CompareTo (other : Rational) =
let this = this.Normalise ()
let other = other.Normalise ()
let cmp =
Rational.Numerator this
* Rational.Denominator other
- Rational.Numerator other
* Rational.Denominator this
cmp.Sign
interface System.IComparable with
member this.CompareTo (other : obj) =
match other with
| :? Rational as other ->
(this :> System.IComparable<Rational>)
.CompareTo other
| _ -> failwith "how?!"
[<RequireQualifiedAccess>]
module Rational =
let inline make a b = Rational.Make a b
let add (r1 : Rational) (r2 : Rational) =
Rational.Make
(Rational.Numerator r1 * Rational.Denominator r2
+ Rational.Numerator r2 * Rational.Denominator r1)
(Rational.Denominator r1 * Rational.Denominator r2)
|> Rational.Normalise
let ofInt (i : int) : Rational = { Num = BigInteger i; Denom = BigInteger.Zero; IsNormal = true }
let times (r1 : Rational) (r2 : Rational) =
Rational.Make
(Rational.Numerator r1 * Rational.Numerator r2)
(Rational.Denominator r1 * Rational.Denominator r2)
|> Rational.Normalise
let subtract (r1 : Rational) (r2 : Rational) : Rational =
Rational.Make
(Rational.Numerator r1 * Rational.Denominator r2
- Rational.Numerator r2 * Rational.Denominator r1)
(Rational.Denominator r1 * Rational.Denominator r2)
|> Rational.Normalise
let reciprocal (r : Rational) : Rational =
Rational.Make (Rational.Denominator r) (Rational.Numerator r)
let divide (r1 : Rational) (r2 : Rational) : Rational =
times r1 (reciprocal r2)

82
RayTracing/Ray.fs
View File

@@ -1,59 +1,73 @@
namespace RayTracing
type Ray<'a> =
{
Origin : Point<'a>
Vector : Vector<'a>
}
type Ray =
private
{
Origin : Point
Vector : UnitVector
}
[<RequireQualifiedAccess>]
module Ray =
let walkAlong<'a> (num : Num<'a>) (ray : Ray<'a>) (magnitude : 'a) : Point<'a> =
let (Point origin) = ray.Origin
let (Vector vector) = ray.Vector |> Vector.unitise num |> Option.get
let make' (origin : Point) (vector : Vector) : Ray option =
match Vector.unitise vector with
| None -> None
| Some v ->
{
Origin = origin
Vector = v
}
|> Some
Array.zip origin vector
|> Array.map (fun (originCoord, directionCoord) -> num.Add originCoord (num.Times directionCoord magnitude))
|> Point
let between<'a> (num : Num<'a>) (p1 : Point<'a>) (p2 : Point<'a>) : Ray<'a> =
let make (origin : Point) (vector : UnitVector) : Ray =
{
Origin = p1
Vector = Point.difference num p2 p1
Origin = origin
Vector = vector
}
/// Given two rays from the same point, what is the angle between them?
let angle<'a> (num : Num<'a>) (r1 : Ray<'a>) (r2 : Ray<'a>) : 'a Radian =
// a.b = |a| |b| cos theta
let v1 = walkAlong num { r1 with Origin = r2.Origin } num.One
let v2 = walkAlong num r2 num.One
let (Radian answer) = num.ArcCos (Vector.dot num (Point.difference num v1 r2.Origin) (Point.difference num v2 r2.Origin))
match num.Compare (num.Double answer) num.Pi with
| Greater ->
num.Subtract num.Pi answer
| _ ->
answer
|> Radian
let walkAlong (ray : Ray) (magnitude : float) : Point =
let (Point origin) = ray.Origin
let (UnitVector (Vector vector)) = ray.Vector
let parallelTo<'a> (p1 : Point<'a>) (ray : Ray<'a>) : Ray<'a> =
Array.init origin.Length (fun i ->
origin.[i] + (vector.[i] * magnitude)
)
|> Point
let between (p1 : Point) (p2 : Point) : Ray option =
make' p1 (Point.difference p2 p1)
let parallelTo (p1 : Point) (ray : Ray) : Ray =
{
Vector = ray.Vector
Origin = p1
}
let liesOn<'a> (num : 'a Num) (point : Point<'a>) (ray : Ray<'a>) : bool =
let liesOn (point : Point) (ray : Ray) : bool =
match point, ray.Origin, ray.Vector with
| Point x, Point y, Vector ray ->
let rec go (state : 'a option) (i : int) =
| Point x, Point y, UnitVector (Vector ray) ->
let rec go (state : float option) (i : int) =
if i >= x.Length then state else
let d = x.[i]
let x = y.[i]
let r = ray.[i]
match state with
| None -> go (Some (num.Divide (num.Subtract d x) r)) (i + 1)
| None -> go (Some ((d - x) / r)) (i + 1)
| Some prevT ->
let t = num.Divide (num.Subtract d x) r
if num.Equal prevT t then go (Some prevT) (i + 1) else None
let t = (d - x) / r
if Float.equal prevT t then go (Some prevT) (i + 1) else None
go None 0
|> Option.isSome
let vector r = r.Vector
let origin r = r.Origin
let flip (r : Ray) =
{
Origin = r.Origin
Vector =
let (UnitVector v) = r.Vector
Vector.scale -1.0 v
|> UnitVector
}

21
RayTracing/Ray.fsi Normal file
View File

@@ -0,0 +1,21 @@
namespace RayTracing
type Ray
[<RequireQualifiedAccess>]
module Ray =
val make' : Point -> Vector -> Ray option
val make : Point -> UnitVector -> Ray
val walkAlong : Ray -> float -> Point
val between : Point -> Point -> Ray option
val parallelTo : Point -> Ray -> Ray
val liesOn : Point -> Ray -> bool
val vector : Ray -> UnitVector
val origin : Ray -> Point
val flip : Ray -> Ray

6
RayTracing/RayTracing.fsproj
View File

@@ -5,12 +5,11 @@
</PropertyGroup>
<ItemGroup>
<Compile Include="Rational.fs" />
<Compile Include="Algebraic.fs" />
<Compile Include="Num.fs" />
<Compile Include="Float.fs" />
<Compile Include="Point.fs" />
<Compile Include="Pixel.fs" />
<Compile Include="Domain.fs" />
<Compile Include="Ray.fsi" />
<Compile Include="Ray.fs" />
<Compile Include="Plane.fs" />
<Compile Include="Sphere.fs" />
@@ -23,7 +22,6 @@
<ItemGroup>
<PackageReference Include="System.IO.Abstractions" Version="13.2.28" />
<PackageReference Include="System.Runtime.Numerics" Version="4.3.0" />
</ItemGroup>
</Project>

18
RayTracing/SampleImages.fs
View File

@@ -34,26 +34,26 @@ module SampleImages =
let random = Random ()
let aspectRatio = 16.0 / 9.0
let camera =
Camera.makeBasic Num.float 4.0 aspectRatio (Point [| 0.0 ; 0.0 ; 0.0 |])
let pixels = 1800
Camera.makeBasic 4.0 aspectRatio (Point [| 0.0 ; 0.0 ; 0.0 |])
let pixels = 200
{
Objects =
[|
Hittable.Sphere (Sphere.make Num.float (SphereStyle.LambertReflection (1.0, { Red = 255uy ; Green = 255uy ; Blue = 0uy }, random)) (Point [| 0.0 ; 0.0 ; 9.0 |]) 1.0)
Hittable.Sphere (Sphere.make Num.float (SphereStyle.PureReflection (1.0, { Red = 0uy ; Green = 255uy ; Blue = 255uy })) (Point [| 1.5 ; 0.5 ; 8.0 |]) 0.5)
Hittable.Sphere (Sphere.make Num.float (SphereStyle.LightSource Colour.Blue) (Point [| -1.5 ; 1.5 ; 8.0 |]) 0.5)
Hittable.Sphere (Sphere.make (SphereStyle.LambertReflection (1.0, { Red = 255uy ; Green = 255uy ; Blue = 0uy }, random)) (Point [| 0.0 ; 0.0 ; 9.0 |]) 1.0)
Hittable.Sphere (Sphere.make (SphereStyle.PureReflection (1.0, { Red = 0uy ; Green = 255uy ; Blue = 255uy })) (Point [| 1.5 ; 0.5 ; 8.0 |]) 0.5)
Hittable.Sphere (Sphere.make (SphereStyle.LightSource Colour.Blue) (Point [| -1.5 ; 1.5 ; 8.0 |]) 0.5)
// Side mirror
Hittable.InfinitePlane (InfinitePlane.make Num.float (InfinitePlaneStyle.PureReflection (1.0, { Colour.White with Green = 240uy })) (Point [| 0.1 ; 0.0 ; 16.0 |]) (Vector [| -2.0 ; 0.0 ; -1.0 |]))
Hittable.InfinitePlane (InfinitePlane.make (InfinitePlaneStyle.PureReflection (1.0, { Colour.White with Green = 240uy })) (Point [| 0.1 ; 0.0 ; 16.0 |]) (Vector [| -2.0 ; 0.0 ; -1.0 |] |> Vector.unitise |> Option.get))
// Floor mirror
Hittable.InfinitePlane (InfinitePlane.make Num.float (InfinitePlaneStyle.PureReflection (0.4, Colour.White)) (Point [| 0.0 ; -1.0 ; 0.0 |]) (Vector [| 0.0 ; 1.0 ; 0.0 |]))
Hittable.InfinitePlane (InfinitePlane.make (InfinitePlaneStyle.PureReflection (0.4, Colour.White)) (Point [| 0.0 ; -1.0 ; 0.0 |]) (Vector [| 0.0 ; 1.0 ; 0.0 |] |> Vector.unitise |> Option.get))
// Back plane
Hittable.InfinitePlane (InfinitePlane.make Num.float (InfinitePlaneStyle.PureReflection (0.6, Colour.White)) (Point [| 0.0 ; 0.0 ; 16.0 |]) (Vector [| 0.0 ; 0.0 ; -1.0 |]))
Hittable.InfinitePlane (InfinitePlane.make (InfinitePlaneStyle.PureReflection (0.6, Colour.White)) (Point [| 0.0 ; 0.0 ; 16.0 |]) (Vector [| 0.0 ; 0.0 ; -1.0 |] |> Vector.unitise |> Option.get))
// Light pad behind us
Hittable.InfinitePlane (InfinitePlane.make Num.float (InfinitePlaneStyle.LightSource Colour.White) (Point [| 0.0 ; 1.0 ; -1.0 |]) (Vector [| 0.0 ; -1.0 ; 1.0 |]))
Hittable.InfinitePlane (InfinitePlane.make (InfinitePlaneStyle.LightSource Colour.White) (Point [| 0.0 ; 1.0 ; -1.0 |]) (Vector [| 0.0 ; -1.0 ; 1.0 |] |> Vector.unitise |> Option.get))
|]
}
|> Scene.render progressIncrement (aspectRatio * (float pixels) |> int) pixels camera

70
RayTracing/Scene.fs
View File

@@ -2,58 +2,55 @@ namespace RayTracing
open System
type Hittable<'a> =
| Sphere of Sphere<'a>
| InfinitePlane of InfinitePlane<'a>
type Hittable =
| Sphere of Sphere
| InfinitePlane of InfinitePlane
[<RequireQualifiedAccess>]
module Hittable =
let hits<'a>
(num : Num<'a>)
(ray : Ray<'a>)
let hits
(ray : Ray)
(incomingColour : Pixel)
(h : Hittable<'a>)
: (Point<'a> * Ray<'a> option * Pixel) option
(h : Hittable)
: (Point * Ray option * Pixel) option
=
match h with
| Sphere s ->
Sphere.intersections num s ray incomingColour
Sphere.intersections s ray incomingColour
|> Array.tryHead
| InfinitePlane plane ->
InfinitePlane.intersections num plane ray incomingColour
InfinitePlane.intersections plane ray incomingColour
|> Array.tryHead
type Scene<'a> =
type Scene =
{
Objects : Hittable<'a> array
Objects : Hittable array
}
[<RequireQualifiedAccess>]
module Scene =
let hitObject<'a>
(num : Num<'a>)
(s : Scene<'a>)
(ray : Ray<'a>)
let hitObject
(s : Scene)
(ray : Ray)
(colour : Pixel)
: (Point<'a> * Ray<'a> option * Pixel) array
: (Point * Ray option * Pixel) array
=
s.Objects
|> Array.choose (Hittable.hits num ray colour)
|> Num.sortInPlaceBy num (fun (a, _, _) -> Vector.normSquared num (Point.difference num a ray.Origin))
|> Array.choose (Hittable.hits ray colour)
|> Float.sortInPlaceBy (fun (a, _, _) -> Vector.normSquared (Point.difference a (Ray.origin ray)))
let internal traceRay<'a>
let internal traceRay
(maxCount : int)
(num : Num<'a>)
(scene : Scene<'a>)
(ray : Ray<'a>)
(scene : Scene)
(ray : Ray)
(colour : Pixel)
: Pixel
=
let rec go (bounces : int) (ray : Ray<'a>) (colour : Pixel) : Pixel =
let rec go (bounces : int) (ray : Ray) (colour : Pixel) : Pixel =
if bounces > maxCount then Colour.Black else
let thingsWeHit = hitObject num scene ray colour
let thingsWeHit = hitObject scene ray colour
match thingsWeHit with
| [||] ->
// Ray goes off into the distance and is never heard from again
@@ -64,20 +61,19 @@ module Scene =
| None ->
colour
| Some outgoingRay ->
go (bounces + 1) { outgoingRay with Vector = Vector.unitise num outgoingRay.Vector |> Option.get } colour
go (bounces + 1) outgoingRay colour
go 0 ray colour
let render<'a>
let render
(progressIncrement : float<progress> -> unit)
(maxWidthCoord : int)
(maxHeightCoord : int)
(camera : Camera<'a>)
(s : Scene<'a>)
(camera : Camera)
(s : Scene)
: float<progress> * Image Async
=
let rand = Random ()
let num = camera.Num
// For each pixel in the output, send a ray from the camera
// in the direction of that pixel.
let rowsIter = 2 * maxHeightCoord + 1
@@ -95,15 +91,17 @@ module Scene =
Array.init camera.SamplesPerPixel (fun _ ->
// TODO make this be deterministic
let pointOnXAxis =
num.DivideInteger (num.Add (num.TimesInteger col camera.ViewportWidth) (num.RandomBetween01 rand)) maxWidthCoord
|> Ray.walkAlong num camera.ViewportXAxis
((float col * camera.ViewportWidth) + (Float.random rand)) / float maxWidthCoord
|> Ray.walkAlong camera.ViewportXAxis
let toWalkUp = Ray.parallelTo pointOnXAxis camera.ViewportYAxis
let endPoint =
num.DivideInteger (num.Add (num.TimesInteger row camera.ViewportHeight) (num.RandomBetween01 rand)) maxHeightCoord
|> Ray.walkAlong num toWalkUp
let ray = Ray.between num camera.View.Origin endPoint
((float row * camera.ViewportHeight) + (Float.random rand)) / float maxHeightCoord
|> Ray.walkAlong toWalkUp
let ray =
Ray.between (Ray.origin camera.View) endPoint
|> Option.get
let result = traceRay 50 num s ray Colour.White
let result = traceRay 50 s ray Colour.White
result
)
|> Pixel.average

140
RayTracing/Sphere.fs
View File

@@ -2,46 +2,43 @@ namespace RayTracing
open System
type Sphere<'a> =
type Sphere =
{
Centre : Point<'a>
Radius : 'a
Centre : Point
Radius : float
/// If an incoming ray has the given colour, and hits the
/// given point (which is guaranteed to be on the surface),
/// what colour ray does it output and in what direction?
Reflection : Ray<'a> -> Pixel -> Point<'a> -> Ray<'a> option * Pixel
Reflection : Ray -> Pixel -> Point -> Ray option * Pixel
}
type SphereStyle<'a> =
type SphereStyle =
/// An emitter of light.
| LightSource of Pixel
/// An absorbing black sphere, with a small light-emitting cap.
| LightSourceCap of Pixel
/// Perfect reflection, as you would see from a smooth flat metal surface.
/// Albedo must be between 0 and 1.
| PureReflection of albedo : 'a * colour : Pixel
| PureReflection of albedo : float * colour : Pixel
/// An ideal matte (diffusely-reflecting) surface: apparent brightness of the
/// surface is the same regardless of the angle of view.
/// Albedo must be between 0 and 1.
| LambertReflection of albedo : 'a * colour : Pixel * Random
| LambertReflection of albedo : float * colour : Pixel * Random
[<RequireQualifiedAccess>]
module Sphere =
let normal<'a> (num : Num<'a>) (centre : Point<'a>) (p : Point<'a>) : Ray<'a> =
{
Origin = p
Vector = Point.difference num p centre
}
let normal (centre : Point) (p : Point) : Ray =
Ray.make' p (Point.difference p centre)
|> Option.get
let reflection<'a>
(num : Num<'a>)
(style : SphereStyle<'a>)
(centre : Point<'a>)
(radius : 'a)
: Ray<'a> -> Pixel -> Point<'a> -> Ray<'a> option * Pixel
let reflection
(style : SphereStyle)
(centre : Point)
(radius : float)
: Ray -> Pixel -> Point -> Ray option * Pixel
=
let normal = normal num centre
let normal = normal centre
fun incomingRay incomingColour strikePoint ->
let normal = normal strikePoint
@@ -52,12 +49,12 @@ module Sphere =
let circleCentreZCoord =
match centre with
| Point v -> Array.head v
let zCoordLowerBound = num.Add circleCentreZCoord (num.Subtract radius (num.DivideInteger radius 4))
let zCoordLowerBound = circleCentreZCoord + (radius - (radius / 4.0))
let strikeZCoord =
match strikePoint with
| Point v -> Array.head v
let colour =
match num.Compare strikeZCoord zCoordLowerBound with
match Float.compare strikeZCoord zCoordLowerBound with
| Greater ->
Pixel.combine colour incomingColour
| _ ->
@@ -66,122 +63,97 @@ module Sphere =
| SphereStyle.LambertReflection (albedo, colour, rand) ->
let outgoing =
{
Origin = strikePoint
Vector =
let (Point centre) = centre
let sphereCentre = Ray.walkAlong num normal num.One
let offset = Vector.randomUnit num rand centre.Length
let target = Ray.walkAlong num { Origin = sphereCentre ; Vector = offset } num.One
Point.difference num target strikePoint
}
let (Point centre) = centre
let sphereCentre = Ray.walkAlong normal 1.0
let offset = UnitVector.random rand centre.Length
let target = Ray.walkAlong (Ray.make sphereCentre offset) 1.0
Point.difference target strikePoint
|> Ray.make' strikePoint
let newColour = Pixel.combine incomingColour colour
Some outgoing, Pixel.darken num newColour albedo
outgoing, Pixel.darken newColour albedo
| SphereStyle.PureReflection (albedo, colour) ->
let plane =
Plane.makeSpannedBy normal incomingRay
|> Plane.orthonormalise num
|> Plane.orthonormalise
let outgoing =
match plane with
| None ->
// Incoming ray is directly along the normal
{
Origin = strikePoint
Vector = incomingRay.Vector |> Vector.scale num (num.Negate num.One)
}
Ray.flip incomingRay
|> Ray.parallelTo strikePoint
|> Some
| Some plane ->
// Incoming ray is (plane1.ray) plane1 + (plane2.ray) plane2
// We want the reflection in the normal, so need (plane1.ray) plane1 - (plane2.ray) plane2
let normalComponent = (Vector.dot num plane.V1 incomingRay.Vector)
let tangentComponent = num.Negate (Vector.dot num plane.V2 incomingRay.Vector)
{
Origin = strikePoint
Vector =
Ray.walkAlong num { Origin = Ray.walkAlong num { Origin = plane.Point ; Vector = plane.V1 } normalComponent ; Vector = plane.V2 } tangentComponent
|> Point.difference num strikePoint
}
let normalComponent = UnitVector.dot plane.V1 (Ray.vector incomingRay)
let tangentComponent = - (UnitVector.dot plane.V2 (Ray.vector incomingRay))
Ray.walkAlong (Ray.make (Ray.walkAlong (Ray.make plane.Point plane.V1) normalComponent) plane.V2) tangentComponent
|> Point.difference strikePoint
|> Ray.make' strikePoint
let newColour = Pixel.combine incomingColour colour
let darkened = Pixel.darken num newColour albedo
Some outgoing, darkened
let darkened = Pixel.darken newColour albedo
outgoing, darkened
let make<'a> (num : Num<'a>) (style : SphereStyle<'a>) (centre : Point<'a>) (radius : 'a) : Sphere<'a> =
let make (style : SphereStyle) (centre : Point) (radius : float) : Sphere =
{
Centre = centre
Radius = radius
Reflection = reflection num style centre radius
Reflection = reflection style centre radius
}
let liesOn<'a> (num : Num<'a>) (point : Point<'a>) (sphere : Sphere<'a>) : bool =
num.Equal (Vector.normSquared num (Point.difference num sphere.Centre point)) (num.Times sphere.Radius sphere.Radius)
let liesOn (point : Point) (sphere : Sphere) : bool =
Float.equal (Vector.normSquared (Point.difference sphere.Centre point)) (sphere.Radius * sphere.Radius)
/// Returns the intersections of this ray with this sphere.
/// The nearest intersection is returned first, if there are multiple.
/// Does not return any intersections which are behind us.
/// If the sphere is made of a material which does not re-emit light, you'll
/// get a None for the outgoing ray.
let intersections<'a>
(num : Num<'a>)
(sphere : Sphere<'a>)
(ray : Ray<'a>)
let intersections
(sphere : Sphere)
(ray : Ray)
(incomingColour : Pixel)
: (Point<'a> * Ray<'a> option * Pixel) array
: (Point * Ray option * Pixel) array
=
// The sphere is all points P such that Point.normSquared (P - sphere.Centre) = sphere.Radius^2
// The ray is all ray.Origin + t ray.Vector for any t.
// So the intersection is all P such that
// Point.normSquared (ray.Origin + t ray.Vector - sphere.Centre) = sphere.Radius^2
// Simplified,
// t^2 Point.normSquared ray.Vector
// + 2 t Vector.dot ray.Vector (ray.Origin - sphere.Centre)
// + Point.normSquared (ray.Origin - sphere.Centre) - sphere.Radius^2
// = 0
// That is:
let difference =
Point.difference num ray.Origin sphere.Centre
Point.difference (Ray.origin ray) sphere.Centre
let vector = ray.Vector |> Vector.unitise num |> Option.get
let a = Vector.normSquared num vector
let b = (UnitVector.dot' (Ray.vector ray) difference) * 2.0
let b =
num.Double (Vector.dot num vector difference)
let c = (Vector.normSquared difference) - (sphere.Radius * sphere.Radius)
let c =
num.Subtract (Vector.normSquared num difference) (num.Times sphere.Radius sphere.Radius)
let discriminant =
num.Subtract (num.Times b b) (num.Double (num.Double (num.Times a c)))
let discriminant = (b * b) - (4.0 * c)
let ts =
match num.Compare discriminant num.Zero with
match Float.compare discriminant 0.0 with
| Comparison.Equal ->
[|
num.Negate (num.Divide b (num.Double a))
- (b / 2.0)
|]
| Comparison.Less -> [||]
| Comparison.Greater ->
let intermediate = num.Sqrt discriminant
let denom = num.Double a
let intermediate = sqrt discriminant
[|
num.Divide (num.Add (num.Negate b) intermediate) denom
num.Divide (num.Add (num.Negate b) (num.Negate intermediate)) denom
(intermediate - b) / 2.0 ;
- (b + intermediate) / 2.0
|]
// Don't return anything that's behind us
|> Array.filter (fun i -> num.Compare i num.Zero = Greater)
|> Array.filter (fun i -> Float.compare i 0.0 = Greater)
ts
|> function
| [||] -> [||]
| [|x|] -> [|x|]
| [|x ; y|] ->
match num.Compare x y with
match Float.compare x y with
| Less -> [|x ; y|]
| Equal -> failwith "Nooo"
| Greater -> [|y ; x|]
| _ -> failwith "Impossible"
|> Array.map (fun pos ->
let strikePoint = Ray.walkAlong num ray pos
let strikePoint = Ray.walkAlong ray pos
let outgoing, colour = sphere.Reflection ray incomingColour strikePoint
strikePoint, outgoing, colour
)