SatcomCoordinates.jl
Documentation for SatcomCoordinates.jl
.
SatcomCoordinates.ForwardOrInverse Type
const ForwardOrInverse{F <: AbstractCRSTransform} = Union{F, InverseTransform{<:Any, <:F}}
Union representing either a forward or reverse transform of F
sourceSatcomCoordinates.WithNumbertype Type
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.
sourceSatcomCoordinates.AER Type
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
SatcomCoordinates.AbstractAffineCRSTransform Type
AbstractAffineCRSTransform{T}
Abstract type representing an affine transform between two CRSs with numbertype T
.
SatcomCoordinates.AbstractCRSTransform Type
AbstractCRSTransform{T}
Abstract type representing a coordinate transform between two CRSs with numbertype T
.
SatcomCoordinates.AbstractFieldValue Type
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):
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 thesvector
fieldraw_properties(::AbstractFieldValue)
: Returns a NamedTuple with the properties of thesvector
fieldproperties_names(::Type{<:AbstractFieldValue{U, CRS}})
: Returns the property names of theCRS
parametric typeDefault constructor taking
N
numbers or a Tuple/Svector, with the values interpreted as quantities of unitU
if not provided directly as quantititesConstructionBase.getproperties(::AbstractFieldValue{U, CRS})
: Returns a NamedTuple with the properties of the field assuming the properties names ofCRS
and returning values with unitU
.Base.getproperty
: extraction of properties directly from the NamedTuple returned byConstructionBase.getproperties
SatcomCoordinates.AbstractGeocentricPosition Type
AbstractGeocentricPosition{T} <: AbstractPosition{T, 3}
Abstract type representing a position in a geocentric coordinate system.
sourceSatcomCoordinates.AbstractLocalPosition Type
AbstractLocalPosition{T} <: AbstractPosition{T, 3}
Abstract type representing a position in a local coordinate system.
sourceSatcomCoordinates.AbstractPointing Type
AbstractPointing{T} <: AbstractSatcomCoordinate{T, 3}
Abstract type representing a pointing direction in 3 dimensions which is backed by fields with shared numbertype
T
.
SatcomCoordinates.AbstractPointingOffset Type
AbstractPointingOffset{T} <: AbstractSatcomCoordinate{T, 2}
Abstract type representing a pointing offset between two pointing directions.
Currently only has two concrete subtypes: UVOffset
and ThetaPhiOffset
.
SatcomCoordinates.AbstractSatcomCoordinate Type
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
SatcomCoordinates.AbstractTopocentricPosition Type
AbstractTopocentricPosition{T} <: AbstractPosition{T, 3}
Abstract type representing a position in a topocentric coordinate system.
sourceSatcomCoordinates.AngularPointing Type
AngularPointing{T}
Abstract type representing a pointing direction identified by two angles in degrees, represented with fields of types Deg{T}
.
SatcomCoordinates.AzEl Type
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
SatcomCoordinates.AzElDistance Type
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.
SatcomCoordinates.AzOverEl Type
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.
SatcomCoordinates.BasicCRSTransform Type
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.
SatcomCoordinates.CRSRotation Type
CRSRotation{T, R <: StaticMatrix{T, 3, 3}} <: AbstractCRSRotation{T}
A type representing a basic rotation of a coordinate system.
sourceSatcomCoordinates.ECEF Type
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 metersy::Met{T}
: Y-coordinate in metersz::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
SatcomCoordinates.ECI Type
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 metersy::Met{T}
: Y-coordinate in metersz::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
SatcomCoordinates.ENU Type
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 metersy::Met{T}
: Y-coordinate in metersz::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
SatcomCoordinates.ElOverAz Type
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
SatcomCoordinates.GeneralizedSpherical Type
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 CRSP
. The name of the property is actually notangle1
but will be the first property name of objects of typeP
angle2::Deg{T}
: Second angle of the referenced angular pointing CRSP
. The name of the property is actually notangle2
but will be the second property name of objects of typeP
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
.
SatcomCoordinates.InverseTransform Type
InverseTransform{T, F <: AbstractCRSTransform{T}} <: AbstractCRSTransform{T}
A type representing an inverse of an AbstractCRSTransform
.
SatcomCoordinates.LLA Type
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
SatcomCoordinates.LocalCartesian Type
struct LocalCartesian{T} <: AbstractLocalPosition{T, 3}
Represents a position in a generic local CRS.
Properties
x::Met{T}
: X-coordinate in metersy::Met{T}
: Y-coordinate in metersz::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
.
SatcomCoordinates.NED Type
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 metersy::Met{T}
: Y-coordinate in metersz::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
SatcomCoordinates.PointingOffset Type
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.
SatcomCoordinates.PointingVersor Type
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 theu
property name.y::T
: The component along the Y axis of the corresponding reference frame. Can also be accessed with thev
property name.z::T
: The component along the Z axis of the corresponding reference frame. Can also be accessed with thew
property name.
See also: UV
, ThetaPhi
, AzOverEl
, ElOverAz
SatcomCoordinates.Spherical Type
const Spherical{T} = GeneralizedSpherical{ThetaPhi, T}
Type representing a position in ISO/Physics spherical coordinates
sourceSatcomCoordinates.ThetaPhi Type
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-calledpolar angle
, representing the angle between theZ
axis and theXY
plane of the reference frame. It is always normalized to fall within the [-90°, 90°] range.φ::Deg{T}
: The so-calledazimuth angle
, representing the angle between they
andx
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
andtp.t
can be used to access theθ
fieldtp.φ
,tp.ϕ
,tp.phi
andtp.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
SatcomCoordinates.UV Type
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
SatcomCoordinates.add_angular_offset Method
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
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
SatcomCoordinates.angle_offset_rotation Method
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
sourceSatcomCoordinates.change_numbertype Function
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
SatcomCoordinates.default_numbertype Method
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
SatcomCoordinates.enforce_numbertype Function
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
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
SatcomCoordinates.get_angular_distance Method
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
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
SatcomCoordinates.get_angular_offset Method
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
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
SatcomCoordinates.has_numbertype Function
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).
SatcomCoordinates.numbertype Function
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
SatcomCoordinates.origin Method
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.
SatcomCoordinates.raw_properties Method
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
SatcomCoordinates.raw_svector Method
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
SatcomCoordinates.rotation Method
rotation(t::AbstractAffineCRSTransform)
This function should return a CRSRotation
object representing the rotation to align the starting CRS to the target CRS.
SatcomCoordinates.wrap_spherical_angles Method
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 thatel ∈ [-90°, 90°]
andaz ∈ [-180°, 180°]
!!!
sourceSatcomCoordinates.wrap_spherical_angles_normalized Method
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 thatel ∈ [-90°, 90°]
andaz ∈ [-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.