How to plot hurricane data on a map

13

As we're in the thick of hurricane season, I'm sure some of my fellow map-makers are wanting to create their own custom hurricane maps. In this blog post, I step you through the process, so you'll know where to find all the pieces, and how to put them together. I demonstrate using SAS software, but the same concepts can be applied to most any mapping software.

Before we get started, here's a picture to get you into the mood for hurricanes... Here's a wicket sky pic my friend Paula took at Ocracoke (one of the barrier islands off the NC coast), as Hurricane Dorian approached. The island was evacuated a short time after she took the photo! (Do you think this picture is beautiful, or scary?)

Hurricane Data

Hurricane Dorian is bearing down on the U.S. East Coast as we speak, therefore I'll use that as the data for my example. I went to the NOAA gis data download page, and found the Shapefile(s) for Dorian's cone of uncertainty, and forecast track.

I downloaded and extracted the zip file. Here's a list of all the files in the download - I'm particularly interested in the polygon (pgn) and line (lin) shapefiles.

I used the following code to import the polygon shapefile for the cone of uncertainty, and rename the X and Y variables to long and lat.

proc mapimport
   datafile='dorian_shapefile/al052019-040_5day_pgn.shp'
   out=cone;
id stormname;
run;
data cone; set cone (rename=(x=long y=lat));
run;

Checking the Data

I then used Proc SGplot to draw the polygon, as a quick sanity-check. (Yep - it's looking good!)

proc sgplot data=cone;
polygon y=lat x=long id=stormname;
run;

Plotting on a Map

It's all well & good to see the polygon for the cone of uncertainty plotted, but unless you have uncanny spatial perception of raw latitude and longitude values, you probably still don't know exactly where (in relation to the U.S. East Coast) Hurricane Dorian is forecast to go. New features in Proc SGmap make it easy to overlay the polygon data on a map, so let's do that!

proc sgmap mapdata=cone noautolegend;
esrimap url="http://services.arcgisonline.com/arcgis/rest/services/NatGeo_World_Map";
choromap / mapid=stormname lineattrs=(color=orange thickness=4px);
run;

The above example draws the cone of uncertainty as an empty polygon, but the choropleth map can also fill in the polygon with color - or more specifically, transparent color. You have to supply some response data for that to happen. I just make up some bogus response data to use.

data cone_attr;
stormname='Dorian';
colorvar='foo';
run;

proc sgmap mapdata=cone maprespdata=cone_attr noautolegend;
esrimap url="http://services.arcgisonline.com/arcgis/rest/services/NatGeo_World_Map";
choromap colorvar / mapid=stormname lineattrs=(color=orange thickness=4px) transparency=.7;
run;

Many of the hurricane forecast maps also show a line and markers in the middle of the cone, therefore let's add that. I import the line (_lin) shapefile, and overlay a series plot of this data on the map to represent the line. And I overlay a scatter plot to represent the markers.

proc mapimport
   datafile='dorian_shapefile/al052019-040_5day_lin.shp'
   out=line;
id stormname;
run;
data line; set line (rename=(x=long y=lat));
run;

proc sgmap mapdata=cone maprespdata=cone_attr plotdata=line noautolegend;
esrimap url="http://services.arcgisonline.com/arcgis/rest/services/NatGeo_World_Map";
choromap colorvar / mapid=stormname lineattrs=(color=orange thickness=4px) transparency=.7;
series x=long y=lat / lineattrs=(color=firebrick thickness=1px) smoothconnect;
scatter x=long y=lat / markerattrs=(color=firebrick size=20px symbol=circlefilled) transparency=.5;
run;

Finishing Touches

I've now got the map just about the way I want it ... except I would rather it also show a little more of the Bahamas and Florida Keys in the background map. By default, the amount of map shown is automatically chosen, based on the polygon/marker/line data being plotted on the map. One way to show more of the background map is to add a scatter plot marker in the area you want to also show on the map, and make it small/invisible.

data corners;
long_corner=-83.0; lat_corner=24.5; output;
run;
data plot_all; set line corners;
run;

proc sgmap mapdata=cone maprespdata=cone_attr plotdata=plot_all noautolegend;
esrimap url="http://services.arcgisonline.com/arcgis/rest/services/NatGeo_World_Map";
choromap colorvar / mapid=stormname lineattrs=(color=orange thickness=4px) transparency=.7;
series x=long y=lat / lineattrs=(color=firebrick thickness=1px) smoothconnect;
scatter x=long y=lat / markerattrs=(color=firebrick size=20px symbol=circlefilled) transparency=.5;
scatter x=long_corner y=lat_corner / markerattrs=(size=1px) transparency=1;
run;

Now that you know a few tricks about working with hurricane data, what ideas do you have for creating your own custom hurricane maps? Feel free to share them in the comments!


Hurricane Dorian: I'm going to Jacksonville.

U.S. Southeast coast: Jacksonville Florida or Jacksonville North Carolina?

Hurricane Dorian: Yes.

Here's a picture (from my friend Raven) of Dorian sliding past Jacksonville. Can you tell which Jacksonville? 🙂

 

Share

About Author

Robert Allison

The Graph Guy!

Robert has worked at SAS for over a quarter century, and his specialty is customizing graphs and maps - adding those little extra touches that help them answer your questions at a glance. His educational background is in Computer Science, and he holds a BS, MS, and PhD from NC State University.

13 Comments

  1. Hi Robert,

    For some reason I am unable to display the smooth connected line inside the cone. I tried running your code on SAS Studio with Base SAS 9.4M6 (SAS OnDemand for Academics) but the line completely is missing.

    What may be the reason?

      • It looks like the SERIES statement in the PROC SGMAP is being ignored for some reason.

        Here is my SAS code based on your code above :

        /******************************************/

        filename odsout '~/Hurricane';
        ODS _ALL_ CLOSE;

        proc mapimport
        datafile='/home/pchegoor/Hurricane/al052019-044A_5day_pgn.shp'
        out=cone;
        id stormname;
        run;
        data cone; set cone (rename=(x=long y=lat));
        run;

        data cone_attr;
        stormname='Dorian';
        colorvar='foo';
        run;

        proc mapimport
        datafile='/home/pchegoor/Hurricane/al052019-044A_5day_lin.shp'
        out=line;
        id stormname;
        run;
        data line; set line (rename=(x=long y=lat));
        run;

        ODS HTML path=odsout body="Hurricane.htm" ;

        ods graphics /
        noscale /* if you don't use this option, the text will be resized */
        imagemap tipmax=2500
        imagefmt=png imagename="Hurricane"
        width=800px height=600px;

        data corners;
        long_corner=-83.0; lat_corner=24.5; output;
        run;
        data plot_all; set line corners;
        run;

        proc sgmap mapdata=cone maprespdata=cone_attr plotdata=plot_all noautolegend;
        esrimap url="http://services.arcgisonline.com/arcgis/rest/services/NatGeo_World_Map";
        choromap colorvar / mapid=stormname lineattrs=(color=orange thickness=4px) transparency=.7;
        series x=long y=lat / lineattrs=(color=firebrick thickness=1px) smoothconnect;
        scatter x=long y=lat / markerattrs=(color=firebrick size=20px symbol=circlefilled) transparency=.5;
        scatter x=long_corner y=lat_corner / markerattrs=(size=1px) transparency=1;
        run;

        ODS HTML CLOSE;
        ODS LISTING;

        /***************************************************/

        The only warnings i see are as below but i don't think they are relateded to missing line:

        117
        118 proc sgmap mapdata=cone maprespdata=cone_attr plotdata=plot_all noautolegend;
        119 esrimap url="http://services.arcgisonline.com/arcgis/rest/services/NatGeo_World_Map";
        120 choromap colorvar / mapid=stormname lineattrs=(color=orange thickness=4px) transparency=.7;
        121 series x=long y=lat / lineattrs=(color=firebrick thickness=1px) smoothconnect;
        122 scatter x=long y=lat / markerattrs=(color=firebrick size=20px symbol=circlefilled) transparency=.5;
        123 scatter x=long_corner y=lat_corner / markerattrs=(size=1px) transparency=1;
        124 run;

        NOTE: Choromap boundaries drawn using LAT and LONG from map data set WORK.CONE.
        NOTE: PROCEDURE SGMAP used (Total process time):
        real time 5.49 seconds
        user cpu time 0.08 seconds
        system cpu time 0.01 seconds
        memory 8405.78k
        OS Memory 27704.00k
        Timestamp 09/04/2019 05:01:18 PM
        Step Count 63 Switch Count 1
        Page Faults 0
        Page Reclaims 2366
        Page Swaps 0
        Voluntary Context Switches 29
        Involuntary Context Switches 0
        Block Input Operations 0
        Block Output Operations 896

        WARNING: Some features may not be displayed on the map because of missing location information in the data.
        WARNING: Some features may not be displayed on the map because of missing location information in the data.
        NOTE: There were 1 observations read from the data set WORK.CONE_ATTR.
        NOTE: There were 9 observations read from the data set WORK.PLOT_ALL.

          • Robert Allison

            Try downloading the latest shapefile, and running this short test, and see if the series line works (modify the path to the shapefiles as needed).
            This works for me, in PC SAS 9.4m6, submitting the code in the DMS interface.

            proc mapimport
            datafile='chegoor/al052019-045_5day_lin.shp'
            out=line;
            id stormname;
            run;
            data line; set line (rename=(x=long y=lat));
            run;

            proc sgmap plotdata=line noautolegend;
            openstreetmap;
            series x=long y=lat / lineattrs=(color=firebrick thickness=1px) smoothconnect;
            scatter x=long y=lat / markerattrs=(color=firebrick size=20px symbol=circlefilled) transparency=.5;
            run;

  2. Here is my TEST Code :

    /****************************************************************/

    filename odsout '~/Hurricane';
    ODS _ALL_ CLOSE;
    proc mapimport
    datafile='~/Hurricane/al052019-045_5day_lin.shp'
    out=line;
    id stormname;
    run;
    data line; set line (rename=(x=long y=lat));
    run;

    ODS HTML path=odsout body="TEST.htm" ;

    ods graphics /
    noscale /* if you don't use this option, the text will be resized */
    imagemap tipmax=2500
    imagefmt=png imagename="TEST"
    width=800px height=600px;

    proc sgmap plotdata=line noautolegend;
    openstreetmap;
    series x=long y=lat / lineattrs=(color=firebrick thickness=1px) smoothconnect;
    scatter x=long y=lat / markerattrs=(color=firebrick size=20px symbol=circlefilled) transparency=.5;
    run;

    ODS HTML CLOSE;
    ODS LISTING;

    /*************************************************************************************/

    Here is my Output , still the line is missing :

Leave A Reply

Back to Top