The GOES-R (Geostationary Operational Environmental Satellite) program images weather phenomena from a set of satellites in geostationary orbits. The GOES-16 and GOES-17 satellites are the first two of four planned GOES-R satellites.
This notebook provides an example of accessing GOES-16 data from blob storage on Azure, including (1) finding the data file corresponding to a date and time, (2) retrieving that file from blob storage, and (3) opening the downloaded file using the xarray library, and (4) rendering the image.
This dataset is stored in the West Europe Azure region, so this notebook will run most efficiently on Azure compute located in the same region. If you are using this data for environmental science applications, consider applying for an AI for Earth grant to support your compute requirements.
This dataset is documented at aka.ms/ai4edata-goes-r.
import os
import tempfile
import numpy as np
import shutil
import urllib
import matplotlib.pyplot as plt
import xarray
from azure.storage.blob import ContainerClient
temp_dir = os.path.join(tempfile.gettempdir(),'goes')
os.makedirs(temp_dir,exist_ok=True)
# Choose a specific product, day, and time to analyze
product = 'ABI-L2-MCMIPF'
syear = '2020'; sday = '002'; shour = '14'
# Storage locations are documented at http://aka.ms/ai4edata-goes16
storage_account_name = 'goeseuwest'
container_name = 'noaa-goes16'
storage_account_url = 'https://' + storage_account_name + '.blob.core.windows.net'
goes_blob_root = storage_account_url + '/' + container_name + '/'
# Create a ContainerClient to enumerate blobs
goes_container_client = ContainerClient(account_url=storage_account_url,
container_name=container_name,
credential=None)
# Download a URL to a temporary file
def download_url(url):
url_as_filename = url.replace('://', '_').replace('/', '_')
destination_filename = os.path.join(temp_dir,url_as_filename)
urllib.request.urlretrieve(url, destination_filename)
return destination_filename
prefix = product + '/' + syear + '/' + sday + '/' + shour + '/'
print('Finding blobs matching prefix: {}'.format(prefix))
generator = goes_container_client.list_blobs(name_starts_with=prefix)
blobs = []
for blob in generator:
blobs.append(blob.name)
print('Found {} images'.format(len(blobs)))
# There will be several scans this hour, we'll take the first
image_index = 0
url = goes_blob_root + blobs[image_index]
print('Processing image:\n{}'.format(url))
Finding blobs matching prefix: ABI-L2-MCMIPF/2020/002/14/ Found 6 images Processing image: https://goeseuwest.blob.core.windows.net/noaa-goes16/ABI-L2-MCMIPF/2020/002/14/OR_ABI-L2-MCMIPF-M6_G16_s20200021400218_e20200021409537_c20200021410039.nc
# GOES-16 MCMIPF files are ~300MB. Not too big to fit in memory, so sometimes it may be
# preferable to download to file first, sometimes it will be better to load straight to
# memory.
download_to_file = True
if download_to_file:
filename = download_url(url)
dataset = xarray.open_dataset(filename)
else:
import netCDF4
import requests
response = requests.get(url)
nc4_ds = netCDF4.Dataset(os.path.basename(url), memory=response.content)
store = xarray.backends.NetCDF4DataStore(nc4_ds)
dataset = xarray.open_dataset(store)
print('Scan starts at: {}'.format(dataset.time_coverage_start))
print('Scan ends at: {}'.format(dataset.time_coverage_end))
Scan starts at: 2020-01-02T14:00:21.8Z Scan ends at: 2020-01-02T14:09:53.7Z
# Bands are documented at:
#
# https://www.ncdc.noaa.gov/data-access/satellite-data/goes-r-series-satellites/glossary
#
# We'll use the red/"veggie"/blue bands with wavelengths 0.64, 0.86, and 0.47, respectively.
#
# This is close enough to RGB for today, but there's a great tutorial on getting closer to
# true color (and doing other fancy rendering tricks with GOES data!) here:
#
# https://unidata.github.io/python-gallery/examples/mapping_GOES16_TrueColor.html
#
r = dataset['CMI_C02'].data; r = np.clip(r, 0, 1)
g = dataset['CMI_C03'].data; g = np.clip(g, 0, 1)
b = dataset['CMI_C01'].data; b = np.clip(b, 0, 1)
# Brighten the image a bit for to look more stylish
gamma = 2.5; r = np.power(r, 1/gamma); g = np.power(g, 1/gamma); b = np.power(b, 1/gamma)
# Create a single RGB image for plotting
rgb = np.dstack((r, g, b))
fig = plt.figure(figsize=(7.5, 7.5), dpi=100)
# This definitely looks slicker with fancy borders on, at the cost of some extra
# imports.
show_fancy_borders = True
if not show_fancy_borders:
plt.imshow(rgb); ax = plt.gca(); ax.axis('off');
else:
import metpy
import cartopy.crs as ccrs
# Pick an arbitrary channel to get the x/y coordinates and projection information
# associated with the scan
dummy_channel = dataset.metpy.parse_cf('CMI_C01')
x = dummy_channel.x; y = dummy_channel.y
ax = fig.add_subplot(1, 1, 1, projection=dummy_channel.metpy.cartopy_crs)
ax.imshow(rgb, origin='upper', extent=(x.min(), x.max(), y.min(), y.max()))
ax.coastlines(resolution='50m', color='black')
ax.add_feature(ccrs.cartopy.feature.BORDERS);
shutil.rmtree(temp_dir)