-
Notifications
You must be signed in to change notification settings - Fork 282
Description
I have LiDAR data that I'd like to convert to radar format (because I want to use pyart package) and then grid it. However, after gridding, I noticed that the ds.rsw values are NaN. Could you please advise on how to retrieve valid values for ds.rsw?
import pyart
import numpy as np
from datetime import datetime
from netCDF4 import Dataset
import warnings
warnings.filterwarnings('ignore')
Load data lidar data
file_path = "D:/lidar/ncfiles/WCS000248_2023-09-23_09-50-36_ppi_40_100m.nc"
data = Dataset(file_path)
sweep_group = data.groups["Sweep_860123"]
time = sweep_group.variables["time"][:]
latitude = data.variables["latitude"][:]
longitude = data.variables["longitude"][:]
altitude = data.variables["altitude"][:]
azimuth = sweep_group.variables["azimuth"][:]
elevation = sweep_group.variables["elevation"][:]
range_ = sweep_group.variables["range"][:]
radial_wind_speed = sweep_group.variables["radial_wind_speed"][:]
rsw = np.array(radial_wind_speed)
Create radar object using pyart
radar = pyart.testing.make_empty_ppi_radar(rsw.shape[1], len(azimuth), 1)
radar.latitude['data'] = np.array([latitude])
radar.longitude['data'] = np.array([longitude])
radar.altitude['data'] = np.array([altitude])
#radar.time['data'] = np.array([(t - time_converted[0]).total_seconds() for t in time_converted])
radar.time = {
'standard_name': 'time',
'long_name': 'time in seconds since volume start',
'calendar': 'gregorian',
'units': 'seconds since 2023-09-23T04:35:05Z',
'comment': 'times are relative to the volume start_time',
'data': np.array([(t - time_converted[0]).total_seconds() for t in time_converted]),
'FillValue':1e+20
}
radar.azimuth['data'] = azimuth
radar.elevation['data'] = elevation
radar.range['data'] = range
radial_wind_speed_dict = {
'long_name': 'radial_wind_speed',
'standard_name': 'radial_wind_speed_of_scatterers_away_from_instrument',
'units': 'm/s',
'sampling_ratio': 1.0,
'_FillValue': -9999 ,
'grid_mapping': 'grid_mapping',
'coordinates': 'time range',
'data': np.ma.masked_invalid(rsw) # Mask invalid data
}
radar.fields = { 'rws': radial_wind_speed_dict}
#success plot this:
import matplotlib.pyplot as plt
from pyart.graph import RadarDisplay
display = RadarDisplay(radar)
fig, ax = plt.subplots(figsize=(10, 8))
display.plot_ppi("rws", sweep=0, ax=ax, cmap="coolwarm")
plt.show()
Now grid
grid_limits = ((10., 4000.), (-4500., 4500.), (-4500., 4500.))
grid_shape = (20, 50, 50)
grid = pyart.map.grid_from_radars([radar], grid_limits=grid_limits, grid_shape=grid_shape)
ds = grid_dv.to_xarray()
print(ds)
print(ds)
<xarray.Dataset> Size: 641kB
Dimensions: (time: 1, z: 20, y: 50, x: 50, nradar: 1)
Coordinates: (12/16)
- time (time) object 8B 2023-09-23 04:35:05
- z (z) float64 160B 10.0 220.0 ... 3.79e+03 4e+03
lat (y, x) float64 20kB 45.02 45.02 ... 45.1 45.1
lon (y, x) float64 20kB 7.603 7.605 ... 7.715 7.717 - y (y) float64 400B -4.5e+03 -4.316e+03 ... 4.5e+03
- x (x) float64 400B -4.5e+03 -4.316e+03 ... 4.5e+03
...
origin_altitude (time) float64 8B nan
radar_altitude (nradar) float64 8B nan
radar_latitude (nradar) float64 8B 45.06
radar_longitude (nradar) float64 8B 7.66
radar_time (nradar) int64 8B 0
radar_name (nradar) <U10 40B 'fake_radar'
Dimensions without coordinates: nradar
Data variables:
rws (time, z, y, x) float64 400kB nan nan ... nan
ROI (time, z, y, x) float32 200kB 500.0 ... 500.0
Attributes:
radar_name: fake_radar
nradar: 1
instrument_name: fake_radar
np.nanmax(ds.rws)
Out[18]: nan
np.nanmin(ds.rws)
Out[19]: nan