Assign GeoCoords#

This notebook demonstrates how to work with geocoordinates (longitude and latitude) instead of radar-centric x and y coordinates.

[1]:

import cmweather  # noqa: F401
import matplotlib.pyplot as plt
import xarray as xr
from open_radar_data import DATASETS

import xradar as xd

Define Functions#

[2]:

def get_geocoords(ds):
    """
    Converts Cartesian coordinates (x, y, z) in a radar dataset to geographic
    coordinates (longitude, latitude, altitude) using CRS transformation.

    Parameters
    ----------
    ds : xarray.Dataset
        Radar dataset with Cartesian coordinates.

    Returns
    -------
    xarray.Dataset
        Dataset with added 'lon', 'lat', and 'alt' coordinates and their attributes.
    """
    from pyproj import CRS, Transformer

    # Convert the dataset to georeferenced coordinates
    ds = ds.xradar.georeference()
    # Define source and target coordinate reference systems (CRS)
    src_crs = ds.xradar.get_crs()
    trg_crs = CRS.from_user_input(4326)  # EPSG:4326 (WGS 84)
    # Create a transformer for coordinate conversion
    transformer = Transformer.from_crs(src_crs, trg_crs)
    # Transform x, y, z coordinates to latitude, longitude, and altitude
    trg_y, trg_x, trg_z = transformer.transform(ds.x, ds.y, ds.z)
    # Assign new coordinates with appropriate attributes
    ds = ds.assign_coords(
        {
            "lon": (ds.x.dims, trg_x, xd.model.get_longitude_attrs()),
            "lat": (ds.y.dims, trg_y, xd.model.get_latitude_attrs()),
            "alt": (ds.z.dims, trg_z, xd.model.get_altitude_attrs()),
        }
    )
    return ds


def fix_sitecoords(ds):
    coords = ["longitude", "latitude", "altitude", "altitude_agl"]
    for coord in coords:
        # Compute median excluding NaN
        data = ds[coord].median(skipna=True).item()
        attrs = ds[coord].attrs if coord in ds else {}
        ds = ds.assign_coords({coord: xr.DataArray(data=data, attrs=attrs)})
    return ds

Load Data#

[3]:

file1 = DATASETS.fetch("cfrad.20080604_002217_000_SPOL_v36_SUR.nc")
file2 = DATASETS.fetch("cfrad.20211011_201557.188_to_20211011_201617.720_DOW8_PPI.nc")

Downloading file 'cfrad.20080604_002217_000_SPOL_v36_SUR.nc' from 'https://github.com/openradar/open-radar-data/raw/main/data/cfrad.20080604_002217_000_SPOL_v36_SUR.nc' to '/home/docs/.cache/open-radar-data'.
Downloading file 'cfrad.20211011_201557.188_to_20211011_201617.720_DOW8_PPI.nc' from 'https://github.com/openradar/open-radar-data/raw/main/data/cfrad.20211011_201557.188_to_20211011_201617.720_DOW8_PPI.nc' to '/home/docs/.cache/open-radar-data'.

Example #2#

[11]:

dtree2 = xd.io.open_cfradial1_datatree(file2)

[12]:

try:
    dtree2 = dtree2.xradar.georeference()
except Exception:
    print("Georeferencing failed!")

Georeferencing failed!

[13]:

print("Longitude", dtree2["sweep_0"]["longitude"].shape)
print("Latitude", dtree2["sweep_0"]["latitude"].shape)

Longitude (731,)
Latitude (731,)

Important Note:

This radar data is from a mobile research radar called Doppler on Wheels (DOW), and its site coordinates (latitude, longitude) often vary slightly during operation, as can be seen from the shape ((731)) of the data in the above cell, while we expect it to be of unity shapes or empty, i.e., (1) or (). As a result, multiple site coordinate values can exist, creating a challenge for assigning consistent x, y, z or lat, lon, and alt coordinates using the current georeferencing system in xradar. To address this, a custom function like fix_sitecoords (defined above) can be created, leveraging the map_over_sweeps function to standardize the site coordinates.

Assign GeoCoords#

Define Functions#

Load Data#

Example #1#

Assign lat, lon, and alt#

Example #2#

Fix Coords#