I have example code on the Weather Radar github demonstrating how to obtain National Weather Service radar images and merge them with basemaps using the GeoTiler library.

For a more friendly, step-by-step walkthrough, let's breakdown the code!

I thought it might be nice to break that code down and go through it step by step!

Libraries

The Weather Radar uses several libraries, which you'll need to download first.

• Geotiler, to construct basemaps.
• datetime and pytz for date time and time zone support.
• requests, for http handling.
• Pillow and numpy for image handling and processing.
• xmltodict, to handle XML parsing.
• Adafruit Blinka, to help run CircuitPython libraries on the RaspberryPi.
• Adafruit RGB display, to drive the TFT display (which in this case is the ili9341).

import time
from datetime import datetime
import pytz

import os
import json
import math
import requests
import logging

from PIL import Image, ImageDraw, ImageFont
from io import BytesIO
import xmltodict
import geotiler
import numpy as np

#Adafruit & CircuitPython libraries
import board
import digitalio
import adafruit_rgb_display.ili9341 as ili9341

#Secrets! (openweather API & lat long coordinates)
from secrets import secrets

Secrets

You'll also need to populate the secrets.py file:

• The NWS api requires a header with User-Agent and contact email to be sent with every request, so these need to be filled in.
• Coordinates (in latitude, longitude format) should be of the point of interest. These are used to get the closest radar station as well as local warnings & alerts.
• If you know it, add the radar station ID as a fallback in case the coordinates don't work.

  secrets = {
      'header' : {
                  'User-Agent': 'DIY radar viewer',
                  'Contact': 'EMAIL_ADDRESS'
                  },
      'coordinates' : (47.168598,-123.559299), #lat, long
      'station': '' #station ID fallback, e.g. klgx
  }

Initialise Display

Here we define some of the pins and settings for the TFT screen, and initialise the display (check out some of the examples on the Adafruit RGB display library).

CURR_DIR = f"{os.path.dirname(__file__)}/"

cs_pin = digitalio.DigitalInOut(board.CE0)
dc_pin = digitalio.DigitalInOut(board.D25)
reset_pin = digitalio.DigitalInOut(board.D24)
BAUDRATE = 24000000
spi = board.SPI()
disp = ili9341.ILI9341(
    spi,
    rotation=270,
    cs=cs_pin,
    dc=dc_pin,
    rst=reset_pin,
    baudrate=BAUDRATE,
)

Now for some functions and settings that we'll run only once during start up.

Get radar Station ID

We'll use the location_to_station function to take the lat long coordinates from the secrets file and make a request to api.weather.gov for that particular point. In return we'll get the assigned radar station ID for that location, the location city, state, timezone, and forecast URL. 

def location_to_station():
    ''' Get the ID of the nearest radar station from lat long coordinates.
        Also nearest population centre, state, time zone, and forecast URLs

        Returns the station code! (must be lowercase)
    '''
    global station
    global forecast_url
    global timeZone

    ######################################
    # Try to get the lat long point file #
    ######################################
    point_url = f"https://api.weather.gov/points/{lat_long[0]},{lat_long[1]}"

    try:
        response = requests.get(point_url, headers=headers, timeout=5)
    except:
        print("Connection Problems getting Lat/Long point data!")
        response = False

    if response:
        point_file = json.load(BytesIO(response.content))
        point_file = point_file['properties']

        station = point_file['radarStation'].lower()
        location_city = point_file['relativeLocation']['properties']['city']
        location_state = point_file['relativeLocation']['properties']['state']
        timeZone = point_file['timeZone']
        forecast_url = point_file['forecast']

        print(f"{location_city} ({location_state}) -- {station.upper()} -- {timeZone}")
    else:
        print(f"Unable to get location or radar ({response})")
        get_station_data(secrets["station"]) #Use a fallback if it doesn't work!

    return station

Bounding Coordinates

Now with a station code, we can use that and the radar layer we're interested in (in this case, 'bohp' which is is the 1 hour precipitation layer) to make a request for the WMS GetCapabilities xml document from the NWS (e.g. this one for Jacksonville, FL). This will give us the SW and NE bounding coordinates that we'll use later to make basemaps and download the actual radar images.

Here's  a page of the OGC compliant layers from the NWS, including alerts, warnings, and layers for the 200+ weather radar stations).

lat_long = secrets['coordinates']
headers = secrets['header']
layer = 'bohp'

location_to_station()

capabilities_url = f'https://opengeo.ncep.noaa.gov:443/geoserver/{station}/{station}_{layer}/wms?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetCapabilities'

#Get the SW and NE coordinates from the WMS GetCapabilities file
minx, miny, maxx, maxy = get_bounding_coordinates(capabilities_url)

def get_bounding_coordinates(url):
    ''' Get and return the bounding coordinates of a WMS layer.
        Param url: The url to the GetCapabilities xml

        Returns minx, miny, maxx, maxy

    '''

    ################################
    # Get the GetCapabilities file #
    ################################
    try:
        response = requests.get(url, headers=headers, timeout=5 )
    except:
        print("Connection Problems: Getting Bounding coordinates")
        response = False

    ################################
    # Get the bounding coordinates #
    ################################
    if response:
        xml_file = BytesIO(response.content)
        capabilties_dict = xmltodict.parse(xml_file.read()) #Parse the XML into a dictionary
        bounding_coordinates = capabilties_dict["WMS_Capabilities"]["Capability"]["Layer"]["EX_GeographicBoundingBox"]

        minx = bounding_coordinates['westBoundLongitude']
        maxx = bounding_coordinates['eastBoundLongitude']
        miny = bounding_coordinates['southBoundLatitude']
        maxy = bounding_coordinates['northBoundLatitude']
        return minx, miny, maxx, maxy
    else:
        print(f"Couldn't get GetCapabilities file ({response})")
        return 0, 0, 0, 0