Brian Añano
Brian Añano

Reputation: 21

How to draw a circle or ellipse on a map (via matplotlib and Cartopy)?

I am trying to plot the track of a tropical cyclone (Typhoon Saola, 2023) using matplotlib and Cartopy. The best track data is from the Japan Meteorological Agency (JMA) and is stored in a CSV file. Plotting the track, including the line and center location points were okay, but I struggled with plotting the gale-force and storm-force wind fields which would require plotting a circle (or an ellipse if the wind field is asymmetrical) on a map. The problem is that the circles/ellipses do not appear on the map.

I tried using matplotlib.patches.Circle and matplotlib.patches.Ellipse for this task, wherein the center is given as lat-lon values of the tropical cyclone center and the radius/semi-major exis/semi-minor axis are provided in nautical miles based on the best track data. But the circle/ellipse are not plotted on the map. Only the tropical cyclone track, the gridlines and the coastlines are visible.

This is my code so far:

# Program for extracting TC track data from a comma-separated values (CSV) file
# esp. JMA best track data

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib.patches as patches
import cartopy.crs as ccrs

def angle(direction):
    if direction == 'N':
        return 0
    elif direction == 'NE':
        return 45
    elif direction == 'E':
        return 90
    elif direction == 'SE':
        return 135
    elif direction == 'S':
        return 180
    elif direction == 'SW':
        return 225
    elif direction == 'W':
        return 270
    elif direction == 'NW':
        return 315

# Importing TC best track data from a CSV file
saola_besttrack = pd.read_csv('C:/Users/Brian/Desktop/saola_track.csv')

# Setting the latitude and longitude values of the PAR region
PAR_longitudes = [120, 135, 135, 115, 115, 120, 120]
PAR_latitudes = [25, 25, 5, 5, 15, 21, 25]

for i in range(len(saola_besttrack)):
    # Taking a subset of the entire data from start to current
    subset = saola_besttrack.iloc[:i+1]
    
    # Taking additional subsets based on wind intensity
    TD_track = subset[subset['Wind (kt)'] <= 33]
    TS_track = subset[(subset['Wind (kt)'] >= 34) & (subset['Wind (kt)'] <= 47)]
    STS_track = subset[(subset['Wind (kt)'] >= 48) & (subset['Wind (kt)'] <= 63)]
    TY_track = subset[(subset['Wind (kt)'] >= 64) & (subset['Wind (kt)'] <= 99)]
    STY_track = subset[subset['Wind (kt)'] >= 100]
    
    # Initializing the map
    fig = plt.figure(figsize=(25., 25.), dpi=250)
    ax = plt.axes(projection=ccrs.PlateCarree())
    ax.set_extent([115, 130, 10, 25], ccrs.PlateCarree())
    
    # Adding gridlines
    gridlines = ax.gridlines(draw_labels=True, xlocs=np.arange(-180, 181, 2),
                             ylocs=np.arange(-90, 91, 2), color='gray',
                             linestyle='--')
    gridlines.top_labels = False
    gridlines.right_labels = False
    gridlines.xlabel_style = {'size': 20}
    gridlines.ylabel_style = {'size': 20}
    
    # Adding the coastlines
    ax.coastlines('10m', edgecolor='black', linewidth=2.5)
    
    # Plotting the PAR boundary
    ax.plot(PAR_longitudes, PAR_latitudes, color='black', linewidth=6,
            transform=ccrs.PlateCarree(), linestyle='-.')
    
    # Plotting the gale-force wind field
    if pd.isna(subset.loc[i, 'Direc. of Major Gale Axis']):
        pass
    elif subset.loc[i, 'Direc. of Major Gale Axis'] == 'symmetric':
        gale_windfield = patches.Circle((subset.loc[i, 'Long.'], subset.loc[i, 'Lat.']),
                                        subset.loc[i, 'Radius of Major Gale Axis (nm)'] * 1.852,
                                        edgecolor='yellow', facecolor='none', linestyle='-',
                                        transform=ccrs.PlateCarree())
        ax.add_patch(gale_windfield)
    else:
        major_axis_angle = angle(subset.loc[i, 'Direc. of Major Gale Axis'])
        gale_windfield = patches.Ellipse((subset.loc[i, 'Long.'], subset.loc[i, 'Lat.']),
                                         2 * subset.loc[i, 'Radius of Major Gale Axis (nm)'] * 1.852,
                                         2 * subset.loc[i, 'Radius of Minor Gale Axis (nm)'] * 1.852,
                                         angle=major_axis_angle, edgecolor='yellow',
                                         facecolor='none', linestyle='-',
                                         transform=ccrs.PlateCarree())
        ax.add_patch(gale_windfield)
    
    # Plotting the storm-force wind field
    if pd.isna(subset.loc[i, 'Direc. of Major Storm Axis']):
        pass
    elif subset.loc[i, 'Direc. of Major Storm Axis'] == 'symmetric':
        storm_windfield = patches.Circle((subset.loc[i, 'Long.'], subset.loc[i, 'Lat.']),
                                         subset.loc[i, 'Radius of Major Storm Axis (nm)'] * 1.852,
                                         edgecolor='red', facecolor='none', linestyle='-',
                                         transform=ccrs.PlateCarree())
        ax.add_patch(storm_windfield)
    else:
        major_axis_angle = angle(subset.loc[i, 'Direc. of Major Storm Axis'])
        storm_windfield = patches.Ellipse((subset.loc[i, 'Long.'], subset.loc[i, 'Lat.']),
                                          2 * subset.loc[i, 'Radius of Major Storm Axis (nm)'] * 1.852,
                                          2 * subset.loc[i, 'Radius of Minor Storm Axis (nm)'] * 1.852,
                                          angle=major_axis_angle, edgecolor='red',
                                          facecolor='none', linestyle='-',
                                          transform=ccrs.PlateCarree())
        ax.add_patch(storm_windfield)
    
    # Plotting the TC track
    ax.plot(subset['Long.'], subset['Lat.'], '-k', transform=ccrs.PlateCarree(),
            linewidth=5)
    ax.plot(TD_track['Long.'], TD_track['Lat.'], 'o', color='blue',
            transform=ccrs.PlateCarree(), markersize=25)
    ax.plot(TS_track['Long.'], TS_track['Lat.'], 'o', color='green',
            transform=ccrs.PlateCarree(), markersize=25)
    ax.plot(STS_track['Long.'], STS_track['Lat.'], 'o', color='orange',
            transform=ccrs.PlateCarree(), markersize=25)
    ax.plot(TY_track['Long.'], TY_track['Lat.'], 'o', color='red',
            transform=ccrs.PlateCarree(), markersize=25)
    ax.plot(STY_track['Long.'], STY_track['Lat.'], 'o', color='purple',
            transform=ccrs.PlateCarree(), markersize=25)
    
    plt.show()

# End of program

And this is the content of the file saola_track.csv:

Year,Month,Day,Hour,Lat.,Long.,Wind (kt),Gust (kt),Direc. of Major Storm Axis,Radius of Major Storm Axis (nm),Radius of Minor Storm Axis (nm),Direc. of Major Gale Axis,Radius of Major Gale Axis (nm),Radius of Minor Gale Axis (nm)

2023,8,22,0,18.8,128.3,0,0,,,,,,

2023,8,22,6,18.6,127.5,0,0,,,,,,

2023,8,22,12,18.5,126.8,0,0,,,,,,

2023,8,22,18,18.5,126.5,0,0,,,,,,

2023,8,23,0,18.5,126.2,0,0,,,,,,

2023,8,23,6,18.6,125.9,0,0,,,,,,

2023,8,23,12,18.8,125.6,0,0,,,,,,

2023,8,23,18,19.5,125.3,0,0,,,,,,

2023,8,24,0,20,125.1,0,0,,,,,,

2023,8,24,6,20.5,124.8,35,0,,,,SE,150,90

2023,8,24,12,20.2,124.5,35,0,,,,SE,150,90

2023,8,24,18,20.1,124.3,40,0,,,,SE,150,90

2023,8,25,0,19.9,124,45,0,,,,SE,150,90

2023,8,25,6,19.8,123.8,50,0,,,,SE,150,90

2023,8,25,12,19.6,123.6,65,0,symmetric,30,30,SE,150,90

2023,8,25,18,19.2,123.4,70,0,symmetric,40,40,SE,150,90

2023,8,26,0,18.5,123.2,75,0,symmetric,40,40,SE,150,90

2023,8,26,6,18,123.2,85,0,symmetric,50,50,SE,150,90

2023,8,26,12,17.6,123,90,0,symmetric,50,50,SE,150,90

2023,8,26,18,17.2,122.9,95,0,symmetric,50,50,SE,150,90

2023,8,27,0,16.8,122.9,95,0,symmetric,50,50,SE,150,90

2023,8,27,6,16.5,123,95,0,symmetric,50,50,SE,150,90

2023,8,27,12,16.2,123.2,95,0,symmetric,50,50,SE,150,90

2023,8,27,18,16,123.7,85,0,symmetric,40,40,SE,150,90

2023,8,28,0,16.8,124.4,80,0,symmetric,40,40,SE,150,90

2023,8,28,6,17.5,124.2,85,0,symmetric,40,40,symmetric,120,120

2023,8,28,12,17.9,124,85,0,symmetric,40,40,symmetric,120,120

2023,8,28,18,18.3,123.8,85,0,symmetric,40,40,symmetric,120,120

2023,8,29,0,18.5,123.5,85,0,symmetric,40,40,symmetric,120,120

2023,8,29,6,18.9,123.1,90,0,symmetric,40,40,symmetric,120,120

2023,8,29,12,19.3,122.7,95,0,symmetric,50,50,NW,150,120

2023,8,29,18,19.9,121.9,100,0,symmetric,50,50,NW,150,120

2023,8,30,0,20.1,121,105,0,symmetric,50,50,NW,150,120

2023,8,30,6,20.4,120.3,105,0,symmetric,50,50,NW,150,120

2023,8,30,12,20.7,119.6,105,0,symmetric,50,50,N,150,120

2023,8,30,18,20.9,118.7,105,0,symmetric,50,50,N,150,120

2023,8,31,0,21.1,118.2,105,0,symmetric,50,50,NE,150,120

2023,8,31,6,21.2,117.8,105,0,symmetric,50,50,NE,150,120

2023,8,31,12,21.5,117.3,100,0,symmetric,50,50,NE,150,120

2023,8,31,18,21.6,116.8,95,0,symmetric,50,50,NE,150,120

2023,9,1,0,21.9,116.2,95,0,symmetric,50,50,NE,150,120

Upvotes: 1

Views: 69

Answers (0)

Related Questions