In [1]:
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.axes_grid1 import make_axes_locatable
import re
import cartopy.crs as ccrs
import wget
import subprocess
# Larger figure size
fig_size = [14, 10]
plt.rcParams['figure.figsize'] = fig_size
Utility functions to dowload IONEX files, get TEC and plot TEC maps. See below for examples of use.
In [2]:
def parse_map(tecmap, exponent = -1):
tecmap = re.split('.*END OF TEC MAP', tecmap)[0]
return np.stack([np.fromstring(l, sep=' ') for l in re.split('.*LAT/LON1/LON2/DLON/H\\n',tecmap)[1:]])*10**exponent
def get_tecmaps(filename):
with open(filename) as f:
ionex = f.read()
return [parse_map(t) for t in ionex.split('START OF TEC MAP')[1:]]
def get_tec(tecmap, lat, lon):
i = round((87.5 - lat)*(tecmap.shape[0]-1)/(2*87.5))
j = round((180 + lon)*(tecmap.shape[1]-1)/360)
return tecmap[i,j]
def ionex_filename(year, day, centre, zipped = True):
return '{}g{:03d}0.{:02d}i{}'.format(centre, day, year % 100, '.Z' if zipped else '')
def ionex_ftp_path(year, day, centre):
return 'ftp://cddis.gsfc.nasa.gov/gnss/products/ionex/{:04d}/{:03d}/{}'.format(year, day, ionex_filename(year, day, centre))
def ionex_local_path(year, day, centre = 'esa', directory = '/tmp', zipped = False):
return directory + '/' + ionex_filename(year, day, centre, zipped)
def download_ionex(year, day, centre = 'esa', output_dir = '/tmp'):
wget.download(ionex_ftp_path(year, day, centre), output_dir)
subprocess.call(['gzip', '-d', ionex_local_path(year, day, centre, output_dir, zipped = True)])
def plot_tec_map(tecmap):
proj = ccrs.PlateCarree()
f, ax = plt.subplots(1, 1, subplot_kw=dict(projection=proj))
ax.coastlines()
h = plt.imshow(tecmap, cmap='viridis', vmin=0, vmax=100, extent = (-180, 180, -87.5, 87.5), transform=proj)
plt.title('VTEC map')
divider = make_axes_locatable(ax)
ax_cb = divider.new_horizontal(size='5%', pad=0.1, axes_class=plt.Axes)
f.add_axes(ax_cb)
cb = plt.colorbar(h, cax=ax_cb)
plt.rc('text', usetex=True)
cb.set_label('TECU ($10^{16} \\mathrm{el}/\\mathrm{m}^2$)')
This downloads all the daily IONEX files for the year 2017.
This takes a while to run, since it downloads one file at a time.
In [3]:
for day in range(1, 366):
download_ionex(2017, day)
Example 1
Plot all the TEC maps for a given day.
In [4]:
for tecmap in get_tecmaps(ionex_local_path(2017, 349)):
plot_tec_map(tecmap)
plt.show()