Skip to content

SatcomCoordinates.jl

Documentation for SatcomCoordinates.jl.

SatcomCoordinates.ForwardOrInverse Type
julia
const ForwardOrInverse{F <: AbstractCRSTransform} = Union{F, InverseTransform{<:Any, <:F}}

Union representing either a forward or reverse transform of F

source
SatcomCoordinates.WithNumbertype Type
julia
const WithNumbertype{T} = Union{AbstractSatcomCoordinate{T}, AbstractCRSTransform{T}}

Union representing the types defined and exported by this package, which always have a numbertype as first parameter.

source
SatcomCoordinates.AER Type
julia
struct AER{T} <: AbstractTopocentricPosition{T}

Represents a position in the Azimuth-Elevation-Range (AER) coordinate system, which is a local coordinate system centered at a point on or above the surface of an Ellipsoid. The Elevation and Azimuth angles are always defined w.r.t. the ENU CRS with the same origin. More specifically:

Properties

  • az::Deg{T}: Azimuth angle, defined angle in the XY (North-East) plane from the +Y (North) direction to the object, positive towards +X (East) direction. It is constrained to be in the range [-180°, 180°]

  • el::Deg{T}: Elevation angle, defined as the angle between the XY plane and the point being described by the AER coordinates, positive towards +Z (Up) direction. It is constrained to be in the range [-90°, 90°]

  • r::Met{T}: Range in meters between the origin of the ENU CRS and the point being described by the AER coordinates.

Basic Constructors

AER(az::ValidAngle, el::ValidAngle, r::ValidDistance)

ValidAngle is either a Real Number, or a subtype of Unitful.Angle. ValidDistance is either a Real Number, or a subtype of Unitful.Length.

If of the provided argument is NaN, the returned AER object will contain NaN for all fields.

Fallback constructors

All subtypes of P <: AbstractSatcomCoordinate can also be constructed using a Tuple or SVector as input, which will be splatted into the standard constructor for P

See also: ENU, NED.

source
SatcomCoordinates.AbstractAffineCRSTransform Type
julia
AbstractAffineCRSTransform{T}

Abstract type representing an affine transform between two CRSs with numbertype T.

source
SatcomCoordinates.AbstractCRSTransform Type
julia
AbstractCRSTransform{T}

Abstract type representing a coordinate transform between two CRSs with numbertype T.

source
SatcomCoordinates.AbstractFieldValue Type
julia
AbstractFieldValue{U, CRS, T}

Abstract type representing the value of a physical field expressed in a specific coordinate reference system CRS, and whose components have an associated unit U and a numbertype T.

Concrete implementations of this subtype are expected to have a single inner field svector which is a SVector{N, T} (N being the number of dimensions of the referenced CRS) and define the BasicTypes.constructor_without_checks method as per other coordinates in order to exploit the convenience functions defined on this abstract type.

An example concrete type representing velocity in a 3D CRS can be implemented as follows (assuming to have Quantity, @u_str and dimension imported from Unitful, and SVector from StaticArrays):

julia
const U = typeof(u"m/s")
const D = dimension(u"m/s")

const V{T} = Quantity{T, D, U}

struct VelocityFieldValue{CRS <: AbstractPosition{<:Any, 3}, T} <: AbstractFieldValue{U, CRS, T}
    svector::SVector{3, T}

    BasicTypes.constructor_without_checks(::Type{VelocityFieldValue{CRS, T}}, sv::SVector{3, T}) where {CRS, T} = new{CRS, T}(sv)
end

fv = VelocityFieldValue{ECEF}(1, 2, 3)

@test fv.x == 1u"m/s"
@test fv.y == 2u"m/s"
@test fv.z == 3u"m/s"

@test raw_svector(fv) == SVector{3, Float64}(1, 2, 3)
@test raw_properties(fv) == (x=1, y=2, z=3)

Concrete subtypes

SatcomCoordinates.jl currently does not implement concrete subtypes of AbstractFieldValue, but only defines the following methods that work if concrete subtypes are impelemented as in the example above:

  • raw_svector(::AbstractFieldValue): Returns the svector field

  • raw_properties(::AbstractFieldValue): Returns a NamedTuple with the properties of the svector field

  • properties_names(::Type{<:AbstractFieldValue{U, CRS}}): Returns the property names of the CRS parametric type

  • Default constructor taking N numbers or a Tuple/Svector, with the values interpreted as quantities of unit U if not provided directly as quantitites

  • ConstructionBase.getproperties(::AbstractFieldValue{U, CRS}): Returns a NamedTuple with the properties of the field assuming the properties names of CRS and returning values with unit U.

  • Base.getproperty: extraction of properties directly from the NamedTuple returned by ConstructionBase.getproperties

source
SatcomCoordinates.AbstractGeocentricPosition Type
julia
AbstractGeocentricPosition{T} <: AbstractPosition{T, 3}

Abstract type representing a position in a geocentric coordinate system.

source
SatcomCoordinates.AbstractLocalPosition Type
julia
AbstractLocalPosition{T} <: AbstractPosition{T, 3}

Abstract type representing a position in a local coordinate system.

source
SatcomCoordinates.AbstractPointing Type
julia
AbstractPointing{T} <: AbstractSatcomCoordinate{T, 3}

Abstract type representing a pointing direction in 3 dimensions which is backed by fields with shared numbertype T.

source
SatcomCoordinates.AbstractPointingOffset Type
julia
AbstractPointingOffset{T} <: AbstractSatcomCoordinate{T, 2}

Abstract type representing a pointing offset between two pointing directions.

Currently only has two concrete subtypes: UVOffset and ThetaPhiOffset.

source
SatcomCoordinates.AbstractSatcomCoordinate Type
julia
AbstractSatcomCoordinate{T, N}

General abstract type identifying a coordinate with N dimensions and an underlying number type T. The number type T is not necessarily the type of the fields in the type instance, but their underlying real type (this is only different for fields whose types are Unitful quantities, where numbertype(::Quantity{T}) where T = T).

Note

The term coordinate is used here in a loose sense, identifying both position in space as well as pointing directions

source
SatcomCoordinates.AbstractTopocentricPosition Type
julia
AbstractTopocentricPosition{T} <: AbstractPosition{T, 3}

Abstract type representing a position in a topocentric coordinate system.

source
SatcomCoordinates.AngularPointing Type
julia
AngularPointing{T}

Abstract type representing a pointing direction identified by two angles in degrees, represented with fields of types Deg{T}.

source
SatcomCoordinates.AzEl Type
julia
AzEl{T} <: AngularPointing{T}

Object specifying a pointing direction in "Elevation/Azimuth" coordinates, defined following the convention used for Azimuth-Elevation-Range (AER) coordinates used by MATLAB and by this package.

This represents the azimuth/elevation definition often used for describing pointing from a user terminal on ground towards a satellite. They are very similar to the ElOverAz definition but with a rotation of the underlying CRS such that elevation is 90° in the direction of the +Z axis and azimuth is computed on the XY plane from the +Y axis towards the +X axis.

Assuming u, v, and w to be direction cosines of the pointing versor ̂p, their relation with the El and Az angles is:

  • u = cos(El) * sin(Az)

  • v = cos(El) * cos(Az)

  • w = sin(El)

Fields

  • az::Deg{T}: The azimuth angle in degrees, constrained to be in the [-180°, 180°] range.

  • el::Deg{T}: The elevation angle in degrees, constrained to be in the [-90°, 90°] range.

Note

The fields of AzEl objects can also be accessed via getproperty using the azimuth and elevation aliases.

See also: ThetaPhi, PointingVersor, UV, ElOverAz, AzOverEl

source
SatcomCoordinates.AzElDistance Type
julia
const AzElDistance{T} = GeneralizedSpherical{AzEl, T}

Type representing a position w.r.t. a local CRS in Azimuth, Elevation and Range coordinates. The difference between AzElDistance and AER is that AER is a always referred to the ENU CRS, while AzElDistance is for a generic local CRS.

source
SatcomCoordinates.AzOverEl Type
julia
AzOverEl{T} <: AngularPointing{T}

Object specifying a pointing direction in "Azimuth over Elevation" coordinates, which specify the elevation and azimuth angles that needs to be fed to an azimuth-over-elevation positioner for pointing to a target towards the pointing direction ̂p.

Following the convention used in most Antenna-related literature, the elevation and azimuth are 0° in the direction of the +Z axis of the reference frame.

Assuming u, v, and w to be direction cosines of the pointing versor ̂p, their relation with the Az and El angles is:

  • u = -sin(Az)

  • v = cos(Az) * sin(El)

  • w = cos(Az) * cos(El)

Note

The equations above are used to represent the "Azimuth over Elevation" coordinates in GRASP. Some textbooks, however, use the opposite convention, meaning that the same equations (with a possible flip in the az/u sign) are used to describe an "Elevation over Azimuth" coordinate system. This is for example the case in the book "Theory and Practice of Modern Antenna Range Measurements" by Clive Parini et al.

Properties

  • az::Deg{T}: The azimuth angle in degrees, constrained to be in the [-180°, 180°] range.

  • el::Deg{T}: The elevation angle in degrees, constrained to be in the [-90°, 90°] range.

Note

The fields of AzOverEl objects can also be accessed via getproperty using the azimuth and elevation aliases.

source
SatcomCoordinates.BasicCRSTransform Type
julia
BasicCRSTransform{T, R <: Union{CRSRotation{T}, Identity}, O <: AbstractPosition{T, 3}} <: AbstractCRSTransform{T}

A type representing a basic transformation (rotation + translation).

Fields:

  • rotation::R: The rotation of the transformation.

  • origin::O: The origin of the transformation.

source
SatcomCoordinates.CRSRotation Type
julia
CRSRotation{T, R <: StaticMatrix{T, 3, 3}} <: AbstractCRSRotation{T}

A type representing a basic rotation of a coordinate system.

source
SatcomCoordinates.ECEF Type
julia
struct ECEF{T} <: AbstractGeocentricPosition{T}

Represents a position in the Earth-Centered, Earth-Fixed (ECEF) coordinate system (or generically for other planets, Ellipsoid-Centered, Ellipsoid-Fixed).

Properties

  • x::Met{T}: X-coordinate in meters

  • y::Met{T}: Y-coordinate in meters

  • z::Met{T}: Z-coordinate in meters

Basic Constructors

ECEF(x::ValidDistance, y::ValidDistance, z::ValidDistance)

ValidDistance is either a Real Number, or a subtype of Unitful.Length.

If of the provided argument is NaN, the returned ECEF object will contain NaN for all fields.

Fallback constructors

All subtypes of P <: AbstractSatcomCoordinate can also be constructed using a Tuple or SVector as input, which will be splatted into the standard constructor for P

See also: ECI, LLA.

source
SatcomCoordinates.ECI Type
julia
struct ECI{T} <: AbstractGeocentricPosition{T}

Represents a position in the Earth-Centered, Inertial (ECI) coordinate system (or generically for other planets, Ellipsoid-Centered, Inertial).

Properties

  • x::Met{T}: X-coordinate in meters

  • y::Met{T}: Y-coordinate in meters

  • z::Met{T}: Z-coordinate in meters

Basic Constructors

ECI(x::ValidDistance, y::ValidDistance, z::ValidDistance)

ValidDistance is either a Real Number, or a subtype of Unitful.Length.

If of the provided argument is NaN, the returned ECI object will contain NaN for all fields.

Fallback constructors

All subtypes of P <: AbstractSatcomCoordinate can also be constructed using a Tuple or SVector as input, which will be splatted into the standard constructor for P

See also: ECEF, LLA.

source
SatcomCoordinates.ENU Type
julia
struct ENU{T} <: AbstractTopocentricPosition{T}

Represents a position in the East-North-Up (ENU) coordinate system, which is a local coordinate system centered at a point on or above the surface of an Ellipsoid. The direction of the ENU axes is uniquely determined by the latitude, longitude, and altitude of the point w.r.t. the referenced ellipsoid.

Properties

  • x::Met{T}: X-coordinate in meters

  • y::Met{T}: Y-coordinate in meters

  • z::Met{T}: Z-coordinate in meters

Basic Constructors

ENU(x::ValidDistance, y::ValidDistance, z::ValidDistance)

ValidDistance is either a Real Number, or a subtype of Unitful.Length.

If of the provided argument is NaN, the returned ENU object will contain NaN for all fields.

Fallback constructors

All subtypes of P <: AbstractSatcomCoordinate can also be constructed using a Tuple or SVector as input, which will be splatted into the standard constructor for P

See also: NED, AER.

source
SatcomCoordinates.ElOverAz Type
julia
ElOverAz{T} <: AngularPointing{T}

Object specifying a pointing direction in "Elevation over Azimuth" coordinates, which specify the azimuth and elevation angles that needs to be fed to an elevation-over-azimuth positioner for pointing to a target towards the pointing direction ̂p.

Following the convention used in most Antenna-related literature, the elevation and azimuth are 0° in the direction of the +Z axis of the reference frame.

Assuming u, v, and w to be direction cosines of the pointing versor ̂p, their relation with the El and Az angles is:

  • u = -sin(Az) * cos(El)

  • v = sin(El)

  • w = cos(Az) * cos(El)

Note

The equations above are used to represent the "Elevation over Azimuth" coordinates in GRASP. Some textbooks, however, use the opposite convention, meaning that the same equations (with sometimes an additional flip in the az/u sign) are used to describe an "Azimuth over Elevation" coordinate system. This is for example the case in the book "Theory and Practice of Modern Antenna Range Measurements" by Clive Parini et al.

Fields

  • az::Deg{T}: The azimuth angle in degrees, constrained to be in the [-180°, 180°] range.

  • el::Deg{T}: The elevation angle in degrees, constrained to be in the [-90°, 90°] range.

Note

The fields of ElOverAz objects can also be accessed via getproperty using the azimuth and elevation aliases.

See also: AzOverEl, ThetaPhi, PointingVersor, UV

source
SatcomCoordinates.GeneralizedSpherical Type
julia
struct GeneralizedSpherical{T, P <: AngularPointing} <: AbstractPosition{T, 3}

Represents a position in a local CRS, defined by an angular pointing direction and a distance.

Note

The parameter P should not be a concrete type of angular pointing but a basic subtype without the numbertype parameter, e.g. AzEl instead of AzEl{Float64}.

Properties

  • angle1::Deg{T}: First angle of the referenced angular pointing CRS P. The name of the property is actually not angle1 but will be the first property name of objects of type P

  • angle2::Deg{T}: Second angle of the referenced angular pointing CRS P. The name of the property is actually not angle2 but will be the second property name of objects of type P

  • r::Met{T}: The distance of the object from the origin of the CRS

Basic Constructors

GeneralizedSpherical(pointing::AngularPointing, r::ValidDistance)
GeneralizedSpherical{P[, T]}(a1::ValidAngle, a2::ValidAngle, r::ValidDistance)

ValidDistance is either a Real Number, or a subtype of Unitful.Length.

If of the provided argument is NaN, the returned PointingAndDistance object will contain NaN for all fields.

See also: Spherical, AzElDistance, LocalCartesian.

source
SatcomCoordinates.InverseTransform Type
julia
InverseTransform{T, F <: AbstractCRSTransform{T}} <: AbstractCRSTransform{T}

A type representing an inverse of an AbstractCRSTransform.

source
SatcomCoordinates.LLA Type
julia
struct LLA{T} <: AbstractGeocentricPosition{T}

Identify a point on or above earth using geodetic coordinates

Properties

  • lat::Deg{T}: Latitude of the point in degrees [-90°, 90°]

  • lon::Deg{T}: Longitude of the point in degrees [-180°, 180°]

  • alt::Met{T}: Altitude of the point above the reference ellipsoid

Basic Constructors

LLA(lat::ValidAngle,lon::ValidAngle,alt::ValidDistance)
LLA(lat::ValidAngle,lon::ValidAngle) # Defaults to 0m altitude

ValidAngle is a either a Real number or a Unitful.Quantity of unit either u"rad" or u"°".

ValidDistance is either a Real Number, or a subtype of Unitful.Length.

If of the provided argument is NaN, the returned LLA object will contain NaN for all fields.

Fallback constructors

All subtypes of P <: AbstractSatcomCoordinate can also be constructed using a Tuple or SVector as input, which will be splatted into the standard constructor for P

See also: ECEF, ECI.

source
SatcomCoordinates.LocalCartesian Type
julia
struct LocalCartesian{T} <: AbstractLocalPosition{T, 3}

Represents a position in a generic local CRS.

Properties

  • x::Met{T}: X-coordinate in meters

  • y::Met{T}: Y-coordinate in meters

  • z::Met{T}: Z-coordinate in meters

Basic Constructors

LocalCartesian(x::ValidDistance, y::ValidDistance, z::ValidDistance)

ValidDistance is either a Real Number, or a subtype of Unitful.Length.

If of the provided argument is NaN, the returned LocalCartesian object will contain NaN for all fields.

Fallback constructors

All subtypes of P <: AbstractSatcomCoordinate can also be constructed using a Tuple or SVector as input, which will be splatted into the standard constructor for P

See also: GeneralizedSpherical.

source
SatcomCoordinates.NED Type
julia
struct NED{T} <: AbstractTopocentricPosition{T}

Represents a position in the North-East-Down (NED) coordinate system, which is a local coordinate system centered at a point on or above the surface of an Ellipsoid. The direction of the NED axes is uniquely determined by the latitude, longitude, and altitude of the point w.r.t. the referenced ellipsoid.

Properties

  • x::Met{T}: X-coordinate in meters

  • y::Met{T}: Y-coordinate in meters

  • z::Met{T}: Z-coordinate in meters

Basic Constructors

NED(x::ValidDistance, y::ValidDistance, z::ValidDistance)

ValidDistance is either a Real Number, or a subtype of Unitful.Length.

If of the provided argument is NaN, the returned NED object will contain NaN for all fields.

Fallback constructors

All subtypes of P <: AbstractSatcomCoordinate can also be constructed using a Tuple or SVector as input, which will be splatted into the standard constructor for P

See also: ENU, AER.

source
SatcomCoordinates.PointingOffset Type
julia
struct PointingOffset{P, T} <: AbstractPointingOffset{T}

Type used to describe an offset between two pointing directions.

Note

The P parameter must be a subtype of Union{UV, AngularPointing} but must be the basic pointing type without specific numbertype (e.g. UV rather than UV{Float64}) as the numbertype is already tracked by the T parameter within the PointingOffset type.

source
SatcomCoordinates.PointingVersor Type
julia
PointingVersor{T} <: AbstractPointing{T}

A unit vector (versor) representing a pointing direction in 3D space. Its components are the x, y, and z components of the unit vector and can also be seen as the u, v, and w direction cosines of the direction identified by the PointingVersor instance.

Properties

  • x::T: The component along the X axis of the corresponding reference frame. Can also be accessed with the u property name.

  • y::T: The component along the Y axis of the corresponding reference frame. Can also be accessed with the v property name.

  • z::T: The component along the Z axis of the corresponding reference frame. Can also be accessed with the w property name.

See also: UV, ThetaPhi, AzOverEl, ElOverAz

source
SatcomCoordinates.Spherical Type
julia
const Spherical{T} = GeneralizedSpherical{ThetaPhi, T}

Type representing a position in ISO/Physics spherical coordinates

source
SatcomCoordinates.ThetaPhi Type
julia
ThetaPhi{T} <: AngularPointing{T}

An object specifying a pointing direction in ThetaPhi coordinates, defined as the θ and φ in the (ISO/Physics definition) spherical coordinates representation

Assuming u, v, and w to be direction cosines of the pointing versor ̂p, their relation with the θ and φ angles is:

  • u = sin(θ) * cos(φ)

  • v = sin(θ) * sin(φ)

  • w = cos(θ)

Properties

  • θ::Deg{T}: The so-called polar angle, representing the angle between the Z axis and the XY plane of the reference frame. It is always normalized to fall within the [-90°, 90°] range.

  • φ::Deg{T}: The so-called azimuth angle, representing the angle between the y and x component of the pointing direction. It is always normalized to fall within the [-180°, 180°] range.

While the field name use the greek letters, the specific fields of an arbitrary ThetaPhi object tp can be accessed with alternative symbols:

  • tp.θ, tp.theta and tp.t can be used to access the θ field

  • tp.φ, tp.ϕ, tp.phi and tp.p can be used to access the φ field

Basic Constructors

ThetaPhi(θ,φ)

The basic constructor takes 2 separate numbers θ, φ and instantiate the object. MarkdownAST.LineBreak()

Provided inputs are intepreted as degrees if not provided using angular quantities from Unitful. If either of the inputs is NaN, the returned ThetaPhi object will contain NaN for both fields.

ThetaPhi(tp)

The ThetaPhi struct can be created using any 2-element Tuple or Vector/StaticVector as input, which will internally call the 2-arguments constructor.

See also: PointingVersor, UV

source
SatcomCoordinates.UV Type
julia
UV{T} <: AbstractPointing{T}

Specify a pointing direction in UV coordinates, which are equivalent to the direction cosines with respect to the X and Y axis of the reference frame. They can also be related to the spherical coordinates (ISO/Physics) spherical coordinates representation by the following equations:

  • u = sin(θ) * cos(φ)

  • v = sin(θ) * sin(φ)

Note

UV coordinates can only be used to represent pointing direction in the half-hemisphere containing the cartesian +Z axis.

Note

The inputs values must also satisfy u^2 + v^2 <= 1 + SatcomCoordinates.UV_CONSTRUCTOR_TOLERANCE[] or an error will be thrown. The SatcomCoordinates.UV_CONSTRUCTOR_TOLERANCE is a Ref{Float64} which defaults to 1e-5 (In case u^2 + v^2 > 1 the inputs are normalized to ensure u^2 + v^2 = 1).

Properties

  • u::T

  • v::T

Basic Constructors

UV{T}(u,v)

The basic constructor takes 2 separate numbers u, v and instantiate the object assuming that u^2 + v^2 <= 1 (condition for valid UV pointing), throwing an error otherwise. MarkdownAST.LineBreak()

If either of the inputs is NaN, the returned UV object will contain NaN for both fields.

UV{T}(uv)

The UV{T} can be created using any 2-element Tuple or StaticVector as input, which will internally call the 2-arguments constructor.

See also: PointingVersor, ThetaPhi

source
SatcomCoordinates.add_angular_offset Method
julia
p = add_angular_offset(p₀::AbstractPointing, offset_angles::Union{ThetaPhi, ThetaPhiOffset})
p = add_angular_offset(p₀::PointingType, θ::ValidAngle, φ::ValidAngle = 0.0)
p = add_angular_offset(output_type, p₀, args...)

Compute the resulting pointing direction p obtained by adding an angular offset expressed as θ and φ angles (following the ISO/Physics convention for spherical coordinates) [deg] to the starting position identified by p₀.

The input starting position p₀ must be any subtype of AbstractPointing The input offset_angles can be provided as an instance of one of the following types:

  • ThetaPhi

  • ThetaPhiOffset

and is converted to ThetaPhiOffset internally, with non-unitful values being interpreted as angles in degrees.

The output is of type output_type if provided or of the same type as p₀ otherwise.

Note

If output_type is UV, the function will throw an error if the final pointing direction is located behind the viewer as the output in UV would be ambiguous. This is not the case for other subtypes of AbstractPointing so an explicit output type (different from UV) should be provided if the target is expected to be behind the viewer.

The offset angles can also be provided separately as 2nd and 3rd argument (optional, defaults to 0.0) to the function using the second method signature. In this case, the inputs are treated as angles in degrees unless explicitly provided using quantitites with ° unit from the Unitful package.

This function performs the inverse operation of get_angular_offset so the following code should return true

julia
using SatcomCoordinates
uv1 = UV(.3,.4)
uv2 = UV(-.2,.5)
offset = get_angular_offset(uv1, uv2)
p = add_angular_offset(uv1, offset)
p  uv2

See also: get_angular_offset, get_angular_distance, ThetaPhi, UV

source
SatcomCoordinates.angle_offset_rotation Method
julia
angle_offset_rotation(θ, φ)

Compute the rotation matrix to find offset points following the procedure in this stackexchnge answer: https://math.stackexchange.com/questions/4343044/rotate-vector-by-a-random-little-amount

source
SatcomCoordinates.change_numbertype Function
julia
change_numbertype(T::Type, x)

Functions that change the underlying numbertype of the provided object x to the first argument T.

It has a fallback default implementation for types defined within this package which calls convert on the provided object x to the type basetype(x){T}.

See also numbertype, enforce_numbertype, has_numbertype, default_numbertype

source
SatcomCoordinates.default_numbertype Method
julia
default_numbertype(args...)

Function that returns the common valid numbertype among the arguments provided as input. It finds the common numbertype via promote_type and either return that (if it's a subtype of AbstractFloat) or Float64 if it's not.

This function is useful to automatically extract from inputs the numbertype to use in constructors.

See also numbertype, enforce_numbertype, has_numbertype, change_numbertype

source
SatcomCoordinates.enforce_numbertype Function
julia
enforce_numbertype(input_type, [default_numbertype]) where {C <: Union{AbstractSatcomCoordinate, AbstractCRSTransform}}

Function that takes as input a type and returns a potentialy more specialized subtype of the input type with the numbertype parameter set if not specified in input_type. Optionally, this function accepts a secon type (or value) as argument and infers the numbertype to set as default (if not alredy present). The default numbertype when the function is called with 1-argument is Float64.

Examples

julia
enforce_numbertype(UV) == UV{Float64} # Provide a default as not present in input type
enforce_numbertype(UV{Float32}) == UV{Float32} # Returns the same input type as it already has a numbertype
enforce_numbertype(UV, Float32) == UV{Float32} # Provide a custom default as not present in input type
enforce_numbertype(UV, 1) == UV{Int64} # Provide a custom default as not present in input type

See also numbertype, enforce_numbertype, has_numbertype, default_numbertype

source
SatcomCoordinates.get_angular_distance Method
julia
get_angular_distance(p₁::AbstractPointing, p₂::AbstractPointing)

Compute the angular distance [°] between the target pointing direction p₂ and the starting pointing direction p₁.

Note

This function's output should be approximately equivalent to the θ (theta) component of get_angular_offset but has a faster implementation. Use this in case the φ (phi) component is not required and speed is important. The following code should evaluate to true

julia
using SatcomCoordinates
uv1 = UV(.3,.4)
uv2 = UV(-.2,.5)
offset = get_angular_offset(uv1, uv2)
Δθ = get_angular_distance(uv1, uv2)
offset.theta  Δθ

See also: add_angular_offset, get_angular_offset

source
SatcomCoordinates.get_angular_offset Method
julia
offset = get_angular_offset(p₁::AbstractPointing, p₂::AbstractPointing)::ThetaPhiOffset

Compute the angular offset required to reach the target pointing direction p₂ from starting pointing direction p₁. The two input pointings can be of any valid AbstractPointing type.

The output is of type ThetaPhiOffset

Note

This function performs the inverse operation of add_angular_offset so the following code should return true

julia
using ReferenceViews
uv1 = UV(.3,.4)
uv2 = UV(-.2,.5)
offset = get_angular_offset(uv1, uv2)
p = add_angular_offset(uv1, offset)
p  uv2

Check out get_angular_distance for a slightly faster implementation in case you only require the angular distance rather than the 2D offset.

See also: add_angular_offset

source
SatcomCoordinates.has_numbertype Function
julia
has_numbertype(T::Type)
has_numbertype(::T)

This function shall return true if the provided type T or object of type T has an associated numbertype.

This function will return false for types defined within this package that do not have the numbertype parameter specified (type T is thus a UnionAll on the numbertype parameter).

source
SatcomCoordinates.numbertype Function
julia
numbertype(T::Type)
numbertype(::T)

This function shall return the underlying numbertype of the provided Type or object.

The concept of numbertype is defined here as the subtype of Real which is used to represent the numerical values in the object's field. It is not directly the type of the fields, mainly as we consider fields of type Unitul.Quantity{T} to have numbertype T.

All the types defined in this package have an assciated parametric numbertype as first parameter.

See also enforce_numbertype, change_numbertype, has_numbertype, default_numbertype

source
SatcomCoordinates.origin Method
julia
origin(t::AbstractAffineCRSTransform)

This function should return an object which is subtype of CartesianPosition and represents the origin of the starting CRS in the target CRS.

source
SatcomCoordinates.raw_properties Method
julia
raw_properties(coords::AbstractSatcomCoordinate)

Generate a NamedTuple starting from the raw SVector holding the coords data and assigning a label to each valid property of the coordinate (as defined by the Base.propertynames(coords) function).

See also raw_svector

source
SatcomCoordinates.raw_svector Method
julia
raw_svector(coord::AbstractSatcomCoordinate)

Extracts the raw SVector storing the data for the provided coordinate. This assumes that coord has a field called svector and just calls getfield(coord, :svector).

Concrete subtypes that do not follow this convention should overload this function.

See also raw_properties

source
SatcomCoordinates.rotation Method
julia
rotation(t::AbstractAffineCRSTransform)

This function should return a CRSRotation object representing the rotation to align the starting CRS to the target CRS.

source
SatcomCoordinates.wrap_spherical_angles Method
julia
az, el = wrap_spherical_angles(az::ValidAngle, el::ValidAngle, ::Type{<:ThetaPhi}) where T <: Deg{<:Real}
θ, φ = wrap_spherical_angles::ValidAngle, φ::ValidAngle, ::Type{<:Union{AzOverEl, ElOverAz}}) where T <: Deg{<:Real}

Function that takes as input two angles representing two orthogonal angular components of spherical coordinates (e.g. θ/φ, el/az, etc.) and returns two angles normalized to a consistent wrapping identifying the full sphere:

  • θ/φ angles are wrapped such that θ ∈ [0°, 180°] and φ ∈ [-180°, 180°]

  • el/az angles are wrapped such that el ∈ [-90°, 90°] and az ∈ [-180°, 180°]

!!!

source
SatcomCoordinates.wrap_spherical_angles_normalized Method
julia
az, el = wrap_spherical_angles_normalized(az::T, el::T, ::Type{<:ThetaPhi}) where T <: Deg{<:Real}
θ, φ = wrap_spherical_angles_normalized::T, φ::T, ::Type{<:Union{AzOverEl, ElOverAz}}) where T <: Deg{<:Real}

Function that takes as input two angles representing two orthogonal angular components of spherical coordinates (e.g. θ/φ, el/az, etc.) and returns two angles normalized to a consistent wrapping identifying the full sphere:

  • θ/φ angles are wrapped such that θ ∈ [0°, 180°] and φ ∈ [-180°, 180°]

  • el/az angles are wrapped such that el ∈ [-90°, 90°] and az ∈ [-180°, 180°]

Note

This function already assumes that the provided input angles are already normalized such that both are in the [-180°, 180°] range. If you want to normalize the inputs automatically use the wrap_first_angle function.

source