BimAnt

EPSG、WKT和PROJ

admin -

在之前的教程中，我们介绍了什么是坐标参考系统 (CRS)、坐标参考系统的组成部分以及投影坐标参考系统和地理坐标参考系统之间的一般差异。 在这个教程中，我们将介绍 CRS 信息的不同存储方式。

1、CRS的格式

有多种格式可用于记录 CRS。 三种常见格式包括：

  • proj.4
  • EPSG
  • Well-known Text(WKT)

通常，我们有一种格式的 CRS 信息，需要将该 CRS 转换为另一种格式，以便在 Python 等工具中使用。 因此，最好熟悉一些可能会遇到的关键格式。

查找 CRS 字符串的最强大的网站之一是 Spatialreference.org。 你可以使用网站上的搜索来查找 EPSG 代码。 一旦找到与感兴趣的 CRS 关联的页面，就可以查看与该 CRS 关联的所有各种格式：EPSG 4326 - WGS84 geographic。

2、PROJ 或 PROJ.4 字符串

PROJ.4 字符串是识别空间或坐标参考系统的一种紧凑方式。 PROJ.4 字符串是 Geopandas 可以接受的格式之一。 但是，请注意许多Python库正在转向更简洁的 EPSG 格式。

使用 PROJ.4 语法，你可以指定完整的参数集，包括定义特定 CRS 的椭圆、基准、投影单位和投影定义。

以下是 proj.4 字符串的示例：

+proj=utm +zone=11 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0

请注意，crs 信息是使用 + 号组合的一串字符和数字构成的。 上述数据的 CRS 为 proj4 格式。 该字符串包含 Python 或其他 GIS 可能需要的所有单个 CRS 元素。 每个元素都用 + 号指定，类似于 .csv 文件由  , 分隔的方式。 在每个 + 之后，我们看到 CRS 元素被定义。 例如 +proj= 和 +datum=。

可以将 proj4 字符串分解为其单独的组件（同样由 + 号分隔），如下所示：

  • +proj=utm：投影是UTM，UTM有几个区域。
  • +zone=11：zone为11，是美国西海岸的zone。
  • datum=WGS84：基准WGS84（基准是指投影中使用的坐标系的0,0参考）
  • +units=m：坐标的单位是米。
  • +ellps=WGS84：数据的椭球体（地球圆度的计算方式）为WGS84

请注意，该区域对于 UTM 投影是唯一的。 并非所有 CRS 都有区域。

另请注意，虽然加利福尼亚位于赤道上方 - 在北半球 - 区域之后没有 N（即11N），在 UTM proj4 规范中明确指定S来表示南半球，但是如果没有 S，那么你可以假设 这是一个北半球投影。

接下来，看看另一个 CRS 定义。

+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0

这是纬度/经度或地理投影。 proj4 字符串的组成部分分解如下。

  • +proj=longlat：数据在地理（经纬度）坐标系中
  • datum=WGS84：基准WGS84（基准是指投影中使用的坐标系的0,0参考）
  • +ellps=WGS84：椭球体（地球圆度的计算方式）为WGS84

请注意，上面没有指定单位。 这是因为这个地理坐标参考系统是纬度和经度，通常以十进制度数记录。

小技巧：每个 proj4 字符串的最后一部分是 +towgs84=0,0,0 。 这是一个转换系数，在需要基准转换时使用。

3、EPSG代码

EPSG 代码是代表 CRS 定义的 4-5 位数字。 首字母缩略词 EPSG 来自现已解散的欧洲石油调查组（European Petroleum Survey Group）。 每个代码都是一个四位五位数字，代表一个特定的 CRS 定义。点击这里查看EPSG代码。

导入在本模块中使用的 worldBoundary 层以探索 CRS。

import os
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import geopandas as gpd
from shapely.geometry import Point
import earthpy as et

# Set working dir & get data
data = et.data.get_data('spatial-vector-lidar')
os.chdir(os.path.join(et.io.HOME, 'earth-analytics'))

# Import world boundary shapefile
worldBound_path = os.path.join("data", "spatial-vector-lidar", "global", 
                               "ne_110m_land", "ne_110m_land.shp")
worldBound = gpd.read_file(worldBound_path)

worldBound.crs

结果如下：

## <Geographic 2D CRS: EPSG:4326>
Name: WGS 84
Axis Info [ellipsoidal]:
- Lat[north]: Geodetic latitude (degree)
- Lon[east]: Geodetic longitude (degree)
Area of Use:
- name: World
- bounds: (-180.0, -90.0, 180.0, 90.0)
Datum: World Geodetic System 1984
- Ellipsoid: WGS 84
- Prime Meridian: Greenwich

请注意，上面返回的 CRS 由两部分组成：

  • “init”告诉 python 将提供 CRS 定义（即 EPSG 代码），并且
  • epsg 代码本身 epsg：4326

4、如何在 Python 中创建 CRS 对象

经常需要为空间对象定义 CRS。 例如，在前面的课程中，我们创建了新的空间点层，并且必须定义点 x,y 位置所在的 CRS。

为此，我们完成了以下步骤：

  • 为单个点 (x,y) 手动创建了一个数组。
  • 把那个 x,y 点变成了一个有形状的点对象
  • 最后将该点对象转换为 pandas GeoDataFrame
# Create a numpy array with x,y location of Boulder
boulder_xy = np.array([[476911.31, 4429455.35]])

# Create shapely point object
boulder_xy_pt = [Point(xy) for xy in boulder_xy]

# Convert to spatial dataframe - geodataframe -- assign the CRS using epsg code
boulder_loc = gpd.GeoDataFrame(boulder_xy_pt,
                               columns=['geometry'],
                               crs={'init': 'epsg:2957'})

# View crs of new spatial points object
boulder_loc.crs

结果如下：

## <Projected CRS: EPSG:2957>
Name: NAD83(CSRS) / UTM zone 13N
Axis Info [cartesian]:
- E[east]: Easting (metre)
- N[north]: Northing (metre)
Area of Use:
- name: Canada - 108°W to 102°W
- bounds: (-108.0, 48.99, -102.0, 84.0)
Coordinate Operation:
- name: UTM zone 13N
- method: Transverse Mercator
Datum: NAD83 Canadian Spatial Reference System
- Ellipsoid: GRS 1980
- Prime Meridian: Greenwic

5、WKT 或Well-known Text

鉴于许多工具，包括 ESRI 的 ArcMap 和 ENVI 使用WKT这种格式，识别这种格式很有用。 Well-known Text (WKT) 是一种紧凑的机器和人类可读的几何对象表示。 它使用方括号 [] 和以逗号 (,) 分隔的元素的组合来定义坐标参考系统 (CRS) 定义的元素。

以下是 WGS84 地理的 WKT 示例：

GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]]

请注意，这里使用全部大写来明确描述元素 - 例如：

  • UNIT
  • DATUM

有时 WKT 结构化 CRS 信息嵌入在元数据文件中——类似于下面的结构：

GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0,
        AUTHORITY["EPSG","8901"]],
    UNIT["degree",0.01745329251994328,
        AUTHORITY["EPSG","9122"]],
    AUTHORITY["EPSG","4326"]]

原文链接：Understand EPSG, WKT and Other CRS Definition Styles

BimAnt翻译整理，转载请标明出处