BUG: Level 2 VAD calculation differs from Level 3 data

[nawendt]

PyART VAD winds speeds appear much too low when reading Level 2 data. When I look at the wind speeds with height from a Level 3 VAD table I get:

array([12., 14., 16., 19., ..., 89., 58., 82.])

These data are in kts. Using PyART on the Level 2 data from IEM and conversion to kts I get:

array([6.09 , 6.31, 6.22, ..., 0.33, 0.12, 0.072])

The altitudes that I am requesting are the same as those present in the Level 3 data and are as follows:

array([487.68, 548.64, 609.6, 670.56, 883.92, 914.4, 1127.76, 1219.2, 1341.12, 1524.0, 1676.4, 1828.8, 2042.16, 2133.6, 2438.4, 2499.36, 2743.2, 3048.0, 3078.48, 3352.8, 3657.6, 3749.04, 3962.4, 4267.2, 4572.0, 4572.0, 4876.8, 5181.6, 5486.4, 5608.32, 5791.2, 6096.0, 6705.6, 7315.2, 7620.0, 9144.0])

The differences become much more exaggerated with height. I have attached a Level 2 (KTLX_*) and Level 3 (sn.*) files from the same time in order to help reproduce the issue.

ktlx_042821_2218.zip

[scollis]

Hey @nawendt we are looking at this now

[scollis]

Z, aliased Z and corrected Z

[zssherman]

@nawendt Do you have the full level 3 vad wind field? We are trying to use that for testing purposes, and trying to read the sn file itself with no luck.

[scollis]

For one, I am concerned about the applicability of VAD here given the convective nature of the scene... We are making a synthetic volume to test

[dopplershift]

Well regardless of the applicability of VAD, if it’s not known why it’s differing from NEXRAD’s results that’s still disconcerting. (or an opportunity to improve the test suite :grin:)

[scollis]

Oh. huge plus one @dopplershift ... @zssherman is going to take the level three winds and make a synthetic radar volume and retrieve winds from that.. To understand the root cause.

[nawendt]

Here is the decoded level 3 data for you in npz format. sn.0029.zip

[zssherman]

Thank you @nawendt !

[zssherman]

So, with the original VAD code, i'm still having issues pinpointing the difference. However, when using Jonathan's VAD code, based on Browning et al, the winds I retrieve in KT are: 20.3240570708959, 20.451149033017053, 20.213930672243716, 19.842817484939356, 16.331657163186538, 15.342313715314917, 12.281523221306688, 12.805500394234128, 13.78001289503375, 13.723269145443385, 20.60172418223752, 17.135698426627375, 13.610807946294253, 13.329828841880865, 15.557544007023141, 16.109913564188926, 16.842267878346817, 17.88028905850967, 17.316152623027715, 23.63088252110223, 14.508157458549444, 15.469422356211929, 19.51005555495896, 17.427313966621313, 23.555126852720747, 23.555126852720747, 30.921264806768647, --, 39.92789815342219, 40.70061709014852, 40.71902491127103, 46.43913105668294, 53.24009637713336, 54.39377835258477, 56.63506496121095, 76.19654012824995

Which are more reasonable, I will continue to try to see what going on what the current VAD code, but back up plan might be to revert to Jonathan's VAD code, i'm still digging into the current VAD code.

[zssherman]

@nawendt In the meantime, I'm adding another vad suite that has more reasonable results. I am still working on a fix for the original vad code as I've narrowed down the problem I believe, but still haven't got a fix working yet. From what I'm seeing is every other sweep in the radar data you provided has masked velocity sweeps, which I believe is throwing the code off, but the fixes I've tried so far don't seem to work, so still working on it. #997

[zssherman]

Let me know if this helps in the meantime. I was using corrected velocity, just trying to get something usable for you for now as I continue trying to fix the problem in the original vad code.

[nawendt]

Thanks for this. So far, things look OK. I'll have a better idea when there are more clouds/weather in the vicinity of the radars I am looking at.

[nawendt]

After looking at a longer time series of VADs I processed from university radars, I was noticing some much higher values than you would expect. I also was using corrected velocity based on dealias_unwrap_phase. Not all values were uniformly high. Some were close, if not a little lower, than nearby WSR-88Ds. This could just be not correcting the velocity properly on my end. I'll keep looking. I attached a file that was particularly off. The requested heights are np.linspace(170,12000,75).

KULM_20210612_111704.gz

[mgrover1]

@nawendt Do you know what method the NWS uses here?

I did some testing with this, and found a similar trend where the Level2 data is lesser in magnitude compared to the Level3 VAD.

I used the following function to calculate the VADs here:

def calculuate_vad(radar, vel_field, percentile=50, method='browning'):
    for idx in range(radar.nsweeps):
        radar_1sweep = radar.extract_sweeps([idx])
        try:
            if method=='michelson':
                vad = pyart.retrieve.vad_michelson(
            radar_1sweep, vel_field, z_want=zlevels
            )

            elif method=='browning':
                vad = pyart.retrieve.vad_browning(
            radar_1sweep, vel_field, z_want=zlevels
            )
            else:
                print('Not a valid method! Use browning or michelson')
            u_allsweeps.append(vad.u_wind)
            v_allsweeps.append(vad.v_wind)
        except ValueError:
            continue
    # Average U and V over all sweeps and compute magnitude and angle
    u_avg = np.nanpercentile(np.array(u_allsweeps), percentile, axis=0)
    v_avg = np.nanpercentile(np.array(v_allsweeps), percentile, axis=0)
    direction = np.rad2deg(np.arctan2(-u_avg, -v_avg)) % 360
    speed = np.sqrt(u_avg**2 + v_avg**2)
    return speed, direction

Using the 90th percentile (percentile = 90) results in the following:

Which is much closer, but still off in some areas.

[nawendt]

Thanks for looking at it. I'm not sure if the exact algorithm, but I'll get into contact with someone who can answer that.

[nawendt]

I was able to get into contact with some folks from the radar operations center. They use an algorithm developed by Rabin and Zrnic (1980). That method gets supplemented by Chrisman and Smith (2009) in order to increase the number of valid wind retrievals.

[mgrover1]

We fixed this with #1428 - it had to do with how we were masking data.