Как я могу генерировать карту плитки из радарных данных Nexrad 3 -го уровня? ⇐ Python

[Anonymous]

Я пытаюсь сгенерировать скользкие плитки карты (например, для использования с листочками или открытыми слоями), используя данные RADAR 3 уровня 3 (например, файлы. NIDS или .NXL3, такие как N0Q, N0B и т. Д.). Моя цель состоит в том, чтобы отобразить эти радарные продукты на интерактивной веб -карте, аналогично тому, как такие услуги, как RadarScope или Radaromega, показывают радиолокационные наложения. Я использовал Pyart и Matplotlib и попытался использовать PNG из этого для создания плиток, но это не сработало. < /P>
код я попробовал < /p>
import os
import subprocess
import numpy as np
import pyart
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
from PIL import Image
import rasterio
from rasterio.transform import from_bounds

# CONFIG
RADAR_FILE = "Level3_IND_N0B_20250710_0228.nids"
PNG_FILE = "radar_kind.png"
TIF_FILE = "radar_kind.tif"
TILE_DIR = "tiles_kind"
ZOOM_RANGE = "5-6"

print("1 ▶ Reading radar file...")
radar = pyart.io.read_nexrad_level3(RADAR_FILE)

# Get radar gate coverage bounds
print("2 ▶ Calculating radar gate bounds...")
lons = radar.gate_longitude['data'].ravel()
lats = radar.gate_latitude['data'].ravel()
valid = np.isfinite(lons) & np.isfinite(lats)
lons = lons[valid]
lats = lats[valid]

min_lon, max_lon = lons.min(), lons.max()
min_lat, max_lat = lats.min(), lats.max()
pad = 0.02
min_lon -= pad
max_lon += pad
min_lat -= pad
max_lat += pad

# Create map with CartoPy using actual lat/lon extent
print("3 ▶ Plotting radar on map using CartoPy...")
fig = plt.figure(figsize=(16, 16), dpi=150)
ax = plt.axes(projection=ccrs.PlateCarree())
ax.set_extent([min_lon, max_lon, min_lat, max_lat], crs=ccrs.PlateCarree())

display = pyart.graph.RadarMapDisplay(radar)
display.plot_ppi_map(
'reflectivity',
sweep=0,
vmin=-32, vmax=64,
cmap="NWSRef",
ax=ax,
projection=ccrs.PlateCarree(),
resolution='50m',
colorbar_flag=False,
title_flag=False
)

# Remove axes and background
ax.axis('off')
ax.set_facecolor('none')
fig.patch.set_alpha(0.0)

plt.savefig(PNG_FILE, bbox_inches='tight', pad_inches=0, transparent=True)
plt.close()

# Convert to GeoTIFF using the exact bounds used in the plot
print("4 ▶ Converting PNG to GeoTIFF with correct bounds...")
img = Image.open(PNG_FILE)
width, height = img.size

transform = from_bounds(min_lon, min_lat, max_lon, max_lat, width, height)

r, g, b, a = img.convert("RGBA").split()
r_np = np.array(r)
g_np = np.array(g)
b_np = np.array(b)
a_np = np.array(a)

with rasterio.open(
TIF_FILE, 'w',
driver='GTiff',
height=height,
width=width,
count=4,
dtype='uint8',
crs='EPSG:4326',
transform=transform
) as dst:
dst.write(r_np, 1)
dst.write(g_np, 2)
dst.write(b_np, 3)
dst.write(a_np, 4)

# Generate tiles using gdal2tiles
print("5 ▶ Generating tiles...")
if os.path.isdir(TILE_DIR):
subprocess.run(["rm", "-rf", TILE_DIR], check=True)

if subprocess.run(["which", "gdal2tiles.py"], capture_output=True).returncode == 0:
cmd = ["gdal2tiles.py", "-p", "raster", "-z", ZOOM_RANGE, "-w", "none", TIF_FILE, TILE_DIR]
else:
cmd = ["python3", "-m", "gdal2tiles", "-p", "raster", "-z", ZOOM_RANGE, "-w", "none", TIF_FILE, TILE_DIR]

subprocess.run(cmd, check=True)
print(f" Tiles generated in ./{TILE_DIR}/")


Подробнее здесь: https://stackoverflow.com/questions/796 ... radar-data