NOAA NHC Wind Speed Probabilities#

Demonstrate the use of geoJSON and shapefile data with PlotGeometry in MetPy’s simplified plotting interface. This example walks through plotting cities, along with 5-day tropical-storm-force wind speed probabilities from NOAA National Hurricane Center.

```import geopandas

from metpy.cbook import get_test_data
from metpy.plots import MapPanel, PanelContainer, PlotGeometry
```

Read in the shapefile file containing the wind probabilities.

```wind_data = geopandas.read_file(get_test_data('nhc_wind_prob_2021082012.zip'))
```

Add the color scheme to the GeoDataFrame. This is the same color scheme used by the National Hurricane Center for their wind speed probability plots.

```wind_data['fill'] = ['none', '#008B00', '#00CD00', '#7FFF00', '#FFFF00', '#FFD700',
'#CD8500', '#FF7F00', '#CD0000', '#8B0000', '#8B008B']
wind_data
```
PERCENTAGE geometry fill
0 <5% MULTIPOLYGON (((-148.82395 16.55476, -148.8239... none
1 5-10% MULTIPOLYGON (((-150.35245 19.02733, -150.3524... #008B00
2 10-20% MULTIPOLYGON (((-150.53228 19.61176, -150.5322... #00CD00
3 20-30% MULTIPOLYGON (((-144.55314 19.88150, -144.5531... #7FFF00
4 30-40% MULTIPOLYGON (((-144.19349 19.97141, -144.1934... #FFFF00
5 40-50% MULTIPOLYGON (((-144.37331 20.10628, -144.3733... #FFD700
6 50-60% MULTIPOLYGON (((-144.23845 20.19619, -144.2384... #CD8500
7 60-70% MULTIPOLYGON (((-142.88977 20.42097, -142.8897... #FF7F00
8 70-80% MULTIPOLYGON (((-144.01367 20.33106, -144.0136... #CD0000
9 80-90% MULTIPOLYGON (((-93.43823 19.47689, -93.43823 ... #8B0000
10 >90% MULTIPOLYGON (((-144.05862 20.46592, -144.0586... #8B008B

Read in the shapefile file containing the cities.

```cities = geopandas.read_file(get_test_data('us_cities.zip'))
```

There are thousands of cities in the United States. We choose a few cities here that we want to display on our plot.

```cities = cities.loc[
((cities['NAME'] == 'Myrtle Beach') & (cities['STATE'] == 'SC'))
| ((cities['NAME'] == 'Hatteras') & (cities['STATE'] == 'NC'))
| ((cities['NAME'] == 'Ocean City') & (cities['STATE'] == 'MD'))
| ((cities['NAME'] == 'New York') & (cities['STATE'] == 'NY'))
| ((cities['NAME'] == 'Nantucket') & (cities['STATE'] == 'MA'))
| ((cities['NAME'] == 'Portland') & (cities['STATE'] == 'ME'))
]
cities
```
GNIS_ID ANSICODE FEATURE FEATURE2 NAME POP_2010 COUNTY COUNTYFIPS STATE STATE_FIPS LATITUDE LONGITUDE PopPlLat PopPlLong ELEV_IN_M ELEV_IN_FT geometry
327 573692.0 582683 Civil County Seat Portland 66194.0 Cumberland 005 ME 23 43.661471 -70.255326 43.661471 -70.255326 19.0 62.0 POINT (-70.25533 43.66147)
4913 616716.0 2378175 Census County Seat Nantucket 7446.0 Nantucket 019 MA 25 41.283456 -70.099461 41.283456 -70.099461 7.0 23.0 POINT (-70.09946 41.28346)
6259 586284.0 2391341 Civil -999 Ocean City 7102.0 Worcester 047 MD 24 38.336503 -75.084906 38.336503 -75.084906 1.0 3.0 POINT (-75.08491 38.33650)
8329 -999.0 2395220 Civil -999 New York 8175133.0 New York 061 NY 36 40.761493 -73.981431 40.663931 -73.938270 32.0 105.0 POINT (-73.98143 40.76149)
14717 1249770.0 2404346 Civil -999 Myrtle Beach 27109.0 Horry 051 SC 45 33.689060 -78.886694 33.689060 -78.886694 8.0 26.0 POINT (-78.88669 33.68906)
25784 1020636.0 2628633 Census -999 Hatteras 504.0 Dare 055 NC 37 35.219342 -75.690161 35.219342 -75.690161 0.0 0.0 POINT (-75.69016 35.21934)

Make sure that both GeoDataFrames have the same coordinate reference system (CRS).

```cities = cities.to_crs(wind_data.crs)
```

We want to find out what the probability of tropical-storm-force winds is for each of the cities we selected above. Geopandas provides a spatial join method, which merges the two GeoDataFrames and can tell us which wind speed probability polygon each of our city points lies within. That information is stored in the ‘PERCENTAGE’ column below.

```cities = geopandas.sjoin(cities, wind_data, how='left', op='within')
cities
```
```/opt/hostedtoolcache/Python/3.11.8/x64/lib/python3.11/site-packages/sphinx_gallery/gen_gallery.py:244: FutureWarning: The `op` parameter is deprecated and will be removed in a future release. Please use the `predicate` parameter instead.
return 0.0, func()
```
GNIS_ID ANSICODE FEATURE FEATURE2 NAME POP_2010 COUNTY COUNTYFIPS STATE STATE_FIPS LATITUDE LONGITUDE PopPlLat PopPlLong ELEV_IN_M ELEV_IN_FT geometry index_right PERCENTAGE fill
327 573692.0 582683 Civil County Seat Portland 66194.0 Cumberland 005 ME 23 43.661471 -70.255326 43.661471 -70.255326 19.0 62.0 POINT (-70.25533 43.66147) 5 40-50% #FFD700
4913 616716.0 2378175 Census County Seat Nantucket 7446.0 Nantucket 019 MA 25 41.283456 -70.099461 41.283456 -70.099461 7.0 23.0 POINT (-70.09946 41.28346) 8 70-80% #CD0000
6259 586284.0 2391341 Civil -999 Ocean City 7102.0 Worcester 047 MD 24 38.336503 -75.084906 38.336503 -75.084906 1.0 3.0 POINT (-75.08491 38.33650) 3 20-30% #7FFF00
8329 -999.0 2395220 Civil -999 New York 8175133.0 New York 061 NY 36 40.761493 -73.981431 40.663931 -73.938270 32.0 105.0 POINT (-73.98143 40.76149) 5 40-50% #FFD700
14717 1249770.0 2404346 Civil -999 Myrtle Beach 27109.0 Horry 051 SC 45 33.689060 -78.886694 33.689060 -78.886694 8.0 26.0 POINT (-78.88669 33.68906) 0 <5% none
25784 1020636.0 2628633 Census -999 Hatteras 504.0 Dare 055 NC 37 35.219342 -75.690161 35.219342 -75.690161 0.0 0.0 POINT (-75.69016 35.21934) 2 10-20% #00CD00

Plot the wind speed probability polygons from the ‘geometry’ column. Use the ‘fill’ column we created above as the fill colors for the polygons, and set the stroke color to ‘none’ for all of the polygons.

```wind_geo = PlotGeometry()
wind_geo.geometry = wind_data['geometry']
wind_geo.fill = wind_data['fill']
wind_geo.stroke = 'none'
```

Plot the cities from the ‘geometry’ column, marked with diamonds (‘D’). Label each point with the name of the city, and it’s probability of tropical-storm-force winds on the line below. Points are set to plot in white and the font color is set to black.

```city_geo = PlotGeometry()
city_geo.geometry = cities['geometry']
city_geo.marker = 'D'
city_geo.labels = cities['NAME'] + '\n(' + cities['PERCENTAGE'] + ')'
city_geo.fill = 'white'
city_geo.label_facecolor = 'black'
```

Add the geometry plots to a panel and container. Finally, we are left with a complete plot of wind speed probabilities, along with some select cities and their specific probabilities.

```panel = MapPanel()
panel.title = 'NHC 5-Day Tropical-Storm-Force Wind Probabilities (Valid 12z Aug 20 2021)'
panel.plots = [wind_geo, city_geo]
panel.area = [-90, -52, 27, 48]
panel.projection = 'mer'
panel.layers = ['lakes', 'land', 'ocean', 'states', 'coastline', 'borders']

pc = PanelContainer()
pc.size = (12, 10)
pc.panels = [panel]
pc.show()
```

Total running time of the script: (0 minutes 12.843 seconds)

Gallery generated by Sphinx-Gallery