Sunday, March 15, 2015

Post 3: Data Collection

                                            Goals and Objectives

This exercise is the first step in the sand mining semester project for University of Wisconsin-Eau Claire's Geographic Information Systems II course. The exercise involves downloading data for the project from several internet sources, including the USGS, USDA, USDOT, as well as county-based files.

In addition to importing the geospatial files for the sand mine project, this exercise will also help familiarize one about manipulating the data in ArcGIS; for example, by joining tables, projecting data into the appropriate coordinate system, building and designing a geodatabase to store the data, and checking all data for accuracy.

The location with which the geospatial data for the exercise relates is Trempealeau County in west-central Wisconsin (fig.1).
 
Figure 1. Location of Trempealeau County in west-central Wisconsin (red).

                                                     Methods


                                               Data Acquisition

The first section of the exercise involved acquiring data from different sources. Six data sets were imported from four websites:
  • USGS:
  1. National land cover data (NLCD) that differentiates the general types of land found in Trempealeau County; e.g. cropland, open water, hay/pasture land, etc. (fig. 2).
  2. Digital elevation model (DEM) which shows the differences in elevation throughout the county (fig. 3).
 


Figure 2. National land cover data map created for Trempealeau County, Wisconsin.


Figure 3.  Digital elevation model created for Trempealeau County, Wisconsin.
Data Source:
http://nationalmap.gov/viewer.html

  • USDA:
     3. National Agricultural Statistics Survey (NASS) is similar to the NCLD data from the USGS in
         that it too differentiates the land use for whatever area of interest (AOI) it covers (fig. 4). How-
         ever, NASS data differentiates land use based on crop land, and is thus more detailed than the
        generalized  NLCD from the USGS site.
     4. Web soil survey displays the soil map units (MU) for Trempealeau County (fig. 5).


Figure 4. National agricultural statistics survey map created for Trempealeau County, Wisconsin.  


Figure 5. Soil index map created for Trempealeau County, Wisconsin. No legend included due to the numerous amount of data.

Data Source:
NASS http://datagateway.nrcs.usda.gov/
Soil Survey http://www.nrcs.usda.gov/wps/portal/nrcs/site/soils/home/

  • USDOT:
       5. The extent of railroad systems in Trempealeau County is shown in this data (fig. 6).

Data Source:
http://www.rita.dot.gov/bts/sites/rita.dot.gov.bts/files/publications/national_transportation_atlas_database/index.html

  • Trempealeau County Land Records (TCLR):
      6. Geospatial data from this website included a geodatabase with many feature classes that are
          specific to Trempealeau County; for example, county boundaries, roads and highways, etc.

Figure 6. Railroad map compiled for Trempealeau County, Wisconsin.

Data Source:
http://www.tremplocounty.com/landrecords/


                                         Data Manipulation

Once data were imported into the appropriate location from the various sources  discussed above, it was manipulated using Python Scripting (PS; see Post 2: Python Scripting for more details). Such manipulation was needed in order to create the maps displayed in figures 2-6 above. Manipulation of the datasets included using PS software included:
  • Projecting data into the appropriate coordinate system (i.e. NAD_1983_HARN_WISCRS_Trempealeau_County_Feet)
  • Extracting raster data to the Trempealeau County boundary
  • Importing raster data into the appropriate geodatabase
In addition to the data manipulation using PS software, several datasets were altered using ArcGIS tools. For instance, clipping national railroad data to the Trempealeau County boundary and joining soil table data from the USDA with other tables previously downloaded with the Trempealeau County geodatabase.

                                             Data Accuracy

An assessment of data accuracy was conducted on the datasets used for the project (table 1). While challenging, it is important to examine all relevant geospatial metadata to ensure that is accurate, precise, and up to date. Some of the information in table 1 is directly dependent upon what is contained in a datasets metadata, e.g. temporal accuracy, scale, and minimum mapping unit. However, several of the fields in table 1 can be calculated indirectly. For instance, planimetric coordinate accuracy (PCA) could be calculated using an equation derived from an Excel graph (fig. 7). Once the PCA is determined, the effective resolution could be determined; being half the PCA.



Figure 7. Graph showing the formula used to derive planimetric coordinate accuracy using known PCA data and the inverse scale data.




Table 1.
 

 Summary dataset accuracy from the five sources used on the project thus far.

                                            

                                                Conclusions

Determining the accuracy of the datasets was one of the most tedious steps of the exercise, but also one of the most important.  By using the most accurate data that is required by a given project we are ensuring that our products are accurate as well. The accuracy portion of the exercise worked well because it forced the user to investigate the metadata from various data sources.

The exercise also did a good job of exposing the user to new software like Python Scripting. While Python Scripting is rather tedious and cumbersome to work with as a beginner, it definitely has potential as far as a data-processing tool. For example, creating a data-processing loop for rastersets was useful because once a user identifies a few parameters in the scripting code, like file location and definition of data, all valid raster data will be processed in succession.

       

  










Post 2: Python Scripting (3 Posts)

Pyscripting 1: March 15, 2015                                                     

                                                Introduction


The purpose of this exercise is to expose the student to Python Scripting (PS) as it pertains to tools found in ArcGIS.

According to ArcGIS Help, PS is a free programming language that can be integrated into a number of computer software programs, like ArcGIS products. Python can be used to create batch-processing of several datasets at one time. For example, in this Python scripting (PS) exercise several raster datasets were processed in succession using PS in conjunction with ArcMap tools, like projection and extract.

In addition to its capability to process numerous datasets at once, PS is a powerful tool because it is a used by numerous users of varying skill levels. As such, PS tutorials and information are found at several websites. In fact, ArcGIS Help pages for Geodatabase toolsets usually have a section that gives an example of PS for each tool (fig. 1).

Following are several websites that are useful for those interested in learning more about Python Scripting:
                              https://wiki.python.org/moin/BeginnersGuide/NonProgrammers





Figure 1. Python Scripting code sample for loading raster sets to a geodatabase. From AcrGIS Help 10.1 (http://resources.arcgis.com/en/help/main/10.1/index.html#//001200000025000000).


                                                    Methods


In order to create the PS for this exercise (fig. 2), several steps were followed:

  • Setting  the environment defined the workspace (i.e. folders) where the data to be processed originated from:
#set environment
env.workspace ="Q:\StudentCoursework\CHupy\g2155.GEOG.337.001\XXX\Ex_5\Tremp_raster2"
arcpy.env.overwriteOutput = True

  • ListRasters returns a list of rasters available in the workspace (defined above) based on name and/or raster type (e.g. tiff, jpeg, etc.):
#Get and print a list all rasters in workspace
listofrasters = arcpy.ListRasters()
for raster in listofrasters:
    print(raster)

  • Before projecting the rasters, they needed to be defined:
#set up output variable
    rasterOut = '{}_Out.tif'.format(raster)
    rasterExtract = '{}_Extract.tif'.format(raster)

  • Rasters were projected into the appropriate coordinate system:
#project raster
    arcpy.ProjectRaster_management(raster, rasterOut, "Q:\StudentCoursework\CHupy\g2155.GEOG.337.001\XXX\Lab\Ex_5\Tremp_raster2\TMP.gdb\Fema\S_BASE_INDEX")

  • Rasters were extracted to the appropriate extents and saved:
#extract rasters
    outExtractByMask = ExtractByMask(rasterOut, "Q:\StudentCoursework\CHupy\g2155.GEOG.337.001XXX\Lab\Ex_5\Tremp_raster2\TMP.gdb\Boundaries\County_Boundary")

  #save extract
    outExtractByMask.save(rasterExtract)

  • The projected and extracted rasters were sent to the appropriate database:
#onvert rasters to TMP geodatabase
    arcpy.RasterToGeodatabase_conversion(outExtractByMask, "Q:\StudentCoursework\CHupy\g2155.GEOG.337.001\XXXLab\Ex_5\Tremp_raster2\TMP.gdb")







Figure 2. Screen-shot of the Python Script created for the exercise.

                                                      Conclusion

Once the rasters were projected, extracted, and saved to the appropriate geodatabase, they were used as part of the semester project for my geography 337 course (see Post 3: 'Data Collection').


 

Pyscripting 2: April 21, 2015

                                                   Introduction


 The Python Scripting assignment in this section correlates to Post 5: Network Analysis. Post 5 deals with the routing of frac sand from mining facilities in Wisconsin to railroad depots for shipment across the country.
  However, not all mining facilities that were geocoded in Post 4 are active. As such, a Python Script (fig. 1) will be written to only select mines that are active and save them as a layer. The Python Script will also create a layer which is comprised of feature classes that are categorized as mines, as opposed to facilities that are classified as processing facilities, for example. Finally, the script will create a feature layer of all the active mine facilities that have no railroads within 1.5 km of them.


Figure 1. Python Script for creating feature classes for active mining facilities at a distance of greater than 1.5 km from railroad depots for network analysis activity in Post 5.

                                                      Methods

  • System modules were imported and the environment was set in the Python Scripting program:
 
#import system modules
import arcpy
from arcpy import env

#set environment
env.workspace= "Q:\StudentCoursework\CHupy\g2155.GEOG.337.001\LUCZAKJN\Lab\ex7\ex7.gdb"

arcpy.env.overwriteOutput = True

  • Variables used in the script were defined:
 
#set variables
all = "all_mines"
active = "active_mines"
status = "status_mines"
nrmines = "mines_norail"
wi = "wi"
rails = "rails_wtm"
finalnr = "mines_norail_final"

  • SQL statements were established to select:
Active Mines:
#write SQL statements to select all active mines
activeSQL = "Status" + " = " + "'Active'"

Mining Facilities:
#write SQL statement to select "mines"
minesSQL = "Facility_T" + " LIKE " + "'%Mine%'"

Mines without a rail depot:
#write SQL statement to select all mines without a rail depot
nrminesSQL = " NOT " + "Facility_T" + " LIKE " + "'%Rail%'"

  • Temporary feature layer is created for each SQL above (i.e. activeSQL, minesSQL, and nrminesSQL:

Active:
#select mines that are active
arcpy.MakeFeatureLayer_management(all, active, activeSQL)

Facility:
#select mines that are active with "mine" facility field
arcpy.MakeFeatureLayer_management(active, status, minesSQL)

No Rails:
#select mines that are active with no rail depot
arcpy.MakeFeatureLayer_management(status, nrmines, nrminesSQL)

  • Mines within 1.5 Km of a rail depot were removed:
 #remove mines within 1.5 miles of a rail line
arcpy.SelectLayerByLocation_management(nrmines, "WITHIN_A_DISTANCE", rails, "1.5 KILOMETER", "REMOVE_FROM_SELECTION")

  • New layers selected above were saved to the ex7 geodatabase:
#save new layers of selected information
arcpy.CopyFeatures_management(nrmines, finalnr)


 

Pyscripting 2: May 13, 2015


Introduction
This python script corresponds to blog Post 6: Raster Analysis of Sand Mine Locations (Blog Post 6 ) The purpose of this python script is to create a weighted index model based on the five risk assessment models generated in the last post (Figure 1). One of the models is given a weighted value and a new risk model was generated (Figure 2).
 
Figure 1. Screenshot of the PyScripter interface.
 
Figure 2. Weighted risk assessment map generated using PyScripter.

 
 
Methods
 
  • System modules were imported into the PyScripter interface, the working environment was set, and the spatial analyst extension was activated:
 
#import system modules
import arcpy
from arcpy import env
 
#set environment
env.workspace = "Q:\StudentCoursework\CHupy\g2155.GEOG.337.001\LUCZAKJN\Lab\Ex_5\Tremp_raster2\NewTMP.gdb"
arcpy.env.overwriteOutput = True
 
#Activate License
arcpy.CheckOutExtension ("spatial")
 
 
  • Variables for the five rasters were specified:
 
#Set Variables
streams   arcpy.Raster("streamsreclassuse1")
park   arcpy.Raster("parkreclassuse1")
school  =  arcpy.Raster("school_reclass1")
prime  =  arcpy.Raster("farmlandreclassuse")
zone  =  arcpy.Raster("zoning_reclassuse1")
 
 
  • Raster math (multiplication) was performed on the variable of interest (i.e. school):
 
#Raster Math
outweight  =  (school*1.5)
 
 
  • The weighted index was created by adding cell values of each risk assessment model
 
#Weighted Index
weighted  =  (outweight + park + streams + prime + zone)
 
  • The resulting index model was saved in PyScripter since rasters file are temporary by default
#Save
weighted.save("weighted_result")
 
print "the script has completed"