Geometrical Analysis

Working with data

Video Locked

Please log in to watch this video

Log In
Chapter Info
Course The Ultimate PostGIS course
Module Working with data
Chapter Geometrical Analysis

Chapter Content

Geometrical analysis forms the backbone of spatial data processing in PostGIS. Whether you're calculating land parcel areas, measuring road lengths, or determining optimal locations, understanding PostGIS geometry analysis functions is essential for any GIS professional. This comprehensive guide covers all the key functions for measuring, analyzing, and manipulating spatial geometries.

Understanding Geometry vs Geography

Before diving into analysis functions, it's crucial to understand the difference between PostGIS geometry and geography types, as this affects how measurements are calculated:

Aspect Geometry Geography
Coordinate System Planar (projected) Spherical (WGS84)
Measurement Units Depends on projection Always meters
Accuracy High for local areas High for global calculations
Performance Faster Slower but more accurate
Best Use Local/regional analysis Global/continental analysis

Setting Up Sample Data

Let's create sample datasets using both geometry and geography types to demonstrate the analysis functions:

-- Enable PostGIS Extension
CREATE EXTENSION IF NOT EXISTS postgis;

-- Geometry table
CREATE TABLE geom_features (
  id SERIAL PRIMARY KEY,
  name TEXT,
  geom GEOMETRY
);

-- Geography table
CREATE TABLE geog_features (
  id SERIAL PRIMARY KEY,
  name TEXT,
  geog GEOGRAPHY
);

-- Insert sample data with real-world coordinates
INSERT INTO geom_features (name, geom) VALUES
('Delhi_Polygon', ST_GeomFromText('POLYGON((77.1 28.6, 77.3 28.6, 77.3 28.8, 77.1 28.8, 77.1 28.6))', 4326)),
('Highway_Line', ST_GeomFromText('LINESTRING(77.2 28.6, 77.25 28.65, 77.3 28.7)', 4326)),
('City_Center', ST_GeomFromText('POINT(77.2 28.65)', 4326)),
('Mumbai_Area', ST_GeomFromText('POLYGON((72.8 19.0, 73.0 19.0, 73.0 19.2, 72.8 19.2, 72.8 19.0))', 4326));

INSERT INTO geog_features (name, geog) VALUES
('Delhi_Geo_Polygon', ST_GeogFromText('POLYGON((77.1 28.6, 77.3 28.6, 77.3 28.8, 77.1 28.8, 77.1 28.6))')),
('Highway_Geo_Line', ST_GeogFromText('LINESTRING(77.2 28.6, 77.25 28.65, 77.3 28.7)')),
('City_Geo_Center', ST_GeogFromText('POINT(77.2 28.65)')),
('Mumbai_Geo_Area', ST_GeogFromText('POLYGON((72.8 19.0, 73.0 19.0, 73.0 19.2, 72.8 19.2, 72.8 19.0))'));

Calculating Area of Polygons

The ST_Area function calculates the area of polygon geometries. For accurate results with geometry types, always project to an appropriate coordinate system.

Area Calculation Examples

-- Geometry (requires projection for accurate results)
SELECT 
  name, 
  ST_Area(geom) AS area_degrees,
  ST_Area(ST_Transform(geom, 32643)) AS area_m2,
  ST_Area(ST_Transform(geom, 32643)) / 10000 AS area_hectares
FROM geom_features 
WHERE GeometryType(geom) = 'POLYGON';

-- Geography (automatically calculates in square meters)
SELECT 
  name, 
  ST_Area(geog) AS area_m2,
  ST_Area(geog) / 10000 AS area_hectares,
  ST_Area(geog) / 1000000 AS area_km2
FROM geog_features 
WHERE GeometryType(geog::geometry) = 'POLYGON';

Real-World Area Analysis

-- Compare areas between different regions
SELECT 
  name,
  ROUND(ST_Area(geog) / 1000000, 2) AS area_km2,
  CASE 
    WHEN ST_Area(geog) / 1000000 > 1000 THEN 'Large'
    WHEN ST_Area(geog) / 1000000 > 100 THEN 'Medium'
    ELSE 'Small'
  END AS size_category
FROM geog_features 
WHERE GeometryType(geog::geometry) = 'POLYGON'
ORDER BY ST_Area(geog) DESC;
Use Case Example Best Practice
Land parcel size Property assessment Use local UTM projection
Forest cover extent Environmental monitoring Use geography for large areas
Building footprints Urban planning Use projected coordinates
Administrative boundaries Census analysis Use geography for accuracy

Measuring Length and Perimeter

Use ST_Length for linear features and ST_Perimeter for polygon boundaries.

Length Calculations

-- Geometry (project for accurate measurements)
SELECT 
  name, 
  ST_Length(geom) AS length_degrees,
  ST_Length(ST_Transform(geom, 32643)) AS length_m,
  ST_Length(ST_Transform(geom, 32643)) / 1000 AS length_km
FROM geom_features 
WHERE GeometryType(geom) = 'LINESTRING';

-- Geography (automatically in meters)
SELECT 
  name, 
  ST_Length(geog) AS length_m,
  ST_Length(geog) / 1000 AS length_km
FROM geog_features 
WHERE GeometryType(geog::geometry) = 'LINESTRING';

Perimeter Calculations

-- Calculate perimeter of polygons
SELECT 
  name,
  ST_Perimeter(geog) AS perimeter_m,
  ST_Perimeter(geog) / 1000 AS perimeter_km,
  ST_Area(geog) / ST_Perimeter(geog) AS compactness_ratio
FROM geog_features 
WHERE GeometryType(geog::geometry) = 'POLYGON';

Advanced Length Analysis

-- Road network analysis
WITH road_segments AS (
  SELECT 
    name,
    ST_Length(geog) AS segment_length,
    ST_StartPoint(geog::geometry) AS start_point,
    ST_EndPoint(geog::geometry) AS end_point
  FROM geog_features 
  WHERE GeometryType(geog::geometry) = 'LINESTRING'
)
SELECT 
  name,
  ROUND(segment_length / 1000, 2) AS length_km,
  ST_AsText(start_point) AS start_coords,
  ST_AsText(end_point) AS end_coords
FROM road_segments;

Finding Centroids and Centers

ST_Centroid computes the geometric center of a feature, useful for labeling and location reference.

-- Calculate centroids
SELECT 
  name, 
  ST_AsText(ST_Centroid(geom)) AS centroid_wkt,
  ST_X(ST_Centroid(geom)) AS centroid_lon,
  ST_Y(ST_Centroid(geom)) AS centroid_lat
FROM geom_features 
WHERE GeometryType(geom) = 'POLYGON';

-- Geography centroids (convert to geometry first)
SELECT 
  name,
  ST_AsText(ST_Centroid(geog::geometry)) AS centroid_wkt,
  ST_X(ST_Centroid(geog::geometry)) AS centroid_lon,
  ST_Y(ST_Centroid(geog::geometry)) AS centroid_lat
FROM geog_features 
WHERE GeometryType(geog::geometry) = 'POLYGON';

Alternative Center Calculations

-- Different types of centers
SELECT 
  name,
  ST_AsText(ST_Centroid(geom)) AS geometric_center,
  ST_AsText(ST_PointOnSurface(geom)) AS point_on_surface,
  ST_AsText(ST_Centroid(ST_Envelope(geom))) AS envelope_center
FROM geom_features 
WHERE GeometryType(geom) = 'POLYGON';
Function Description Use Case
ST_Centroid Geometric center (may be outside polygon) Mathematical analysis
ST_PointOnSurface Point guaranteed to be inside polygon Label placement
ST_Centroid(ST_Envelope()) Center of bounding box Quick approximation

Bounding Boxes and Envelopes

ST_Envelope returns the minimum bounding rectangle, useful for spatial indexing and viewport calculations.

-- Calculate bounding boxes
SELECT 
  name,
  ST_AsText(ST_Envelope(geom)) AS bbox_wkt,
  ST_XMin(geom) AS min_lon,
  ST_YMin(geom) AS min_lat,
  ST_XMax(geom) AS max_lon,
  ST_YMax(geom) AS max_lat,
  (ST_XMax(geom) - ST_XMin(geom)) * (ST_YMax(geom) - ST_YMin(geom)) AS bbox_area
FROM geom_features;

-- Extent of entire dataset
SELECT 
  ST_AsText(ST_Envelope(ST_Collect(geom))) AS dataset_extent,
  ST_XMin(ST_Collect(geom)) AS dataset_min_lon,
  ST_YMin(ST_Collect(geom)) AS dataset_min_lat,
  ST_XMax(ST_Collect(geom)) AS dataset_max_lon,
  ST_YMax(ST_Collect(geom)) AS dataset_max_lat
FROM geom_features;

Convex Hull Analysis

ST_ConvexHull creates the smallest convex polygon that encloses geometries.

-- Individual convex hulls
SELECT 
  name,
  ST_AsText(ST_ConvexHull(geom)) AS convex_hull,
  ST_Area(ST_ConvexHull(geom)) AS hull_area,
  ST_Area(geom) / ST_Area(ST_ConvexHull(geom)) AS convexity_ratio
FROM geom_features 
WHERE GeometryType(geom) = 'POLYGON';

-- Convex hull of all features combined
SELECT 
  ST_AsText(ST_ConvexHull(ST_Collect(geom))) AS overall_convex_hull,
  ST_Area(ST_ConvexHull(ST_Collect(geom))) AS total_hull_area
FROM geom_features;

Directional Analysis

ST_Azimuth determines the bearing between two points, useful for navigation and directional analysis.

-- Calculate bearings between points
SELECT 
  'Delhi to Mumbai' AS route,
  ST_Azimuth(
    ST_Point(77.2, 28.65), 
    ST_Point(72.85, 19.1)
  ) AS bearing_radians,
  ST_Azimuth(
    ST_Point(77.2, 28.65), 
    ST_Point(72.85, 19.1)
  ) * 180 / pi() AS bearing_degrees,
  CASE 
    WHEN ST_Azimuth(ST_Point(77.2, 28.65), ST_Point(72.85, 19.1)) * 180 / pi() BETWEEN 0 AND 90 THEN 'Northeast'
    WHEN ST_Azimuth(ST_Point(77.2, 28.65), ST_Point(72.85, 19.1)) * 180 / pi() BETWEEN 90 AND 180 THEN 'Southeast'
    WHEN ST_Azimuth(ST_Point(77.2, 28.65), ST_Point(72.85, 19.1)) * 180 / pi() BETWEEN 180 AND 270 THEN 'Southwest'
    ELSE 'Northwest'
  END AS direction;

-- Azimuth analysis for line segments
SELECT 
  name,
  ST_Azimuth(ST_StartPoint(geom), ST_EndPoint(geom)) * 180 / pi() AS line_bearing
FROM geom_features 
WHERE GeometryType(geom) = 'LINESTRING';

Buffer Analysis

ST_Buffer creates buffer zones around features, essential for proximity analysis.

-- Create buffers around features
SELECT 
  name,
  ST_AsText(ST_Buffer(geom, 0.01)) AS buffer_01_degrees,
  ST_Area(ST_Buffer(geom, 0.01)) AS buffer_area_degrees
FROM geom_features 
WHERE GeometryType(geom) = 'POINT';

-- Geography buffers (in meters)
SELECT 
  name,
  ST_Area(ST_Buffer(geog, 1000)) AS buffer_1km_area_m2,
  ST_Area(ST_Buffer(geog, 5000)) AS buffer_5km_area_m2,
  ST_Area(ST_Buffer(geog, 10000)) AS buffer_10km_area_m2
FROM geog_features 
WHERE GeometryType(geog::geometry) = 'POINT';

Multi-Ring Buffer Analysis

-- Create multiple buffer rings
WITH buffer_rings AS (
  SELECT 
    name,
    generate_series(1000, 5000, 1000) AS buffer_distance
  FROM geog_features 
  WHERE GeometryType(geog::geometry) = 'POINT'
)
SELECT 
  br.name,
  br.buffer_distance,
  ST_Area(ST_Buffer(gf.geog, br.buffer_distance)) / 1000000 AS buffer_area_km2
FROM buffer_rings br
JOIN geog_features gf ON br.name = gf.name
ORDER BY br.name, br.buffer_distance;

Geometry Simplification

ST_Simplify reduces geometry complexity while preserving overall shape.

-- Simplify geometries for different zoom levels
SELECT 
  name,
  ST_NPoints(geom) AS original_points,
  ST_NPoints(ST_Simplify(geom, 0.001)) AS simplified_001_points,
  ST_NPoints(ST_Simplify(geom, 0.01)) AS simplified_01_points,
  ST_NPoints(ST_Simplify(geom, 0.1)) AS simplified_1_points,
  ST_AsText(ST_Simplify(geom, 0.01)) AS simplified_geometry
FROM geom_features 
WHERE GeometryType(geom) IN ('POLYGON', 'LINESTRING');

3D Geometry Analysis

PostGIS supports 3D analysis for geometries with Z coordinates.

-- Create 3D geometries
INSERT INTO geom_features (name, geom) VALUES
('3D_Line', ST_GeomFromText('LINESTRING Z(77.2 28.6 100, 77.25 28.65 150, 77.3 28.7 200)', 4326)),
('3D_Point', ST_GeomFromText('POINT Z(77.2 28.65 250)', 4326));

-- 3D measurements
SELECT 
  name,
  ST_Length(geom) AS 2d_length,
  ST_3DLength(geom) AS 3d_length,
  ST_3DLength(geom) - ST_Length(geom) AS elevation_difference,
  ST_Z(ST_StartPoint(geom)) AS start_elevation,
  ST_Z(ST_EndPoint(geom)) AS end_elevation
FROM geom_features 
WHERE GeometryType(geom) = 'LINESTRING' AND ST_Is3D(geom);

Advanced Geometric Analysis

Shape Complexity Metrics

-- Calculate shape complexity metrics
SELECT 
  name,
  ST_Area(geog) AS area_m2,
  ST_Perimeter(geog) AS perimeter_m,
  -- Compactness ratio (1 = perfect circle)
  (4 * pi() * ST_Area(geog)) / POWER(ST_Perimeter(geog), 2) AS compactness,
  -- Shape index
  ST_Perimeter(geog) / (2 * sqrt(pi() * ST_Area(geog))) AS shape_index,
  -- Convexity ratio
  ST_Area(geog) / ST_Area(ST_ConvexHull(geog::geometry)::geography) AS convexity
FROM geog_features 
WHERE GeometryType(geog::geometry) = 'POLYGON';

Spatial Statistics

-- Dataset-wide spatial statistics
SELECT 
  COUNT(*) AS total_features,
  AVG(ST_Area(geog)) AS avg_area_m2,
  STDDEV(ST_Area(geog)) AS stddev_area,
  MIN(ST_Area(geog)) AS min_area_m2,
  MAX(ST_Area(geog)) AS max_area_m2,
  SUM(ST_Area(geog)) / 1000000 AS total_area_km2
FROM geog_features 
WHERE GeometryType(geog::geometry) = 'POLYGON';

Performance Optimization

Optimization Description When to Use
Spatial Indexes GIST indexes on geometry columns Always for spatial queries
Appropriate Projections Use local projections for measurements Geometry-based calculations
Geography Type Use for global/large-scale analysis Cross-continental data
Geometry Simplification Reduce complexity for visualization Web mapping applications
Bounding Box Filters Pre-filter with ST_Intersects(bbox) Large dataset queries

Real-World Applications

Urban Planning

-- Analyze urban density and accessibility
SELECT 
  name,
  ST_Area(geog) / 1000000 AS area_km2,
  -- Calculate building density (assuming point features are buildings)
  (SELECT COUNT(*) FROM geog_features gf2 
   WHERE GeometryType(gf2.geog::geometry) = 'POINT' 
   AND ST_Within(gf2.geog, gf1.geog)) AS building_count,
  -- Buildings per square kilometer
  (SELECT COUNT(*) FROM geog_features gf2 
   WHERE GeometryType(gf2.geog::geometry) = 'POINT' 
   AND ST_Within(gf2.geog, gf1.geog)) / (ST_Area(gf1.geog) / 1000000) AS buildings_per_km2
FROM geog_features gf1
WHERE GeometryType(gf1.geog::geometry) = 'POLYGON';

Environmental Analysis

-- Environmental impact assessment
SELECT 
  name,
  ST_Area(geog) / 10000 AS area_hectares,
  ST_Perimeter(geog) / 1000 AS perimeter_km,
  -- Edge-to-area ratio (habitat fragmentation indicator)
  (ST_Perimeter(geog) / 1000) / (ST_Area(geog) / 10000) AS edge_area_ratio
FROM geog_features 
WHERE GeometryType(geog::geometry) = 'POLYGON'
ORDER BY edge_area_ratio DESC;

Best Practices

  1. Choose appropriate coordinate systems: Use local projections for accurate measurements
  2. Validate geometries: Use ST_IsValid() before analysis
  3. Consider scale: Use geography for global analysis, geometry for local
  4. Index spatial columns: Create GIST indexes for better performance
  5. Simplify when appropriate: Reduce complexity for visualization
  6. Handle edge cases: Check for NULL geometries and empty results
  7. Document assumptions: Clearly state coordinate systems and units used