pyart icon indicating copy to clipboard operation
pyart copied to clipboard

User defined sweep in det_sys_phase

Open joshua-wx opened this issue 4 years ago • 3 comments
trafficstars

Hi,

Our radar volumes in Australia now use a top-down scanning strategy, so the lowest sweep is last. When reading odimh5 files in pyart, it preserves this order. Most functions in pyart are not impacted (as far as I'm aware), expect for det_sys_phase, which defaults to the first sweep. I've made some changes to this function and the internal function it calls to permit user defined sweep. Let me know if this can be merged into pyart! See attached

def det_sys_phase(radar, ncp_lev=0.4, rhohv_lev=0.6, ncp_field=None, rhv_field=None, phidp_field=None, sweep=0): """ Determine the system phase.

Parameters
----------
radar : Radar
    Radar object for which to determine the system phase.
ncp_lev : float, optional
    Miminum normal coherent power level. Regions below this value will
    not be included in the phase calculation.
rhohv_lev : float, optional
    Miminum copolar coefficient level. Regions below this value will not
    be included in the phase calculation.
ncp_field, rhv_field, phidp_field : str, optional
    Field names within the radar object which represent the normal
    coherent power, the copolar coefficient, and the differential phase
    shift. A value of None for any of these parameters will use the
    default field name as defined in the Py-ART configuration file.
sweep : int, optional
    Define the sweep index from the radar object to use for detecting
    the system phase. Default to the first sweep.
Returns
-------
sys_phase : float or None
    Estimate of the system phase. None is not estimate can be made.

"""
# parse the field parameters
if ncp_field is None:
    ncp_field = get_field_name('normalized_coherent_power')
if rhv_field is None:
    rhv_field = get_field_name('cross_correlation_ratio')
if phidp_field is None:
    phidp_field = get_field_name('differential_phase')

ncp = radar.fields[ncp_field]['data'][:, 30:]
rhv = radar.fields[rhv_field]['data'][:, 30:]
phidp = radar.fields[phidp_field]['data'][:, 30:]
first_ray_idx = radar.sweep_start_ray_index['data'][sweep]
last_ray_idx = radar.sweep_end_ray_index['data'][sweep]
return _det_sys_phase(ncp, rhv, phidp, first_ray_idx, last_ray_idx, ncp_lev,
                      rhohv_lev)

def _det_sys_phase(ncp, rhv, phidp, first_ray_idx, last_ray_idx, ncp_lev=0.4, rhv_lev=0.6): """ Determine the system phase, see :py:func:det_sys_phase. """ good = False phases = [] for radial in range(first_ray_idx, last_ray_idx + 1): meteo = np.logical_and(ncp[radial, :] > ncp_lev, rhv[radial, :] > rhv_lev) mpts = np.where(meteo) if len(mpts[0]) > 25: good = True msmth_phidp = smooth_and_trim(phidp[radial, mpts[0]], 9) phases.append(msmth_phidp[0:25].min()) if not good: return None return np.median(phases)

Cheers, Joshua

joshua-wx avatar Jun 09 '21 06:06 joshua-wx

Hi @joshua-wx, yeah user defined sweep should work. If you could a PR would be much appreciated! Thank you!

zssherman avatar Jun 09 '21 18:06 zssherman

Hey @joshua-wx - any interest in adding this to Py-ART? Following up here :)

mgrover1 avatar Nov 22 '22 15:11 mgrover1