**Home**»

**Software**»

**PostGIS: Spatially enabled Relational Database Sytem**»

**Analysis of SSURGO Data in PostGIS: An Overview**» Profile water storage as calculated from SSURGO

# Profile water storage as calculated from SSURGO

**Overview**

This example illustrates a component-percent weighted mean water storage, based on total profile storage values. The `comppct_r` column is used as a weight to adjust the influence of the profile water storage values for each component: larger weights (larger components) result in more influence. When performing weighted means, use care in the selection of an appropriate weighting parameter. A selection of a weighting variable not direclty related to the property which is being averaged can result in an effect known as Simpson's paradox.

Some background on Soil Water Holding Capacity and Irrigation Management.

The general formula for computing profile water storage is defined as follows. Note that the `awc_r` column in the SSURGO database (chorizon table) is pre-adjusted to compensate for coarse fragments (see relevant SSURGO metadata below).

*Calculate the water storage within a given component*

**component_whc** = profile_sum_awc = sum( (hzdepb_r - hzdept_r) * awc_r )

*Calculate the weighted mean value of the profile water storage, for a given map unit*

weighted mean = sum(weights * x) / sum(weights)

**mapunit_whc** = sum( comppct_r * component_whc ) / sum( comppct_r )

Note: weights corresponding to non-soil components must be filtered out

From the SSURGO 2.1 table definitions:

**Column Label: AWC - Representative Value (awc_r)**

The amount of water that an increment of soil depth, inclusive of fragments, can store that is available to plants. AWC is expressed as a volume fraction, and is commonly estimated as the difference between the water contents at 1/10 or 1/3 bar (field capacity) and 15 bars

(permanent wilting point) tension and adjusted for salinity, and fragments.

1. Reduce to 1:1 relationship with **cokey** by aggrigating AWC horizon data in the `chorizon` table, and linking to the `component` table.

```
SELECT mukey, compname, comppct_r, a.* FROM component
JOIN
(
SELECT cokey, sum( (hzdepb_r - hzdept_r) * awc_r) AS component_whc, sum((hzdepb_r - hzdept_r)) AS depth
FROM chorizon WHERE areasymbol = 'ca113'
GROUP BY cokey
) AS a
ON component.cokey = a.cokey
WHERE component.areasymbol = 'ca113'
ORDER BY mukey ;
```

Query Results

mukey | compname | comppct_r | cokey | component_whc | depth

--------+---------------------+-----------+----------------+---------------+-------

459225 | Millsholm | 30 | 459225:624008 | 6.82 | 58

459225 | Dibble | 55 | 459225:624007 | 12.36 | 86

459226 | Millsholm | 40 | 459226:624013 | 6.82 | 58

459226 | Dibble | 45 | 459226:624012 | 12.36 | 86

459227 | Millsholm | 30 | 459227:624019 | 6.82 | 58

459227 | Dibble | 60 | 459227:624018 | 12.36 | 99

...

2. Reduce to 1:1 relationship with **mukey** by aggreigating based on *component percent*:

```
SELECT mukey,
-- note that weights from non-soil components must be removed
-- otherwise, weighted mean values will be too low
SUM(comppct_r * component_whc) / SUM(comppct_r) AS comppct_weighted_whc_cm
FROM component
JOIN
(
SELECT cokey, sum( (hzdepb_r - hzdept_r) * awc_r) AS component_whc,
sum((hzdepb_r - hzdept_r)) AS depth
FROM chorizon
WHERE areasymbol = 'ca113'
GROUP BY cokey
) AS a
USING (cokey)
WHERE component.areasymbol = 'ca113'
-- filter out components that are missing soils data
AND a.component_whc IS NOT NULL
GROUP BY mukey ;
```

mukey | comppct_weighted_whc_cm

--------+-------------------------

459225 | 10

459226 | 10

459227 | 11

...

Link to polygons and setup access permissions:

```
-- create the new table with both geometry and attributes
CREATE TABLE yolo_whc AS
SELECT ogc_fid, wkb_geometry AS wkb_geometry, b.mukey, b.comppct_weighted_whc_cm
FROM mapunit_poly
-- use LEFT JOIN to include non-soil polygons in the result set
-- alternatively use JOIN to ignore non-soil polygons
LEFT JOIN
(
SELECT mukey,
-- note that weights from non-soil components must be removed
-- otherwise, weighted mean values will be too low
SUM(comppct_r * component_whc) / SUM(comppct_r) AS comppct_weighted_whc_cm
FROM component
JOIN
(
SELECT cokey, sum( (hzdepb_r - hzdept_r) * awc_r) AS component_whc,
sum((hzdepb_r - hzdept_r)) AS depth
FROM chorizon
WHERE areasymbol = 'ca113'
GROUP BY cokey
) AS a
USING (cokey)
WHERE component.areasymbol = 'ca113'
-- filter out components that are missing soils data
AND a.component_whc IS NOT NULL
GROUP BY mukey
) AS b
-- JOIN constraint
USING (mukey)
-- optional constraint to limit geometry search in mapunit_poly table
WHERE mapunit_poly.areasymbol = 'ca113' ;
```

Create indexes and register the new geometry:

```
-- create attribute and spatial index:
CREATE UNIQUE INDEX yolo_whc_idx ON yolo_whc (ogc_fid) ;
CREATE INDEX yolo_whc_spatial_idx ON yolo_whc USING gist (wkb_geometry gist_geometry_ops);
-- register in geometry_columns table:
INSERT INTO geometry_columns VALUES ('','public','yolo_whc','wkb_geometry',2,9001,'POLYGON');
-- optional:
-- fix permissions
-- GRANT SELECT on table yolo_whc to [user] ;
-- cleanup
vacuum analyze yolo_whc ;
```

4. Dump to shapefile with `ogr2ogr`:

` ogr2ogr ywhc.shp PG:"dbname=ssurgo_combined user=xxxx password=xxxx" yolo_whc`

## Software

- General Purpose Programming with Scripting Languages
- LaTeX Tips and Tricks
- PostGIS: Spatially enabled Relational Database Sytem
- PROJ: forward and reverse geographic projections
- GDAL and OGR: geodata conversion and re-projection tools
- R: advanced statistical package
- GRASS GIS: raster, vector, and imagery analysis
- Generic Mapping Tools: high quality map production