Wednesday, May 15, 2013

Tutorial: Shapefile; A look at the DBF file.

This is the 2nd tutorial on using shapefiles in GrADS.  If this is your first time using shapefiles in GrADS, it is recommended you check out this tutorial.  This tutorial will dig a little deeper and look at the 'q dbf' command to look at the different shapefile properties.  For this tutorial, we will take a break from the usual focus on weather data, and instead focus on politics!  What we will do is draw the United States colored by which way each state voted in the 2013 presidential election.

Blue: Democrat, Red: Republican

So, to start you first need to download the required shapefile: Shapefile
Once you have this file put in your shapefile folder, we can begin.  Now we do need to open a file to get the environmental scaling, but that's all we need it for.  I'm going to use a simple GFS file from the GrADS Data Server.

    'sdfopen http://nomads.ncep.noaa.gov:9090/dods/gfs_hd/gfs_hd20130512/gfs_hd_00z'
    'set lat 22 50'
    'set lon -125 -68'
    'set clevs 10000000'
    'd pressfc/100'


Now, we have our file open and our map set up.  Now we are ready to plot our shapefiles.  However, we need to specify the color of each state.  This is where the 'q dfb' command.  The 'q dbf' takes a look a the database (.dbf) file that came in the shapefile zip folder, which is essentially the metadata.  So for this example, we would query the database file of the state shapefile.

    'q dbf Shapefiles/s_06se12.dbf'
   say result

What you will see is a list of information separated by commas.  What this information tells you, is the state abbreviation, shapefile id number, lat/lon, state name.  So, all we need to do, is extract this information and compare it to our previous knowledge of who won the election.  The hardest part of this is dealing with the fact that each value is separated by a comma instead of a space; we can't just use the 'subwrd()' command.  So, since we know that arrays can be set with strings as well as numbers, we can simply set a list of different colors, corresponding to each state abbreviation.

The first few numbers might look like:
    *2=Republican
    *4=Democrat

    col.AL=2
    col.AR=2
    col.AZ=2
    col.CA=4
    ....
    ....


Once, you have your list set, we simply need to loop through each shape within the shapefile, and determine its state abbreviation, then set the shapefile fill color to the color corresponding to the list above.

    i=2
    check=1
    while(check=1)
      'q dbf Shapefiles/s_06se12.dbf'
      line=sublin(result,i)
      if(line='' | line='(NULL)');check=0;else;

        com_count=1
        vars=1
        while(vars<=3)
          com_check=1
          len=1
          c=com_count
          while(com_check=1)
            comma=substr(line,c,1)
            if(comma=',');var.vars=substr(line,com_count,len-1);com_check=0;endif
            if(comma='' & vars=11);var.vars=substr(line,com_count,len-1);com_check=0;endif
            len=len+1
            c=c+1
         endwhile
         com_count=com_count+len-1
         vars=vars+1
       endwhile
       state=subwrd(var.3,1)
       'set shpopts 'col.state
       'draw shp Shapefiles/s_06se12.shp 'i-2
      endif
     i=i+1
   endwhile


This block of code may look complicated, but mostly it is just looking for commas within the string, and setting variables to the values in-between each comma.  You will also see the option "i-2" added on at the end of the 'drawshp' command.  This is telling GrADS to only plot the shapefile that has the id number i-2.  Since I starts at 2 (the first line of data in the shapefile), and since the first ID number is 0, we need to offset the pointer by 2.  So what this loop does, is loop through all of the shapes in the shapefile, and set the color to the list prescribed at the beginning of the script. 

If you run this script in GrADS, you will see a the shapes start to show up colored on the map.  So, that does it for this tutorial, hopefully it helped you work a little bit with the 'q dbf' command. 

Download Example Script Here

1 comment: