1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
|
#!/usr/bin/env python3
from __future__ import annotations
import typing
import argparse
import matplotlib.pyplot as mpl
import numpy as np
import pymap3d as pm
p = argparse.ArgumentParser()
p.add_argument("alt_m", help="altitude [meters]", type=float, default=0.0, nargs="?")
args = p.parse_args()
lat, lon = np.meshgrid(np.arange(-90, 90, 0.1), np.arange(-180, 180, 0.2))
x, y, z = pm.geodetic2ecef(lat, lon, args.alt_m)
def panel(ax, val, name: str, cmap: str | None = None):
hi = ax.pcolormesh(lon, lat, val, cmap=cmap)
ax.set_title(name)
fg.colorbar(hi, ax=ax).set_label(name + " [m]")
ax.set_xlabel("longitude [deg]")
fg = mpl.figure(figsize=(16, 5))
axs: typing.Any = fg.subplots(1, 3, sharey=True)
fg.suptitle("geodetic2ecef")
panel(axs[0], x, "x", "bwr")
panel(axs[1], y, "y", "bwr")
panel(axs[2], z, "z", "bwr")
axs[0].set_ylabel("latitude [deg]")
fg = mpl.figure(figsize=(16, 5))
axs = fg.subplots(1, 3, sharey=True)
fg.suptitle(r"|$\nabla$ geodetic2ecef|")
panel(axs[0], np.hypot(*np.gradient(x)), r"|$\nabla$ x|")
panel(axs[1], np.hypot(*np.gradient(y)), r"|$\nabla$ y|")
panel(axs[2], np.hypot(*np.gradient(z)), r"|$\nabla$ z|")
axs[0].set_ylabel("latitude [deg]")
mpl.show()
|
