minor fixes + basic directions

Author: andreapasquale94

src/Core/axes.jl

@@ -17,10 +17,10 @@ Add a new axes node to `frames`.
     [`add_axes_fixedoffset!`](@ref) and [`add_axes_root!`](@ref).
 """
 function add_axes!(
-    frames::FrameSystem{O,N}, name::Symbol, id::Int,
-    funs::FrameAxesFunctions{O,N}=FrameAxesFunctions{O,N}(),
+    frames::FrameSystem{O,T}, name::Symbol, id::Int,
+    funs::FrameAxesFunctions{O,T}=FrameAxesFunctions{O,T}(),
     parentid=nothing
-) where {O,N<:Number}
+) where {O,T<:Number}
 
     if has_axes(frames, id)
         # Check if a set of axes with the same ID is already registered within 
@@ -64,7 +64,7 @@ function add_axes!(
 
 
     # Create point
-    node = FrameAxesNode{O,N}(name, id, parentid, funs)
+    node = FrameAxesNode{O,T}(name, id, parentid, funs)
 
     # Insert new point in the graph
     add_axes!(frames, node)
@@ -86,10 +86,10 @@ represented by `dcm`, a Direction Cosine Matrix (DCM).
 See also [`add_axes!`](@ref).
 """
 function add_axes_fixedoffset!(
-    frames::FrameSystem{O,N}, name::Symbol, id::Int, parent, dcm::DCM{N}
-) where {O,N}
+    frames::FrameSystem{O,T}, name::Symbol, id::Int, parent, dcm::DCM{T}
+) where {O,T}
 
-    funs = FrameAxesFunctions{O,N}(t -> Rotation{O}(dcm))
+    funs = FrameAxesFunctions{O,T}(t -> Rotation{O}(dcm))
     add_axes!(frames, name, id, funs, axes_id(frames, parent))
 end
 
@@ -107,9 +107,9 @@ despite the rotation depends on time).
 See also [`add_axes!`](@ref).
 """
 function add_axes_projected!(
-    frames::FrameSystem{O,N}, name::Symbol, id::Int, parent, fun::Function
-) where {O,N}
-    funs = FrameAxesFunctions{O,N}(t -> Rotation{O}(fun(t)))
+    frames::FrameSystem{O,T}, name::Symbol, id::Int, parent, fun::Function
+) where {O,T}
+    funs = FrameAxesFunctions{O,T}(t -> Rotation{O}(fun(t)))
     add_axes!(frames, name, id, funs, axes_id(frames, parent))
 end
 
@@ -134,17 +134,17 @@ If `δfun`, `δ²fun` or `δ³fun` are not provided, they are computed via autom
     function does not perform any checks on the output types. 
 """
 function add_axes_rotating!(
-    frames::FrameSystem{O,N}, name::Symbol, id::Int, parent, fun,
+    frames::FrameSystem{O,T}, name::Symbol, id::Int, parent, fun,
     δfun=nothing, δ²fun=nothing, δ³fun=nothing,
-) where {O,N}
+) where {O,T}
 
     for (order, fcn) in enumerate([δfun, δ²fun, δ³fun])
         if (O < order + 1 && !isnothing(fcn))
             @warn "ignoring $fcn, frame system order is less than $(order+1)"
         end
     end
 
-    funs = FrameAxesFunctions{O,N}(
+    funs = FrameAxesFunctions{O,T}(
         t -> Rotation{O}(fun(t)),
 
         # First derivative 

src/Core/directions.jl

@@ -1,32 +1,4 @@
 
-"""
-    add_direction!(frames, name::Symbol, axesid, funs)
-
-Add a new direction node to `frames`.
-
-### Inputs 
-- `frames` -- Target frame system 
-- `name` -- Direction name, must be unique within `frames` 
-- `axesid` -- ID of the axes the direction is expressed in
-- `funs` -- `DirectionFunctions` object storing the functions to compute the direction and, 
-            eventually, its time derivatives. It must match the type and order of `frames`.
-"""
-function add_direction!(
-    frames::FrameSystem{O,N}, name::Symbol, axesid::Int, funs::DirectionFunctions{O,N}
-) where {O,N<:Number}
-    if name in directions(frames)
-        throw(
-            ArgumentError(
-                "A direction with name $name is already registered in the input frame system.",
-            ),
-        )
-    end
-
-    dir = Direction{O,N}(name, length(directions(frames)) + 1, axesid, funs)
-    push!(directions(frames), Pair(name, dir))
-    nothing
-end
-
 """
     add_direction!(frames, name::Symbol, axes, fun, δfun=nothing, δ²fun=nothing, δ³fun=nothing)
 
@@ -47,7 +19,7 @@ If `δfun`, `δ²fun` or `δ³fun` are not provided, they are computed via autom
     function does not perform any checks on the output types. 
 """
 function add_direction!(
-    frames::FrameSystem{O,N}, name::Symbol, ax, fun::Function,
+    frames::FrameSystem{O,N}, name::Symbol, axes, fun::Function,
     δfun=nothing, δ²fun=nothing, δ³fun=nothing
 ) where {O,N}
 
@@ -58,7 +30,7 @@ function add_direction!(
     end
 
     funs = DirectionFunctions{O,N}(
-        t -> Translation{O}fun(t),
+        t -> Translation{O}(fun(t)),
 
         # First derivative
         if isnothing(δfun)
@@ -86,18 +58,22 @@ function add_direction!(
                 if isnothing(δ²fun)
                     (
                         if isnothing(δfun)
-                            t -> Translation{O}(vcat(fun(t), D¹(fun, t), D²(fun, t), D³(fun, t)))
+                            t -> Direction{O}(vcat(fun(t), D¹(fun, t), D²(fun, t), D³(fun, t)))
                         else
-                            t -> Translation{O}(vcat(δfun(t), D²(fun, t), D³(fun, t)))
+                            t -> Direction{O}(vcat(δfun(t), D²(fun, t), D³(fun, t)))
                         end
                     )
                 else
-                    t -> Translation{O}(vcat(δ²fun(t), D³(fun, t)))
+                    t -> Direction{O}(vcat(δ²fun(t), D³(fun, t)))
                 end
             )
         else
             t -> Translation{O}(δ³fun(t))
         end,
     )
-    return add_direction!(frames, name, axes_id(frames, ax), funs)
+
+    axid = axes_id(frames, axes)
+    dir = DirectionDefinition{O,N}(name, length(directions(frames)) + 1, axid, funs)
+    push!(directions(frames), Pair(name, dir))
+    nothing
 end

src/Core/graph.jl

@@ -8,7 +8,7 @@ end
     FrameSystem{O, T, S}
 
 A `FrameSystem` instance manages a collection of user-defined `FramePointNode`, 
-`FrameAxesNode` and `Direction` objects, enabling computation of arbitrary transformations 
+`FrameAxesNode` and `DirectionDefinition` objects, enabling computation of arbitrary transformations 
 between them. It is created by specifying the maximum transformation order `O`, the outputs 
 datatype `T` and an `AbstractTimeScale` instance `S`.
 
@@ -28,7 +28,7 @@ The parameter `S` can be dropped, in case the default (`BarycentricDynamicalTime
 struct FrameSystem{O,T<:Number,S<:AbstractTimeScale}
     points::AliasGraph{PointsGraph{O,T},Dict{Symbol,Int}}
     axes::AliasGraph{AxesGraph{O,T},Dict{Symbol,Int}}
-    dir::Dict{Symbol,Direction{O,T}}
+    dir::Dict{Symbol,DirectionDefinition{O,T}}
 end
 
 function FrameSystem{O,T,S}() where {O,T,S}

src/Core/nodes.jl

@@ -243,7 +243,7 @@ const AxesGraph{O,T} = MappedNodeGraph{FrameAxesNode{O,T},SimpleGraph{Int}}
 const DirectionFunctions{O,T} = FramePointFunctions{O,T}
 
 """
-    Direction{O, T}
+    DirectionDefinition{O, T}
 
 Define a new direction.
 
@@ -252,18 +252,15 @@ Define a new direction.
 - `id` -- direction ID
 - `f` -- `DirectionFunctions` container 
 """
-struct Direction{O,T}
+struct DirectionDefinition{O,T}
     name::Symbol
     id::Int
+    axesid::Int
 
     # internals
     f::DirectionFunctions{O,T}
 end
 
-function Direction{O,T}(name::Symbol, id::Int, dfun::DirectionFunctions{O,T}) where {O,T}
-    return Direction{O,T}(name, id, dfun)
-end
-
-function Base.show(io::IO, d::Direction{O,T}) where {O,T}
-    return println(io, "Direction{$O, $T}(name=$(d.name), id=$(d.id))")
+function Base.show(io::IO, d::DirectionDefinition{O,T}) where {O,T}
+    return println(io, "DirectionDefinition{$O, $T}(name=$(d.name), id=$(d.id), axesid=$(d.axesid))")
 end

src/Core/points.jl

@@ -138,17 +138,17 @@ If `δfun`, `δ²fun` or `δ³fun` are not provided, they are computed with auto
     dimensions. 
 """
 function add_point_dynamical!(
-    frames::FrameSystem{O,N}, name::Symbol, id::Int, parent, ax, fun,
+    frames::FrameSystem{O,T}, name::Symbol, id::Int, parent, ax, fun,
     δfun=nothing, δ²fun=nothing, δ³fun=nothing
-) where {O,N}
+) where {O,T}
 
     for (order, fcn) in enumerate([δfun, δ²fun, δ³fun])
         if (O < order + 1 && !isnothing(fcn))
             @warn "ignoring $fcn, frame system order is less than $(order+1)"
         end
     end
 
-    funs = FramePointFunctions{O,N}(
+    funs = FramePointFunctions{O,T}(
         t -> Translation{O}(fun(t)),
 
         # First derivative

src/Core/transform.jl

@@ -153,7 +153,7 @@ for (order, axfun, _axfun, pfun, _pfun, _pfwd, _pbwd, dfun) in zip(
         an observing point, in a given set of axes, at the desired time `t` expressed in 
         seconds since `J2000`. 
         """
-        function ($pfun)(fr::FrameSystem{O,N}, from, to, ax, t::Number) where {O,N}
+        function ($pfun)(fr::FrameSystem{O,T}, from, to, ax, t::Number) where {O,T}
             if O < $order
                 throw(
                     ErrorException(
@@ -167,7 +167,7 @@ for (order, axfun, _axfun, pfun, _pfun, _pfwd, _pbwd, dfun) in zip(
             toid = point_id(fr, to)
             axid = axes_id(fr, ax)
 
-            fromid == toid && return @SVector zeros(N, 3 * $order)
+            fromid == toid && return @SVector zeros(T, 3 * $order)
 
             # Check to ensure that the two points are registerd
             for id in (fromid, toid)
@@ -266,4 +266,66 @@ for (order, axfun, _axfun, pfun, _pfun, _pfwd, _pbwd, dfun) in zip(
         end
 
     end
+
+    # --------------------------------------------------------------------------------------
+    # Directions transformations
+    # --------------------------------------------------------------------------------------
+
+    @eval begin
+
+        """
+            $($dfun)(frames::FrameSystem, name::Symbol, axes, ep::Epoch) 
+
+        Compute the direction vector `name` of order $(3*$order) at epoch `ep` expressed in 
+        the `axes` frame.
+
+        Requires a frame system of order ≥ $($order).
+        """
+        @inline function ($dfun)(
+            frames::FrameSystem{<:Any,<:Any,S}, name::Symbol, axes, ep::Epoch{S}
+        ) where {S}
+            return $(dfun)(frames, name, axes, j2000s(ep))
+        end
+
+        """
+            $($dfun)(frames::FrameSystem, name::Symbol, axes, t::Number) 
+
+        Compute the direction vector `name` of order $(3*$order) at epoch `t`, where `t` is 
+        expressed in seconds since `J2000`.
+
+        Requires a frame system of order ≥ $($order).
+        """
+        function ($dfun)(frames::FrameSystem{O,N}, name::Symbol, axes, t::Number) where {O,N}
+            if O < $order
+                throw(
+                    ErrorException(
+                        "Insufficient frame system order: " *
+                        "transformation requires at least order $($order).",
+                    ),
+                )
+            end
+
+            if !has_direction(frames, name)
+                throw(
+                    ErrorException(
+                        "No direction with name $(name) registered in the frame system."
+                    )
+                )
+            end
+
+            node = directions(frames)[name]
+            stv = Translation{$order}(node.f[$order](t))
+
+            thisaxid = node.axesid
+            axid = axes_id(frames, axes)
+            if thisaxid != axid
+                stv = ($axfun)(frames, thisaxid, axid, t) * stv
+            end
+
+            D = 3 * $order
+            return @views SVector(stv)[1:D]
+        end
+
+    end
+
 end

src/Definitions/directions.jl

@@ -0,0 +1,90 @@
+
+"""
+    add_direction_position!(frames, name::Symbol, origin, target, axes)
+
+Add a direction based on the position vector from `origin` to `target` in the specified `axes`.
+"""
+function add_direction_position!(
+    frames::FrameSystem{O,T}, name::Symbol, from, to, axes
+) where {O,T}
+
+    fromid, toid = point_id(frames, from), point_id(frames, to)
+    axid = axes_id(frames, axes)
+
+    for id in (fromid, toid)
+        !has_point(frames, id) && throw(
+            ArgumentError(
+                "no points with id $id registered in the given frame system."
+            )
+        )
+    end
+
+    !has_axes(frames, axid) && throw(
+        ArgumentError(
+            "no axes with id $axes registered in the given frame system"
+        )
+    )
+
+    fun = t -> unitvec(vector3(frames, fromid, toid, axid, t))
+    return add_direction!(frames, name, axid, fun)
+end
+
+"""
+    add_direction_velocity!(frames, name::Symbol, origin, target, axes)
+
+Add a direction based on the velocity vector from `origin` to `target` in the specified `axes`.
+"""
+function add_direction_velocity!(
+    frames::FrameSystem{O,T}, name::Symbol, from, to, axes
+) where {O,T}
+
+    fid, tid = point_id(frames, from), point_id(frames, to)
+    axid = axes_id(frames, axes)
+
+    for id in (fid, tid)
+        !has_point(frames, id) && throw(
+            ArgumentError(
+                "no points with id $id registered in the given frame system."
+            )
+        )
+    end
+
+    !has_axes(frames, axid) && throw(
+        ArgumentError(
+            "no axes with name $axes registered in the given frame system"
+        )
+    )
+
+    fun = t -> @views(unitvec(vector6(frames, fid, tid, axid, t)[4:end]))
+    return add_direction!(frames, name, axid, fun)
+end
+
+"""
+    add_direction_orthogonal!(frames, name::Symbol, dir1, dir2)
+
+Add a direction as the cross product between two existing directions (i.e. `dir1` and `dir2`).
+"""
+function add_direction_orthogonal!(
+    frames::FrameSystem{O,T}, name::Symbol, dir1::Symbol, dir2::Symbol, axes
+) where {O,T}
+
+    for d in (dir1, dir2)
+        if !(d in keys(directions(frames)))
+            throw(
+                ArgumentError(
+                    "no direction with name $d registered in the given frame system"
+                )
+            )
+        end
+    end
+
+    axid = axes_id(frames, axes)
+    !has_axes(frames, axid) && throw(
+        ArgumentError(
+            "no axes with name $axes registered in the given frame system"
+        )
+    )
+
+    fun = t -> unitvec(cross3(direction3(frames, dir1, axid, t), direction3(frames, dir2, axid, t)))
+    return add_direction!(frames, name, axid, fun)
+end