From 2ae60d036b643cac2bce9f1d6d3808d0f11f2b6e Mon Sep 17 00:00:00 2001 From: Graham McIntire Date: Mon, 30 Mar 2026 08:36:21 -0500 Subject: [PATCH] Add SRTM auto-download on tile miss from AWS Terrain Tiles When a local .hgt tile is missing, download it from the public AWS S3 skadi bucket, decompress with zlib, and write to the tiles directory before retrying the lookup. Falls back to Open-Meteo/OpenTopo APIs if the download fails. --- config/dev.exs | 3 + config/runtime.exs | 1 + config/test.exs | 1 + lib/microwaveprop/terrain/elevation_client.ex | 17 ++ lib/microwaveprop/terrain/srtm.ex | 175 +++++++++++++++++ .../terrain/elevation_client_test.exs | 50 +++++ test/microwaveprop/terrain/srtm_test.exs | 180 ++++++++++++++++++ 7 files changed, 427 insertions(+) create mode 100644 lib/microwaveprop/terrain/srtm.ex create mode 100644 test/microwaveprop/terrain/srtm_test.exs diff --git a/config/dev.exs b/config/dev.exs index 73379e92..93533b11 100644 --- a/config/dev.exs +++ b/config/dev.exs @@ -73,6 +73,9 @@ config :microwaveprop, MicrowavepropWeb.Endpoint, # Enable dev routes for dashboard and mailbox config :microwaveprop, dev_routes: true +# Use local SRTM1 tiles for elevation lookups instead of the Open-Meteo API +config :microwaveprop, srtm_tiles_dir: Path.expand("~/srtm/tiles") + # Initialize plugs at runtime for faster development compilation config :phoenix, :plug_init_mode, :runtime diff --git a/config/runtime.exs b/config/runtime.exs index 9db82ec3..82577ee9 100644 --- a/config/runtime.exs +++ b/config/runtime.exs @@ -66,6 +66,7 @@ if config_env() == :prod do secret_key_base: secret_key_base config :microwaveprop, :dns_cluster_query, System.get_env("DNS_CLUSTER_QUERY") + config :microwaveprop, srtm_tiles_dir: "/srtm" # ## SSL Support # diff --git a/config/test.exs b/config/test.exs index 13f5f531..4d13c776 100644 --- a/config/test.exs +++ b/config/test.exs @@ -34,6 +34,7 @@ config :microwaveprop, iem_req_options: [plug: {Req.Test, Microwaveprop.Weather. # Route HTTP requests through Req.Test stubs config :microwaveprop, solar_req_options: [plug: {Req.Test, Microwaveprop.Weather.SolarClient}] +config :microwaveprop, srtm_req_options: [plug: {Req.Test, Microwaveprop.Terrain.Srtm}, retry: false] # Initialize plugs at runtime for faster test compilation config :phoenix, :plug_init_mode, :runtime diff --git a/lib/microwaveprop/terrain/elevation_client.ex b/lib/microwaveprop/terrain/elevation_client.ex index 8b038682..ab669629 100644 --- a/lib/microwaveprop/terrain/elevation_client.ex +++ b/lib/microwaveprop/terrain/elevation_client.ex @@ -1,11 +1,26 @@ defmodule Microwaveprop.Terrain.ElevationClient do @moduledoc false + alias Microwaveprop.Terrain.Srtm + @batch_size 100 @spec fetch_elevation_profile(float(), float(), float(), float(), pos_integer()) :: {:ok, list(map())} | {:error, String.t()} def fetch_elevation_profile(lat1, lon1, lat2, lon2, n \\ 64) do + case srtm_tiles_dir() do + nil -> + fetch_elevation_profile_api(lat1, lon1, lat2, lon2, n) + + tiles_dir -> + case Srtm.fetch_elevation_profile(lat1, lon1, lat2, lon2, tiles_dir, n) do + {:ok, _profile} = ok -> ok + {:error, _} -> fetch_elevation_profile_api(lat1, lon1, lat2, lon2, n) + end + end + end + + defp fetch_elevation_profile_api(lat1, lon1, lat2, lon2, n) do pts = sample_path(lat1, lon1, lat2, lon2, n) dist_km = haversine_km(lat1, lon1, lat2, lon2) @@ -25,6 +40,8 @@ defmodule Microwaveprop.Terrain.ElevationClient do end end + defp srtm_tiles_dir, do: Application.get_env(:microwaveprop, :srtm_tiles_dir) + def sample_path(lat1, lon1, lat2, lon2, n) do for i <- 0..n do f = i / n diff --git a/lib/microwaveprop/terrain/srtm.ex b/lib/microwaveprop/terrain/srtm.ex new file mode 100644 index 00000000..d9cd86c0 --- /dev/null +++ b/lib/microwaveprop/terrain/srtm.ex @@ -0,0 +1,175 @@ +defmodule Microwaveprop.Terrain.Srtm do + @moduledoc false + + @samples 3601 + @void -32_768 + @base_url "https://elevation-tiles-prod.s3.amazonaws.com/skadi" + + @spec tile_filename(float(), float()) :: String.t() + def tile_filename(lat, lon) do + lat_floor = floor(lat) + lon_floor = floor(lon) + + lat_prefix = if lat_floor >= 0, do: "N", else: "S" + lon_prefix = if lon_floor >= 0, do: "E", else: "W" + + lat_str = lat_floor |> abs() |> Integer.to_string() |> String.pad_leading(2, "0") + lon_str = lon_floor |> abs() |> Integer.to_string() |> String.pad_leading(3, "0") + + "#{lat_prefix}#{lat_str}#{lon_prefix}#{lon_str}.hgt" + end + + @spec download_tile(float(), float(), String.t()) :: {:ok, String.t()} | {:error, term()} + def download_tile(lat, lon, tiles_dir) do + filename = tile_filename(lat, lon) + lat_dir = String.slice(filename, 0, 3) + url = "#{@base_url}/#{lat_dir}/#{filename}.gz" + path = Path.join(tiles_dir, filename) + + case Req.get(url, req_options()) do + {:ok, %{status: 200, body: body}} -> + decompressed = :zlib.gunzip(body) + File.write!(path, decompressed) + {:ok, path} + + {:ok, %{status: 404}} -> + {:error, :not_available} + + {:ok, %{status: status}} -> + {:error, "SRTM download HTTP #{status}"} + + {:error, reason} -> + {:error, "SRTM download error: #{inspect(reason)}"} + end + end + + @spec lookup(float(), float(), String.t()) :: + {:ok, integer()} | {:error, :no_tile} | {:error, :void} + def lookup(lat, lon, tiles_dir) do + path = Path.join(tiles_dir, tile_filename(lat, lon)) + + case :file.open(path, [:read, :binary, :raw]) do + {:ok, fd} -> + read_elevation(fd, lat, lon) + + {:error, :enoent} -> + if File.dir?(tiles_dir) do + case download_tile(lat, lon, tiles_dir) do + {:ok, _path} -> + case :file.open(path, [:read, :binary, :raw]) do + {:ok, fd} -> read_elevation(fd, lat, lon) + {:error, _} -> {:error, :no_tile} + end + + {:error, _} -> + {:error, :no_tile} + end + else + {:error, :no_tile} + end + end + end + + @spec fetch_elevation_profile(float(), float(), float(), float(), String.t(), pos_integer()) :: + {:ok, list(map())} | {:error, term()} + def fetch_elevation_profile(lat1, lon1, lat2, lon2, tiles_dir, n \\ 64) do + pts = sample_path(lat1, lon1, lat2, lon2, n) + dist_km = haversine_km(lat1, lon1, lat2, lon2) + + results = + Enum.reduce_while(pts, {:ok, []}, fn pt, {:ok, acc} -> + case lookup(pt.lat, pt.lon, tiles_dir) do + {:ok, elev} -> + entry = %{ + lat: pt.lat, + lon: pt.lon, + d: pt.d, + elev: elev, + dist_km: pt.d * dist_km + } + + {:cont, {:ok, acc ++ [entry]}} + + {:error, reason} -> + {:halt, {:error, reason}} + end + end) + + results + end + + defp read_elevation(fd, lat, lon) do + row = round((floor(lat) + 1 - lat) * (@samples - 1)) + col = round((lon - floor(lon)) * (@samples - 1)) + offset = (row * @samples + col) * 2 + + result = + case :file.pread(fd, offset, 2) do + {:ok, <>} when elev == @void -> + {:error, :void} + + {:ok, <>} -> + {:ok, elev} + + _ -> + {:error, :void} + end + + :file.close(fd) + result + end + + defp req_options do + defaults = [ + compressed: false, + decode_body: false, + retry: &retry?/2, + max_retries: 3, + retry_delay: &retry_delay/1 + ] + + overrides = Application.get_env(:microwaveprop, :srtm_req_options, []) + Keyword.merge(defaults, overrides) + end + + defp retry?(_request, response) do + case response do + %Req.Response{status: status} when status in [429, 500, 502, 503, 504] -> true + %{__exception__: true} -> true + _ -> false + end + end + + defp retry_delay(n) do + base = Integer.pow(2, n) * 1_000 + jitter = :rand.uniform(1_000) + base + jitter + end + + defp sample_path(lat1, lon1, lat2, lon2, n) do + for i <- 0..n do + f = i / n + + %{ + lat: lat1 + f * (lat2 - lat1), + lon: lon1 + f * (lon2 - lon1), + d: f + } + end + end + + defp haversine_km(lat1, lon1, lat2, lon2) do + dlat = deg_to_rad(lat2 - lat1) + dlon = deg_to_rad(lon2 - lon1) + rlat1 = deg_to_rad(lat1) + rlat2 = deg_to_rad(lat2) + + a = + :math.sin(dlat / 2) ** 2 + + :math.cos(rlat1) * :math.cos(rlat2) * :math.sin(dlon / 2) ** 2 + + 2 * 6371.0 * :math.asin(:math.sqrt(a)) + end + + defp deg_to_rad(deg), do: deg * :math.pi() / 180 +end diff --git a/test/microwaveprop/terrain/elevation_client_test.exs b/test/microwaveprop/terrain/elevation_client_test.exs index 37f56918..3c87ad83 100644 --- a/test/microwaveprop/terrain/elevation_client_test.exs +++ b/test/microwaveprop/terrain/elevation_client_test.exs @@ -111,3 +111,53 @@ defmodule Microwaveprop.Terrain.ElevationClientTest do end end end + +defmodule Microwaveprop.Terrain.ElevationClientSrtmTest do + use ExUnit.Case, async: false + + alias Microwaveprop.Terrain.ElevationClient + + @tiles_dir Path.expand("~/srtm/tiles") + @has_tiles File.exists?(Path.join(@tiles_dir, "N32W097.hgt")) + + setup do + prev = Application.get_env(:microwaveprop, :srtm_tiles_dir) + + on_exit(fn -> + if prev, + do: Application.put_env(:microwaveprop, :srtm_tiles_dir, prev), + else: Application.delete_env(:microwaveprop, :srtm_tiles_dir) + end) + + :ok + end + + describe "SRTM integration" do + @tag :srtm + @tag skip: !@has_tiles && "SRTM tiles not available" + test "uses local SRTM tiles when configured" do + Application.put_env(:microwaveprop, :srtm_tiles_dir, @tiles_dir) + + assert {:ok, profile} = + ElevationClient.fetch_elevation_profile(32.78, -96.8, 32.5, -96.5, 4) + + assert length(profile) == 5 + assert hd(profile).elev > 0 + end + + @tag :srtm + test "falls back to API when SRTM tiles missing" do + Application.put_env(:microwaveprop, :srtm_tiles_dir, "/nonexistent/tiles") + + Req.Test.stub(ElevationClient, fn conn -> + Req.Test.json(conn, %{"elevation" => [200.0, 250.0, 180.0]}) + end) + + assert {:ok, profile} = + ElevationClient.fetch_elevation_profile(32.9, -97.0, 30.3, -97.7, 2) + + assert length(profile) == 3 + assert hd(profile).elev == 200.0 + end + end +end diff --git a/test/microwaveprop/terrain/srtm_test.exs b/test/microwaveprop/terrain/srtm_test.exs new file mode 100644 index 00000000..680a98f6 --- /dev/null +++ b/test/microwaveprop/terrain/srtm_test.exs @@ -0,0 +1,180 @@ +defmodule Microwaveprop.Terrain.SrtmTest do + use ExUnit.Case, async: true + + alias Microwaveprop.Terrain.Srtm + + @tiles_dir Path.expand("~/srtm/tiles") + @has_tiles File.exists?(Path.join(@tiles_dir, "N32W097.hgt")) + + # A full SRTM tile: 3601 * 3601 samples, each 2 bytes big-endian signed + @samples 3601 + @test_elevation 150 + @synthetic_hgt :binary.copy(<<@test_elevation::signed-big-integer-size(16)>>, @samples * @samples) + + describe "tile_filename/2" do + test "northern hemisphere, western hemisphere" do + assert Srtm.tile_filename(32.78, -96.8) == "N32W097.hgt" + end + + test "southern hemisphere, eastern hemisphere" do + assert Srtm.tile_filename(-0.5, 10.3) == "S01E010.hgt" + end + + test "equator, prime meridian boundary" do + assert Srtm.tile_filename(0.5, -0.5) == "N00W001.hgt" + end + + test "exact integer coordinates" do + assert Srtm.tile_filename(32.0, -97.0) == "N32W097.hgt" + end + + test "negative lon floor wrapping" do + assert Srtm.tile_filename(35.2, -106.7) == "N35W107.hgt" + end + end + + describe "download_tile/3" do + test "downloads and decompresses tile on success" do + tmp_dir = Path.join(System.tmp_dir!(), "srtm_test_#{System.unique_integer([:positive])}") + File.mkdir_p!(tmp_dir) + + on_exit(fn -> File.rm_rf!(tmp_dir) end) + + gzipped = :zlib.gzip(@synthetic_hgt) + + Req.Test.stub(Srtm, fn conn -> + Plug.Conn.send_resp(conn, 200, gzipped) + end) + + assert {:ok, path} = Srtm.download_tile(32.0, -97.0, tmp_dir) + assert path == Path.join(tmp_dir, "N32W097.hgt") + assert File.exists?(path) + assert File.read!(path) == @synthetic_hgt + end + + test "returns {:error, :not_available} on 404" do + tmp_dir = Path.join(System.tmp_dir!(), "srtm_test_#{System.unique_integer([:positive])}") + File.mkdir_p!(tmp_dir) + + on_exit(fn -> File.rm_rf!(tmp_dir) end) + + Req.Test.stub(Srtm, fn conn -> + Plug.Conn.send_resp(conn, 404, "Not Found") + end) + + assert {:error, :not_available} = Srtm.download_tile(80.0, 0.0, tmp_dir) + end + + test "returns error message on other HTTP status" do + tmp_dir = Path.join(System.tmp_dir!(), "srtm_test_#{System.unique_integer([:positive])}") + File.mkdir_p!(tmp_dir) + + on_exit(fn -> File.rm_rf!(tmp_dir) end) + + Req.Test.stub(Srtm, fn conn -> + Plug.Conn.send_resp(conn, 500, "Internal Server Error") + end) + + assert {:error, "SRTM download HTTP 500"} = Srtm.download_tile(32.0, -97.0, tmp_dir) + end + end + + describe "lookup/3 with auto-download" do + test "auto-downloads missing tile and returns elevation" do + tmp_dir = Path.join(System.tmp_dir!(), "srtm_test_#{System.unique_integer([:positive])}") + File.mkdir_p!(tmp_dir) + + on_exit(fn -> File.rm_rf!(tmp_dir) end) + + gzipped = :zlib.gzip(@synthetic_hgt) + + Req.Test.stub(Srtm, fn conn -> + Plug.Conn.send_resp(conn, 200, gzipped) + end) + + assert {:ok, @test_elevation} = Srtm.lookup(32.0, -97.0, tmp_dir) + end + + test "returns {:error, :no_tile} when download fails" do + tmp_dir = Path.join(System.tmp_dir!(), "srtm_test_#{System.unique_integer([:positive])}") + File.mkdir_p!(tmp_dir) + + on_exit(fn -> File.rm_rf!(tmp_dir) end) + + Req.Test.stub(Srtm, fn conn -> + Plug.Conn.send_resp(conn, 404, "Not Found") + end) + + assert {:error, :no_tile} = Srtm.lookup(80.0, 0.0, tmp_dir) + end + end + + describe "lookup/3" do + @describetag :srtm + + @tag skip: !@has_tiles && "SRTM tiles not available at #{@tiles_dir}" + test "returns elevation for Dallas area" do + assert {:ok, elev} = Srtm.lookup(32.78, -96.8, @tiles_dir) + # Dallas is roughly 120-200m elevation + assert is_number(elev) + assert elev > 100 + assert elev < 300 + end + + @tag skip: !@has_tiles && "SRTM tiles not available at #{@tiles_dir}" + test "returns {:error, :no_tile} for missing tile" do + Req.Test.stub(Srtm, fn conn -> + Plug.Conn.send_resp(conn, 404, "Not Found") + end) + + assert {:error, :no_tile} = Srtm.lookup(80.0, 0.0, @tiles_dir) + end + end + + describe "fetch_elevation_profile/6" do + @describetag :srtm + + @tag skip: !@has_tiles && "SRTM tiles not available at #{@tiles_dir}" + test "returns profile with correct number of points" do + # Dallas to a point within the same tile + assert {:ok, profile} = + Srtm.fetch_elevation_profile(32.78, -96.8, 32.5, -96.5, @tiles_dir, 10) + + assert length(profile) == 11 + + first = hd(profile) + assert Map.has_key?(first, :lat) + assert Map.has_key?(first, :lon) + assert Map.has_key?(first, :d) + assert Map.has_key?(first, :elev) + assert Map.has_key?(first, :dist_km) + end + + @tag skip: !@has_tiles && "SRTM tiles not available at #{@tiles_dir}" + test "all elevations are positive for Dallas area path" do + assert {:ok, profile} = + Srtm.fetch_elevation_profile(32.78, -96.8, 32.5, -96.5, @tiles_dir, 8) + + Enum.each(profile, fn pt -> + assert pt.elev > 0, "Expected positive elevation, got #{pt.elev}" + end) + end + + @tag skip: !@has_tiles && "SRTM tiles not available at #{@tiles_dir}" + test "first and last dist_km are correct" do + assert {:ok, profile} = + Srtm.fetch_elevation_profile(32.78, -96.8, 32.5, -96.5, @tiles_dir, 4) + + first = hd(profile) + last = List.last(profile) + + assert_in_delta first.dist_km, 0.0, 0.01 + assert last.dist_km > 0 + end + + test "returns error when tiles dir does not exist" do + assert {:error, _} = + Srtm.fetch_elevation_profile(32.78, -96.8, 32.5, -96.5, "/nonexistent", 4) + end + end +end