Unidata / MetPy

MetPy is a collection of tools in Python for reading, visualizing and performing calculations with weather data.
https://unidata.github.io/MetPy/
BSD 3-Clause "New" or "Revised" License
1.24k stars 413 forks source link

Level2File only calculates slant range #1025

Open ytao579 opened 5 years ago

ytao579 commented 5 years ago

https://github.com/Unidata/MetPy/blob/2e7deea7798176bcbadf0860a2cbe6c392f8954d/examples/formats/NEXRAD_Level_2_File.py#L52

Hi, I was looking at this coord calculation in Metpy and comparing how pyart calculate the coordinates (x, y = display._get_x_y(sweep, True, True)), seems very different.

In Metpy example, the coords are calculated as (NEXRAD_Level_2_File):

_xlocs = ref_range * np.sin(np.deg2rad(az[:, np.newaxis]))
_ylocs_ = ref_range * np.cos(np.deg2rad(az[:, np.newaxis]))

However, in Pyart, it's a much more complex calculation: in its antenna_to_cartesian method (called by get_gate_x_y_z in radar.py), it is doing:

theta_e = elevations * np.pi / 180.0    # elevation angle in radians.
theta_a = azimuths * np.pi / 180.0      # azimuth angle in radians.
R = 6371.0 * 1000.0 * 4.0 / 3.0     # effective radius of earth in meters.
r = ranges * 1000.0                 # distances to gates in meters.

z = (r ** 2 + R ** 2 + 2.0 * r * R * np.sin(theta_e)) ** 0.5 - R
s = R * np.arcsin(r * np.cos(theta_e) / (R + z))  # arc length in m.
x = s * np.sin(theta_a)
y = s * np.cos(theta_a)
return x, y, z

And I was using both Metpy and Pyart example to plot a NEXRAD data file, the result are very different, can some one please help me to understand why they calculate the coordinates differently and which one is right?

Thanks, ytao

dopplershift commented 5 years ago

I've answered the question regarding the difference from PyART on StackOverflow.

I'm going to leave this open as a reminder that we should probably do better (when we clean up the rest of the NEXRAD reading code) and do the proper ground range calculation.