Using Siphon to get NEXRAD Level 3 data from a TDS

First, point to the top-level thredds radar server accessor to find what datasets are available.



In [1]:

    
from siphon.radarserver import RadarServer, get_radarserver_datasets
ds = get_radarserver_datasets('http://thredds.ucar.edu/thredds/')
print(list(ds.keys()))









    



['NEXRAD Level III Radar for Case Study CCS039', 'TDWR Level III Radar from IDD', 'NEXRAD Level II Radar for Case Study CCS039', 'NEXRAD Level II Radar from IDD', 'NEXRAD Level III Radar from IDD']

Now create an instance of RadarServer to point to the appropriate radar server access URL. This is pulled from the catalog reference url.



In [2]:

    
url = ds['NEXRAD Level III Radar from IDD'].follow().catalog_url
rs = RadarServer(url)

Look at the variables available in this dataset



In [3]:

    
rs.variables









    Out[3]:





{'DAA',
 'DHR',
 'DOD',
 'DPA',
 'DPR',
 'DSD',
 'DSP',
 'DTA',
 'DU3',
 'DU6',
 'DVL',
 'EET',
 'HHC',
 'N0C',
 'N0H',
 'N0K',
 'N0M',
 'N0Q',
 'N0R',
 'N0S',
 'N0U',
 'N0V',
 'N0X',
 'N0Z',
 'N1C',
 'N1H',
 'N1K',
 'N1M',
 'N1P',
 'N1Q',
 'N1S',
 'N1U',
 'N1X',
 'N2C',
 'N2H',
 'N2K',
 'N2M',
 'N2Q',
 'N2S',
 'N2U',
 'N2X',
 'N3C',
 'N3H',
 'N3K',
 'N3M',
 'N3Q',
 'N3S',
 'N3U',
 'N3X',
 'NAC',
 'NAH',
 'NAK',
 'NAM',
 'NAQ',
 'NAU',
 'NAX',
 'NBC',
 'NBH',
 'NBK',
 'NBM',
 'NBQ',
 'NBU',
 'NBX',
 'NCR',
 'NET',
 'NMD',
 'NST',
 'NTP',
 'NVL',
 'NVW',
 'OHA',
 'PTA'}

Create a new query object to help request the data. Using the chaining methods, ask for data from radar FTG (Denver) for now for the product N0Q, which is reflectivity data for the lowest tilt. We see that when the query is represented as a string, it shows the encoded URL.



In [4]:

    
from datetime import datetime
query = rs.query()
query.stations('FTG').time(datetime.utcnow()).variables('N0Q')









    Out[4]:





var=N0Q&time=2015-07-15T22%3A50%3A49.691861&stn=FTG

We can use the RadarServer instance to check our query, to make sure we have required parameters and that we have chosen valid station(s) and variable(s)



In [5]:

    
rs.validate_query(query)









    Out[5]:





True

Make the request, which returns an instance of TDSCatalog. This handles parsing the catalog



In [6]:

    
catalog = rs.get_catalog(query)

We can look at the datasets on the catalog to see what data we found by the query. We find one NIDS file in the return



In [7]:

    
catalog.datasets









    Out[7]:





{'Level3_FTG_N0Q_20150715_2247.nids': <siphon.catalog.Dataset at 0x10b831e48>}

We can pull that dataset out of the dictionary and look at the available access URLs. We see URLs for OPeNDAP, CDMRemote, and HTTPServer (direct download).



In [8]:

    
ds = list(catalog.datasets.values())[0]
ds.access_urls









    Out[8]:





{'CdmRemote': 'http://thredds.ucar.edu/thredds/cdmremote/nexrad/level3/IDD/N0Q/FTG/20150715/Level3_FTG_N0Q_20150715_2247.nids',
 'HTTPServer': 'http://thredds.ucar.edu/thredds/fileServer/nexrad/level3/IDD/N0Q/FTG/20150715/Level3_FTG_N0Q_20150715_2247.nids',
 'OPENDAP': 'http://thredds.ucar.edu/thredds/dodsC/nexrad/level3/IDD/N0Q/FTG/20150715/Level3_FTG_N0Q_20150715_2247.nids'}

We'll use the CDMRemote reader in Siphon and pass it the appropriate access URL.



In [9]:

    
from siphon.cdmr import Dataset
data = Dataset(ds.access_urls['CdmRemote'])

The CDMRemote reader provides an interface that is almost identical to the usual python NetCDF interface. We pull out the variables we need for azimuth and range, as well as the data itself.



In [10]:

    
rng = data.variables['gate'][:] / 1000.
az = data.variables['azimuth'][:]
ref = data.variables['BaseReflectivityDR'][:]

Then convert the polar coordinates to Cartesian



In [11]:

    
import numpy as np
x = rng * np.sin(np.deg2rad(az))[:, None]
y = rng * np.cos(np.deg2rad(az))[:, None]
ref = np.ma.array(ref, mask=np.isnan(ref))

Finally, we plot them up using matplotlib.



In [12]:

    
import matplotlib.pyplot as plt
fig, ax = plt.subplots(1, 1, figsize=(9, 8))
ax.pcolormesh(x, y, ref)
ax.set_aspect('equal', 'datalim')
ax.set_xlim(-460, 460)
ax.set_ylim(-460, 460)









    














    











    Out[12]:





(-460, 460)