Last active
December 1, 2023 20:13
-
-
Save mrange/c19170f3b3892e6a31521600980f8dd0 to your computer and use it in GitHub Desktop.
Raymarcher in F#
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| #nowarn "9" | |
| open System | |
| open System.Diagnostics | |
| open System.Numerics | |
| open SixLabors.ImageSharp | |
| open SixLabors.ImageSharp.Processing | |
| open SixLabors.ImageSharp.PixelFormats | |
| open SixLabors.ImageSharp.Advanced | |
| open System.Threading.Tasks | |
| open System.Runtime.Intrinsics | |
| open System.Runtime.Intrinsics.X86 | |
| open FSharp.NativeInterop | |
| // Turn off Tiered Compilation | |
| // Set-Item env:COMPlus_TieredCompilation -value 0 | |
| module RayMarcher = | |
| type V1 = float32 | |
| let inline v1v v : V1 = float32 v | |
| type Mat = Matrix3x2 | |
| let inline matv m11 m12 m13 m21 m22 m23 : Mat = | |
| Mat (v1v m11, v1v m12, v1v m21, v1v m22, v1v m13, v1v m23) | |
| type V3 = Vector3 | |
| let inline v3v x y z : V3 = V3 (v1v x, v1v y, v1v z) | |
| let inline v3normalize v : V3 = V3.Normalize v | |
| let inline v3cross x y: V3 = V3.Cross (x, y) | |
| type V8 = Vector256<float32> | |
| let inline v8add (a : V8) (b : V8) : V8 = Avx.Add (a, b) | |
| let inline v8sub (a : V8) (b : V8) : V8 = Avx.Subtract (a, b) | |
| let inline v8mul (a : V8) (b : V8) : V8 = Avx.Multiply (a, b) | |
| let inline v8man (a : V8) (b : V8) : V8 = Avx.Max (a, b) | |
| let inline v8min (a : V8) (b : V8) : V8 = Avx.Min (a, b) | |
| let inline v8sqrt (a : V8) : V8 = Avx.Sqrt (a) | |
| let inline v8isqrt(a : V8) : V8 = Avx.ReciprocalSqrt (a) | |
| let inline v8inv (a : V8) : V8 = Avx.Reciprocal (a) | |
| let inline v8mask (a : V8) : int = Avx.MoveMask (a) | |
| let inline v8le (a : V8) (b : V8) : int = v8mask (Avx.Compare (a, b, FloatComparisonMode.LessThanOrEqualOrderedNonSignaling)) | |
| let v8v1 v : V8 = | |
| let sa = NativePtr.stackalloc<float32> 1 | |
| NativePtr.set sa 0 v | |
| Avx.BroadcastScalarToVector256 sa | |
| let _0 = v8v1 0.F | |
| let _1 = v8v1 1.F | |
| let _2 = v8v1 2.F | |
| let _l2 = v8v1 2.F | |
| module FastMath = | |
| let inline padd a b = | |
| let a0, a1 = a | |
| let b0, b1 = b | |
| v8add a0 b0, v8add a1 b1 | |
| let inline psub a b = | |
| let a0, a1 = a | |
| let b0, b1 = b | |
| v8sub a0 b0, v8sub a1 b1 | |
| let inline pmul a b = | |
| let a0, a1 = a | |
| let b0, b1 = b | |
| v8mul a0 b0, v8mul a1 b1 | |
| let inline ple a b = | |
| let a0, a1 = a | |
| let b0, b1 = b | |
| ((v8le a0 b0) <<< 8) ||| (v8le a1 b1) | |
| let inline l2 v = | |
| let x0, y0, z0, x1, y1, z1 = v | |
| let x0 = v8mul x0 x0 | |
| let y0 = v8mul y0 y0 | |
| let z0 = v8mul z0 z0 | |
| let x1 = v8mul x1 x1 | |
| let y1 = v8mul y1 y1 | |
| let z1 = v8mul z1 z1 | |
| v8add (v8add x0 y0) z0, v8add (v8add x1 y1) z1 | |
| let inline l8 v = | |
| let x0, y0, z0, x1, y1, z1 = v | |
| let x0 = v8mul x0 x0 | |
| let y0 = v8mul y0 y0 | |
| let z0 = v8mul z0 z0 | |
| let x0 = v8mul x0 x0 | |
| let y0 = v8mul y0 y0 | |
| let z0 = v8mul z0 z0 | |
| let x0 = v8mul x0 x0 | |
| let y0 = v8mul y0 y0 | |
| let z0 = v8mul z0 z0 | |
| let x1 = v8mul x1 x1 | |
| let y1 = v8mul y1 y1 | |
| let z1 = v8mul z1 z1 | |
| let x1 = v8mul x1 x1 | |
| let y1 = v8mul y1 y1 | |
| let z1 = v8mul z1 z1 | |
| let x1 = v8mul x1 x1 | |
| let y1 = v8mul y1 y1 | |
| let z1 = v8mul z1 z1 | |
| v8add (v8add x0 y0) z0, v8add (v8add x1 y1) z1 | |
| let inline l1 v = | |
| let l20, l21 = l2 v | |
| v8sqrt l20, v8sqrt l21 | |
| let inline scale s v = | |
| let x0, y0, z0, x1, y1, z1 = v | |
| v8mul s x0, v8mul s y0, v8mul s z0, v8mul s x1, v8mul s y1, v8mul s z1 | |
| let inline normal v = | |
| let x0, y0, z0, x1, y1, z1 = v | |
| let l20, l21 = l2 v | |
| let i0 = v8isqrt l20 | |
| let i1 = v8isqrt l21 | |
| v8mul i0 x0, v8mul i0 y0, v8mul i0 z0, v8mul i1 x1, v8mul i1 y1, v8mul i1 z1 | |
| let inline sphere v r = | |
| let l10, l11 = l1 v | |
| v8sub l10 r, v8sub l11 r | |
| let inline roundbox v r = | |
| let l80, l81 = l8 v | |
| let l0 = v8sqrt (v8sqrt (v8sqrt l80)) | |
| let l1 = v8sqrt (v8sqrt (v8sqrt l81)) | |
| v8sub l0 r, v8sub l1 r | |
| let inline translate v t = | |
| let tx, ty, tz = t | |
| let x0, y0, z0, x1, y1, z1 = v | |
| v8add x0 tx, v8add y0 ty, v8add z0 tz, v8add x1 tx, v8add y1 ty, v8add z1 tz | |
| let inline ray o d t = | |
| let ox0, oy0, oz0, ox1, oy1, oz1 = o | |
| let dx0, dy0, dz0, dx1, dy1, dz1 = d | |
| let t0, t1 = t | |
| v8add ox0 (v8mul dx0 t0), v8add oy0 (v8mul dy0 t0), v8add oz0 (v8mul dz0 t0), v8add ox1 (v8mul dx1 t1), v8add oy1 (v8mul dy1 t1), v8add oz1 (v8mul dz1 t1) | |
| let inline df r p = | |
| sphere p r | |
| type [<Struct>] V = | |
| // x0 y0 z0 x1 y1 z1 | |
| | V of V8*V8*V8*V8*V8*V8 | |
| static member inline New x0 y0 z0 x1 y1 z1 = V (x0, y0, z0, x1, y1, z1) | |
| static member Zero = V.New _0 _0 _0 _0 _0 _0 | |
| let inline vto v = | |
| let x0, y0, z0, x1, y1, z1 = v | |
| V (x0, y0, z0, x1, y1, z1) | |
| let inline vfrom (V (x0, y0, z0, x1, y1, z1)) = | |
| x0, y0, z0, x1, y1, z1 | |
| let dfOnlyLocals r p = | |
| let p = vfrom p | |
| let d0, d1 = FastMath.df r p | |
| struct (d0, d1) | |
| let df v = dfOnlyLocals _1 v | |
| let rec rayMarch l2 ro rd = | |
| let ro = vfrom ro | |
| let rd = vfrom rd | |
| let l2 = l2, l2 | |
| let mutable i = 50 | |
| let mutable t0 = _0 | |
| let mutable t1 = _0 | |
| let mutable notIntersected = true | |
| let mutable notOutOfBounds = true | |
| while notIntersected && notOutOfBounds do | |
| let t = t0, t1 | |
| let p = FastMath.ray ro rd t | |
| let struct (d0, d1) = df (vto p) | |
| let d = (d0, d1) | |
| let d2 = FastMath.pmul d d // One way to take abs of d | |
| if (FastMath.ple d2 l2) = 0xFFFF then | |
| notIntersected <- false | |
| else | |
| let t = FastMath.padd t d | |
| t0 <- fst t | |
| t1 <- snd t | |
| let i = i - 1 | |
| notOutOfBounds <- i > 0 | |
| type Scene = | |
| { | |
| ViewTransform : Mat | |
| Origin : V3 | |
| Width : int | |
| Height : int | |
| W : V3 | |
| U : V3 | |
| V : V3 | |
| } | |
| static member New width height origin lookAt up : Scene = | |
| let w = v1v width | |
| let h = v1v height | |
| let ratio = w/h | |
| let xscale = 2.F*ratio/h | |
| let xoffset = -1.F*ratio | |
| let yscale = 2.F/h | |
| let yoffset = -1.F | |
| let m = matv xscale 0 xoffset 0 yscale yoffset | |
| let w = v3normalize (lookAt - origin) | |
| let u = v3normalize (v3cross w up) | |
| let v = v3normalize (v3cross w u) | |
| { | |
| ViewTransform = m | |
| Origin = origin | |
| Width = width | |
| Height = height | |
| W = w | |
| U = u | |
| V = v | |
| } | |
| member this.Pixels (pxs : Rgba32 array) x y = | |
| let w = this.Width | |
| let yo = y*w | |
| let xo = x | |
| for yy = 0 to 3 do | |
| for xx = 0 to 3 do | |
| pxs.[yo + xo + w*yy + xx] <- Rgba32 (10uy, 20uy, 30uy) | |
| type [<Struct>] Measure = | |
| | Measure of int64*int*int*int | |
| override x.ToString () = | |
| let (Measure (ms, cc0, cc1, cc2)) = x | |
| sprintf "{ Measured, %d ms, (%d, %d, %d) cc }" ms cc0 cc1 cc2 | |
| static member (-) (Measure (ams, acc0, acc1, acc2), Measure (bms, bcc0, bcc1, bcc2)) = | |
| Measure (ams - bms, acc0 - bcc0, acc1 - bcc1, acc2 - bcc2) | |
| let measure = | |
| let sw = Stopwatch () | |
| sw.Start () | |
| fun () -> | |
| let cc i = GC.CollectionCount i | |
| Measure (sw.ElapsedMilliseconds, cc 0, cc 1, cc 2) | |
| let sequentialFor b e a = | |
| for i = b to e-1 do | |
| a i | |
| let parallelFor b e a = | |
| Parallel.For (b, e, Action<int> a) |> ignore | |
| open RayMarcher | |
| let render w h = | |
| let s = Scene.New w h (v3v 0 0 -3) (v3v 0 0 0) (v3v 0 1 0) | |
| let w4 = w >>> 2 | |
| let h4 = h >>> 2 | |
| let pxs = Array.zeroCreate<Rgba32> (w*h) | |
| let row y4 = | |
| let y = y4 <<< 2 | |
| for x4 = 0 to w4 - 1 do | |
| let x = x4 <<< 2 | |
| s.Pixels pxs x y | |
| sequentialFor 0 h4 row | |
| pxs | |
| [<EntryPoint>] | |
| let main argv = | |
| let w = 512*2 | |
| let h = 384*2 | |
| if w % 4 <> 0 then failwith "width expected to be divisible by 4" | |
| if h % 4 <> 0 then failwith "height expected to be divisible by 4" | |
| let before = measure () | |
| let pixels = render w h | |
| let after = measure () | |
| printfn "A %dx%d image cost: %s" w h (string (after - before)) | |
| use img = new Image<Rgba32>(w, h) | |
| let ps = img.GetPixelSpan() | |
| for i = 0 to (ps.Length - 1) do | |
| ps.[i] <- pixels.[i] | |
| img.Save "avx.png" | |
| 0 |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| #nowarn "9" | |
| open System | |
| open System.Diagnostics | |
| open System.Runtime.Intrinsics | |
| open System.Runtime.Intrinsics.X86 | |
| open FSharp.NativeInterop | |
| // Turn off Tiered Compilation | |
| // Set-Item env:COMPlus_TieredCompilation -value 0 | |
| type [<Struct>] Measure = | |
| | Measure of int64*int*int*int | |
| override x.ToString () = | |
| let (Measure (ms, cc0, cc1, cc2)) = x | |
| sprintf "{ Measured, %d ms, (%d, %d, %d) cc }" ms cc0 cc1 cc2 | |
| static member (-) (Measure (ams, acc0, acc1, acc2), Measure (bms, bcc0, bcc1, bcc2)) = | |
| Measure (ams - bms, acc0 - bcc0, acc1 - bcc1, acc2 - bcc2) | |
| let measure = | |
| let sw = Stopwatch () | |
| sw.Start () | |
| fun () -> | |
| let cc i = GC.CollectionCount i | |
| Measure (sw.ElapsedMilliseconds, cc 0, cc 1, cc 2) | |
| let time a = | |
| let before = measure () | |
| let r = a () | |
| let after = measure () | |
| r, (after - before) | |
| type VF = Vector256<float32> | |
| let inline ( +++ ) (a : VF) (b : VF) : VF = Avx.Add (a, b) | |
| let inline ( --- ) (a : VF) (b : VF) : VF = Avx.Subtract (a, b) | |
| let inline ( *** ) (a : VF) (b : VF) : VF = Avx.Multiply (a, b) | |
| let vf_ v = | |
| let sa = NativePtr.stackalloc<float32> 1 | |
| NativePtr.set sa 0 v | |
| Avx.BroadcastScalarToVector256 sa | |
| type [<Struct>] V3 = | |
| | V3 of VF*VF*VF | |
| static member Zero : V3 = V3 (VF.Zero, VF.Zero, VF.Zero) | |
| static member inline (+) (V3 (lx, ly, lz), V3 (rx, ry, rz)) : V3 = | |
| V3 (lx +++ rx, ly +++ ry, lz +++ rz) | |
| static member inline (-) (V3 (lx, ly, lz), V3 (rx, ry, rz)) : V3 = | |
| V3 (lx --- rx, ly --- ry, lz --- rz) | |
| static member inline (*) (V3 (lx, ly, lz), V3 (rx, ry, rz)) : V3 = | |
| V3 (lx *** rx, ly *** ry, lz *** rz) | |
| module TupleTest = | |
| type [<Struct>] V = | |
| // x0 y0 z0 x1 y1 z1 | |
| | V of VF*VF*VF*VF*VF*VF | |
| static member Zero : V = V (VF.Zero, VF.Zero, VF.Zero, VF.Zero, VF.Zero, VF.Zero) | |
| let inline l2 xyz = | |
| let x0, y0, z0, x1, y1, z1 = xyz | |
| x0 *** x0 +++ y0 *** y0 +++ z0 *** z0, x1 *** x1 +++ y1 *** y1 +++ z1 *** z1 | |
| let inline l1 xyz = | |
| let l20, l21 = l2 xyz | |
| Avx.Sqrt l20, Avx.Sqrt l21 | |
| let inline sphere xyz r = | |
| let l10, l11 = l1 xyz | |
| // l10 --- r, l11 --- r | |
| l10, l11 | |
| let inline mind (ld0 : VF, rd0 : VF) (ld1 : VF, rd1 : VF) : VF*VF = | |
| Avx.Min (ld0, rd0), Avx.Min (ld1, rd1) | |
| let _r0 = vf_ 1.F | |
| let _r1 = vf_ 2.F | |
| let inline df (V (x0, y0, z0, x1, y1, z1)) = | |
| let r0 = _r0 | |
| let r1 = _r1 | |
| // Two spheres translated | |
| let xyz = x0, y0, z0, x1, y1, z1 | |
| let s0 = sphere xyz r0 | |
| let s1 = sphere xyz r1 | |
| mind s0 s1 | |
| let test () = | |
| let z = V.Zero | |
| let n = 10000000 | |
| let mutable acc0 = VF.Zero | |
| let mutable acc1 = VF.Zero | |
| for x = 1 to n do | |
| let d0, d1 = df z | |
| acc0 <- acc0 +++ d0 | |
| acc1 <- acc1 +++ d1 | |
| acc0, acc1 | |
| [<EntryPoint>] | |
| let main argv = | |
| for x = 0 to 100 do | |
| let _, m = time TupleTest.test | |
| printfn "... it cost: %s" (string m) | |
| 0 |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| open System | |
| open System.Diagnostics | |
| open System.Numerics | |
| open SixLabors.ImageSharp | |
| open SixLabors.ImageSharp.Processing | |
| open SixLabors.ImageSharp.PixelFormats | |
| open SixLabors.ImageSharp.Advanced | |
| open System.Threading.Tasks | |
| // Turn off Tiered Compilation | |
| // Set-Item env:COMPlus_TieredCompilation -value 0 | |
| module Raymarcher = | |
| type V1 = float32 | |
| type V2 = Vector2 | |
| type V3 = Vector3 | |
| type Mat = Matrix3x2 | |
| let inline v1 x = float32 x | |
| let inline clamp_1 (x : V1) (y : V1) (z : V1) = if x < y then y elif x > z then z else x | |
| let inline lerp_1 (x: V1) (y : V1) (z : V1) = z*(y - x) + x | |
| let inline mod_1 (x : V1) (y : V1) = x - y*floor (x/y) | |
| (* | |
| let inline mod_1 (x : V1) (y : V1) = | |
| let m = x % y | |
| let s = Math.Sign x + Math.Sign y | |
| if s = 0 then m + y else m | |
| *) | |
| let inline positive_1 x = if x > 0.F then x else 0.F | |
| let inline v2 x y = V2 (v1 x, v1 y) | |
| let inline abs_2 (x : V2) = V2.Abs x | |
| let inline length_2 (x : V2) = x.Length () | |
| let inline length8_2 (x : V2) = | |
| let x2 = x*x | |
| let x4 = x2*x2 | |
| let x8 = x4*x4 | |
| let d8 = x8.X + x8.Y | |
| Math.Pow (float d8, 1./8.) |> v1 | |
| let inline cmax_2 (x : V2) = max x.X x.Y | |
| let inline cmin_2 (x : V2) = min x.X x.Y | |
| let inline v3 x y z = V3 (v1 x, v1 y, v1 z) | |
| let inline normalize_3 v = V3.Normalize v | |
| let inline dot_3 x y = V3.Dot (x, y) | |
| let inline cross_3 x y = V3.Cross (x, y) | |
| let inline reflect_3 i n = i - 2.F*n*dot_3 n i | |
| let inline clamp_3 (x : V3) (y : V3) (z : V3) = V3.Clamp (x, y, z) | |
| let inline lerp_3 (x : V3) (y : V3) (z : V1) = V3.Lerp (x, y, z) | |
| let inline mod_3 (x : V3) (y : V3) = v3 (mod_1 x.X y.X) (mod_1 x.Y y.Y) (mod_1 x.Z y.Z) | |
| let inline cmax_3 (x : V3) = max x.X (max x.Y x.Z) | |
| let inline cmin_3 (x : V3) = min x.X (min x.Y x.Z) | |
| let inline length_3 (x : V3) = x.Length () | |
| let inline length8_3 (x : V3) = | |
| let x2 = x*x | |
| let x4 = x2*x2 | |
| let x8 = x4*x4 | |
| let d8 = x8.X + x8.Y + x8.Z | |
| Math.Pow (float d8, 1./8.) |> v1 | |
| let inline mat m11 m12 m13 m21 m22 m23 = Mat (v1 m11, v1 m12, v1 m21, v1 m22, v1 m13, v1 m23) | |
| let pi = v1 System.Math.PI | |
| let tau = 2.0F*pi | |
| let inline rgb r g b = v3 (v1 r/255.F) (v1 g/255.F) (v1 b/255.F) | |
| type [<Struct>] Material = | |
| { | |
| DistanceToSurface : float32 | |
| Color : V3 | |
| Reflecting : float32 | |
| } | |
| static member inline New d c r : Material = { DistanceToSurface = d; Color = c; Reflecting = r } | |
| static member Empty : Material = Material.New 0.F (v3 0 0 0) 0.F | |
| type DistanceFunction = V3 -> Material | |
| type SkyFunction = V3 -> V3 | |
| type PostFunction = V1 -> V1 -> V3 -> V3 | |
| type PostFunction_ = OptimizedClosures.FSharpFunc<V1, V1, V3, V3> | |
| type [<Struct>] Intersection = I of float32*V3*int*Material | |
| module Metrics = | |
| let now : unit -> int64 = | |
| let sw = Stopwatch () | |
| sw.Start () | |
| fun () -> sw.ElapsedMilliseconds | |
| #if METRICS | |
| let mutable distanceEstimations = 0L | |
| let mutable rayMarches = 0L | |
| let mutable shadowMarches = 0L | |
| let mutable outOfBounds = 0L | |
| let mutable insideBounds = 0L | |
| #endif | |
| open Metrics | |
| module Details = | |
| let maxBounces = 5 | |
| let maxMarches = 90 | |
| let maxShadowMarches = 36 | |
| let maxDistance = 16.F | |
| let initialStep = 0.01F | |
| let epsilon = 0.0001F | |
| let minShadow = 0.25F | |
| let tolerance = 0.0001F | |
| let epsilonX = v3 epsilon 0 0 | |
| let epsilonY = v3 0 epsilon 0 | |
| let epsilonZ = v3 0 0 epsilon | |
| // F# compiler generates tail. attribute which is great in principle | |
| // but causes severe performance degradation here (about 3x slower) | |
| let inline suppressTailCall (m : Material) : Material = | |
| if m.DistanceToSurface = 0.F then m else m | |
| let inline suppressTailCall_ (v : V3) : V3 = | |
| if v.X = 0.F then v else v | |
| let inline isnanf f = System.Single.IsNaN f | |
| let inline boxDistance_3 p b = | |
| let d = (V3.Abs p) - b | |
| length_3 (V3.Max (d, V3.Zero)) + cmax_3 (V3.Min (d, V3.Zero)) | |
| let computeNormal (df : DistanceFunction) (pos : V3) : V3 = | |
| let inline c off = (df (pos + off)).DistanceToSurface | |
| let vmax = v3 (c epsilonX) (c epsilonY) (c epsilonZ) | |
| let vmin = v3 (c -epsilonX) (c -epsilonY) (c -epsilonZ) | |
| normalize_3 (vmax - vmin) | |
| module Loops = | |
| let rec rayMarch | |
| (df : DistanceFunction) | |
| (ray : V3 ) | |
| (pos : V3 ) | |
| (travelled: V1 ) | |
| (iter : int ) : Intersection = | |
| if travelled < maxDistance then | |
| let m = df pos | |
| let d = m.DistanceToSurface | |
| if iter <= 0 then | |
| #if METRICS | |
| Interlocked.Add (&rayMarches, int64 maxMarches) |> ignore | |
| #endif | |
| I (travelled, pos, maxMarches, m) | |
| elif d < tolerance then | |
| #if METRICS | |
| Interlocked.Add (&rayMarches, int64 (maxMarches - iter)) |> ignore | |
| #endif | |
| I (travelled, pos, maxMarches - iter, m) | |
| else | |
| rayMarch df ray (pos + d*ray) (travelled + d) (iter - 1) | |
| else | |
| #if METRICS | |
| Interlocked.Add (&rayMarches, int64 (maxMarches - iter)) |> ignore | |
| #endif | |
| I (nanf, pos, maxMarches - iter, Material.Empty) | |
| let rec softShadow | |
| (df : DistanceFunction) | |
| (pos : V3 ) | |
| (lightDir : V3 ) | |
| (minStep : V1 ) | |
| (maxt : V1 ) | |
| (k : V1 ) | |
| (shadow : V1 ) | |
| (travelled: V1 ) | |
| (iter : int ) : V1 = | |
| if iter > 0 && travelled < maxt && minShadow < shadow then | |
| let m = df pos | |
| let d = m.DistanceToSurface | |
| let shadow = min shadow (k*d/travelled) | |
| let step = max minStep d | |
| softShadow df (pos + step*lightDir) lightDir minStep maxt k shadow (travelled + step) (iter - 1) | |
| else | |
| #if METRICS | |
| Interlocked.Add (&shadowMarches, int64 (maxShadowMarches - iter)) |> ignore | |
| #endif | |
| clamp_1 shadow minShadow 1.F | |
| let raymarch | |
| (df : DistanceFunction) | |
| (ray : V3 ) | |
| (origin : V3 ) | |
| (initial : V1 ) : Intersection = | |
| Loops.rayMarch df ray (origin + initial*ray) initial maxMarches | |
| let softShadow | |
| (df : DistanceFunction) | |
| (pos : V3 ) | |
| (lightDir : V3 ) | |
| (mint : V1 ) | |
| (maxt : V1 ) | |
| (k : V1 ) : V1 = | |
| Loops.softShadow df (pos + mint*lightDir) lightDir (0.5F*mint) maxt k 1.F mint maxShadowMarches | |
| module Loops2 = | |
| let rec pixel | |
| (df : DistanceFunction) | |
| (skyf : SkyFunction ) | |
| (light : V3 ) | |
| (ray : V3 ) | |
| (origin : V3 ) | |
| (limit : V1 ) | |
| (bounce : int ) : V3 = | |
| let (I (dist, pos, iter, mat)) = raymarch df ray origin initialStep | |
| if isnanf dist then | |
| skyf ray |> suppressTailCall_ | |
| else | |
| let nor = computeNormal df pos | |
| let lv = light - pos | |
| let ll = length_3 lv | |
| let ld = (1.F / ll)*lv | |
| let dif = positive_1 (dot_3 ld nor) | |
| let occ = (v1 (maxMarches - iter)) / v1 maxMarches | |
| let sha = softShadow df pos ld initialStep ll 24.F | |
| let tot = lerp_1 0.1F 1.F (occ*dif*sha) | |
| let col = mat.Color*tot | |
| let fres= pown (1.F - abs (dot_3 nor ray)) 4 | |
| let refl= lerp_1 mat.Reflecting 1.F fres | |
| if bounce > 0 && refl > limit then | |
| let ref = reflect_3 ray nor | |
| let rcol = pixel df skyf light ref pos (limit/refl) (bounce - 1) | |
| (1.F - refl)*col + refl*rcol | |
| else | |
| col | |
| open Details | |
| module Primitives = | |
| let inline plane (c: V3) (refl : V1) (n : V3) (o : V1) : DistanceFunction = | |
| let n = normalize_3 n | |
| fun p -> | |
| let dist = dot_3 n p + o | |
| Material.New dist c refl | |
| let inline sphere (c: V3) (refl : V1) (r : V1) : DistanceFunction = | |
| fun p -> | |
| let dist = length_3 p - r | |
| Material.New dist c refl | |
| let inline softCube (c: V3) (refl : V1) (r : V1) : DistanceFunction = | |
| fun p -> | |
| let dist = length8_3 p - r | |
| Material.New dist c refl | |
| let inline box (c: V3) (refl : V1) (b : V3) : DistanceFunction = | |
| fun p -> | |
| let dist = boxDistance_3 p b | |
| Material.New dist c refl | |
| let inline boundingBox (limit : V1) (c: V3) (refl : V1) (b : V3) (df : DistanceFunction) : DistanceFunction = | |
| fun p -> | |
| let dist = boxDistance_3 p b | |
| if dist > limit then | |
| #if METRICS | |
| Interlocked.Increment &outOfBounds |> ignore | |
| #endif | |
| Material.New dist c refl | |
| else | |
| #if METRICS | |
| Interlocked.Increment &insideBounds |> ignore | |
| #endif | |
| let m2 = df p | |
| if dist > m2.DistanceToSurface then | |
| Material.New dist c refl | |
| else | |
| m2 | |
| module Transformations = | |
| let inline rotateAroundX (a : V1) (df : DistanceFunction): DistanceFunction = | |
| let c = cos a | |
| let s = sin a | |
| fun p -> | |
| let pp = v3 p.X (c*p.Y + s*p.Z) (-s*p.Y + c*p.Z) | |
| df pp |> suppressTailCall | |
| let inline rotateAroundY (a : V1) (df : DistanceFunction): DistanceFunction = | |
| let c = cos a | |
| let s = sin a | |
| fun p -> | |
| let pp = v3 (c*p.X + s*p.Z) p.Y (-s*p.X + c*p.Z) | |
| df pp |> suppressTailCall | |
| let inline rotateAroundZ (a : V1) (df : DistanceFunction): DistanceFunction = | |
| let c = cos a | |
| let s = sin a | |
| fun p -> | |
| let pp = v3 (c*p.X + s*p.Y) (-s*p.X + c*p.Y) p.Z | |
| df pp |> suppressTailCall | |
| let inline scale (s : V1) (df : DistanceFunction): DistanceFunction = | |
| fun p -> | |
| let m = df (p/s) | |
| { m with DistanceToSurface = m.DistanceToSurface*s } | |
| let inline translate (t : V3) (df : DistanceFunction): DistanceFunction = | |
| fun p -> | |
| df (p - t) |> suppressTailCall | |
| let inline unionWith (from : DistanceFunction) (s : DistanceFunction): DistanceFunction = | |
| fun p -> | |
| let m1 = from p | |
| let m2 = s p | |
| if m1.DistanceToSurface < m2.DistanceToSurface then | |
| m1 | |
| else | |
| m2 | |
| let inline intersectWith (from : DistanceFunction) (s : DistanceFunction): DistanceFunction = | |
| fun p -> | |
| let m1 = from p | |
| let m2 = s p | |
| if m1.DistanceToSurface > m2.DistanceToSurface then | |
| m1 | |
| else | |
| m2 | |
| let inline subtractFrom (from : DistanceFunction) (s : DistanceFunction): DistanceFunction = | |
| fun p -> | |
| let m1 = from p | |
| let m2 = s p | |
| let d2 = -m2.DistanceToSurface | |
| if m1.DistanceToSurface > d2 then | |
| m1 | |
| else | |
| { m2 with DistanceToSurface = d2 } | |
| let inline subtract (s: DistanceFunction) (from : DistanceFunction): DistanceFunction = | |
| subtractFrom from s | |
| let inline softUnionWith power (from : DistanceFunction) (s : DistanceFunction): DistanceFunction = | |
| fun p -> | |
| let m1 = from p | |
| let m2 = s p | |
| let sum = exp (-power*m1.DistanceToSurface) + exp (-power*m2.DistanceToSurface) | |
| let d = -log sum / power | |
| if abs (d - m1.DistanceToSurface) < abs (d - m2.DistanceToSurface) then | |
| { m1 with DistanceToSurface = d } | |
| else | |
| { m2 with DistanceToSurface = d } | |
| let inline repeat (rep : V3) (df : DistanceFunction): DistanceFunction = | |
| let rep2 = 0.5F*rep | |
| fun p -> | |
| let pp = mod_3 p rep - rep2 | |
| df pp |> suppressTailCall | |
| let inline distort (d : V3 -> V1) (df : DistanceFunction): DistanceFunction = | |
| fun p -> | |
| let m = df p | |
| let dd = d p | |
| { m with DistanceToSurface = m.DistanceToSurface + dd } | |
| type Scene = | |
| { | |
| DistanceFunction : DistanceFunction | |
| SkyFunction : SkyFunction | |
| PostFunction : PostFunction_ | |
| Light : V3 | |
| ViewTransform : Mat | |
| Origin : V3 | |
| W : V3 | |
| U : V3 | |
| V : V3 | |
| } | |
| static member New df skyf postf width height origin lookAt up light : Scene = | |
| let width = v1 width | |
| let height = v1 height | |
| let ratio = width/height | |
| let xscale = 2.F*ratio/width | |
| let xoffset = -1.F*ratio | |
| let yscale = 2.F/height | |
| let yoffset = -1.F | |
| let m = mat xscale 0 xoffset 0 yscale yoffset | |
| let w = normalize_3 (lookAt - origin) | |
| let u = normalize_3 (cross_3 w up) | |
| let v = normalize_3 (cross_3 w u) | |
| #if METRICS | |
| let df : DistanceFunction = | |
| fun p -> | |
| Interlocked.Increment &distanceEstimations |> ignore | |
| df p |> suppressTailCall | |
| #endif | |
| { | |
| DistanceFunction = df | |
| SkyFunction = skyf | |
| PostFunction = PostFunction_.Adapt postf | |
| Light = light | |
| ViewTransform = m | |
| Origin = origin | |
| W = w | |
| U = u | |
| V = v | |
| } | |
| member t.Pixel (x : int) (y : int) : V3 = | |
| let v = Vector2.Transform (Vector2 (v1 x, v1 y), t.ViewTransform) | |
| let ray = normalize_3 (v.X*t.U + v.Y*t.V + 2.5F*t.W) | |
| let col = Loops2.pixel t.DistanceFunction t.SkyFunction t.Light ray t.Origin 0.01F maxBounces | |
| t.PostFunction.Invoke (v.X, v.Y, col) |> suppressTailCall_ | |
| module MyScene = | |
| open Raymarcher | |
| open Raymarcher.Primitives | |
| open Raymarcher.Transformations | |
| let white = v3 1 1 1 | |
| let gray = v3 0.5 0.5 0.5 | |
| let sunPos = v3 1 3 -3 | |
| let sunDir = normalize_3 sunPos | |
| let sunColor = v3 1 0.75 0.5 | |
| let bigBox = | |
| softCube white 0.F 4.F | |
| |> translate (v3 0 -4 0) | |
| let holeyBox = | |
| sphere gray 0.F 1.2F | |
| |> subtractFrom (softCube gray 0.F 1.F) | |
| |> unionWith (sphere white 0.75F 1.1F) | |
| let post x y col = | |
| let col = V3.Clamp (col, V3.Zero, V3.One) | |
| let col = v3 (sqrt col.X) (sqrt col.Y) (sqrt col.Z) | |
| let col = 0.6F*col + 0.4F*(col*col*(v3 3 3 3 - 2.F*col)) | |
| let sat = dot_3 col (v3 0.33 0.33 0.33) | |
| let col = lerp_3 col (v3 sat sat sat) -0.4F | |
| col | |
| let sky (ray : V3) : V3 = | |
| let tf = dot_3 ray (v3 0 1 0) | |
| let sf = positive_1 (dot_3 ray sunDir) | |
| let s1 = sunColor*pown sf 80 | |
| let s2 = sunColor*pown sf 20 | |
| let t1 = 0.5F*(tf + 1.F)*(v3 0.20 0.3F 0.1) + v3 0 0 0.25F | |
| t1 + s1 + s2 | |
| let scene w h = | |
| let df = | |
| (sphere white 0.75F 0.5F |> translate (v3 -2 -0.5 0)) | |
| |> unionWith holeyBox | |
| |> translate (v3 0 1 0) | |
| |> scale 0.75F | |
| |> unionWith bigBox | |
| Scene.New df sky post w h (v3 -1 1 -3) (v3 0 0 0) (v3 0 1 0) sunPos | |
| type [<Struct>] Measure = | |
| | Measure of int64*int*int*int | |
| override x.ToString () = | |
| let (Measure (ms, cc0, cc1, cc2)) = x | |
| sprintf "{ Measured, %d ms, (%d, %d, %d) cc }" ms cc0 cc1 cc2 | |
| static member (-) (Measure (ams, acc0, acc1, acc2), Measure (bms, bcc0, bcc1, bcc2)) = | |
| Measure (ams - bms, acc0 - bcc0, acc1 - bcc1, acc2 - bcc2) | |
| let measure = | |
| let sw = Stopwatch () | |
| sw.Start () | |
| fun () -> | |
| let cc i = GC.CollectionCount i | |
| Measure (sw.ElapsedMilliseconds, cc 0, cc 1, cc 2) | |
| let sequentialFor b e a = | |
| for i = b to e-1 do | |
| a i | |
| let parallelFor b e a = | |
| Parallel.For (b, e, Action<int> a) |> ignore | |
| let render w h = | |
| let pxs = Array.zeroCreate<Rgba32> (w*h) | |
| let s = MyScene.scene w h | |
| let row y = | |
| let yo = y*w | |
| for x = 0 to w - 1 do | |
| let p = s.Pixel x y | |
| pxs.[yo + x] <- Rgba32 p | |
| parallelFor 0 h row | |
| pxs | |
| [<EntryPoint>] | |
| let main argv = | |
| let w = 512*2 | |
| let h = 384*2 | |
| let before = measure () | |
| let pixels = render w h | |
| let after = measure () | |
| printfn "A %dx%d image cost: %s" w h (string (after - before)) | |
| use img = new Image<Rgba32>(w, h) | |
| let ps = img.GetPixelSpan() | |
| for i = 0 to (ps.Length - 1) do | |
| ps.[i] <- pixels.[i] | |
| img.Save "reference.png" | |
| 0 |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment