diff --git a/.editorconfig b/.editorconfig new file mode 100644 index 0000000..575e700 --- /dev/null +++ b/.editorconfig @@ -0,0 +1,13 @@ +root = true + +[*.fs] +fsharp_space_before_uppercase_invocation=true +fsharp_space_before_member=true +fsharp_space_before_colon=true +fsharp_multiline_block_brackets_on_same_column=true +fsharp_newline_between_type_definition_and_members=true +fsharp_keep_if_then_in_same_line=true +fsharp_align_function_signature_to_indentation=true +fsharp_alternative_long_member_definitions=true +fsharp_disable_elmish_syntax=true +fsharp_multi_line_lambda_closing_newline=true diff --git a/RayTracing.App/Program.fs b/RayTracing.App/Program.fs index 4c8aaa8..4ba8122 100644 --- a/RayTracing.App/Program.fs +++ b/RayTracing.App/Program.fs @@ -7,23 +7,28 @@ open Spectre.Console module Program = type ProgressTask with + member this.Increment (prog : float) = this.Increment (prog / 1.0) - let go (ctx : ProgressContext) = - let fs = FileSystem() + let go (sample : SampleImages) (ctx : ProgressContext) = + let fs = FileSystem () + let output = fs.Path.GetTempFileName () |> fun s -> fs.Path.ChangeExtension (s, ".ppm") |> fs.FileInfo.FromFileName let task = ctx.AddTask "[green]Generating image[/]" - let maxProgress, image = SampleImages.gradient task.Increment + let maxProgress, image = SampleImages.get sample task.Increment task.MaxValue <- maxProgress / 1.0 let image = image |> Async.RunSynchronously let outputTask = ctx.AddTask "[green]Writing image[/]" - let maxProgress, writer = ImageOutput.toPpm outputTask.Increment image output + + let maxProgress, writer = + ImageOutput.toPpm outputTask.Increment image output + outputTask.MaxValue <- maxProgress / 1.0 writer |> Async.RunSynchronously @@ -31,18 +36,28 @@ module Program = printfn "%s" output.FullName [] - let main (_ : string []) : int = + let main (argv : string []) : int = + let sample = + argv + |> Array.exactlyOne + |> function + | "spheres" -> SampleImages.Spheres + | "gradient" -> SampleImages.Gradient + | s -> failwithf "Unrecognised arg: %s" s + let prog = - AnsiConsole.Progress() - .Columns( - TaskDescriptionColumn(), - ProgressBarColumn(), - PercentageColumn(), - RemainingTimeColumn(), - SpinnerColumn() + AnsiConsole + .Progress() + .Columns ( + TaskDescriptionColumn (), + ProgressBarColumn (), + PercentageColumn (), + RemainingTimeColumn (), + SpinnerColumn () ) + prog.HideCompleted <- false prog.AutoClear <- false - prog.Start go + prog.Start (go sample) 0 diff --git a/TestRayTracing/PpmOutputExample.txt b/RayTracing.Test/PpmOutputExample.txt similarity index 100% rename from TestRayTracing/PpmOutputExample.txt rename to RayTracing.Test/PpmOutputExample.txt diff --git a/TestRayTracing/TestRayTracing.fsproj b/RayTracing.Test/RayTracing.Test.fsproj similarity index 70% rename from TestRayTracing/TestRayTracing.fsproj rename to RayTracing.Test/RayTracing.Test.fsproj index 0392368..39725bd 100644 --- a/TestRayTracing/TestRayTracing.fsproj +++ b/RayTracing.Test/RayTracing.Test.fsproj @@ -7,10 +7,17 @@ + + + + + + + diff --git a/RayTracing.Test/TestPixel.fs b/RayTracing.Test/TestPixel.fs new file mode 100644 index 0000000..5f87da7 --- /dev/null +++ b/RayTracing.Test/TestPixel.fs @@ -0,0 +1,108 @@ +namespace RayTracing.Test + +open RayTracing +open FsCheck +open FsUnitTyped +open NUnit.Framework +open System + +[] +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 + + Check.QuickThrowOnFailure property + + [] + let ``Average of a few pixels, case 1`` () = + let pixels = + [| + (0uy, 234uy, 0uy) + (0uy, 212uy, 0uy); (0uy, 59uy, 0uy); (0uy, 225uy, 0uy); (0uy, 132uy, 0uy); + (0uy, 69uy, 0uy); (0uy, 207uy, 0uy); (0uy, 212uy, 0uy); (0uy, 30uy, 0uy); + (0uy, 0uy, 0uy); (0uy, 179uy, 0uy); (0uy, 234uy, 0uy); (0uy, 54uy, 0uy); + (0uy, 43uy, 0uy) + |] + |> Array.map (fun (r, g, b) -> { Red = r ; Green = g ; Blue = b }) + + let avg = Pixel.average pixels + + avg.Green |> shouldEqual (pixels |> Seq.map (fun i -> float i.Green) |> Seq.average |> Math.Round |> byte) + avg.Red |> shouldEqual (pixels |> Seq.map (fun i -> float i.Red) |> Seq.average |> Math.Round |> byte) + avg.Blue |> shouldEqual (pixels |> Seq.map (fun i -> float i.Blue) |> Seq.average |> Math.Round |> byte) + + [] + let ``Average of a few pixels, case 2`` () = + let pixels = + [| + (0uy, 0uy, 136uy) + (0uy, 0uy, 90uy); (0uy, 0uy, 109uy); (0uy, 0uy, 204uy); (0uy, 0uy, 209uy); + (0uy, 0uy, 31uy); (0uy, 0uy, 244uy); (0uy, 0uy, 67uy); (0uy, 0uy, 139uy); + (0uy, 0uy, 161uy); (0uy, 0uy, 179uy); (0uy, 0uy, 173uy); (0uy, 0uy, 100uy); + (0uy, 0uy, 109uy); (0uy, 0uy, 122uy); (0uy, 0uy, 27uy); (0uy, 0uy, 249uy); + (0uy, 0uy, 54uy) + |] + |> Array.map (fun (r, g, b) -> { Red = r ; Green = g ; Blue = b }) + + let avg = Pixel.average pixels + + avg.Green |> shouldEqual (pixels |> Seq.map (fun i -> float i.Green) |> Seq.average |> Math.Round |> byte) + avg.Red |> shouldEqual (pixels |> Seq.map (fun i -> float i.Red) |> Seq.average |> Math.Round |> byte) + avg.Blue |> shouldEqual (pixels |> Seq.map (fun i -> float i.Blue) |> Seq.average |> Math.Round |> byte) + + [] + let ``Average of a few pixels, case 3`` () = + let pixels = + [| + (0uy, 0uy, 0uy) + (0uy, 0uy, 123uy) + |] + |> Array.map (fun (r, g, b) -> { Red = r ; Green = g ; Blue = b }) + + let avg = Pixel.average pixels + + avg.Green |> shouldEqual (pixels |> Seq.map (fun i -> float i.Green) |> Seq.average |> Math.Round |> byte) + avg.Red |> shouldEqual (pixels |> Seq.map (fun i -> float i.Red) |> Seq.average |> Math.Round |> byte) + avg.Blue |> shouldEqual (pixels |> Seq.map (fun i -> float i.Blue) |> Seq.average |> Math.Round |> byte) + + [] + let ``Average of a few pixels`` () = + let property (fst : byte * byte * byte) (values : (byte * byte * byte) list) : bool = + let values = fst :: values + let pixels = + values + |> List.map (fun (a, b, c) -> { Pixel.Red = a ; Green = b ; Blue = c }) + + let avg = Pixel.average 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) + && avg.Blue = (pixels |> List.map (fun i -> float i.Blue) |> List.average |> Math.Round |> byte) + + Check.QuickThrowOnFailure property + + [] + let ``Combine pixels with white`` () = + let property (r : byte) (g : byte) (b : byte) : bool = + let original = { Red = r ; Green = g ; Blue = b } + let combined = + original + |> Pixel.combine Colour.White + combined = original + + Check.QuickThrowOnFailure property + + [] + let ``Combine pixels with black`` () = + let property (r : byte) (g : byte) (b : byte) : bool = + let original = { Red = r ; Green = g ; Blue = b } + original + |> Pixel.combine Colour.Black + |> (=) Colour.Black + + Check.QuickThrowOnFailure property diff --git a/RayTracing.Test/TestPlane.fs b/RayTracing.Test/TestPlane.fs new file mode 100644 index 0000000..4e075c2 --- /dev/null +++ b/RayTracing.Test/TestPlane.fs @@ -0,0 +1,21 @@ +namespace RayTracing.Test + +open NUnit.Framework +open FsCheck +open RayTracing + +[] +module TestPlane = + + [] + let ``Orthogonalise does make orthogonal vectors`` () = + let property (p : Plane) : bool = + let orth = Plane.orthonormalise Num.float 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 + + property + |> Prop.forAll (Arb.fromGen TestUtils.planeGen) + |> Check.QuickThrowOnFailure \ No newline at end of file diff --git a/RayTracing.Test/TestPpmOutput.fs b/RayTracing.Test/TestPpmOutput.fs new file mode 100644 index 0000000..d1b7d76 --- /dev/null +++ b/RayTracing.Test/TestPpmOutput.fs @@ -0,0 +1,51 @@ +namespace RayTracing.Test + +open RayTracing +open NUnit.Framework +open FsUnitTyped +open System.IO.Abstractions.TestingHelpers + +[] +module TestRayTracing = + + [] + let ``Wikipedia example of PPM output`` () = + let fs = MockFileSystem () + + let expected = + TestUtils.getEmbeddedResource "PpmOutputExample.txt" + + let image = + [| + [| + { Red = 255uy; Blue = 0uy; Green = 0uy } + { Red = 0uy; Blue = 0uy; Green = 255uy } + { Red = 0uy; Blue = 255uy; Green = 0uy } + |] + [| + { + Red = 255uy + Blue = 0uy + Green = 255uy + } + { + Red = 255uy + Blue = 255uy + Green = 255uy + } + { Red = 0uy; Blue = 0uy; Green = 0uy } + |] + |] + |> Image + + let outputFile = + fs.Path.GetTempFileName () + |> fs.FileInfo.FromFileName + + let _, writer = + ImageOutput.toPpm ignore image outputFile + + writer |> Async.RunSynchronously + + fs.File.ReadAllText outputFile.FullName + |> shouldEqual expected diff --git a/RayTracing.Test/TestRational.fs b/RayTracing.Test/TestRational.fs new file mode 100644 index 0000000..27bb765 --- /dev/null +++ b/RayTracing.Test/TestRational.fs @@ -0,0 +1,34 @@ +namespace RayTracing.Test + +open NUnit.Framework +open FsCheck +open RayTracing + +[] +module TestRational = + + [] + 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 + + [] + 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 diff --git a/RayTracing.Test/TestRay.fs b/RayTracing.Test/TestRay.fs new file mode 100644 index 0000000..e3f8f75 --- /dev/null +++ b/RayTracing.Test/TestRay.fs @@ -0,0 +1,74 @@ +namespace RayTracing.Test + +open RayTracing +open NUnit.Framework +open FsCheck + +[] +module TestRay = + + [] + 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 + ) + : bool + = + let origin1 = [| originX; originY; originZ |] |> Point + + 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 output2 = + Ray.walkAlong num ray2 magnitude + + let actual = + Point.difference num output output2 + + let expected = + Point.difference num origin1 origin2 + + Vector.equal num actual expected + + let gen : Gen = + Arb.generate + |> Gen.map NormalFloat.op_Explicit + + let gen = + Gen.zip (Gen.three (Gen.three gen)) gen + + property Num.float + |> Prop.forAll (Arb.fromGen gen) + |> Check.QuickThrowOnFailure + + [] + let ``walkAlong walks the right distance`` () = + let property (ray : Ray, 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) + + property + |> Prop.forAll (Arb.fromGen (Gen.zip TestUtils.rayGen (Arb.generate |> Gen.map NormalFloat.op_Explicit))) + |> Check.QuickThrowOnFailure + + [] + let ``walkAlong stays on the ray`` () = + let property (ray : Ray, distance : float) = + let walked = Ray.walkAlong Num.float ray distance + Ray.liesOn Num.float walked ray + + property + |> Prop.forAll (Arb.fromGen (Gen.zip TestUtils.rayGen (Arb.generate |> Gen.map NormalFloat.op_Explicit))) + |> Check.QuickThrowOnFailure diff --git a/RayTracing.Test/TestSphere.fs b/RayTracing.Test/TestSphere.fs new file mode 100644 index 0000000..5a2a834 --- /dev/null +++ b/RayTracing.Test/TestSphere.fs @@ -0,0 +1,45 @@ +namespace RayTracing.Test + +open RayTracing +open NUnit.Framework +open FsCheck + +[] +module TestSphere = + + [] + let ``Point at distance r from centre lies on sphere`` () = + 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 gen : Gen * float * Point> = + gen { + let! centre = TestUtils.pointGen + let! radius = Arb.generate |> Gen.map NormalFloat.op_Explicit + let! theta = + Arb.generate + |> Gen.map NormalFloat.op_Explicit + |> Gen.map Radian + let! phi = + Arb.generate + |> 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 + |] + |> Point + |> Point.difference Num.float centre + |> fun (Vector v) -> Point v + return centre, radius, surfacePoint + } + + property + |> Prop.forAll (Arb.fromGen gen) + |> Check.QuickThrowOnFailure \ No newline at end of file diff --git a/RayTracing.Test/TestSphereIntersection.fs b/RayTracing.Test/TestSphereIntersection.fs new file mode 100644 index 0000000..25c59b2 --- /dev/null +++ b/RayTracing.Test/TestSphereIntersection.fs @@ -0,0 +1,62 @@ +namespace RayTracing.Test + +open NUnit.Framework +open FsCheck +open FsUnitTyped +open RayTracing + +[] +module TestSphereIntersection = + + + let sphere : Gen> = + gen { + let! origin = TestUtils.pointGen + let! radius = Arb.generate + return Sphere.make Num.float (SphereStyle.LightSource Colour.White) origin radius.Get + } + + [] + let ``Intersection of sphere and ray does lie on both`` () = + let property (ray : Ray, sphere : Sphere) : bool = + let intersections = Sphere.intersections Num.float 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 + rayOk && sphereOk + ) + && + intersections + |> Array.map (fun (intersection, _, _) -> Vector.normSquared Num.float (Point.difference Num.float ray.Origin intersection)) + |> Seq.pairwise + |> Seq.forall (fun (i, j) -> Num.float.Compare i j = Less) + + property + |> Prop.forAll (Arb.fromGen (Gen.zip TestUtils.rayGen sphere)) + |> Check.QuickThrowOnFailure + + [] + 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 + + let intersections = Sphere.intersections Num.float sphere ray Colour.White + + intersections + |> Array.map (fun (intersection, _, _) -> Vector.normSquared Num.float (Point.difference Num.float ray.Origin intersection)) + |> Seq.pairwise + |> Seq.forall (fun (i, j) -> Num.float.Compare i j = Less) + |> shouldEqual true + + intersections + |> Seq.forall (fun (p, _, _) -> Ray.liesOn Num.float p ray) + |> shouldEqual true + + intersections + |> Seq.forall (fun (p, _, _) -> Sphere.liesOn Num.float p sphere) + |> shouldEqual true diff --git a/RayTracing.Test/TestUtils.fs b/RayTracing.Test/TestUtils.fs new file mode 100644 index 0000000..2148f5a --- /dev/null +++ b/RayTracing.Test/TestUtils.fs @@ -0,0 +1,76 @@ +namespace RayTracing.Test + +open System.Numerics +open RayTracing +open System.IO +open System.Reflection +open FsCheck + +[] +module TestUtils = + + type Dummy = + class + end + + let getEmbeddedResource (filename : string) : string = + let filename = + Assembly + .GetAssembly(typeof) + .GetManifestResourceNames () + |> Seq.filter (fun s -> s.EndsWith filename) + |> Seq.exactlyOne + + use stream = + Assembly + .GetAssembly(typeof) + .GetManifestResourceStream filename + + use reader = new StreamReader (stream) + reader.ReadToEnd().Replace ("\r\n", "\n") + + let rationalGen : Gen = + 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 = + [ + 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 |> Gen.map NormalFloat.op_Explicit + + let pointGen = + Gen.three Arb.generate + |> Gen.map (fun (i, j, k) -> Point [| i.Get ; j.Get ; k.Get |]) + + let vectorGen = + Gen.three Arb.generate + |> Gen.map (fun (i, j, k) -> Vector [| i.Get ; j.Get ; k.Get |]) + + let rayGen : Gen> = + gen { + let! origin = pointGen + let! direction = vectorGen + return { Origin = origin ; Vector = direction } + } + + let planeGen = + gen { + let! origin = pointGen + let! v1 = vectorGen + let! v2 = vectorGen + return Plane.makeSpannedBy { Origin = origin ; Vector = v1 } { Origin = origin ; Vector = v2 } + } diff --git a/RayTracing.sln b/RayTracing.sln index 98ba0d9..7817ba4 100644 --- a/RayTracing.sln +++ b/RayTracing.sln @@ -4,18 +4,14 @@ VisualStudioVersion = 16.6.30114.105 MinimumVisualStudioVersion = 10.0.40219.1 Project("{F2A71F9B-5D33-465A-A702-920D77279786}") = "RayTracing", "RayTracing\RayTracing.fsproj", "{AE2CD8DA-FFA5-49BE-B5C9-1A96A3928325}" EndProject -Project("{F2A71F9B-5D33-465A-A702-920D77279786}") = "TestRayTracing", "TestRayTracing\TestRayTracing.fsproj", "{3C554C00-650A-4DEA-B37F-F0831D21CF37}" +Project("{F2A71F9B-5D33-465A-A702-920D77279786}") = "RayTracing.Test", "RayTracing.Test\RayTracing.Test.fsproj", "{3C554C00-650A-4DEA-B37F-F0831D21CF37}" EndProject Project("{F2A71F9B-5D33-465A-A702-920D77279786}") = "RayTracing.App", "RayTracing.App\RayTracing.App.fsproj", "{94792D42-3D22-4FAE-A15B-B89AAF6DB732}" EndProject Global GlobalSection(SolutionConfigurationPlatforms) = preSolution Debug|Any CPU = Debug|Any CPU - Debug|x64 = Debug|x64 - Debug|x86 = Debug|x86 Release|Any CPU = Release|Any CPU - Release|x64 = Release|x64 - Release|x86 = Release|x86 EndGlobalSection GlobalSection(SolutionProperties) = preSolution HideSolutionNode = FALSE @@ -23,29 +19,15 @@ Global GlobalSection(ProjectConfigurationPlatforms) = postSolution {3C554C00-650A-4DEA-B37F-F0831D21CF37}.Debug|Any CPU.ActiveCfg = Debug|Any CPU {3C554C00-650A-4DEA-B37F-F0831D21CF37}.Debug|Any CPU.Build.0 = Debug|Any CPU - {3C554C00-650A-4DEA-B37F-F0831D21CF37}.Debug|x64.ActiveCfg = Debug|Any CPU - {3C554C00-650A-4DEA-B37F-F0831D21CF37}.Debug|x64.Build.0 = Debug|Any CPU - {3C554C00-650A-4DEA-B37F-F0831D21CF37}.Debug|x86.ActiveCfg = Debug|Any CPU - {3C554C00-650A-4DEA-B37F-F0831D21CF37}.Debug|x86.Build.0 = Debug|Any CPU {3C554C00-650A-4DEA-B37F-F0831D21CF37}.Release|Any CPU.ActiveCfg = Release|Any CPU {3C554C00-650A-4DEA-B37F-F0831D21CF37}.Release|Any CPU.Build.0 = Release|Any CPU - {3C554C00-650A-4DEA-B37F-F0831D21CF37}.Release|x64.ActiveCfg = Release|Any CPU - {3C554C00-650A-4DEA-B37F-F0831D21CF37}.Release|x64.Build.0 = Release|Any CPU - {3C554C00-650A-4DEA-B37F-F0831D21CF37}.Release|x86.ActiveCfg = Release|Any CPU - {3C554C00-650A-4DEA-B37F-F0831D21CF37}.Release|x86.Build.0 = Release|Any CPU {94792D42-3D22-4FAE-A15B-B89AAF6DB732}.Debug|Any CPU.ActiveCfg = Debug|Any CPU {94792D42-3D22-4FAE-A15B-B89AAF6DB732}.Debug|Any CPU.Build.0 = Debug|Any CPU - {94792D42-3D22-4FAE-A15B-B89AAF6DB732}.Debug|x64.ActiveCfg = Debug|Any CPU - {94792D42-3D22-4FAE-A15B-B89AAF6DB732}.Debug|x64.Build.0 = Debug|Any CPU - {94792D42-3D22-4FAE-A15B-B89AAF6DB732}.Debug|x86.ActiveCfg = Debug|Any CPU - {94792D42-3D22-4FAE-A15B-B89AAF6DB732}.Debug|x86.Build.0 = Debug|Any CPU {94792D42-3D22-4FAE-A15B-B89AAF6DB732}.Release|Any CPU.ActiveCfg = Release|Any CPU {94792D42-3D22-4FAE-A15B-B89AAF6DB732}.Release|Any CPU.Build.0 = Release|Any CPU - {94792D42-3D22-4FAE-A15B-B89AAF6DB732}.Release|x64.ActiveCfg = Release|Any CPU - {94792D42-3D22-4FAE-A15B-B89AAF6DB732}.Release|x64.Build.0 = Release|Any CPU - {94792D42-3D22-4FAE-A15B-B89AAF6DB732}.Release|x86.ActiveCfg = Release|Any CPU - {94792D42-3D22-4FAE-A15B-B89AAF6DB732}.Release|x86.Build.0 = Release|Any CPU {AE2CD8DA-FFA5-49BE-B5C9-1A96A3928325}.Debug|Any CPU.ActiveCfg = Debug|Any CPU {AE2CD8DA-FFA5-49BE-B5C9-1A96A3928325}.Debug|Any CPU.Build.0 = Debug|Any CPU + {AE2CD8DA-FFA5-49BE-B5C9-1A96A3928325}.Release|Any CPU.ActiveCfg = Release|Any CPU + {AE2CD8DA-FFA5-49BE-B5C9-1A96A3928325}.Release|Any CPU.Build.0 = Release|Any CPU EndGlobalSection EndGlobal diff --git a/RayTracing/Algebraic.fs b/RayTracing/Algebraic.fs new file mode 100644 index 0000000..9479f0b --- /dev/null +++ b/RayTracing/Algebraic.fs @@ -0,0 +1,32 @@ +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" diff --git a/RayTracing/Camera.fs b/RayTracing/Camera.fs new file mode 100644 index 0000000..ba794d9 --- /dev/null +++ b/RayTracing/Camera.fs @@ -0,0 +1,66 @@ +namespace RayTracing + +type Camera<'a> = + { + Num : Num<'a> + /// How tall is our viewport? + ViewportHeight : 'a + /// How wide is our viewport? + ViewportWidth : 'a + /// In which direction is the camera pointing? + View : Ray<'a> + /// 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> + /// 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> + /// How far away from the camera is the imaginary plane + /// onto which we're collecting the pixels of the result? + FocalLength : 'a + /// How many samples will we take per pixel, for anti-aliasing? + SamplesPerPixel : int + } + +[] +module Camera = + let makeBasic<'a> + (n : Num<'a>) + (focalLength : 'a) + (aspectRatio : 'a) + (origin : Point<'a>) + : Camera<'a> + = + let height = n.Double n.One + let view = + { + Origin = origin + Vector = Vector [| n.Zero ; n.Zero ; n.One |] + } + let xAxis = + { + Origin = origin + Vector = Vector [| n.One ; n.Zero ; n.Zero |] + } + let yAxis = + { + Origin = origin + Vector = Vector [| n.Zero ; n.One ; n.Zero |] + } + + { + Num = n + FocalLength = focalLength + ViewportHeight = height + ViewportWidth = n.Times aspectRatio height + View = view + ViewportXAxis = + Ray.parallelTo (Ray.walkAlong n view focalLength) xAxis + ViewportYAxis = + Ray.parallelTo (Ray.walkAlong n view focalLength) yAxis + SamplesPerPixel = 10 + } \ No newline at end of file diff --git a/RayTracing/Domain.fs b/RayTracing/Domain.fs index 4e88e0b..318f4fa 100644 --- a/RayTracing/Domain.fs +++ b/RayTracing/Domain.fs @@ -1,24 +1,66 @@ namespace RayTracing +open System +open RayTracing + [] type progress -[] -type Pixel = - { - Red : byte - Green : byte - Blue : byte - } - -type Image = - | Image of Pixel [] [] +type Image = Image of Pixel [] [] [] module Image = - let rowCount (Image i) : int = - i.Length + let rowCount (Image i) : int = i.Length - let colCount (Image i) : int = - i.[0].Length + let colCount (Image i) : int = i.[0].Length +[] +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 diff --git a/RayTracing/ImageOutput.fs b/RayTracing/ImageOutput.fs index 62a2c04..ce0e322 100644 --- a/RayTracing/ImageOutput.fs +++ b/RayTracing/ImageOutput.fs @@ -18,8 +18,9 @@ module ImageOutput = (file : IFileInfo) : float * Async = - (float (Image.rowCount image)) * 1.0, + (float (Image.rowCount image) + 1.0) * 1.0, async { + progressIncrement 1.0 use outputStream = file.OpenWrite () use writer = new StreamWriter (outputStream) writer.Write "P3\n" @@ -29,13 +30,14 @@ module ImageOutput = match image with | Image arr -> let writeRow (row : Pixel []) = - for pixel in 0..row.Length - 2 do + for pixel in 0 .. row.Length - 2 do writer.Write (PixelOutput.toPpm row.[pixel]) writer.Write " " + writer.Write (PixelOutput.toPpm row.[row.Length - 1]) progressIncrement 1.0 - for row in 0..arr.Length - 2 do + for row in 0 .. arr.Length - 2 do writeRow arr.[row] writer.Write "\n" diff --git a/RayTracing/InfinitePlane.fs b/RayTracing/InfinitePlane.fs new file mode 100644 index 0000000..352b679 --- /dev/null +++ b/RayTracing/InfinitePlane.fs @@ -0,0 +1,117 @@ +namespace RayTracing + +open System + +type InfinitePlaneStyle<'a> = + /// 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 + /// 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 + +type InfinitePlane<'a> = + { + Normal : Vector<'a> + Point : Point<'a> + /// 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 + } + +[] +module InfinitePlane = + + /// Returns the intersections of this ray with this plane. + /// The nearest intersection is returned first, if there are multiple. + /// 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>) + (incomingColour : Pixel) + : (Point<'a> * Ray<'a> 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 [||] + 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 + | Greater -> + let strikePoint = Ray.walkAlong num { Origin = ray.Origin ; Vector = rayVec } 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 + = + fun incomingRay incomingColour strikePoint -> + match style with + | InfinitePlaneStyle.LightSource colour -> + None, Pixel.combine incomingColour colour + + | 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 newColour = Pixel.combine incomingColour colour + Some outgoing, Pixel.darken num newColour albedo + + | InfinitePlaneStyle.PureReflection (albedo, colour) -> + let plane = + Plane.makeSpannedBy { Origin = strikePoint ; Vector = normal } incomingRay + |> Plane.orthonormalise num + 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) + } + | 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 newColour = Pixel.combine incomingColour colour + let darkened = Pixel.darken num newColour albedo + Some outgoing, darkened + + + let make<'a> (num : Num<'a>) (style : InfinitePlaneStyle<'a>) (pointOnPlane : Point<'a>) (normal : Vector<'a>) : InfinitePlane<'a> = + { + Point = pointOnPlane + Normal = normal + Reflection = reflection num style pointOnPlane normal + } + diff --git a/RayTracing/Num.fs b/RayTracing/Num.fs new file mode 100644 index 0000000..64b4e94 --- /dev/null +++ b/RayTracing/Num.fs @@ -0,0 +1,109 @@ +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 + +[] +module Num = + let float : Num = + 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 = + { + 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 + +[] +module Radian = + let add<'a> (n : Num<'a>) (Radian r1) (Radian r2) = n.Add r1 r2 |> Radian \ No newline at end of file diff --git a/RayTracing/Pixel.fs b/RayTracing/Pixel.fs new file mode 100644 index 0000000..a6d3c08 --- /dev/null +++ b/RayTracing/Pixel.fs @@ -0,0 +1,79 @@ +namespace RayTracing + +open System + +[] +type Pixel = + { + Red : byte + Green : byte + Blue : byte + } + +[] +module Colour = + let Black = + { + Red = 0uy + Green = 0uy + Blue = 0uy + } + let White = + { + Red = 255uy + Green = 255uy + Blue = 255uy + } + let Red = + { + Red = 255uy + Green = 0uy + Blue = 0uy + } + let Green = + { + Red = 0uy + Green = 255uy + Blue = 0uy + } + let Blue = + { + Red = 0uy + Green = 0uy + Blue = 255uy + } + +[] +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 + { + Red = byte (Math.Round (r / count)) + Green = byte (Math.Round (g / count)) + Blue = byte (Math.Round (b / count)) + } + + let combine (p1 : Pixel) (p2 : Pixel) : Pixel = + { + Red = (int p1.Red * int p2.Red) / 255 |> byte + Green = (int p1.Green * int p2.Green) / 255 |> byte + Blue = (int p1.Blue * int p2.Blue) / 255 |> byte + } + + let darken<'a> (num : Num<'a>) (p : Pixel) (albedo : 'a) : 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 + } \ No newline at end of file diff --git a/RayTracing/Plane.fs b/RayTracing/Plane.fs new file mode 100644 index 0000000..8d8d03a --- /dev/null +++ b/RayTracing/Plane.fs @@ -0,0 +1,47 @@ +namespace RayTracing + +/// A plane spanned by two rays from a common origin. +type 'a Plane = + private + { + V1 : 'a Vector + V2 : 'a Vector + Point : 'a Point + } + +type 'a OrthonormalPlane = + { + V1 : 'a Vector + V2 : 'a Vector + Point : 'a Point + } + + +[] +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 vec2 = + Vector.difference num plane.V2 (Vector.scale num coeff v1) + |> Vector.unitise num + match vec2 with + | None -> None + | Some v2 -> + { + V1 = v1 + V2 = v2 + Point = plane.Point + } + |> Some + + let makeSpannedBy<'a> (r1 : 'a Ray) (r2 : 'a Ray) : 'a Plane = + { + V1 = r1.Vector + V2 = r2.Vector + Point = r1.Origin + } + + let basis<'a> (plane : 'a OrthonormalPlane) : 'a Ray * 'a Ray = + { Origin = plane.Point ; Vector = plane.V1 }, { Origin = plane.Point ; Vector = plane.V2 } diff --git a/RayTracing/Point.fs b/RayTracing/Point.fs new file mode 100644 index 0000000..97e768d --- /dev/null +++ b/RayTracing/Point.fs @@ -0,0 +1,44 @@ +namespace RayTracing + +/// An n-dimensional point. +[] +type Point<'a> = Point of 'a array + +[] +type Vector<'a> = Vector of 'a array + +[] +module Point = + let difference<'a> (num : Num<'a>) (p1 : Point<'a>) (p2 : Point<'a>) : Vector<'a> = + match p1, p2 with + | Point p1, Point p2 -> + Array.zip p1 p2 + |> Array.map (fun (a, b) -> num.Subtract a b) + |> Vector + + let sum<'a> (num : Num<'a>) (p1 : Point<'a>) (p2 : Point<'a>) : Point<'a> = + 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 + diff --git a/RayTracing/Rational.fs b/RayTracing/Rational.fs new file mode 100644 index 0000000..209d713 --- /dev/null +++ b/RayTracing/Rational.fs @@ -0,0 +1,124 @@ +namespace RayTracing + +open System.Numerics + +[] +[] +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 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) + .CompareTo other + | _ -> failwith "how?!" + +[] +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) \ No newline at end of file diff --git a/RayTracing/Ray.fs b/RayTracing/Ray.fs new file mode 100644 index 0000000..4035c21 --- /dev/null +++ b/RayTracing/Ray.fs @@ -0,0 +1,59 @@ +namespace RayTracing + +type Ray<'a> = + { + Origin : Point<'a> + Vector : Vector<'a> + } + +[] +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 + + 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> = + { + Origin = p1 + Vector = Point.difference num p2 p1 + } + + /// 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 parallelTo<'a> (p1 : Point<'a>) (ray : Ray<'a>) : Ray<'a> = + { + Vector = ray.Vector + Origin = p1 + } + + let liesOn<'a> (num : 'a Num) (point : Point<'a>) (ray : Ray<'a>) : bool = + match point, ray.Origin, ray.Vector with + | Point x, Point y, Vector ray -> + let rec go (state : 'a 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) + | Some prevT -> + let t = num.Divide (num.Subtract d x) r + if num.Equal prevT t then go (Some prevT) (i + 1) else None + + go None 0 + |> Option.isSome diff --git a/RayTracing/RayTracing.fsproj b/RayTracing/RayTracing.fsproj index d2114e7..33cbd48 100644 --- a/RayTracing/RayTracing.fsproj +++ b/RayTracing/RayTracing.fsproj @@ -5,13 +5,25 @@ + + + + + + + + + + + + diff --git a/RayTracing/SampleImages.fs b/RayTracing/SampleImages.fs index 0ea4e88..f93fa7e 100644 --- a/RayTracing/SampleImages.fs +++ b/RayTracing/SampleImages.fs @@ -1,5 +1,11 @@ namespace RayTracing +open System + +type SampleImages = + | Gradient + | Spheres + [] module SampleImages = @@ -13,10 +19,47 @@ module SampleImages = 256.0, async { - return Array.init 256 (fun height -> - let output = Array.init 256 (pixelAt height) - progressIncrement 1.0 - output - ) - |> Image - } \ No newline at end of file + return + Array.init + 256 + (fun height -> + let output = Array.init 256 (pixelAt height) + progressIncrement 1.0 + output + ) + |> Image + } + + let spheres (progressIncrement : float -> unit) : float * Image Async = + 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 + { + 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) + + // 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 |])) + + // 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 |])) + + // 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 |])) + + // 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 |])) + |] + } + |> Scene.render progressIncrement (aspectRatio * (float pixels) |> int) pixels camera + + let get (s : SampleImages) : (float -> unit) -> float * Image Async = + match s with + | Gradient -> gradient + | Spheres -> spheres + diff --git a/RayTracing/Scene.fs b/RayTracing/Scene.fs new file mode 100644 index 0000000..fdf727b --- /dev/null +++ b/RayTracing/Scene.fs @@ -0,0 +1,116 @@ +namespace RayTracing + +open System + +type Hittable<'a> = + | Sphere of Sphere<'a> + | InfinitePlane of InfinitePlane<'a> + +[] +module Hittable = + let hits<'a> + (num : Num<'a>) + (ray : Ray<'a>) + (incomingColour : Pixel) + (h : Hittable<'a>) + : (Point<'a> * Ray<'a> option * Pixel) option + = + match h with + | Sphere s -> + Sphere.intersections num s ray incomingColour + |> Array.tryHead + | InfinitePlane plane -> + InfinitePlane.intersections num plane ray incomingColour + |> Array.tryHead + +type Scene<'a> = + { + Objects : Hittable<'a> array + } + +[] +module Scene = + + let hitObject<'a> + (num : Num<'a>) + (s : Scene<'a>) + (ray : Ray<'a>) + (colour : Pixel) + : (Point<'a> * Ray<'a> 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)) + + let internal traceRay<'a> + (maxCount : int) + (num : Num<'a>) + (scene : Scene<'a>) + (ray : Ray<'a>) + (colour : Pixel) + : Pixel + = + let rec go (bounces : int) (ray : Ray<'a>) (colour : Pixel) : Pixel = + if bounces > maxCount then Colour.Black else + + let thingsWeHit = hitObject num scene ray colour + match thingsWeHit with + | [||] -> + // Ray goes off into the distance and is never heard from again + Colour.Black + | arr -> + let _strikePoint, outgoingRay, colour = arr.[0] + match outgoingRay with + | None -> + colour + | Some outgoingRay -> + go (bounces + 1) { outgoingRay with Vector = Vector.unitise num outgoingRay.Vector |> Option.get } colour + + go 0 ray colour + + let render<'a> + (progressIncrement : float -> unit) + (maxWidthCoord : int) + (maxHeightCoord : int) + (camera : Camera<'a>) + (s : Scene<'a>) + : float * 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 + let colsIter = 2 * maxWidthCoord + 1 + 1.0 * float (rowsIter * colsIter), async { + return + Array.init rowsIter (fun row -> + let row = row - maxHeightCoord + Array.Parallel.init colsIter (fun col -> + let col = col - maxWidthCoord + // Where does this pixel correspond to, on the imaginary canvas? + // For the early prototype, we'll just take the upper right quadrant + // from the camera. + let ret = + 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 + 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 + + let result = traceRay 50 num s ray Colour.White + result + ) + |> Pixel.average + progressIncrement 1.0 + ret + ) + ) + |> Array.rev + |> Image + } \ No newline at end of file diff --git a/RayTracing/Sphere.fs b/RayTracing/Sphere.fs new file mode 100644 index 0000000..a9d733c --- /dev/null +++ b/RayTracing/Sphere.fs @@ -0,0 +1,187 @@ +namespace RayTracing + +open System + +type Sphere<'a> = + { + Centre : Point<'a> + Radius : 'a + /// 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 + } + +type SphereStyle<'a> = + /// 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 + /// 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 + +[] +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 reflection<'a> + (num : Num<'a>) + (style : SphereStyle<'a>) + (centre : Point<'a>) + (radius : 'a) + : Ray<'a> -> Pixel -> Point<'a> -> Ray<'a> option * Pixel + = + let normal = normal num centre + fun incomingRay incomingColour strikePoint -> + let normal = normal strikePoint + + match style with + | SphereStyle.LightSource colour -> + None, Pixel.combine incomingColour colour + | SphereStyle.LightSourceCap colour -> + let circleCentreZCoord = + match centre with + | Point v -> Array.head v + let zCoordLowerBound = num.Add circleCentreZCoord (num.Subtract radius (num.DivideInteger radius 4)) + let strikeZCoord = + match strikePoint with + | Point v -> Array.head v + let colour = + match num.Compare strikeZCoord zCoordLowerBound with + | Greater -> + Pixel.combine colour incomingColour + | _ -> + Colour.Black + None, colour + + | 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 newColour = Pixel.combine incomingColour colour + Some outgoing, Pixel.darken num newColour albedo + + | SphereStyle.PureReflection (albedo, colour) -> + let plane = + Plane.makeSpannedBy normal incomingRay + |> Plane.orthonormalise num + 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) + } + | 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 newColour = Pixel.combine incomingColour colour + let darkened = Pixel.darken num newColour albedo + Some outgoing, darkened + + let make<'a> (num : Num<'a>) (style : SphereStyle<'a>) (centre : Point<'a>) (radius : 'a) : Sphere<'a> = + { + Centre = centre + Radius = radius + Reflection = reflection num 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) + + /// 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>) + (incomingColour : Pixel) + : (Point<'a> * Ray<'a> 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 + + let vector = ray.Vector |> Vector.unitise num |> Option.get + let a = Vector.normSquared num vector + + let b = + num.Double (Vector.dot num vector difference) + + 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 ts = + match num.Compare discriminant num.Zero with + | Comparison.Equal -> + [| + num.Negate (num.Divide b (num.Double a)) + |] + | Comparison.Less -> [||] + | Comparison.Greater -> + let intermediate = num.Sqrt discriminant + let denom = num.Double a + [| + num.Divide (num.Add (num.Negate b) intermediate) denom + num.Divide (num.Add (num.Negate b) (num.Negate intermediate)) denom + |] + // Don't return anything that's behind us + |> Array.filter (fun i -> num.Compare i num.Zero = Greater) + ts + |> function + | [||] -> [||] + | [|x|] -> [|x|] + | [|x ; y|] -> + match num.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 outgoing, colour = sphere.Reflection ray incomingColour strikePoint + strikePoint, outgoing, colour + ) diff --git a/TestRayTracing/TestPpmOutput.fs b/TestRayTracing/TestPpmOutput.fs deleted file mode 100644 index 3b496f8..0000000 --- a/TestRayTracing/TestPpmOutput.fs +++ /dev/null @@ -1,36 +0,0 @@ -namespace TestRayTracing - -open RayTracing -open NUnit.Framework -open FsUnitTyped -open System.IO.Abstractions.TestingHelpers - -[] -module TestRayTracing = - - [] - let ``Wikipedia example of PPM output`` () = - let fs = MockFileSystem () - let expected = TestUtils.getEmbeddedResource "PpmOutputExample.txt" - - let image = - [| - [| - { Red = 255uy ; Blue = 0uy ; Green = 0uy } - { Red = 0uy ; Blue = 0uy ; Green = 255uy } - { Red = 0uy ; Blue = 255uy ; Green = 0uy } - |] - [| - { Red = 255uy ; Blue = 0uy ; Green = 255uy } - { Red = 255uy ; Blue = 255uy ; Green = 255uy } - { Red = 0uy ; Blue = 0uy ; Green = 0uy } - |] - |] - |> Image - - let outputFile = fs.Path.GetTempFileName () |> fs.FileInfo.FromFileName - let _, writer = ImageOutput.toPpm ignore image outputFile - writer |> Async.RunSynchronously - - fs.File.ReadAllText outputFile.FullName - |> shouldEqual expected \ No newline at end of file diff --git a/TestRayTracing/TestUtils.fs b/TestRayTracing/TestUtils.fs deleted file mode 100644 index d14747b..0000000 --- a/TestRayTracing/TestUtils.fs +++ /dev/null @@ -1,20 +0,0 @@ -namespace TestRayTracing - -open System.IO -open System.Reflection - -[] -module TestUtils = - - type Dummy = class end - - let getEmbeddedResource (filename : string) : string = - let filename = - Assembly.GetAssembly(typeof).GetManifestResourceNames () - |> Seq.filter (fun s -> s.EndsWith filename) - |> Seq.exactlyOne - use stream = - Assembly.GetAssembly(typeof).GetManifestResourceStream filename - - use reader = new StreamReader (stream) - reader.ReadToEnd().Replace("\r\n", "\n") \ No newline at end of file