I tried reading the data in the Fortran code and it did not make a perceivable difference. But note that the comparisons are dependent on the number of things in the scene. The original julia version runs faster than the original Fortran version if the number of “things” is 10.
The above version that I posted (the one that mimicked the Fortran version by using only one type of thing with a label, instead of different types), and is thus faster, is almost as fast as the nim version with 10 objects:
Nim:
nim c -d:danger -d:lto -d:intpow -d:quake --passC:"-march=nativ
e" -r nim10.nim
Hint: used config file '/etc/nim/nim.cfg' [Conf]
Hint: operation successful (300 lines compiled; 0.043 sec total; 5.129MiB peakmem; Dangerous Release B
uild) [SuccessX]
Hint: /home/leandro/Drive/Work/JuliaPlay/nim10 [Exec]
Completed in 154.819779 ms
(and the times becomes similar if the intpow and quake compiler options are removed.
Julia:
leandro@pitico:~/Drive/Work/JuliaPlay% julia RayTracer_LMSingle10.jl
178.229 ms (40 allocations: 979.20 KiB)
Fortran:
leandro@pitico:~/Drive/Work/JuliaPlay% ./ray10.exe
Render time: 253.07 ms.
The original Julia version from the site runs in 229ms for 10 objects here. Thus, faster than the Fortran version and 50% slower than the optimized nim version.
Stressing that the only thing I did there was to increse the number of objects in the nim and Fortran versions, like this:
diff RayTracer.nim nim10.nim
153,154c153,161
< initSphere(Vector(x: -1.0, y: 0.5, z: 1.5), 0.5, ShinySurface),
< initSphere(Vector(x: 0.0, y: 1.0, z: -0.25), 1.0, ShinySurface)
---
> initPlane(Vector(x: 1.0, y: 1.0, z: 0.0), 0.0, CheckerBoardSurface),
> initPlane(Vector(x: 2.0, y: 1.0, z: 0.0), 0.0, CheckerBoardSurface),
> initPlane(Vector(x: 3.0, y: 1.0, z: 0.0), 0.0, CheckerBoardSurface),
> initPlane(Vector(x: 4.0, y: 1.0, z: 0.0), 0.0, CheckerBoardSurface),
> initSphere(Vector(x: 0.0, y: 0.5, z: 1.5), 0.5, ShinySurface),
> initSphere(Vector(x: -1.0, y: 1.0, z: -0.25), 1.0, ShinySurface),
> initSphere(Vector(x: -2.0, y: 1.0, z: -0.25), 1.0, ShinySurface),
> initSphere(Vector(x: -3.0, y: 1.0, z: -0.25), 1.0, ShinySurface),
> initSphere(Vector(x: -4.0, y: 1.0, z: -0.25), 1.0, ShinySurface)
diff main.f95 main_10.f95
538c538
< allocate(scene%things(3))
---
> allocate(scene%things(10))
543,544c543,551
< scene%things(2) = CreateSphere(TVector( 0.0, 1.0, -0.25), 1.0d+0, SHINY_SURFACE)
< scene%things(3) = CreateSphere(TVector(-1.0, 0.5, 1.5 ), 0.5d+0, SHINY_SURFACE)
---
> scene%things(2) = CreatePlane (TVector( 1.0, 1.0, 0.0 ), 0.0d+0, CHECKERBOARD_SURFACE)
> scene%things(3) = CreatePlane (TVector( 2.0, 1.0, 0.0 ), 0.0d+0, CHECKERBOARD_SURFACE)
> scene%things(4) = CreatePlane (TVector( 3.0, 1.0, 0.0 ), 0.0d+0, CHECKERBOARD_SURFACE)
> scene%things(5) = CreatePlane (TVector( 4.0, 1.0, 0.0 ), 0.0d+0, CHECKERBOARD_SURFACE)
> scene%things(6) = CreateSphere(TVector( 0.0, 1.0, -0.25), 1.0d+0, SHINY_SURFACE)
> scene%things(7) = CreateSphere(TVector(-1.0, 0.5, 1.5 ), 0.5d+0, SHINY_SURFACE)
> scene%things(8) = CreateSphere(TVector(-2.0, 0.5, 1.5 ), 0.5d+0, SHINY_SURFACE)
> scene%things(9) = CreateSphere(TVector(-3.0, 0.5, 1.5 ), 0.5d+0, SHINY_SURFACE)
> scene%things(10) = CreateSphere(TVector(-4.0, 0.5, 1.5 ), 0.5d+0, SHINY_SURFACE)
565c572
The Julia version I am running here is:
Code
using StaticArrays
using LinearAlgebra
@enum Surface begin
Shiny = 1
CheckerBoard = 2
end
struct Vec3{T} <: FieldVector{3,T}
x::T
y::T
z::T
end
dot(a::Vec3, b::Vec3) = a.x * b.x + a.y * b.y + a.z * b.z
struct Color{T} <: FieldVector{3,T}
r::T
g::T
b::T
end
dot(a::Color, b::Color) = a.x * b.x + a.y * b.y + a.z * b.z
struct RGBColor <: FieldVector{4,UInt8}
b::UInt8
g::UInt8
r::UInt8
a::UInt8
end
struct Image
width::Int
height::Int
data::Vector{RGBColor}
end
function Image(w::Int, h::Int)
data = Vector{RGBColor}(undef, w * h);
Image(w, h, data)
end
struct Camera{T}
pos::Vec3{T}
forward::Vec3{T}
right::Vec3{T}
up::Vec3{T}
end
function Camera{T}(pos::Vec3{T}, lookAt::Vec3{T}) where T
_down = Vec3{T}(0.0, -1.0, 0.0);
_forward = lookAt - pos
forward = normalize(_forward)
right = normalize(cross(forward, _down)) * 1.5
up = normalize(cross(forward, right)) * 1.5
Camera{T}(pos, forward, right, up)
end
struct Ray{T}
start::Vec3{T}
dir::Vec3{T}
end
struct Thing{T}
type::Int # 1 plane, 2 sphere
vec::Vec3{T} # normal or center
r::T # offset or radius
surface::Surface
end
struct Intersection{T}
thing::Thing{T}
ray::Ray{T}
dist::T
end
struct Light{T}
pos::Vec3{T}
color::Color{T}
end
struct SurfaceProperties{T}
diffuse::Color{T}
specular::Color{T}
reflect::T
roughness::T
end
const WORD = UInt16
const DWORD = UInt32
const LONG = Int32
struct BITMAPINFOHEADER
biSize::DWORD
biWidth::LONG
biHeight::LONG
biPlanes::WORD
biBitCount::WORD
biCompression::DWORD
biSizeImage::DWORD
biXPelsPerMeter::LONG
biYPelsPerMeter::LONG
biClrUsed::DWORD
biClrImportant::DWORD
end
function BITMAPINFOHEADER(w::LONG, h::LONG)
size = sizeof(BITMAPINFOHEADER)
sizeImage = w * h * 4
BITMAPINFOHEADER(size, w, -h, 1, 32, 0, sizeImage, 0, 0, 0, 0)
end
struct BITMAPFILEHEADER
bfType::WORD
bfSize::DWORD
bfReserved::DWORD
bfOffBits::DWORD
end
function BITMAPFILEHEADER(imageSize::DWORD)
bfType = 0x4D42 # 'B' + ('M' << 8)
offBits = 54
size = convert(DWORD, offBits + imageSize);
BITMAPFILEHEADER(bfType, size, 0, offBits)
end
struct Scene{T}
things::Vector{Thing{T}}
lights::Vector{Light{T}}
camera::Camera{T}
function Scene(T)
things = [ Thing{T}(1,Vec3{T}(0.0, 1.0, 0.0), 0.0, CheckerBoard), # plane
Thing{T}(1,Vec3{T}(1.0, 1.0, 0.0), 0.0, CheckerBoard), # plane
Thing{T}(1,Vec3{T}(2.0, 1.0, 0.0), 0.0, CheckerBoard), # plane
Thing{T}(1,Vec3{T}(3.0, 1.0, 0.0), 0.0, CheckerBoard), # plane
Thing{T}(1,Vec3{T}(4.0, 1.0, 0.0), 0.0, CheckerBoard), # plane
Thing{T}(2,Vec3{T}(0.0, 1.0, -0.25), 1.0, Shiny), # sphere
Thing{T}(2,Vec3{T}(-1.0, 0.5, 1.5), 0.5*0.5, Shiny),
Thing{T}(2,Vec3{T}(-2.0, 0.5, 1.5), 0.5*0.5, Shiny),
Thing{T}(2,Vec3{T}(-3.0, 0.5, 1.5), 0.5*0.5, Shiny),
Thing{T}(2,Vec3{T}(-4.0, 0.5, 1.5), 0.5*0.5, Shiny)
]
lights = [
Light{T}(Vec3{T}(-2.0, 2.5, 0.0), Color{T}(0.49, 0.07, 0.07))
Light{T}(Vec3{T}(1.5, 2.5, 1.5), Color{T}(0.07, 0.07, 0.49))
Light{T}(Vec3{T}(1.5, 2.5, -1.5), Color{T}(0.07, 0.49, 0.071))
Light{T}(Vec3{T}(0.0, 3.5, 0.0), Color{T}(0.21, 0.21, 0.35))
]
camera = Camera{T}(Vec3{T}(3.0, 2.0, 4.0), Vec3{T}(-1.0, 0.5, 0.0));
new{T}(things, lights, camera)
end
end
const maxDepth = 5
White(T) = Color{T}(1,1,1)
Grey(T) = Color{T}(0.5,0.5,0.5)
Black(T) = Color{T}(0,0,0)
Background(T) = Black(T)
Defaultcolor(T) = Black(T)
function GetNormal(t::Thing, pos::Vec3)
t.type == 1 && return t.vec
t.type == 2 && return normalize(pos-t.vec)
end
function ToRgbColor(c::Color)
b = floor(UInt8, clamp(c.b, 0.0, 1.0) * 255.0)
g = floor(UInt8, clamp(c.g, 0.0, 1.0) * 255.0)
r = floor(UInt8, clamp(c.r, 0.0, 1.0) * 255.0)
a = convert(UInt8, 255)
return RGBColor(b, g, r, a)
end
function GetPoint(camera::Camera, x::Int, y::Int, screenWidth::Int, screenHeight::Int)
recenterX = (x - (screenWidth / 2.0)) / 2.0 / screenWidth
recenterY = -(y - (screenHeight / 2.0)) / 2.0 / screenHeight
return normalize(camera.forward + ((camera.right * recenterX) + (camera.up * recenterY)))
end
@inline function GetInteresection(t::Thing{T}, ray::Ray{T}) where T
if t.type == 1
denom = dot(t.vec, ray.dir)
if denom <= 0
dist = (dot(t.vec, ray.start) + t.r) / (-denom)
return Intersection(t, ray, dist)
end
elseif t.type == 2
eo = (t.vec - ray.start)
v = dot(eo, ray.dir)
dist = zero(T)
if (v >= 0)
disc = t.r - (dot(eo, eo) - v * v)
dist = (disc >= 0) ? v - sqrt(disc) : dist
end
if dist > 0
return Intersection(t, ray, dist)
end
end
nothing
end
function GetSurfaceProperties(s::Surface, pos::Vec3{T}) where T
if s == Shiny
return SurfaceProperties(White(T), Grey(T), T(0.7), T(250.0))
else
xz = floor(pos.x) + floor(pos.z)
condition = mod(xz, 2) != 0
color = condition ? White(T) : Black(T)
reflect = condition ? 0.1 : 0.7
SurfaceProperties(color, White(T), T(reflect), T(250.0))
end
end
function Save(image::Image, fileName::String)
w = convert(LONG, image.width)
h = convert(LONG, image.height)
infoHeader = BITMAPINFOHEADER(w, h)
offBits = convert(DWORD, 54)
size = convert(DWORD, offBits + infoHeader.biSizeImage)
fileHeader = BITMAPFILEHEADER(infoHeader.biSizeImage)
f = open(fileName, "w")
write(f, fileHeader.bfType)
write(f, fileHeader.bfSize)
write(f, fileHeader.bfReserved)
write(f, fileHeader.bfOffBits)
write(f, infoHeader.biSize)
write(f, infoHeader.biWidth)
write(f, infoHeader.biHeight)
write(f, infoHeader.biPlanes)
write(f, infoHeader.biBitCount)
write(f, infoHeader.biCompression)
write(f, infoHeader.biSizeImage)
write(f, infoHeader.biXPelsPerMeter)
write(f, infoHeader.biYPelsPerMeter)
write(f, infoHeader.biClrUsed)
write(f, infoHeader.biClrImportant)
for c in image.data
write(f, c.b)
write(f, c.g)
write(f, c.r)
write(f, c.a)
end
close(f)
end
function GetClosestIntersection(scene::Scene{T}, ray::Ray{T}) where T
closest = floatmax(T)
closestInter = nothing
for i in eachindex(scene.things)
inter = GetInteresection(scene.things[i], ray)
if (inter != nothing && inter.dist < closest)
closestInter = inter
closest = inter.dist
end
end
return closestInter;
end
function TraceRay(scene::Scene{T}, ray::Ray, depth::Int) where T
isect = GetClosestIntersection(scene, ray);
isect != nothing ? Shade(scene, isect, depth) : Background(T)
end
function Shade(scene::Scene{T}, isect::Intersection, depth::Int) where T
d = isect.ray.dir;
pos = (d * isect.dist) + isect.ray.start;
normal = GetNormal(isect.thing, pos);
reflectDir = normalize(d .- (2*(normal * dot(normal,d))))
surface = GetSurfaceProperties(isect.thing.surface, pos);
naturalColor = Background(T) + GetNaturalColor(scene, surface, pos, normal, reflectDir);
reflectedColor = (depth >= maxDepth) ? Grey(T) : GetReflectionColor(scene, surface, pos, reflectDir, depth);
return naturalColor + reflectedColor;
end
function GetReflectionColor(scene::Scene, surface::SurfaceProperties, pos::Vec3, reflectDir::Vec3, depth::Int)
ray = Ray(pos, reflectDir)
color = TraceRay(scene, ray, depth + 1)
return color*surface.reflect
end
function GetNaturalColor(scene::Scene{T}, surface::SurfaceProperties, pos::Vec3, normv::Vec3, reflectDir::Vec3) where T
result = Black(T)
for light in scene.lights
ldist = light.pos - pos
livec = normalize(ldist)
ray = Ray(pos, livec)
neatIsect = GetClosestIntersection(scene, ray)
isInShadow = neatIsect != nothing && neatIsect.dist < norm(ldist)
if (!isInShadow)
illum = dot(livec,normv)
specular = dot(livec,reflectDir)
lcolor = (illum > 0) ? light.color*illum : Defaultcolor(T)
scolor = (specular > 0) ? light.color*(specular^surface.roughness) : Defaultcolor(T)
result = result + (lcolor .* surface.diffuse) + (scolor .* surface.specular)
end
end
return result
end
function Render(scene::Scene{T}) where T
camera = scene.camera
image = Image(500, 500)
w = image.width - 1
h = image.height - 1
for y in 0:h, x in 0:w
pt = GetPoint(camera, x, y, w, h)
ray = Ray{T}(scene.camera.pos, pt)
color = TraceRay(scene, ray, 0)
image.data[y * image.width + x + 1] = ToRgbColor(color)
end
return image
end
scene = Scene(Float64)
image = Render(scene)
Save(image, "julia-ray64.bmp")
using BenchmarkTools
@btime begin
scene = Scene(Float64)
image = Render($scene)
end