1、什么是坐标参考系
要定义某物的位置，您通常会使用坐标系。 该系统由位于 2（或更多）维空间内的 X 和 Y 值组成。您使用具有 X、Y（有时还有 Z 轴）的坐标系来定义对象在空间中的位置。
虽然上面的坐标系是二维的，但你生活在一个恰好是“圆”的三维地球上。 要定义物体在圆形地球上的位置，您需要一个适应地球形状的坐标系。 当您在纸上或平面计算机屏幕上制作地图时，您从 3 维空间（地球）移动到 2 维空间（您的计算机屏幕或一张纸）。 CRS 的组件定义了如何“扁平化”存在于 3-D 地球空间中的数据。 CRS 还定义了坐标系本身。
CRS 定义了圆形地球上的一个位置与扁平的二维坐标系上的同一位置之间的转换。
2、CRS 的组成部分
坐标参考系统由几个关键部分组成：
坐标系：数据叠加的 X、Y 网格以及定义点在空间中的位置的方式。
水平和垂直单位：用于沿 x、y（和 z）轴定义网格的单位。
基准面：地球形状的模型版本，它定义了用于在空间中放置坐标系的原点。 您将在下面进一步解释这一点。
投影信息：用于使圆形表面（例如地球）上的物体变平的数学方程式，以便您可以在平坦表面（例如您的计算机屏幕或纸质地图）上查看它们。
3、为什么 CRS 很重要
了解您的数据使用的坐标系非常重要 - 特别是当您处理存储在不同坐标系中的不同数据时。 如果您有来自同一位置的数据存储在不同的坐标参考系中，它们将不会在任何 GIS 或其他程序中排列，除非您有像 ArcGIS 或 QGIS 这样支持动态投影的程序。 即使您使用支持动态投影的工具工作，您也会希望所有数据都在同一投影中以执行分析和处理任务。
（1）坐标系和单位
您可以使用 x 和 y 值定义空间位置，例如绘图位置，类似于上图中显示的笛卡尔坐标系。
例如，下图在地理坐标参考系统中显示了世界上所有的大陆。 单位是度，坐标系本身是纬度和经度，原点是赤道与地球中央子午线 (0,0) 相交的位置。
接下来，您将通过探索一些数据来了解有关 CRS 的更多信息。
import os
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from matplotlib.ticker import ScalarFormatter
import seaborn as sns
import geopandas as gpd
from shapely.geometry import Point
import earthpy as et
# Adjust plot font sizes
sns.set(font_scale=1.5)
sns.set_style("white")
# Set working dir & get data
data = et.data.get_data('spatial-vector-lidar')
os.chdir(os.path.join(et.io.HOME, 'earth-analytics'))首先，使用 geopandas 加载一个 shapefile。
# 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)（2）数据绘制
# Plot worldBound data using geopandas
fig, ax = plt.subplots(figsize=(10, 5))
worldBound.plot(color='darkgrey',
ax=ax)
# Set the x and y axis labels
ax.set(xlabel="Longitude (Degrees)",
ylabel="Latitude (Degrees)",
title="Global Map - Geographic Coordinate System - WGS84 Datum\n Units: Degrees - Latitude / Longitude")
# Add the x y graticules
ax.set_axisbelow(True)
ax.yaxis.grid(color='gray',
linestyle='dashed')
ax.xaxis.grid(color='gray',
linestyle='dashed')
地理坐标系 WGS84 基准中的全球地图。
（3）创建空间点对象
接下来，将三个坐标位置添加到您的地图。 请注意，单位是十进制度数（纬度、经度）：
科罗拉多州博尔德：40.0274，-105.2519
挪威奥斯陆：59.9500、10.7500
西班牙马略卡岛：39.6167、2.9833
要在空间上绘制这些点，您将
创建点位置的 numpy 数组
使用 for 循环填充 shapely Point 对象
# Create numpy array of x,y point locations
add_points = np.array([[-105.2519, 40.0274],
[ 10.75 , 59.95 ],
[ 2.9833, 39.6167]])
# Turn points into list of x,y shapely points
city_locations = [Point(xy) for xy in add_points]
city_locations
[<shapely.geometry.point.Point at 0x7fb0b90558b0>,
<shapely.geometry.point.Point at 0x7fb0b9055910>,
<shapely.geometry.point.Point at 0x7fb0b9055970>]
# Create geodataframe using the points
city_locations = gpd.GeoDataFrame(city_locations,
columns=['geometry'],
crs=worldBound.crs)# Create geodataframe using the points
city_locations = gpd.GeoDataFrame(city_locations,
columns=['geometry'],
crs=worldBound.crs)
city_locations.head(3)最后你可以在你的世界地图上绘制点.
# Plot point locations
fig, ax = plt.subplots(figsize=(12, 8))
worldBound.plot(figsize=(10, 5), color='k',
ax=ax)
# Add city locations
city_locations.plot(ax=ax,
color='springgreen',
marker='*',
markersize=45)
# Setup x y axes with labels and add graticules
ax.set(xlabel="Longitude (Degrees)", ylabel="Latitude (Degrees)",
title="Global Map - Geographic Coordinate System - WGS84 Datum\n Units: Degrees - Latitude / Longitude")
ax.set_axisbelow(True)
ax.yaxis.grid(color='gray', linestyle='dashed')
ax.xaxis.grid(color='gray', linestyle='dashed')
地理坐标参考系统中的全球地图，点位置覆盖在顶部。
接下来，导入可以投影到特定坐标参考系中的适当经纬网。
# Import graticule & world bounding box shapefile data
graticule_path = os.path.join("data", "spatial-vector-lidar", "global",
"ne_110m_graticules_all", "ne_110m_graticules_15.shp")
graticule = gpd.read_file(graticule_path)
bbox_path = os.path.join("data", "spatial-vector-lidar", "global",
"ne_110m_graticules_all", "ne_110m_wgs84_bounding_box.shp")
bbox = gpd.read_file(bbox_path)
# Create map axis object
fig, ax = plt.subplots(1, 1, figsize=(15, 8))
# Add bounding box and graticule layers
bbox.plot(ax=ax, alpha=.1, color='grey')
graticule.plot(ax=ax, color='lightgrey')
worldBound.plot(ax=ax, color='black')
# Add points to plot
city_locations.plot(ax=ax,
markersize=60,
color='springgreen',
marker='*')
# Add title and axes labels
ax.set(title="World Map - Geographic Coordinate Reference System (long/lat degrees)",
xlabel="X Coordinates (meters)",
ylabel="Y Coordinates (meters)");
带有经纬网的地理坐标参考系统中的全球地图。
4、地理 CRS优缺点
当您需要定位地球上的位置时，采用十进制度数的地理坐标系非常有用。 但是，纬度和经度位置并不是使用统一的测量单位来定位的。 因此，地理 CRS 不是测量距离的理想选择。 这就是开发其他预计 CRS 的原因。
地理坐标系使用角度定位纬度和经度位置。 因此每条纬线南北移动的间距是不均匀的。
（1）Projected CRS-Robinson
您可以在另一个 CRS - Robinson 中查看上面的相同数据。 Robinson是一个预计的 CRS。 请注意，地图上的国家/地区边界与您在 CRS 中创建的地图相比具有不同的形状：地理纬度/经度 WGS84。
下面您首先将数据重新投影到 robinson 项目 (+proj=robin) 中。 然后你再次绘制数据。
# Reproject the data
worldBound_robin = worldBound.to_crs('+proj=robin')
graticule_robin = graticule.to_crs('+proj=robin')
# Plot the data
fig, ax = plt.subplots(figsize=(12, 8))
worldBound_robin.plot(ax=ax,
color='k')
graticule_robin.plot(ax=ax, color='lightgrey')
ax.set(title="World Map: Robinson Coordinate Reference System",
xlabel="X Coordinates (meters)",
ylabel="Y Coordinates (meters)")
for axis in [ax.xaxis, ax.yaxis]:
formatter = ScalarFormatter()
formatter.set_scientific(False)
axis.set_major_formatter(formatter)
如果将上面使用的相同纬度/经度坐标位置添加到地图中，会发生什么情况？ 请记住，地图上的数据位于CRS - Robinson
# Plot the data
fig, ax = plt.subplots(1, 1, figsize=(12, 8))
worldBound_robin.plot(ax=ax,
color='k')
graticule_robin.plot(ax=ax,
color='lightgrey')
city_locations.plot(ax=ax,
marker='*',
color='springgreen',
markersize=100)
ax.set(title="World Map: Robinson Coordinate Reference System",
xlabel="X Coordinates (meters)",
ylabel="Y Coordinates (meters)")
for axis in [ax.xaxis, ax.yaxis]:
formatter = ScalarFormatter()
formatter.set_scientific(False)
axis.set_major_formatter(formatter)
plt.axis('equal'); 
如果您将地理坐标参考系统 (WGS84) 中的点叠加在使用 Robinson 投影的地图之上，请注意它们没有正确对齐。
请注意，当您尝试将以度为单位的纬度/经度坐标添加到不同 CRS 中的地图时，这些点不在正确的位置。 您需要先将这些点转换为其他数据所在的同一 CRS。将数据集从一个 CRS 转换为另一个 CRS 的过程通常称为重新投影。
在 python 中，您使用 .to_crs 方法重新投影数据。
# Reproject point locations to the Robinson projection
city_locations_robin = city_locations.to_crs(worldBound_robin.crs)
fig, ax = plt.subplots(1, 1, figsize=(12, 8))
worldBound_robin.plot(ax=ax,
cmap='Greys')
ax.set(title="World map (robinson)",
xlabel="X Coordinates (meters)",
ylabel="Y Coordinates (meters)")
city_locations_robin.plot(ax=ax, markersize=100, color='springgreen')
for axis in [ax.xaxis, ax.yaxis]:
formatter = ScalarFormatter()
formatter.set_scientific(False)
axis.set_major_formatter(formatter)
plt.axis('equal');
如果您将您的点重新投影到 Robinson CRS，那么它们会在它们记录的空间位置中很好地绘制在底图的顶部。
（2）地图比较
上面的两个图在视觉上看起来不同，并且使用了不同的坐标系。 并排查看地图上呈现的实际经纬度线或纬度线。
# Reproject graticules and bounding box to robinson
graticule_robinson = graticule.to_crs('+proj=robin')
bbox_robinson = bbox.to_crs('+proj=robin')
# Setup plot with 2 "rows" one for each map and one column
fig, (ax0, ax1) = plt.subplots(2, 1, figsize=(13, 12))
# First plot
bbox.plot(ax=ax0,
alpha=.1,
color='grey')
graticule.plot(ax=ax0,
color='lightgrey')
worldBound.plot(ax=ax0,
color='k')
city_locations.plot(ax=ax0,
markersize=100,
color='springgreen')
ax0.set(title="World Map - Geographic (long/lat degrees)")
# Second plot
bbox_robinson.plot(ax=ax1,
alpha=.1,
color='grey')
graticule_robinson.plot(ax=ax1,
color='lightgrey')
worldBound_robin.plot(ax=ax1,
color='k')
city_locations_robin.plot(ax=ax1,
markersize=100,
color='springgreen')
ax1.set(title="World Map Projected - Robinson (Meters)")
for axis in [ax1.xaxis, ax1.yaxis]:
formatter = ScalarFormatter()
formatter.set_scientific(False)
axis.set_major_formatter(formatter)
请注意，这两张地图上的 x 轴和 y 轴单位完全不同。 这就是为什么地理 WGS84 CRS 中的点在绘制在另一个 CRS（如 Robinson）的地图上时不能正确排列的原因之一。 但是，如果您重新投影数据，那么它们将正确排列。
（3）地理 CRS vs. 投影CRS
上述地图提供了两种主要坐标系的示例：
地理坐标系：跨越整个地球的坐标系（例如纬度/经度）。
投影坐标系：局部坐标系以最小化特定区域的视觉失真（例如 Robinson、UTM、State Plane）
1万+
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
