Degrees

NEWS.md

@@ -2,6 +2,12 @@
 
 We started keeping track of changes in the `NEWS.md` file after version 0.1.0.
 
+# PlantRayTracer 1.0.1 (2026-06-11)
+
+All angles are now in hexadecimal degrees, both for public API and internal calculations.
+
+This affects `DirectionalSource` and `FixedSource` but also all internal calculations.
+
 # PlantRayTracer 1.0.0 (2026-06-10)
 
 No actual changes. We move to version 1.0.0 because:

Project.toml

@@ -1,7 +1,7 @@
 name = "PlantRayTracer"
 uuid = "78485975-e4aa-407d-b3bb-ded5a4265d05"
 authors = ["Alejandro Morales Sierra <alejandro.moralessierra@wur.nl> and contributors"]
-version = "1.0.0"
+version = "1.0.1"
 
 [deps]
 Atomix = "a9b6321e-bd34-4604-b9c9-b65b8de01458"
@@ -20,7 +20,7 @@ Unrolled = "9602ed7d-8fef-5bc8-8597-8f21381861e8"
 Atomix = "1.1.2"
 CoordinateTransformations = "0.6.3"
 LinearAlgebra = "1.11"
-PlantGeomPrimitives = "1.0.0"
+PlantGeomPrimitives = "1.0.1"
 PrecompileTools = "1.3.3"
 Random = "1.11"
 Rotations = "1.5.1"

src/Materials/Lambertian.jl

@@ -54,8 +54,8 @@ Calculate the wavelength-weighted reflectance/transmittance of each type of inte
 This function actually performs sampling of angles (required for Phong calculations)
 =#
 function calculate_interaction(material::Lambertian, power, ray, intersection, rng)
-    Φ = 2π * rand(rng)
-    θ = acos(sqrt(rand(rng)))
+    Φ = 360.0 * rand(rng)
+    θ = acosd(sqrt(rand(rng)))
     mode, coef = choose_outcome(material, power, rng)
     return (mode = mode, coef = coef, θ = θ, Φ = Φ)
 end

src/Materials/Phong.jl

@@ -58,13 +58,13 @@ Calculate the wavelength-weighted reflectance/transmittance of each type of inte
 This function actually performs sampling of angles (required for Phong calculations)
 =#
 function calculate_interaction(material::Phong, power, ray, intersection, rng)
-    Φ = 2 * π * rand(rng)
+    Φ = 360.0 * rand(rng)
     ρs, θ = calc_specular(material, intersection.axes, rng, ray.dir, Φ)
     mode, coef = choose_outcome(material, ρs, power, rng)
     if mode == ρs
         return (mode = mode, coef = coef, θ = θ, Φ = Φ)
     else
-        θ = acos(sqrt(rand(rng)))
+        θ = acosd(sqrt(rand(rng)))
         return (mode = mode, coef = coef, θ = θ, Φ = Φ)
     end
 end
@@ -115,7 +115,7 @@ function calc_specular(m::Phong{nw}, axes, rng, idir, Φ) where {nw}
         # Sample an angle from the Phong lobe
         θ = sample_phong(m, axes, rng, idir, Φ)
         # Estimator of specular reflectance (Lafortune & Willems, 1994)
-        coef = (m.n + 2) / (m.n + 1) * cos(θ)
+        coef = (m.n + 2) / (m.n + 1) * cosd(θ)
         ρs = m.ρsmax .* coef
         ρs = clamp.(ρs, ρ0, m.ρsmax)
         return (ρs = ρs, θ = θ)
@@ -156,9 +156,9 @@ function sample_phong(m::Phong, axes, rng, idir, Φ)
     e2 = normalize(e1 × axes.n)
     n = normalize(e2 × e1)
     # Calculate zenith angle with respect to reflection direction
-    α = acos(rand(rng)^(1 / (m.n + 1)))
+    α = acosd(rand(rng)^(1 / (m.n + 1)))
     # From polar to Cartessian
     dir = polar_to_cartesian((e1 = e1, e2 = e2, n = n), α, Φ)
     # Calculate θ
-    θ = acos(dir ⋅ axes.n)
+    θ = acosd(dir ⋅ axes.n)
 end

src/Sources/Directional.jl

@@ -2,34 +2,34 @@
 # DirectionalSource
 
 """
-    DirectionalSource(box::AABB; θ, Φ, α = π, alpha_soil = 0, beta_soil = π, radiosity, nrays)
-    DirectionalSource(mesh::Mesh; θ, Φ, α = π, alpha_soil = 0, beta_soil = π, radiosity, nrays)
+    DirectionalSource(box::AABB; θ, Φ, α = 180, alpha_soil = 0, beta_soil = 180, radiosity, nrays)
+    DirectionalSource(mesh::Mesh; θ, Φ, α = 180, alpha_soil = 0, beta_soil = 180, radiosity, nrays)
 
 Create a Directional source (including geometry and angle components) by providing an axis-aligned
 bounding box (`box`) or an `Mesh` object as well as the zenith (`θ`), azimuth (`Φ`), azimuth of
 X axis (`α`), slope inclination (`alpha_soil`) and azimuth of slope normal (`beta_soil`) angles,
 the radiosity of the source projected on the horizontal plane and the number of rays to be
-generated. Directional sources may generate incorrect results in the absence of a grid cloner
-that extends the mesh. This is because the rays are generated from the upper face of the mesh's
-bounding box. See VPL documentation for details on light sources.
+generated. All angles are in degrees. Directional sources may generate incorrect results in the
+absence of a grid cloner that extends the mesh. This is because the rays are generated from the
+upper face of the mesh's bounding box. See VPL documentation for details on light sources.
 
 ## Examples
 ```jldoctest
 julia> using PlantGeomPrimitives;
 
 julia> mesh = Ellipse();
 
-julia> source = DirectionalSource(mesh, θ = 0.0, Φ = 0.0, α = π, radiosity = 1.0, nrays = 1_000);
+julia> source = DirectionalSource(mesh, θ = 0.0, Φ = 0.0, α = 180.0, radiosity = 1.0, nrays = 1_000);
 ```
 """
-function DirectionalSource(box::AABB; θ, Φ, α =  π, alpha_soil = 0.0, beta_soil = π, radiosity, nrays)
+function DirectionalSource(box::AABB; θ, Φ, α = 180.0, alpha_soil = 0.0, beta_soil = 180.0, radiosity, nrays)
     dir_geom = create_directional(box)
     # Radiosity is projected onto horizontal plane as we sample from the top of the bounding box
     # The code below ensures that we get the right irradiance onto the mesh
     power = radiosity .* area(dir_geom)
     out = Source(dir_geom, FixedSource(θ, Φ, α, alpha_soil, beta_soil), power./nrays, nrays)
 end
-function DirectionalSource(mesh::PGP.Mesh; θ, Φ, α = π, alpha_soil = 0.0, beta_soil = π, radiosity, nrays)
+function DirectionalSource(mesh::PGP.Mesh; θ, Φ, α = 180.0, alpha_soil = 0.0, beta_soil = 180.0, radiosity, nrays)
     box = AABB(mesh)
     DirectionalSource(box, θ = θ, Φ = Φ, α = α, alpha_soil = alpha_soil, beta_soil = beta_soil, radiosity = radiosity, nrays = nrays)
 end

src/Sources/SourceAngle.jl

@@ -4,11 +4,12 @@
 
 """
     FixedSource(dir)
-    FixedSource(θ, Φ, α = π, alpha_soil = 0, beta_soil = π)
+    FixedSource(θ, Φ, α = 180, alpha_soil = 0, beta_soil = 180)
 
 Create a fixed irradiance source by given a vector with the direction of the
 rays (dir) or zenith (θ), azimuth (Φ), azimuth of X axis (α), slope inclination
-(alpha_soil) and azimuth of slope normal (beta_soil) angles.
+(alpha_soil) and azimuth of slope normal (beta_soil) angles. All angles are in
+degrees.
 
 ## Examples
 ```jldoctest
@@ -22,7 +23,7 @@ julia> source_dir = FixedSource(PG.Vec(0.0, 0.0, -1.0));
 struct FixedSource{FT} <: SourceAngle
     dir::PGP.Vec{FT}
 end
-function FixedSource(θ, Φ::Real, α = π, alpha_soil = 0.0, beta_soil = π)
+function FixedSource(θ, Φ::Real, α = 180.0, alpha_soil = 0.0, beta_soil = 180.0)
     FixedSource(rotate_coordinates(θ, Φ, α, alpha_soil, beta_soil).z)
 end
 generate_direction(a::FixedSource, rng) = a.dir
@@ -54,7 +55,7 @@ end
 LambertianSource(axes) = LambertianSource(axes...)
 
 function generate_direction(a::LambertianSource{FT}, rng) where {FT}
-    Φ = FT(2) * FT(π) * rand(rng, FT)
-    θ = acos(sqrt(rand(rng, FT)))
+    Φ = FT(360) * rand(rng, FT)
+    θ = acosd(sqrt(rand(rng, FT)))
     polar_to_cartesian((e1 = a.x, e2 = a.y, n = a.z), θ, Φ)
 end

src/utils.jl

@@ -36,12 +36,12 @@ rho(vals...) = SA.SVector(vals)
 ###############################################################################
 
 # Calculate direction unit vector on Cartesian system axes from the
-# angles θ and  Φ
+# angles θ and Φ (in degrees)
 function polar_to_cartesian(axes, θ, Φ)
     e1, e2, n = axes
-    dir1 = @. n * cos(θ)
-    dir2 = @. e1 * cos(Φ) * sin(θ)
-    dir3 = @. e2 * sin(Φ) * sin(θ)
+    dir1 = @. n * cosd(θ)
+    dir2 = @. e1 * cosd(Φ) * sind(θ)
+    dir3 = @. e2 * sind(Φ) * sind(θ)
     dir = @. dir1 + dir2 + dir3
     PGP.Vec(normalize(dir))
 end
@@ -63,22 +63,22 @@ rotatez(x) = CT.LinearMap(Rotations.RotZ(x))
  Given angles θ (solar zenith), Φ (solar azimuth), α (azimuth of X axis),
  alpha_soil (slope inclination) and beta_soil (azimuth of slope normal),
  calculate the flipped coordinate system of the plane (i.e. a coordinate system
- for rays generated from the sky towards the Earth)
- By default we assume:
- X axis points towards South (azimuth of pi radians)
- Y axis points towards East (azimuth of pi/2 radians)
+ for rays generated from the sky towards the Earth).
+ All angles are in degrees. By default we assume:
+ X axis points towards South (azimuth of 180 degrees)
+ Y axis points towards East (azimuth of 90 degrees)
  Z axis points towards zenith (away from ground)
  Φ clockwise looking against Z - East is positive
- θ counterclockwise - 0 = Zenith, Sunrise = pi/2
- α is the geocentric azimuth of the X axis (e.g. α = pi/2 means X points East)
- alpha_soil is the inclination of the slope (0 = horizontal, pi/2 = vertical)
- beta_soil is the geocentric azimuth of the slope normal (e.g. pi = south-facing slope)
+ θ counterclockwise - 0 = Zenith, Sunrise = 90
+ α is the geocentric azimuth of the X axis (e.g. α = 90 means X points East)
+ alpha_soil is the inclination of the slope (0 = horizontal, 90 = vertical)
+ beta_soil is the geocentric azimuth of the slope normal (e.g. 180 = south-facing slope)
  PGP = PlantGeomPrimitives.jl package (just generates unit vectors of a Cartesian system)
 =#
-function rotate_coordinates(θ, Φ, α = π, alpha_soil = 0.0, beta_soil = π)
-    # To deal with Irrational and different precision levels
+function rotate_coordinates(θ, Φ, α = 180.0, alpha_soil = 0.0, beta_soil = 180.0)
+    # To deal with different precision levels
     FT = promote_type(typeof(θ), typeof(Φ), typeof(α), typeof(alpha_soil), typeof(beta_soil))
-    rot = rotatez(α - beta_soil) ∘ rotatey(-alpha_soil) ∘ rotatez(beta_soil - FT(π) - Φ) ∘ rotatey(-θ)
+    rot = rotatez(deg2rad(α - beta_soil)) ∘ rotatey(deg2rad(-alpha_soil)) ∘ rotatez(deg2rad(beta_soil - FT(180) - Φ)) ∘ rotatey(deg2rad(-θ))
     (x = .-rot(PGP.X(FT)), y = .-rot(PGP.Y(FT)), z = .-rot(PGP.Z(FT)))
 end
 

test/elements_raytracer.jl

@@ -237,22 +237,22 @@ let
     @test RT.generate_direction(source, rng) == source.dir
 
     # FixedSource: row orientation (α parameter)
-    # α=π is the default (rows North-South); sun on horizon from North → ray along +X
-    source = RT.FixedSource(π/2, 0.0, π)
+    # α=180 is the default (rows North-South); sun on horizon from North → ray along +X
+    source = RT.FixedSource(90.0, 0.0, 180.0)
     @test source.dir ≈ PGP.Vec(1, 0, 0)
-    # α=π/2 means rows East-West; same sun position shifts ray to -Y
-    source = RT.FixedSource(π/2, 0.0, π/2)
+    # α=90 means rows East-West; same sun position shifts ray to -Y
+    source = RT.FixedSource(90.0, 0.0, 90.0)
     @test source.dir ≈ PGP.Vec(0, -1, 0)
 
     # FixedSource: tilted soil surface (alpha_soil, beta_soil)
     # Flat surface (alpha_soil=0) reproduces the default result
-    source = RT.FixedSource(0.0, 0.0, π, 0.0, π)
+    source = RT.FixedSource(0.0, 0.0, 180.0, 0.0, 180.0)
     @test source.dir ≈ PGP.Vec(0, 0, -1.0)
     # South-facing 45° slope, overhead sun: ray enters at 45° to the slope normal
-    source = RT.FixedSource(0.0, 0.0, π, π/4, π)
+    source = RT.FixedSource(0.0, 0.0, 180.0, 45.0, 180.0)
     @test source.dir ≈ PGP.Vec(1/sqrt(2), 0, -1/sqrt(2))
     # Sun perpendicular to that slope: ray points straight into the slope normal
-    source = RT.FixedSource(π/4, π, π, π/4, π)
+    source = RT.FixedSource(45.0, 180.0, 180.0, 45.0, 180.0)
     @test source.dir ≈ PGP.Vec(0, 0, -1.0)
 
     # Lambertian source
@@ -301,16 +301,16 @@ let
     @test dsource.angle.dir ≈ PGP.Vec(0, 0, -1.0)
 
     # DirectionalSource: row orientation (α)
-    dsource = RT.DirectionalSource(gbox, θ = π/2, Φ = 0.0, α = π/2, radiosity = 1.0, nrays = 1_000)
+    dsource = RT.DirectionalSource(gbox, θ = 90.0, Φ = 0.0, α = 90.0, radiosity = 1.0, nrays = 1_000)
     @test dsource.angle.dir ≈ PGP.Vec(0, -1, 0)
 
     # DirectionalSource: tilted surface (alpha_soil, beta_soil)
-    dsource = RT.DirectionalSource(gbox, θ = 0.0, Φ = 0.0, alpha_soil = π/4, beta_soil = π, radiosity = 1.0, nrays = 1_000)
+    dsource = RT.DirectionalSource(gbox, θ = 0.0, Φ = 0.0, alpha_soil = 45.0, beta_soil = 180.0, radiosity = 1.0, nrays = 1_000)
     @test dsource.angle.dir ≈ PGP.Vec(1/sqrt(2), 0, -1/sqrt(2))
 
     # DirectionalSource from Mesh with all new kwargs
     dmesh = PGP.Rectangle()
-    dsource = RT.DirectionalSource(dmesh, θ = 0.0, Φ = 0.0, α = π, alpha_soil = 0.0, beta_soil = π, radiosity = 1.0, nrays = 1_000)
+    dsource = RT.DirectionalSource(dmesh, θ = 0.0, Φ = 0.0, α = 180.0, alpha_soil = 0.0, beta_soil = 180.0, radiosity = 1.0, nrays = 1_000)
     @test dsource isa RT.Source
     @test dsource.angle.dir ≈ PGP.Vec(0, 0, -1.0)
 

test/raytracer_tiles.jl

@@ -105,10 +105,10 @@ let
     test_results(irradiance, sum(radiosity) * 0.7)
 
     # Source is at an angle
-    θ = π/3
-    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cos(θ), nrays = nrays, θ = θ)
+    θ = 60.0
+    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cosd(θ), nrays = nrays, θ = θ)
     pow, irradiance = get_power(graph, N, nw)
-    test_results(irradiance, sum(radiosity) * 0.7 * cos(θ))
+    test_results(irradiance, sum(radiosity) * 0.7 * cosd(θ))
 
     # Tile is also at an angle (because of cloner only a fraction is exposed!)
     N = 1
@@ -119,11 +119,11 @@ let
     graph = PG.Graph(axiom = axiom)
     scene = PGP.Mesh(graph)
     #render(scene)
-    θ = π/3
-    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cos(θ), nrays = nrays, θ = θ)
+    θ = 60.0
+    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cosd(θ), nrays = nrays, θ = θ)
     pow, irradiance = get_power(graph, N, nw)
     d = 0.5
-    sin_alpha = sin(pi/3)
+    sin_alpha = sind(60)
     h = d*sin_alpha
     dp = d - h*h/d
     test_results(irradiance, sum(radiosity) * dp/d)
@@ -159,10 +159,10 @@ let
     test_results(irradiance, sum(radiosity) * 0.7)
 
     # Source is at an angle
-    θ = π/3
-    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cos(θ), nrays = nrays, θ = θ)
+    θ = 60.0
+    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cosd(θ), nrays = nrays, θ = θ)
     pow, irradiance = get_power(graph, N, nw)
-    test_results(irradiance, sum(radiosity) * 0.7 * cos(θ))
+    test_results(irradiance, sum(radiosity) * 0.7 * cosd(θ))
 
     ##################### One tile with multiple triangles and multiple wavelength ####################
     N = 10
@@ -196,10 +196,10 @@ let
     test_results(irradiance, sum(radiosity) * 0.7)
 
     # Source is at an angle
-    θ = π/3
-    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cos(θ), nrays = nrays, θ = θ)
+    θ = 60.0
+    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cosd(θ), nrays = nrays, θ = θ)
     pow, irradiance = get_power(graph, N, nw)
-    test_results(irradiance, sum(radiosity) * 0.7 * cos(θ))
+    test_results(irradiance, sum(radiosity) * 0.7 * cosd(θ))
 
     ##################### Many tiles with multiple triangles and single wavelength ####################
     N = 10
@@ -233,10 +233,10 @@ let
     test_results(irradiance, sum(radiosity) * 0.7)
 
     # Source is at an angle
-    θ = π/3
-    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cos(θ), nrays = nrays, θ = θ)
+    θ = 60.0
+    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cosd(θ), nrays = nrays, θ = θ)
     pow, irradiance = get_power(graph, N, nw)
-    test_results(irradiance, sum(radiosity) * 0.7 * cos(θ))
+    test_results(irradiance, sum(radiosity) * 0.7 * cosd(θ))
 
 
     ##################### Many tiles with multiple triangles and multiple wavelengths ####################
@@ -271,10 +271,10 @@ let
     test_results(irradiance, sum(radiosity) * 0.7)
 
     # Source is at an angle
-    θ = π/3
-    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cos(θ), nrays = nrays, θ = θ)
+    θ = 60.0
+    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cosd(θ), nrays = nrays, θ = θ)
     pow, irradiance = get_power(graph, N, nw)
-    test_results(irradiance, sum(radiosity) * 0.7 * cos(θ))
+    test_results(irradiance, sum(radiosity) * 0.7 * cosd(θ))
 
 
     ##################### Two graphs - many tiles with multiple triangles and multiple wavelengths ####################
@@ -319,12 +319,12 @@ let
     test_results(irradiance, sum(radiosity) * 0.7)
 
     # Source is at an angle
-    θ = π/3
-    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cos(θ), nrays = nrays, θ = θ)
+    θ = 60.0
+    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cosd(θ), nrays = nrays, θ = θ)
     pow, irradiance = get_power(graph1, N, nw)
-    test_results(irradiance, sum(radiosity) * 0.7 * cos(θ))
+    test_results(irradiance, sum(radiosity) * 0.7 * cosd(θ))
     pow, irradiance = get_power(graph2, N, nw)
-    test_results(irradiance, sum(radiosity) * 0.7 * cos(θ))
+    test_results(irradiance, sum(radiosity) * 0.7 * cosd(θ))
 
 
     ##################### Two graphs - many tiles with multiple triangles and multiple wavelengths - use add! ####################
@@ -376,11 +376,11 @@ let
     test_results(irradiance, sum(radiosity) * 0.7)
 
     # Source is at an angle
-    θ = π/3
-    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cos(θ), nrays = nrays, θ = θ)
+    θ = 60.0
+    ray_trace!(scene, psettings, PRT.BVH, radiosity = radiosity*cosd(θ), nrays = nrays, θ = θ)
     pow, irradiance = get_power(graph1, N, nw)
-    test_results(irradiance, sum(radiosity) * 0.7 * cos(θ))
+    test_results(irradiance, sum(radiosity) * 0.7 * cosd(θ))
     pow = sum(mat.power)
     irradiance = pow / PGP.area(scene_extra)
-    test_results(irradiance, sum(radiosity) * 0.7 * cos(θ))
+    test_results(irradiance, sum(radiosity) * 0.7 * cosd(θ))
 end