Geospatial Analysis I: map algebra, neighborhood operations

1
Geospatial Analysis Imap algebra, neighborhood
operations

• Geospatial Analysis and Modeling
• Lecture notes
• Helena Mitasova, NCSU MEAS

2
Outline

• summary statistics global and zonal operations
• neighborhood (focal) operations
• map algebra and local operations
• expressions, operators, functions and variables
• basic calculations, integer and floating point
  data
• "if" conditions, handling NULLs, creating masks
• special operators
• patching, mosaicking and overlay

3
Summary statistics

• Continuous raster data stored numbers are values
  quantifying the phenomenon
• univariate statistics min, max, mean, standard
  deviation,
• histogram requires discretization into bins

4
Summary statistics

• Continuous raster data stored numbers are values
  quantifying the phenomenon
• univariate statistics min, max, mean, standard
  deviation,
• histogram requires discretization into bins
• Discrete raster data stored numbers can be
  values (quantitative data) or category numbers
  (qualitative data)
• univariate statistics can be applied only to
  quantitative attributes
• mode (most frequent cat), diversity (number of
  different cats) apply to category data
• Can be applied as global, zonal, or focal
  operations

5
Global statistics continuous data
Mean and standard deviation for elevation raster
maps m res. N mean stdv 1m 525000
120.769 6.619 10m 5250 120.935 6.600
Same source (lidar elevation), different
interpolation and resolution except 30m SRTM
6
Global statistics continuous data
Mean and standard deviation for elevation raster
maps m res. N mean stdv 1m 525000
120.769 6.619 10m 5250 120.935 6.600
Histograms
Number of cells
1m 10m
elevationm
Same source (lidar elevation), different
interpolation and resolution except 30m SRTM
7
Global statistics continuous data
Histograms
Mean and standard deviation for elevation raster
maps m res. N mean stdv 1m 525000
120.769 6.619 10m 5250 120.935
6.600 30m 624 120.284 6.981
1m 10m 30m
Same source (lidar elevation), different
interpolation and resolution except 30m SRTM
8
Global statistics continuous data
Histograms
Mean and standard deviation for elevation raster
maps m res. N mean stdv 1m 525000
120.769 6.619 10m 5250 120.935
6.600 30m 624 120.284 6.698 SR30m 624
121.991 6.424 500m 1.7 116.845 2.436
1m 10m 30m 30m SRTM
Same source (lidar elevation), different
interpolation and resolution except 30m SRTM
9
Global statistics discrete objects

• Works for soil erodibility factor Does not
  apply to soil ID
• Mean 0.248 stdv 0.028 Mean
  29379.2 stdv 349.302
• but not for erodibility classes

10
Zonal statistics continuous data

• Zonal summary statistics for continuous data
• Raster 1 rasterized polygon (zipcodes)

11
Zonal statistics continuous data

• Zonal summary statistics for continuous data
• Raster 1 rasterized polygon (zipcodes)
• Raster 2 continuous field (elevation)

12
Zonal statistics continuous data

• Zonal summary statistics for continuous data
• Raster 1 rasterized polygon (zipcodes)
• Raster 2 continuous field (elevation)
• Resulting raster 3 mean elevation value for
  each zipcode

zonal statistics applicable to continuous data
min, max, mean, stdv, median
13
Zonal statistics category data

• Zonal statistics for category data
• Raster 1 rasterized polygon (zipcodes)
• Raster 2 category map (land use 1996 with 24
  classes)

HI developed LI developed South. pine Bott.
hardwoods Mixed hardwoods Cultivated Herbaceous
14
Zonal statistics category data

• Zonal statistics for category data
• Raster 1 rasterized polygon (zipcodes)
• Raster 2 category map (land use 1996 with 24
  classes)
• Resulting raster 3 most frequent land use for
  each zipcode
• Guess the categories in the result below

HI developed LI developed South. pine Bott.
hardwoods Mixed hardwoods Cultivated Herbaceous
15
Neighborhood operations

• Focal operations value at a grid cell is
  function of its neighborhood values.
• Grid cell neighborhood moving window, square or
  circular
• Same handling of continuous and discrete data as
  for global and zonal operations

16
Neighborhood operations

• Functions
• mean, stdv,max, min, median, sum, variance
• mode, diversity, interspersion most frequent
  cell value, number of different cell values
• filters
• fluxes
• Implemented as modules or in map algebra
• Image processing applications

17
Neighborhood operations

• Continuous data
• smoothing
• SRTM DEM
• focal mean operator
• 5x5 window

18
Neighborhood operations

• Continuous data
• smoothing
• SRTM DEM
• focal mean operator
• 5x5 window

19
Neighborhood operations

• Discrete data
• land use 1996
• focal diversity
• 7x7 window

number of land use classes in the cell 7x7
neighborhood
20
Neighborhood operations

• Discrete data
• land use 1996
• focal diversity
• 7x7 window

number of land use classes in the cell 7x7
neighborhood
21
Raster map algebra

• Computes a new raster map using expression built
  by applying logical and/or arithmetic operators,
  mathematical functions to existing raster maps
  representing variables
• newmap expression (map1, map2, ...map3)
• Local Expression is applied on per-cell basis
• Focal moving window calculations apply
  expression to the given cell and its neighboring
  cells

22
Raster map algebra

• Each software has its own syntax and set of
  operators and functions, examples
• Logical ?
• Arithm. Operators ?
• Math. Functions ?
• Special operators ?

23
Raster map algebra

• Each software has its own syntax and set of
  operators and functions, examples
• Logical operators and functions
• less than, equal,
• and, or, not, ...
• if(x), if(x,a,b)
• Arithm. Operators
• , -, , /
• Math. Functions
• exp(x,y),sin(x),log(x)
• min(x1,x2,..), max(), median()

24
Map algebra int and fp rules

• Integer and floating point matters for arithmetic
  expressions
• F(map1_int, map2_int) -gt map_int
• F(map1_int, map2_fp) -gt map_fp
• Example compute ndvi index from landsat integer
  maps
• ndvi(tm4-tm3)/(tm4tm3)
• ndvifloat(tm4-tm3)/float(tm4tm3)
• ndvi1.(tm4-tm3)/(tm4tm3)

25
Map algebra int and fp rules
tm3 1-255

• ndvi(tm4-tm3)/(tm4tm3)

?
tm4 1-255
26
Map algebra int and fp rules
tm3 1-255

• ndvi(tm4-tm3)/(tm4tm3)

0
tm4 1-255
27
Map algebra int and fp rules

• tm3 1-255

ndvi1.(tm4-tm3)/(tm4tm3)
tm4 1-255
28
Map algebra int and fp rules

• tm3 1-255

ndvi1.(tm4-tm3)/(tm4tm3)
tm4 1-255
29
Map algebra if conditions

• logical expressions apply to both continuous
  and discrete (category) data and their
  combination
• use for complex reclassification, masking and
  overlays
• Example find all forested areas with elevation gt
  120m

highforestif(landclass965 elevation gt 120.,
elevation, 0)
30
Map algebra if conditions

• logical expressions apply to both continuous
  and discrete (category) data and their
  combination
• use for complex reclassification, masking and
  overlays
• Example find all forested areas with elevation gt
  120m

highforestif(landclass965 elevation gt 120.,
elevation, 0)
31
Map algebra null data handling

• Raster maps can include NULL data (no-data), each
  software has its own rules how to handle them in
  map algebra operations
• General rule If a cell is null in at least one
  map (variable) then it is null in the resulting
  map
• If statements can test for null and/or assign
  a cell null value based on the if condition
• Special operators can be implemented to extend
  the rules applied to nulls

32
Map algebra null data handling

• Compute average elevation from 30m SRTM and one
  tile of 6m NCFlood DEM,
• elev_avg(elev_srtm_30melevlid_D782_6m)/2.

NULL
?
zero
102m
Use if statement to limit the computation to
no-NULL areas
33
Map algebra null data handling

• Compute average elevation from 30m SRTM and one
  tile of 6m NCFlood DEM, cells that do not have
  NCFlood tile values are NULL
• elev_avg(elev_srtm_30melevlid_D782_6m)/2.

NULL
NULL
zero
102m
51m
Use if statement to limit the computation to
no-NULL areas
34
Map algebra relative coordinates

• Focal operations - moving window - operations
• filters or simulations slower than C modules
  but more flexible,
• edge cell problem
• Operations on raster map subset, creating new
  maps
• Replacing values in a subregion of a raster map
• Generating tilted plane
• Modeling dynamic processes
• Fluxes and diffusion

35
Map algebra relative coordinates

• Example for applying operation to subset of a
  given region
• For an area where
• x lt 637033. and y lt 225552. and
• x gt 633985. and y gt222504.
• use lidar-based DEM elevlid_D783_6m
• everywhere else within the define region
• use SRTM-based elev_srtm_30m

36
Map algebra relative coordinates

• Examples generating tilted plane (e.g.
  geological fault)
• zaxbyc,
• zxy100 plane dipping to NW starting at 100m
• planerow()col()98

Enough map algebra for now, more comes along later
37
Raster maps patching and overlay

• Merging several neighboring raster maps into a
  single raster
• Filling-in nulls in base raster map with values
  from additional raster maps
• Order of maps matters
• Issues resolution, extension that is not
  aligned, no-data slivers due to projection

38
Raster maps patching

• Merging several neighboring raster maps into a
  single raster
• example patch DEM tiles A, B, C, D into map E

A B C D
39
Raster maps patching and overlay
Filling-in nulls combine roads, lakes,
building footprint and elevation map into a
single raster
40
Reclassification and rescaling

• Reclassification / recoding transformation rules
  used to convert between raster value types and
  classes
• general form
• change interval or list of values ltwi, wjgt to a
  new value vk or interval of values ltvk,vlgt
• Rescaling
• applies to continuous data or values associated
  with discrete data
• wnew kwold for each cell
• histogram equalized values are distributed
  according to a cumulative histogram

41
Reclassification examples

• Reclassification examples
• How is this handled in DBMS? When do you want to
  do this?

INT zipcodes to INT cities 27511 27513 1 Cary
42
Reclassification examples

• Reclassification examples

Zipcodes to cities 27511 27513 1 Cary
FP aspect to INT cardinal directions
0- 45 4 East 45-135 1 North 135-225 2
West 225-315 3 South 315-360 4 East
43
Rescaling examples

• Slope equal interval
  

44
Rescaling examples

• Slope equal interval
  hist. equalized

45
Summary and references

• summary statistics (global, zonal, focal)
• Neteler ch. 5, Bolstad ch. 10
• map algebra
• Neteler ch. 5.2 , Bolstad ch. 10
• patching, rescaling, reclassification