Drainage Network and Basin Analysis

Southern California

Introduction

Using topographic, drainage network, and land use/land cover (LULC) data for our study area in southern California, we analyzed land use and landcover in relation to the drainage network and basins.

Part 1:  Summarize Land Use/Land Cover
    a.    within the entire watershed above Fillmore
    b.    within 250 meters of channels upstream of Fillmore
    c.    describe and explain differences between these two summaries

Part 2: Map all 5th order drainage basins - identify the basins with:
    a.    the highest discharge
    b.    the highest mean slope value
    c.    the highest total erosion potential
    d.    the greatest percentage of forest land
    e. discuss possible implications of these results - in terms of the patterns of erosion you might expect
         and/or in terms of refinements to a model for erosion potential

Part 3: Identify the areas downstream from all industrial sites
            (both industrial sites and downstream areas within 100 meters of a channel)

ACQUIRING THE DATA:

Data Source:  USGS Land use/Land cover (LULC) data
Record Descriptions

Geographic Information Retrieval and Analysis System (GIRAS)
     The GIRAS format involves a standard character fixed length record (usually
     ASCII-coded, 80 character card images). Each record may consist of 1 to 16 data-element
     fields. Each data-element field may be one of three different types: 1) 16 bit binary integer;
     2) 32 bit binary integer or 3) a string of text characters. The GIRAS file structure is
     comprised of a map header, section header, arc records subfile, coordinate subfile, polygon
     records subfile, file of arcs by polygon (FAP) subfile, text subfile and an associated data
     subfile.

Initially we used Don's "hawked" coverage la_lulc located in the  group1 folder.  Then we tried to change the projection from Albers to UTM using  "project define".

However, this didn't work out because the coverage still appeared to be in Albers when we looked at it in A/P. So we decided to download the data from the USGS homepage.  Then we unzipped it, ran
Arc: fold -80 land_use > la80,
Arc:  girasarc2 la80 lulc utm -119 0 0

Here is a display of our land use land cover data downloaded in UTM:
 
 

Data Manipulation

Ran a dissolve to get rid of the neat lines
Arc: dissolve lulc lulc2 lulc-id poly

1.a.     Summarizing land use within the watershed above Fillmore, CA.

To summarize landuse within the entire watershed above Fillmore, we first had to identify the area that drains into Fillmore.  This was done in gridspace, and subsequently vectorized.  We then clipped the landuse cover to this watershed, and summarized landuse.

We used the flow direction grid from project 2 and the watershed command to interactively create the watershed grid.  The point for Fillmore
Grid: fm_ws = watershed ( flowdirgrid, selectpoint(flowdirgrid, *))
Define the points <9 to END>
Running... 100%
Delineating watershed...
 
 

Then we coverted the Fillmore watershed grid to polygon coverage to use for clipping the landuse coverage
Grid: fm_ws_poly = gridpoly(fm_ws)

The landuse cover was clipped to the watershed using the clip command:
Arc: clip lulc2 fm_ws_poly fmlu
 
 

The minimum polygon size for the landuse classification is 4 ha, with minimum width of 200 meters.
Land use classifications are listed here
 
 

Fillmore land use

Closeup of Fillmore land use

Summary of Fillmore Land Use

1.b.    Summarize LULC within 250 meters of channels upstream of Fillmore

We switched to “grid space” to complete this part of the project in order to have greater flexibility in the generation and selection of data than that provided by the zonal operations in A/I.  By constructing grids for the Fillmore land use and streams vector coverages we had the ability to easily select for values from overlaying grid cells within the same geographic space using grid math and conditional arguments.  Our process is as follows:

1.  Construct the grids from coverages:
Grid: fmlugrid = polygrid(fmlu, lulc-id, #, #, 50)
 Converting polygons from FMLU to grid FMLUGRID

Grid: sbufgrid = polygrid(fmstream_buf, inside, #, #, 50)
 Converting polygons from FMSTREAM_BUF to grid SBUFGRID

2.  Next we set the values of cells that were “inside polygons” within the stream vector coverage to null.  These cells originally had a value of 1.

Grid: sbufgrid2 = con(sbufgrid == 100,100)

3.  We used the CON statement again to construct a grid that contained the land use values from that fell within the stream buffer grid.  Since the stream buffer grid had only values of 100 or null, the null values were ignored in the Con statement and only land use values that existed in overlapping stream buffer cells equal to 100 were selected.

GRID:  fmlu_s_grid =  con(sbufgrid2, fmlugrid)

The following diagrams illustrate the process:
 

Land use grid	Stream buffer grid	Land use within stream buffer grid

The following chart summarized the land use area in square meters by land use category.

1.c.    Comparison of differences in land use

The land use within the stream buffers is relatively similar that of the overall study area.  There is no indication that the land use occuring within the buffer zone is disproportionately residential, industrial, or otherwise developed.  In fact, shrub and brush rangeland and evergreen forest make up the largest land use classes by percent area withiin the buffers.  Our stream coverage may have included more channels than really exist in the landscape as we used a relatively low flow accumulation threshold to determine where channels were located.  This may have demphasized the percent of development within the buffers.
 

2.a.      Map all 5th order drainage basins - identify the basins with the highest discharge

Our method to map the 5th order drainage basins was to select the end points of the 5th order streams and generate basins based on these points.  To do this, we selected only those nodes of 5th order streams that intersected 6th order streams.  We also considered the ocean shoreline and the border of the study area coverages as proxy 6th order streams.  We constructed the watersheds and conducted analysis of land use in gridspace.  The process is as follows:

1.    The stream cover was appended with the shoreline and border arcs of those resective coverages.  We selected for 5th order streams and 6th order streams and placed each in respective coverages using the PUT command in ARC/EDIT.

2.    We converted the nodes of the fifth order streams into points
Arc: nodepoint 5thorder 5points
Arc: build 5points points

3.    This allowed us to select the points that overlapped with the 6th order streams
Arcplot: resel 5points points overlap 6thorder arcs
27 of 1179 selected

4.    These selected points were made into a new coverage using the create command.
Arcplot: create 5pts select 5points points

5.    There were several places were additional points were needed.  These areas were on the edges of the extent, where
Arcedit:  ec 5pts
Arcedit: ef points
Arcedit: add                                     {Interactively added points at places where the 5thorder streams went off the extent of the map}
Arcedit: save

6.    Then we converted these points into a grid:
Arc: pointgrid 5pts 5ptsgrid
 Converting points from 5pts to grid 5ptsgrid
Cell Size (square cell): 200
Convert the Entire Coverage(Y/N)?: y
Enter background value (NODATA | ZERO): nodata

2.b.  Discharge
1.    Before weighting the source points with flow accumulation data and creating the watershed, they must all have the same value.  We used the following commands to set the values of these cells to 1.
GRID:  docell
:: if( 5ptsgrid > 0) 5pts-b = 1
:: end

2.    We then weighted the points used to create the watershed using the flow accumulation grid.
GRID: 5ptsdis = 5ptsgrid * flowacc
GRID: 5dis-b = 5dis / 100000               {This command was used to make numbers easier to visualize}

3.    The watersheds were created with the following command:
GRID: 5dis-b-ws = watershed( flowdirgrid, 5dis-b )

2.c.    Mean slope by watershed

To determine the mean slope for each watershed we used the zonalmean command:
GRID: meanslope = zonalmean( slopegrid, 5thbasin )
 
 

2.d.    Potential erosion by watershed

We used a simple formula to estimate erosion intensity (slope * discharge)
GRID: erosion = meanslope * 5dis-b-ws
 
 

2.e.    Percent forest by watershed

1.     To get a landuse grid for the entire region we used the following command:
Usage:  (*) POLYGRID(<cover>,{item},{lookup_table},{weight_table},{cellsize})
Grid: lugrid = polygrid( lulc2, lulc-id, #, #, 200)
 Converting polygons from LULC2 to grid LUGRID
Number of Rows    = 573
Number of Columns = 930

2.        We selected forested areas from the landuse grid using the con command.  Forested landuse classification are: 41, 42, 43, 61.  We set each forested cell equal to one.
GRID: forest = con( lugrid == 41or lugrid == 42 or lugrid == 43 or lugrid == 61, 1)

3.    Then we summed the number forest cells in each basin:
GRID: sumforest = zonalsum( 5thbasin, forest )

4.    Each basin was multiplied by the area of a cell to yield the total forest area per basin
GRID: forestarea = sumforest * 40000

5.    The zonal area command was used to determine the total area of each basin:
GRID: basinarea = zonalarea( 5thbasin )

6.    The percent forest for each basin was determined using grid math.  The float command was used to create a floating point output grid.
GRID: percentforest = float( forestarea / basinarea )
 
 

2.f.  Discussion of potential erosion
 

Mean Slope	Mean Annual Discharge	Proportion Forested

Potential Erosion determined by mean slope x discharge / % forest

Although this discussion of potential erosion is constrained by a lack of data concerning many of the important processes and structures that influence erosion (e.g. soils), some of the more important factors regulating erosion potential include mean slope, mean discharge, and vegetative cover.  We used a simple, evenly weighted formula to determine potential erosion (mean slope in degrees x mean annual discharge in cubic meters / percent forest cover).  Increasing slope and discharge increase erosion proportionally, while increasing forest cover reduces the potential for erosion by increased soil holding capabilities and increased infiltration.  This method may have disproportionally represented areas with very high or very low forest cover.  Further analysis could be improved by using a weighted method for forest cover.

The area with the highest potential erosion is in the southwest of the project area where discharge is highest and forest cover is lowest.  Mean slope in this area is also relatively high.  This area could be targeted for erosion control measures, especially if human land use is responsible for reduced forest cover.

3.a  Identify the areas downstream from all industrial sites

Industrial sites are identified with the classifications 13 (Industrial) and 15 (Industrial and Commercial) in the landuse cover.  If we select the industrial polygons, and buffer them by 100 meters, we can select stream nodes that overlap with buffered industrial areas.  Then using the TRACE command we can create arcs that flow downstream from the selected nodes.  We buffered these selected streams, converted them to gridspace and analysed land use.

TRACE <BOTH | UPSTREAM | DOWNSTREAM> <cover> <out_selection_file> {ITEM flow_direction_item | FILE flow_direction_file} {* | selection_file} {end_selection_file} {previous_arc_file}
This command finds the arcs and nodes upstream and/or downstream from a set of specified locations and writes them to out_selection_file.  The selection_file specifies that the origin node(s) of the trace are to come from a selection file. This selection file must contain node features for the coverage specified for TRACE. A pathname may be given.

1.    The following process was used to identify industrial areas:
Arcedit: ec lulc2
The edit coverage is now /EXPORT/DATA1/GEO595A/STUDENTS/GROUP4/PROJECT3/LULC2
WARNING the Map extent is not defined
Defaulting the map extent to the BND of /EXPORT/DATA1/GEO595A/STUDENTS/GROUP4/PROJECT3/LULC2
Arcedit: ef polys
5234 element(s) for edit feature POLYS
Arcedit: drawe polys
Arcedit: draw
Arcedit: resel
0 element(s) now selected
Arcedit: asel lulc-id = 13   {This selects industrial areas)
229 element(s) now selected
Arcedit: drawsel
Arcedit: asel lulc-id = 15   {This selects industrial and commercial areas}
239 element(s) now selected
Arcedit: put ind
Creating /EXPORT/DATA1/GEO595A/STUDENTS/GROUP4/PROJECT3/IND2
Copying the polygon(s) into /EXPORT/DATA1/GEO595A/STUDENTS/GROUP4/PROJECT3/IND2.
239 polygon(s) copied
901 arc(s) copied
239 label(s) copied

The cover containing industrial polygons (ind) was then built for polygons using the build command

2.    The industrial polygons were then buffered by 100 meters
Usage: BUFFER <in_cover> <out_cover> {buffer_item} {buffer_table}
              {buffer_distance} {fuzzy_tolerance} {LINE | POLY | POINT | NODE}
Arc: buffer ind ind_buf # # 100

3.    Then we reselected streams that overlapped with the buffered industrial polygons
Arcplot: resel str arcs overlap ind_buf polys
STR arcs : 1152 of 35580 selected.

4.    This selected set was put into a new coverage using the create command.  This was done because selecting only the nodes that overlapped with the industrial polygons missed many streams that passed through these areas, but did not contain node points within the polygons.  By creating a new coverage with these arcs we were able to select the nodes from these arcs to be used in the tracing process.
Arcplot: create indstr select str arcs
Arc: build indstr node

5.    Then we reselected the nodes from all streams that overlapped with the streams that overlapped industrial areas.
Arcplot: resel str nodes overlap indstr nodes
STR nodes : 1545 of 35986 selected.

6.    These nodes were written to a selection file to be used in the TRACE command:
Arcplot: writeselect strnode.sel str nodes
1 selection set(s) written to strnode.sel

7.  The stream coverage needed to be cleaned before running the TRACE command.  This may be due to duplicate nodes which overlap, but are not connected.  This disrupted the tracing since the streams were not continuous without being cleaned.
Arc: clean str # # # line
Cleaning /EXPORT/DATA1/GEO595A/STUDENTS/GROUP4/PROJECT3/STR
 Sorting...
 Intersecting...
 Assembling arcs...

8.    The TRACE command was used to identify the areas downstream from the previous selected nodes of streams that passed through industrial areas.
Arcplot: trace downstream str downstr.sel # strnode.sel
str nodes : 2966 of 35986 traced.
str arcs : 2914 of 35580 traced.

9.   A new coverage was created from the selected set of downstream arcs.
Arcplot: create downstr select str arcs

10.  These downstream arcs were buffer by 100 meters
Arc: buffer downstr downstrbuf # # 100

11.  The buffered streams were converted to gridspace to allow manipulation using the zonal functions
Arc: polygrid downstrbuf dstrgrid inside
Converting polygons from downstrbuf to grid dstrgrid
Cell Size (square cell): 200
Convert the Entire Coverage? (Y/N): y
Number of Rows    = 573
Number of Columns = 901

12. The interior polygons created by buffering were then removed with a con command
Grid: dstrgrid2 = con( dstrgrid == 100, 100)

13.  The gridded, buffered streams were overlapped with the landuse grid to analyse land use within 100 meters of streams downstream from industrial areas.
Grid: dstrlu = con( dstrgrid, lugrid )

Residential and Rangeland uses are the most prevalent land uses for this area.  Croplands, pasture and other industrial uses are also common.  This information could be used to determine what land uses would be most affected by industrial pollutants.