PostGIS 教程：从入门到项目应用

第一部分：环境搭建与基础操作

1.1 环境准备

Docker 快速部署

# docker-compose.yml
version: '3.8'
services:
  postgis:
    image: postgis/postgis:15-3.3
    environment:
      POSTGRES_DB: gisdb
      POSTGRES_USER: gisuser
      POSTGRES_PASSWORD: gispassword
    ports:
      - "5432:5432"
    volumes:
      - postgis_data:/var/lib/postgresql/data

volumes:
  postgis_data:

本地安装验证

-- 检查扩展是否可用
SELECT name, default_version, installed_version 
FROM pg_available_extensions 
WHERE name LIKE 'postgis%';

-- 启用 PostGIS 扩展
CREATE EXTENSION IF NOT EXISTS postgis;
CREATE EXTENSION IF NOT EXISTS postgis_topology;
CREATE EXTENSION IF NOT EXISTS postgis_raster;

-- 验证安装
SELECT PostGIS_Version();

1.2 实战项目：城市地理信息系统

第二部分：数据准备与导入

2.1 创建项目数据库结构

-- 创建项目数据库
CREATE DATABASE city_gis;
\c city_gis;

-- 启用扩展
CREATE EXTENSION postgis;

-- 创建城市表
CREATE TABLE cities (
    id SERIAL PRIMARY KEY,
    name VARCHAR(100) NOT NULL,
    population INTEGER,
    area_sq_km DECIMAL(10,2),
    elevation INTEGER,
    geom GEOMETRY(Point, 4326),
    created_at TIMESTAMP DEFAULT NOW()
);

-- 创建行政区表
CREATE TABLE districts (
    id SERIAL PRIMARY KEY,
    city_id INTEGER REFERENCES cities(id),
    name VARCHAR(100) NOT NULL,
    population_density DECIMAL(10,2),
    geom GEOMETRY(Polygon, 4326)
);

-- 创建道路表
CREATE TABLE roads (
    id SERIAL PRIMARY KEY,
    name VARCHAR(100),
    type VARCHAR(20), -- highway, primary, secondary, etc.
    length_km DECIMAL(8,2),
    geom GEOMETRY(LineString, 4326)
);

-- 创建兴趣点表
CREATE TABLE pois (
    id SERIAL PRIMARY KEY,
    name VARCHAR(100) NOT NULL,
    category VARCHAR(50), -- restaurant, hospital, school, etc.
    address TEXT,
    geom GEOMETRY(Point, 4326)
);

2.2 数据导入实战

从 CSV 导入点数据

-- 准备示例数据 cities.csv
-- name,population,area_sq_km,elevation,latitude,longitude
-- Beijing,21540000,16410,44,39.9042,116.4074
-- Shanghai,24280000,6340,4,31.2304,121.4737

-- 创建临时表
CREATE TEMP TABLE temp_cities (
    name VARCHAR(100),
    population INTEGER,
    area_sq_km DECIMAL(10,2),
    elevation INTEGER,
    latitude DECIMAL(9,6),
    longitude DECIMAL(9,6)
);

-- 导入 CSV
\copy temp_cities FROM 'cities.csv' WITH CSV HEADER;

-- 转换并插入数据
INSERT INTO cities (name, population, area_sq_km, elevation, geom)
SELECT 
    name,
    population,
    area_sq_km,
    elevation,
    ST_SetSRID(ST_MakePoint(longitude, latitude), 4326)
FROM temp_cities;

-- 验证数据
SELECT name, ST_AsText(geom) FROM cities;

从 GeoJSON 导入

-- 假设有 districts.geojson 文件
CREATE OR REPLACE FUNCTION import_geojson_districts()
RETURNS void AS $$
DECLARE
    geojson_text TEXT;
    feature JSON;
BEGIN
    -- 读取 GeoJSON 文件（实际应用中可能从文件或API读取）
    geojson_text := '{
        "type": "FeatureCollection",
        "features": [
            {
                "type": "Feature",
                "properties": {"name": "海淀区", "density": 125.5},
                "geometry": {
                    "type": "Polygon",
                    "coordinates": [[[116.2,39.9],[116.3,39.9],[116.3,40.0],[116.2,40.0],[116.2,39.9]]]
                }
            }
        ]
    }';
    
    -- 解析并插入数据
    INSERT INTO districts (name, population_density, geom)
    SELECT 
        feature->'properties'->>'name',
        (feature->'properties'->>'density')::DECIMAL,
        ST_GeomFromGeoJSON(feature->'geometry')
    FROM json_array_elements(geojson_text::json->'features') AS feature;
    
END;
$$ LANGUAGE plpgsql;

-- 执行导入
SELECT import_geojson_districts();

第三部分：空间索引与性能优化

3.1 创建空间索引

-- 为所有空间列创建索引
CREATE INDEX idx_cities_geom ON cities USING GIST (geom);
CREATE INDEX idx_districts_geom ON districts USING GIST (geom);
CREATE INDEX idx_roads_geom ON roads USING GIST (geom);
CREATE INDEX idx_pois_geom ON pois USING GIST (geom);

-- 创建复合索引
CREATE INDEX idx_pois_category_geom ON pois (category, geom);
CREATE INDEX idx_cities_population_geom ON cities (population, geom);

3.2 性能优化实战

-- 查看索引使用情况
EXPLAIN (ANALYZE, BUFFERS) 
SELECT name FROM cities 
WHERE ST_DWithin(geom, ST_MakePoint(116.4, 39.9), 1.0);

-- 空间查询性能测试函数
CREATE OR REPLACE FUNCTION benchmark_spatial_query()
RETURNS TABLE(query_name TEXT, execution_time INTERVAL, rows_returned BIGINT) AS $$
BEGIN
    RETURN QUERY EXECUTE '
        WITH queries AS (
            SELECT ''邻近城市查询'' as name, clock_timestamp() as start_time,
            (SELECT COUNT(*) FROM cities WHERE ST_DWithin(geom, ST_MakePoint(116.4, 39.9), 2.0)) as row_count
            UNION ALL
            SELECT ''区域包含查询'', clock_timestamp(),
            (SELECT COUNT(*) FROM pois WHERE ST_Within(geom, (SELECT geom FROM districts WHERE name = ''海淀区'')))
        )
        SELECT 
            name,
            clock_timestamp() - start_time as duration,
            row_count
        FROM queries';
END;
$$ LANGUAGE plpgsql;

-- 运行性能测试
SELECT * FROM benchmark_spatial_query();

第四部分：实战查询与分析

4.1 基础空间查询

-- 1. 查找最近的城市
SELECT name, 
       ST_Distance(geom, ST_MakePoint(116.4, 39.9)) as distance_deg
FROM cities 
ORDER BY geom <-> ST_MakePoint(116.4, 39.9)
LIMIT 5;

-- 2. 查找指定半径内的兴趣点
SELECT p.name, p.category,
       ST_Distance(p.geom, c.geom) as distance_from_city
FROM pois p
JOIN cities c ON c.name = 'Beijing'
WHERE ST_DWithin(p.geom, c.geom, 0.5)  -- 约50公里半径
ORDER BY distance_from_city;

-- 3. 计算行政区的面积和周长
SELECT 
    name,
    ST_Area(geom::geography) / 1000000 as area_sq_km,  -- 转换为平方公里
    ST_Perimeter(geom::geography) / 1000 as perimeter_km  -- 转换为公里
FROM districts;

4.2 高级空间分析

-- 1. 缓冲区分析：查找道路附近的兴趣点
SELECT DISTINCT p.name, p.category, r.name as road_name,
       ST_Distance(p.geom::geography, r.geom::geography) as distance_meters
FROM pois p
CROSS JOIN LATERAL (
    SELECT name, geom
    FROM roads 
    WHERE ST_DWithin(p.geom, geom, 0.01)  -- 约1公里
    ORDER BY ST_Distance(p.geom, geom)
    LIMIT 1
) r
WHERE ST_Distance(p.geom::geography, r.geom::geography) < 1000;  -- 1公里内

-- 2. 空间聚合：统计每个行政区的兴趣点数量
SELECT 
    d.name,
    COUNT(p.id) as poi_count,
    ST_Area(d.geom::geography) / 1000000 as area_sq_km,
    ROUND(COUNT(p.id) / (ST_Area(d.geom::geography) / 1000000), 2) as density_per_sq_km
FROM districts d
LEFT JOIN pois p ON ST_Within(p.geom, d.geom)
GROUP BY d.id, d.name, d.geom
ORDER BY density_per_sq_km DESC;

-- 3. 路径分析：查找两点间的最短路径（简化版）
WITH start_point AS (
    SELECT geom FROM pois WHERE name = '起点' LIMIT 1
),
end_point AS (
    SELECT geom FROM pois WHERE name = '终点' LIMIT 1
),
nearest_roads AS (
    SELECT r.id, r.name, r.geom
    FROM roads r, start_point s, end_point e
    WHERE ST_DWithin(r.geom, s.geom, 0.02)  -- 约2公里
       OR ST_DWithin(r.geom, e.geom, 0.02)
)
SELECT ST_Length(ST_ShortestLine(s.geom, e.geom)::geography) as straight_line_distance
FROM start_point s, end_point e;

4.3 地理编码实战

-- 创建地理编码函数
CREATE OR REPLACE FUNCTION geocode_address(
    address_text TEXT,
    city_name TEXT DEFAULT NULL
)
RETURNS TABLE(
    matched_address TEXT,
    latitude DECIMAL(9,6),
    longitude DECIMAL(9,6),
    confidence_score INTEGER
) AS $$
BEGIN
    RETURN QUERY 
    SELECT 
        p.name as matched_address,
        ST_Y(p.geom) as latitude,
        ST_X(p.geom) as longitude,
        CASE 
            WHEN p.name ILIKE '%' || address_text || '%' THEN 100
            WHEN p.address ILIKE '%' || address_text || '%' THEN 80
            ELSE 50
        END as confidence_score
    FROM pois p
    WHERE (p.name ILIKE '%' || address_text || '%' 
           OR p.address ILIKE '%' || address_text || '%')
      AND (city_name IS NULL OR EXISTS (
          SELECT 1 FROM cities c 
          WHERE c.name = city_name AND ST_DWithin(p.geom, c.geom, 1.0)
      ))
    ORDER BY confidence_score DESC
    LIMIT 1;
    
    -- 如果没有找到精确匹配，返回城市中心
    IF NOT FOUND THEN
        RETURN QUERY 
        SELECT 
            city_name as matched_address,
            ST_Y(c.geom) as latitude,
            ST_X(c.geom) as longitude,
            10 as confidence_score
        FROM cities c
        WHERE c.name = city_name
        LIMIT 1;
    END IF;
END;
$$ LANGUAGE plpgsql;

-- 使用地理编码
SELECT * FROM geocode_address('天安门', 'Beijing');

第五部分：Web 应用集成

5.1 RESTful API 设计

-- 创建空间查询 API 函数
CREATE OR REPLACE FUNCTION api_nearby_pois(
    lat DECIMAL,
    lng DECIMAL,
    radius_km DECIMAL DEFAULT 10,
    categories TEXT[] DEFAULT NULL,
    limit_count INTEGER DEFAULT 50
)
RETURNS TABLE(
    id INTEGER,
    name TEXT,
    category TEXT,
    address TEXT,
    latitude DECIMAL,
    longitude DECIMAL,
    distance_km DECIMAL
) AS $$
DECLARE
    search_point GEOMETRY := ST_SetSRID(ST_MakePoint(lng, lat), 4326);
BEGIN
    RETURN QUERY
    SELECT 
        p.id,
        p.name,
        p.category,
        p.address,
        ST_Y(p.geom) as latitude,
        ST_X(p.geom) as longitude,
        ST_Distance(p.geom::geography, search_point::geography) / 1000 as distance_km
    FROM pois p
    WHERE ST_DWithin(p.geom::geography, search_point::geography, radius_km * 1000)
      AND (categories IS NULL OR p.category = ANY(categories))
    ORDER BY p.geom <-> search_point
    LIMIT limit_count;
END;
$$ LANGUAGE plpgsql;

-- 测试 API 函数
SELECT * FROM api_nearby_pois(39.9042, 116.4074, 5.0, ARRAY['restaurant', 'hotel'], 10);

5.2 GeoJSON 输出

-- 生成 GeoJSON 格式数据
CREATE OR REPLACE FUNCTION api_geojson_pois(
    bounds_box BOX2D DEFAULT NULL,
    category_filter TEXT DEFAULT NULL
)
RETURNS JSON AS $$
DECLARE
    result JSON;
BEGIN
    SELECT json_build_object(
        'type', 'FeatureCollection',
        'features', json_agg(
            json_build_object(
                'type', 'Feature',
                'id', p.id,
                'geometry', ST_AsGeoJSON(p.geom)::json,
                'properties', json_build_object(
                    'name', p.name,
                    'category', p.category,
                    'address', p.address
                )
            )
        )
    ) INTO result
    FROM pois p
    WHERE (bounds_box IS NULL OR p.geom && bounds_box)
      AND (category_filter IS NULL OR p.category = category_filter);
    
    RETURN COALESCE(result, '{"type":"FeatureCollection","features":[]}'::json);
END;
$$ LANGUAGE plpgsql;

-- 使用示例
SELECT api_geojson_pois(ST_MakeBox2D(ST_Point(116.3, 39.9), ST_Point(116.5, 40.1)), 'restaurant');

第六部分：高级应用场景

6.1 实时位置追踪系统

-- 创建位置追踪表
CREATE TABLE vehicle_locations (
    id SERIAL PRIMARY KEY,
    vehicle_id INTEGER NOT NULL,
    geom GEOMETRY(Point, 4326),
    speed_kmh DECIMAL(5,2),
    timestamp TIMESTAMP DEFAULT NOW(),
    created_at TIMESTAMP DEFAULT NOW()
);

CREATE INDEX idx_vehicle_locations_geom ON vehicle_locations USING GIST (geom);
CREATE INDEX idx_vehicle_locations_timestamp ON vehicle_locations (timestamp);
CREATE INDEX idx_vehicle_locations_vehicle ON vehicle_locations (vehicle_id, timestamp);

-- 查找附近车辆
CREATE OR REPLACE FUNCTION find_nearby_vehicles(
    center_lat DECIMAL,
    center_lng DECIMAL,
    radius_km DECIMAL DEFAULT 5.0,
    max_age_minutes INTEGER DEFAULT 10
)
RETURNS TABLE(
    vehicle_id INTEGER,
    latitude DECIMAL,
    longitude DECIMAL,
    speed_kmh DECIMAL,
    last_update INTERVAL,
    distance_km DECIMAL
) AS $$
BEGIN
    RETURN QUERY
    SELECT 
        vl.vehicle_id,
        ST_Y(vl.geom) as latitude,
        ST_X(vl.geom) as longitude,
        vl.speed_kmh,
        NOW() - vl.timestamp as last_update,
        ST_Distance(
            vl.geom::geography, 
            ST_SetSRID(ST_MakePoint(center_lng, center_lat), 4326)::geography
        ) / 1000 as distance_km
    FROM vehicle_locations vl
    WHERE vl.timestamp > NOW() - (max_age_minutes || ' minutes')::INTERVAL
      AND ST_DWithin(
          vl.geom::geography, 
          ST_SetSRID(ST_MakePoint(center_lng, center_lat), 4326)::geography,
          radius_km * 1000
      )
    ORDER BY vl.timestamp DESC, distance_km;
END;
$$ LANGUAGE plpgsql;

6.2 空间数据质量监控

-- 数据质量检查函数
CREATE OR REPLACE FUNCTION check_spatial_data_quality()
RETURNS TABLE(
    table_name TEXT,
    check_type TEXT,
    issue_count INTEGER,
    details TEXT
) AS $$
BEGIN
    -- 检查无效几何
    RETURN QUERY
    SELECT 
        'cities' as table_name,
        '无效几何' as check_type,
        COUNT(*) as issue_count,
        string_agg(name || ': ' || ST_IsValidReason(geom), '; ') as details
    FROM cities 
    WHERE NOT ST_IsValid(geom)
    UNION ALL
    
    -- 检查空几何
    SELECT 
        'districts' as table_name,
        '空几何' as check_type,
        COUNT(*) as issue_count,
        string_agg(name, ', ') as details
    FROM districts 
    WHERE geom IS NULL OR ST_IsEmpty(geom)
    UNION ALL
    
    -- 检查坐标范围
    SELECT 
        'pois' as table_name,
        '坐标超出范围' as check_type,
        COUNT(*) as issue_count,
        '经纬度超出合理范围' as details
    FROM pois 
    WHERE ST_X(geom) < -180 OR ST_X(geom) > 180 
       OR ST_Y(geom) < -90 OR ST_Y(geom) > 90;
END;
$$ LANGUAGE plpgsql;

-- 运行数据质量检查
SELECT * FROM check_spatial_data_quality();

第七部分：部署与维护

7.1 自动化维护任务

-- 创建维护函数
CREATE OR REPLACE FUNCTION perform_spatial_maintenance()
RETURNS TABLE(task TEXT, result TEXT) AS $$
BEGIN
    -- 更新空间索引统计信息
    RETURN QUERY SELECT '分析表' as task, '完成' as result FROM (SELECT ANALYZE cities) t;
    RETURN QUERY SELECT '分析表' as task, '完成' as result FROM (SELECT ANALYZE districts) t;
    RETURN QUERY SELECT '分析表' as task, '完成' as result FROM (SELECT ANALYZE roads) t;
    RETURN QUERY SELECT '分析表' as task, '完成' as result FROM (SELECT ANALYZE pois) t;
    
    -- 修复无效几何（如果存在）
    UPDATE cities SET geom = ST_MakeValid(geom) WHERE NOT ST_IsValid(geom);
    IF FOUND THEN
        RETURN QUERY SELECT '修复几何' as task, '修复了无效几何' as result;
    END IF;
    
    -- 清理旧的位置数据
    DELETE FROM vehicle_locations WHERE timestamp < NOW() - INTERVAL '30 days';
    IF FOUND THEN
        RETURN QUERY SELECT '清理数据' as task, '清理了旧位置数据' as result;
    END IF;
END;
$$ LANGUAGE plpgsql;

-- 设置定时任务（需要 pg_cron 扩展）
-- SELECT cron.schedule('0 2 * * *', 'SELECT perform_spatial_maintenance()');

7.2 监控和日志

-- 创建查询日志表
CREATE TABLE spatial_query_log (
    id SERIAL PRIMARY KEY,
    query_type TEXT,
    parameters JSONB,
    execution_time INTERVAL,
    rows_returned INTEGER,
    executed_at TIMESTAMP DEFAULT NOW()
);

-- 包装查询函数
CREATE OR REPLACE FUNCTION logged_spatial_query(
    function_name TEXT,
    parameters JSONB
)
RETURNS JSONB AS $$
DECLARE
    start_time TIMESTAMP := clock_timestamp();
    result JSONB;
    row_count INTEGER;
BEGIN
    EXECUTE 'SELECT ' || function_name || '($1)' 
    INTO result 
    USING parameters;
    
    GET DIAGNOSTICS row_count = ROW_COUNT;
    
    INSERT INTO spatial_query_log 
    (query_type, parameters, execution_time, rows_returned)
    VALUES 
    (function_name, parameters, clock_timestamp() - start_time, row_count);
    
    RETURN result;
END;
$$ LANGUAGE plpgsql;

总结

1. 环境搭建 - Docker 和本地安装
2. 数据管理 - 表设计、数据导入、质量监控
3. 空间分析 - 查询、缓冲区、聚合分析
4. 性能优化 - 索引、查询优化
5. 应用集成 - API 设计、GeoJSON 输出
6. 高级应用 - 实时追踪、地理编码
7. 运维管理 - 自动化维护、监控