Skip to content

add simple projected coordinates #43

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
aplavin wants to merge 3 commits into
base: master
Choose a base branch
from
+61 −3
Open

add simple projected coordinates #43

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
7f5297d
fix isapprox
aplavin Nov 10, 2022
17fdff6
add simple projected coordinates
aplavin Nov 10, 2022
56fbeae
remove _center_angle
aplavin Feb 7, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions src/SkyCoords.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,7 @@ export AbstractSkyCoords,
GalCoords,
FK5Coords,
CartesianCoords,
ProjectedCoords,
separation,
position_angle,
offset,
Expand All @@ -18,6 +19,7 @@ export AbstractSkyCoords,

include("types.jl")
include("cartesian.jl")
include("projected.jl")

# -----------------------------------------------------------------------------
# Helper functions: Create rotation matrix about a given axis (x, y, z)
Expand Down
30 changes: 30 additions & 0 deletions src/projected.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,30 @@
abstract type AbstractProjectedCoords <: AbstractSkyCoords end

""" Projected coordinates with origin and offset.
Approximation for small FoV. """
struct ProjectedCoords{TC,T} <: AbstractProjectedCoords
origin::TC
offset::SVector{2,T}
end

origin(c::ProjectedCoords) = c.origin
lon(c::AbstractProjectedCoords) = lon(origin(c)) + c.offset[1] / cos(lat(origin(c)))
lat(c::AbstractProjectedCoords) = lat(origin(c)) + c.offset[2]

Base.convert(::Type{T}, c::T) where {T<:AbstractProjectedCoords} = c
Copy link
Member

@giordano giordano Dec 25, 2022

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is this method needed? Isn't it already covered by https://github.com/JuliaLang/julia/blob/162ee48e1c34b2a2cd797395353f19a7aca21aa2/base/Base.jl#L64?

Copy link
Member Author

@aplavin aplavin Feb 7, 2023

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The Base method is less specific and won't be called when converting skycoords. Same motivation as

SkyCoords.jl/src/types.jl

Line 76 in 6c8849b

Base.convert(::Type{T}, c::T) where {T<:AbstractSkyCoords} = c
.

Base.convert(::Type{TCto}, c::AbstractProjectedCoords) where {TCto<:AbstractSkyCoords} =
convert(TCto,
constructorof(typeof(origin(c)))(
lon(c), lat(c)
)
)

Base.isapprox(a::ProjectedCoords, b::ProjectedCoords; kwargs...) = isapprox(origin(a), origin(b); kwargs...) && isapprox(a.offset, b.offset; kwargs...)

function project(origin::AbstractSkyCoords, c::AbstractSkyCoords)
cc = convert(typeof(origin), c)
Δlon = rem2pi(lon(cc) - lon(origin), RoundNearest)
offset = SVector(Δlon * cos(lat(origin)), lat(cc) - lat(origin))
ProjectedCoords(origin, offset)
end

15 changes: 12 additions & 3 deletions src/types.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -82,6 +82,15 @@ function Base.convert(::Type{T}, c::S) where {T<:AbstractSkyCoords,S<:AbstractSk
end

Base.:(==)(a::T, b::T) where {T<:AbstractSkyCoords} = lon(a) == lon(b) && lat(a) == lat(b)
Base.isapprox(a::ICRSCoords, b::ICRSCoords; kwargs...) = isapprox(SVector(lon(a), lat(a)), SVector(lon(b), lat(b)); kwargs...)
Base.isapprox(a::GalCoords, b::GalCoords; kwargs...) = isapprox(SVector(lon(a), lat(a)), SVector(lon(b), lat(b)); kwargs...)
Base.isapprox(a::FK5Coords{e}, b::FK5Coords{e}; kwargs...) where {e} = isapprox(SVector(lon(a), lat(a)), SVector(lon(b), lat(b)); kwargs...)
Base.isapprox(a::ICRSCoords, b::ICRSCoords; kwargs...) = isapprox(
SVector(lon(a), lat(a)),
SVector(lon(a) + rem2pi(lon(b) - lon(a), RoundNearest), lat(b));
kwargs...)
Base.isapprox(a::GalCoords, b::GalCoords; kwargs...) = isapprox(
SVector(lon(a), lat(a)),
SVector(lon(a) + rem2pi(lon(b) - lon(a), RoundNearest), lat(b));
kwargs...)
Base.isapprox(a::FK5Coords{e}, b::FK5Coords{e}; kwargs...) where {e} = isapprox(
SVector(lon(a), lat(a)),
SVector(lon(a) + rem2pi(lon(b) - lon(a), RoundNearest), lat(b));
kwargs...)
17 changes: 17 additions & 0 deletions test/runtests.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@ using ConstructionBase: setproperties
using DelimitedFiles
using LinearAlgebra: normalize
using SkyCoords
using SkyCoords: project, origin
using StableRNGs
using Statistics
using Test
Expand All @@ -15,6 +16,19 @@ rad2arcsec(r) = 3600 * rad2deg(r)
# tests against astropy.coordinates
include("astropy.jl")

@testset "projected coords" begin
c0 = ICRSCoords(0.1, -0.2)
c1 = ICRSCoords(0.1 + 1e-5, -0.2 + 3e-5)
cp = project(c0, c1)::ProjectedCoords
@test origin(cp) == c0
@test cp.offset[1] ≈ 0.98 * 1e-5 rtol=1e-4
@test cp.offset[2] ≈ 3e-5
@test convert(ICRSCoords, cp) ≈ c1
@test convert(GalCoords, cp) ≈ convert(GalCoords, c1)
@test cp == cp
@test cp ≈ cp
end

# Test separation between coordinates and conversion with mixed floating types.
@testset "Separation" begin
c1 = ICRSCoords(ℯ, pi / 2)
Expand Down Expand Up @@ -177,6 +191,9 @@ end
@test !(c1 ≈ c4)
@test c1 ≈ c2 rtol=1e-3
@test c1 ≈ c4 rtol=1e-3

@test ICRSCoords(eps(), 1) ≈ ICRSCoords(0, 1)
@test ICRSCoords(eps(), 1) ≈ ICRSCoords(-eps(), 1)
end

@test_broken (!(ICRSCoords(1, 2) ≈ FK5Coords{2000}(1, 2)); true)
Expand Down