Skip to content
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.

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

Make the SatelliteToolbox be an extension instead of a direct dependency. #70

Merged
merged 3 commits into from
Apr 13, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
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
10 changes: 5 additions & 5 deletions .github/workflows/ci.yml
Original file line number Diff line number Diff line change
Expand Up @@ -24,20 +24,20 @@ jobs:
fail-fast: false
matrix:
version:
- '1.8' # Replace this with the minimum Julia version that your package supports. E.g. if your package requires Julia 1.5 or higher, change this to '1.5'.
- '1' # Replace this with the minimum Julia version that your package supports. E.g. if your package requires Julia 1.5 or higher, change this to '1.5'.
#- '1' # Leave this line unchanged. '1' will automatically expand to the latest stable 1.x release of Julia.
#- 'nightly'
os:
#- ubuntu-latest
- ubuntu-20.04
- ubuntu-latest
#- ubuntu-20.04
#- macos-latest
arch:
- x64
include:
# Linux
- name: Linux - Compile only
os: ubuntu-20.04
#os: ubuntu-latest
#os: ubuntu-20.04
os: ubuntu-latest
run_in_pr : true

steps:
Expand Down
12 changes: 10 additions & 2 deletions Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -8,19 +8,27 @@ GMT = "5752ebe1-31b9-557e-87aa-f909b540aa54"
Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7"
Dates = "ade2ca70-3891-5945-98fb-dc099432e06a"
Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a"

[weakdeps]
SatelliteToolboxTle = "7ff27aeb-5fff-4337-a9ee-a9fe6b7ed35e"
SatelliteToolboxPropagators = "c2b69894-ea78-4e2b-9ba6-cedbbc3d14d7"
SatelliteToolboxTransformations = "6b019ec1-7a1e-4f04-96c7-a9db1ca5514d"
PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a"

[extensions]
RemoteSSatTbxExt = ["SatelliteToolboxTle", "SatelliteToolboxPropagators", "SatelliteToolboxTransformations"]

[extras]
Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
SatelliteToolboxTle = "7ff27aeb-5fff-4337-a9ee-a9fe6b7ed35e"
SatelliteToolboxPropagators = "c2b69894-ea78-4e2b-9ba6-cedbbc3d14d7"
SatelliteToolboxTransformations = "6b019ec1-7a1e-4f04-96c7-a9db1ca5514d"

[targets]
test = ["Test"]

[compat]
julia = "1.6"
julia = "1.9"
GMT = "1.9"
PrecompileTools = "1.0"
SatelliteToolboxTle = "1"
Expand Down
105 changes: 105 additions & 0 deletions ext/RemoteSSatTbxExt/RemoteSSatTbxExt.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,105 @@
module RemoteSSatTbxExt
using RemoteS, GMT, Dates
using SatelliteToolboxTle, SatelliteToolboxPropagators, SatelliteToolboxTransformations
using PrecompileTools

function RemoteS.sat_tracks_ext(; geocentric::Bool=false, tiles::Bool=false, position::Bool=false, kwargs...)
(position && tiles) && error("Cannot require tiles and a single position. Makes no sense.")
d = GMT.KW(kwargs)

start = ((val = find_in_dict(d, [:start])[1]) === nothing) ? now(Dates.UTC) : RemoteS.getitDTime(val)

(tiles) && (sat_name = get_sat_name(d))
(tiles) && (start = round(start, Dates.Minute(5))) # This is for MODIS only

if ((val = find_in_dict(d, [:duration])[1]) !== nothing)
if (isa(val, String)) # Accept duration in D(ays), h(ours), m(inutes) or s(econds)
if (endswith(val,"D")) dur = Day(parse(Int, val[1:end-1]))
elseif (endswith(val,"h")) dur = Hour(parse(Int, val[1:end-1]))
elseif (endswith(val,"m")) dur = Minute(parse(Int, val[1:end-1]))
elseif (endswith(val,"s")) dur = Second(parse(Int, val[1:end-1]))
else error("Only 'D', 'h', 'm' or 's' are accepted in duration")
end
elseif (isa(val, Real)) # Assume duration was given in minutes
dur = Minute(trunc(Int, val))
end
stop = start + dur
else
if ((val = find_in_dict(d, [:stop])[1]) !== nothing)
stop = RemoteS.getitDTime(val)
else
stop = start + Minute(100) # Default is ~Terra period
end
end
if ((val = find_in_dict(d, [:step :inc :dt])[1]) !== nothing) # Steps are in seconds
if (isa(val, String))
if (endswith(val,"m")) dt = parse(Int, val[1:end-1]) * 60
elseif (endswith(val,"s")) dt = parse(Int, val[1:end-1])
else error("Only 's' or 'm' are accepted in increment")
end
else
dt = trunc(Int, val)
end
else
dt = 30
end

if ((val_tle = find_in_dict(d, [:tle_obj])[1]) !== nothing) tle = val_tle # Some other fun already got it
else tle = RemoteS.loadTLE(d)
end
orbp = SatelliteToolboxPropagators.Propagators.init(Val(:SGP4), tle)

epoch_jd = orbp.sgp4d.epoch
startmfe = (datetime2julian(DateTime(start)) - epoch_jd) * 24 * 3600
stopmfe = (datetime2julian(DateTime(stop)) - epoch_jd) * 24 * 3600
(tiles) && (dt = 60) # Arbitrary choice that works well for MODIS but may need revision for others
t = startmfe:dt:stopmfe

(position) && (t = [t[1]]) # Single position. Doing it here wastes work above but code is way cleaner

out = Matrix{Float64}(undef, length(t), 4)
#r = SatelliteToolboxPropagators.Propagators.propagate!.(orbp, t)[1] # Doesn't work. Why?

for n = 1:GMT.numel(t)
r = SatelliteToolboxPropagators.Propagators.propagate!(orbp, t[n])[1]
jd = epoch_jd + t[n] / (24 * 3600)
tt = SatelliteToolboxTransformations.r_eci_to_ecef(SatelliteToolboxTransformations.TEME(), SatelliteToolboxTransformations.PEF(), jd) * r
out[n,1], out[n,2], out[n,3], out[n, 4] = tt[1], tt[2], tt[3], jd
end

if (tiles)
out = mapproject(out, E=true, I=true)
return make_sat_tiles(out[1].data, SCENE_HALFW[sat_name], sat_name)
end

if (geocentric)
D = mat2ds(out, geom=UInt32(2), colnames=["X", "Y", "Z", "JulianDay"])
else
D = mapproject(out, E=true, I=true)
D.colnames, D.proj4, D.geom = ["lon", "lat", "alt", "JulianDay"], GMT.prj4WGS84, UInt32(2)
end
D
end

# --------------------------------------------------------------------------------------------
function RemoteS.loadTLE_ext(d::Dict)
# Load a TLE or use a default one. In a function because it's used by two functions
if ((val = find_in_dict(d, [:tle :TLE])[1]) !== nothing)
if (isa(val, String)) tle = SatelliteToolboxTle.read_tle(val)
elseif (isa(val, Vector{String}) && length(val) == 2)
tle = SatelliteToolboxTle.read_tle(val[1], val[2])
else
error("BAD input TLE data")
end
else
tle = SatelliteToolboxTle.read_tle("C:\\v\\AQUA.tle")
end
return tle
end

@setup_workload begin
RemoteS.sat_tracks_ext(tle=["1 27424U 02022A 21245.83760660 .00000135 00000-0 39999-4 0 9997"; "2 27424 98.2123 186.0654 0002229 67.6025 313.3829 14.57107527 28342"], duration=100, geocentric=true)
RemoteS.sat_tracks_ext(position=true, tle=["1 27424U 02022A 21245.83760660 .00000135 00000-0 39999-4 0 9997"; "2 27424 98.2123 186.0654 0002229 67.6025 313.3829 14.57107527 28342"]);
end

end
3 changes: 0 additions & 3 deletions src/RemoteS.jl
Original file line number Diff line number Diff line change
@@ -1,7 +1,6 @@
module RemoteS

using GMT, Printf, Statistics, Dates#, Requires
using SatelliteToolboxTle, SatelliteToolboxPropagators, SatelliteToolboxTransformations
using PrecompileTools

const SCENE_HALFW = Dict("AQUA" => 1163479, "TERRA" => 1163479, "LANDSAT8" => 92500) # half widths
Expand All @@ -22,8 +21,6 @@ include("utils.jl")
include("sat_tracks.jl")

@setup_workload begin
sat_tracks(tle=["1 27424U 02022A 21245.83760660 .00000135 00000-0 39999-4 0 9997"; "2 27424 98.2123 186.0654 0002229 67.6025 313.3829 14.57107527 28342"], duration=100, geocentric=true)
sat_tracks(position=true, tle=["1 27424U 02022A 21245.83760660 .00000135 00000-0 39999-4 0 9997"; "2 27424 98.2123 186.0654 0002229 67.6025 313.3829 14.57107527 28342"]);
truecolor(mat2img(rand(UInt16,128,128)), mat2img(rand(UInt16,128,128)), mat2img(rand(UInt16,128,128)));
get_MODIS_scene_name(datetime2julian(DateTime("2020-09-20")), "A");
reportbands(mat2img(rand(UInt16, 4,4,3), names=["Band 1", "Band 2", "Band 3"]), 3)[1];
Expand Down
110 changes: 15 additions & 95 deletions src/sat_tracks.jl
Original file line number Diff line number Diff line change
Expand Up @@ -37,84 +37,21 @@ and the orbit track can be visualized with
imshow(orb, proj=:Robinson, region=:global, coast=true)
"""
function sat_tracks(; geocentric::Bool=false, tiles::Bool=false, position::Bool=false, kwargs...)
# ...
(position && tiles) && error("Cannot require tiles and a single position. Makes no sense.")
d = KW(kwargs)

start = ((val = find_in_dict(d, [:start])[1]) === nothing) ? now(Dates.UTC) : getitDTime(val)

(tiles) && (sat_name = get_sat_name(d))
(tiles) && (start = round(start, Dates.Minute(5))) # This is for MODIS only

if ((val = find_in_dict(d, [:duration])[1]) !== nothing)
if (isa(val, String)) # Accept duration in D(ays), h(ours), m(inutes) or s(econds)
if (endswith(val,"D")) dur = Day(parse(Int, val[1:end-1]))
elseif (endswith(val,"h")) dur = Hour(parse(Int, val[1:end-1]))
elseif (endswith(val,"m")) dur = Minute(parse(Int, val[1:end-1]))
elseif (endswith(val,"s")) dur = Second(parse(Int, val[1:end-1]))
else error("Only 'D', 'h', 'm' or 's' are accepted in duration")
end
elseif (isa(val, Real)) # Assume duration was given in minutes
dur = Minute(trunc(Int, val))
end
stop = start + dur
else
if ((val = find_in_dict(d, [:stop])[1]) !== nothing)
stop = getitDTime(val)
else
stop = start + Minute(100) # Default is ~Terra period
end
end
if ((val = find_in_dict(d, [:step :inc :dt])[1]) !== nothing) # Steps are in seconds
if (isa(val, String))
if (endswith(val,"m")) dt = parse(Int, val[1:end-1]) * 60
elseif (endswith(val,"s")) dt = parse(Int, val[1:end-1])
else error("Only 's' or 'm' are accepted in increment")
end
else
dt = trunc(Int, val)
end
else
dt = 30
end

if ((val_tle = find_in_dict(d, [:tle_obj])[1]) !== nothing) tle = val_tle # Some other fun already got it
else tle = loadTLE(d)
end
orbp = SatelliteToolboxPropagators.Propagators.init(Val(:SGP4), tle)

epoch_jd = orbp.sgp4d.epoch
startmfe = (datetime2julian(DateTime(start)) - epoch_jd) * 24 * 3600
stopmfe = (datetime2julian(DateTime(stop)) - epoch_jd) * 24 * 3600
(tiles) && (dt = 60) # Arbitrary choice that works well for MODIS but may need revision for others
t = startmfe:dt:stopmfe

(position) && (t = [t[1]]) # Single position. Doing it here wastes work above but code is way cleaner

out = Matrix{Float64}(undef, length(t), 4)
#r = SatelliteToolboxPropagators.Propagators.propagate!.(orbp, t)[1] # Doesn't work. Why?

for n = 1:length(t)
r = SatelliteToolboxPropagators.Propagators.propagate!(orbp, t[n])[1]
jd = epoch_jd + t[n] / (24 * 3600)
tt = SatelliteToolboxTransformations.r_eci_to_ecef(SatelliteToolboxTransformations.TEME(), SatelliteToolboxTransformations.PEF(), jd) * r
out[n,1], out[n,2], out[n,3], out[n, 4] = tt[1], tt[2], tt[3], jd
end
!(isdefined(Main, :SatelliteToolboxPropagators)) &&
(printstyled("Satellite Toolbox Propagators not loaded. Load them with:\n\n"; color=:blue); printstyled("using SatelliteToolboxTle, SatelliteToolboxPropagators, SatelliteToolboxTransformations"; color=:yellow); return nothing)
sat_tracks_ext(; geocentric=geocentric, tiles=tiles, position=position, kwargs...)
end

if (tiles)
out = mapproject(out, E=true, I=true)
return make_sat_tiles(out[1].data, SCENE_HALFW[sat_name], sat_name)
end
function sat_tracks_ext end

if (geocentric)
D = mat2ds(out, geom=UInt32(2), colnames=["X", "Y", "Z", "JulianDay"])
else
D = mapproject(out, E=true, I=true)
D.colnames, D.proj4, D.geom = ["lon", "lat", "alt", "JulianDay"], GMT.prj4WGS84, UInt32(2)
end
D
function loadTLE(d::Dict)
!(isdefined(Main, :SatelliteToolboxTle)) &&
(printstyled("Satellite Toolbox Propagators not loaded. Load them with:\n\n"; color=:blue); printstyled("using SatelliteToolboxTle, SatelliteToolboxPropagators, SatelliteToolboxTransformations"; color=:yellow); return nothing)
loadTLE_ext(d)
end

function loadTLE_ext end

# --------------------------------------------------------------------------------------------
function getitDTime(val)
if isa(val, String) || isa(val, Tuple) ret::DateTime = DateTime(val)
Expand All @@ -124,23 +61,6 @@ function getitDTime(val)
return ret
end

# --------------------------------------------------------------------------------------------
function loadTLE(d::Dict)
# Load a TLE or use a default one. In a function because it's used by two functions
if ((val = find_in_dict(d, [:tle :TLE])[1]) !== nothing)
if (isa(val, String)) tle = SatelliteToolboxTle.read_tle(val)
elseif (isa(val, Vector{String}) && length(val) == 2)
tle = SatelliteToolboxTle.read_tle(val[1], val[2])
else
error("BAD input TLE data")
end
else
#tle = SatelliteToolbox.read_tle("C:\\v\\Landsat8.tle")
tle = SatelliteToolboxTle.read_tle("C:\\v\\AQUA.tle")
end
tle
end

# --------------------------------------------------------------------------------------------
function get_sat_name(d::Dict)::String
# Get the satellite name from kwargs (encoded in 'd'). Used by at least two functions
Expand Down Expand Up @@ -189,7 +109,7 @@ function sat_scenes(track, sat_name::String)
ll14 = geod(track[k,1:2], [azim[k]+90, azim[k]-90], halfwidth)[1]
ll23 = geod(track[k+5,1:2], [azim[k+5]+90, azim[k+5]-90], halfwidth)[1]
sc = get_MODIS_scene_name(track[k,4], sat_name)
D[n+=1] = GMTdataset([ll14[1:1,:]; ll23[1:1,:]; ll23[2:2,:]; ll14[2:2,:]; ll14[1:1,:]], Float64[], Float64[], Dict{String, String}(), ["lon","lat"], String[], sc, String[], "", "", 0, GMT.Gdal.wkbPolygon)
D[n+=1] = GMTdataset([ll14[1:1,:]; ll23[1:1,:]; ll23[2:2,:]; ll14[2:2,:]; ll14[1:1,:]], Float64[], Float64[], GMT.DictSvS(), ["lon","lat"], String[], sc, String[], "", "", 0, GMT.Gdal.wkbPolygon)
end
D[1].proj4 = GMT.prj4WGS84
D
Expand Down Expand Up @@ -237,7 +157,7 @@ function clip_orbits(track, bb::Vector{<:Real})
end_seg = [inds[2:2:length(inds)]; length(azim)+1]
D = Vector{GMTdataset}(undef, length(begin_seg))
colnames, prj4 = isa(track, Vector) ? (track[1].colnames, track[1].proj4) : (track.colnames, track.proj4)
for k = 1:length(begin_seg)
for k = 1:GMT.numel(begin_seg)
D[k] = GMTdataset(isa(segments, Vector) ? segments[1][begin_seg[k]:end_seg[k], :] : segments[begin_seg[k]:end_seg[k], :], Float64[], Float64[], GMT.DictSvS(), colnames, String[], "", String[], prj4, "", 0, GMT.Gdal.wkbLineString)
end
D[1].proj4 = GMT.prj4WGS84
Expand Down Expand Up @@ -312,7 +232,7 @@ function findscenes(lon::Real, lat::Real; kwargs...)

scenes = Vector{String}(undef, 0)
dists = mapproject(D, G=(lon,lat))
for k = 1:length(dists)
for k = 1:GMT.numel(dists)
d, ind = findmin(view(dists[k].data, :,5))
if (d <= SCENE_HALFW[sat])
jd = datetime2julian(floor(julian2datetime(dists[k][ind,4]), Dates.Minute(5)))
Expand All @@ -329,7 +249,7 @@ function day_night_orbits(D; day::Bool=false, night::Bool=false)
(!day && !night) && return D # No selection requested

pass = zeros(Bool, length(D))
for k = 1:length(D)
for k = 1:GMT.numel(D)
lon, lat = D[k][1,1:2]
jd = D[k][1,4]
raise, set = solar(I=@sprintf("%.4f/%.4f+d%s", lon, lat, string(julian2datetime(jd))), C=true)[5:6]
Expand Down
Loading