第5章 Python与地理信息系统

本章主要学习Python处理矢量数据，包含以下内容：

距离测量
坐标转换
矢量数据重投影
Shapefile 文件编辑
海量数据过滤
专题地图创建
非GIS数据类型转换
地理化编码

本篇博文记录第5.1距离测量

5.1 距离测量

地理学第一定律：“任何事物都相关，相近的事物关联更紧密”。因此距离测量显得格外重要。
地球是不规则椭圆状，计算距离时有3种地球模型可供选择：

- 平面

地球被当做没有曲率的平面，只需处理直线即可，可以利用地图投影把地球模型压扁转换成二维平面，利用十进制度数来表示，十进制度数系统会被当做包含x和y的笛卡尔直角坐标系，测量距离则就可以很方便的使用欧式几何的勾股定理了。

但是求得距离的误差比较大。

- 球体
使用一个完美的球体来表达地球的特征，但是地球实际上是一个不规则的椭球体

- 椭球体

存在多种椭球体地球模型基准，其中一套基准就是一组地球形状定义的集合，也被称为大地测量坐标系统。和其他地理坐标系统类似，基准也可以对特定区域进行优化。最常用的基准是基于全球的***WGS84基准***（会优化，2004年是最新的版本），比正圆形球体模型更符合实际。

补充：大地水准面（大地水平面模型）：是指与平均海水面重合并延伸到大陆内部的水准面。

5.1.1勾股定理（也叫欧氏距离公式）——平面地球模型

例如两个点之间，计算距离

a2+b2=c2

import math
# 两个投影点，经密西西比横轴墨卡托投影
x1 = 456456.23
y1 = 1279721.064
x2 = 576628.34
y2 = 1071740.33x_dist = x1 - x2
y_dist = y1 - y2dist_sq = x_dist**2 + y_dist**2dist = math.sqrt(dist_sq)print(dist)

240202.668047278
大该是240202m，也就是240.2km

使用角度进行测量，必须将角度转换为弧，才能在坐标系统之间计算曲面上的距离。


import math
#经纬度坐标
x1 = -90.21
y1 = 32.31
x2 = -88.95
y2 = 30.43x_dist = math.radians(x1 - x2)
y_dist = math.radians(y1 - y2)
dist_squre = x_dist**2 + y_dist**2
dis_rad = math.sqrt(dist_squre)
# 6371251.46m 地球半径
print("弧度是：%f"% dis_rad)
print("距离是：%f"% (dis_rad*6371251.46))

弧度是：0.039500
距离是：251664.466877，大概是251km，比第一次结果大了11km
因此不同测量方式（不用的坐标系统）和地球模型的结果可能存在很大差异

5.1.2 半正矢量公式——球体地球模型

大圆距离（也叫球面距离）：指球体表面上两点之间最短距离（根据两点的经度和纬度来确定大圆上两点之间距离的计算方法）。

球面上任意两点A和B以及球心都可以确定唯一的大圆，这个大圆被称为黎曼圆，而在大圆上连接这两点的较短的一条弧的长度就是大圆距离。

若这两点和球心正好都在球的直径上，则过这三点可以有无数大圆，但两点之间的
弧长都相等，因此该大圆能够将球体等分。

最佳的十进制度数距离测量方法：半正矢量公式。
半正矢公式：根据经度纬度计算两点间距离

import math
x1 = -90.212452861859035
y1 = 32.316272202663704
x2 = -88.952170968942525
y2 = 30.438559624660321
# Convert angle x from degrees to radians
x_dist = math.radians(x1-x2)
y_dist = math.radians(y1-y2)
y1_rad = math.radians(y1)
y2_rad = math.radians(y2)
a = math.sin(y_dist/2)**2 + math.sin(x_dist/2)**2 * math.cos(y1_rad) * math.cos(y2_rad)
c = 2 * math.asin(math.sqrt(a))
# 6371251.46m 地球半径
# 弧长 = 弧度*半径
distance = c * 6371 # kilometers
print(distance)

240.6359762909508
半正矢公式计算结果是240.6km，与使用优化过的重投影方法获得的240.2km相比，它们之间的误差小于0.5km。
半正矢公式是最常用的距离计算公式，因为它相对来说代码量较少并且大部分情况下足够精确，它的精度范围能够控制在1m以内

目前遇到的大部分坐标点都是未经投影变换的十进制度数格式。所以，从距离
测量的角度来看，有以下几个可供选择的方案：

重投影到一个精确的笛卡儿坐标系统中并计算距离；
根据你的精度要求只使用半正矢公式；
使用精度更高的Vincenty公式。

5.1.3Vincenty——椭球体地球模型

Vincenty公式是基于椭球体地球模型的,如果你使用的是一个本地化的椭球体模型，
那么它的计算精度可以远远小于1m。
可以自定义半长轴和扁平率来匹配任何权威的椭球体模型。


import math
distance = None
x1 = -90.212452861859035
y1 = 32.316272202663704
x2 = -88.952170968942525
y2 = 30.438559624660321
#定义椭球体参数，例如NAD83
a = 6378137 #半长轴
f = 1/298.257222101 #逆扁平度
b = abs((f*a)-a) #半短轴
L   =  math.radians(x2-x1)
U1 = math.atan((1 - f) * math.tan(math.radians(y1)))
U2 = math.atan((1 - f) * math.tan(math.radians(y2)))
sinU1 = math.sin(U1)
cosU1 = math.cos(U1)
sinU2 = math.sin(U2)
cosU2 = math.cos(U2)
lam = L
for i in range(100):sinLam = math.sin(lam)cosLam = math.cos(lam)sinSigma = math.sqrt((cosU2 * sinLam) ** 2 +(cosU1 * sinU2 - sinU1 * cosU2 * cosLam) ** 2)if (sinSigma == 0):# 重合点distance = 0breakcosSigma = sinU1 * sinU2 + cosU1 * cosU2 * cosLamsigma = math.atan2(sinSigma, cosSigma)sinAlpha = cosU1 * cosU2 * sinLam / sinSigmacosSqAlpha = 1 - sinAlpha ** 2cos2SigmaM = cosSigma - 2 * sinU1 * sinU2 / cosSqAlphaif math.isnan(cos2SigmaM):cos2SigmaM = 0#   赤道线C = f / 16 * cosSqAlpha * (4 + f * (4 - 3 * cosSqAlpha))LP = lamlam = L + (1 - C) * f * sinAlpha *\(sigma + C * sinSigma * (cos2SigmaM + C * cosSigma *(-1+2 * cos2SigmaM * cos2SigmaM)))if not abs(lam - LP) > 1e-12:breakuSq = cosSqAlpha * (a ** 2 - b ** 2) / b ** 2A = 1 + uSq / 16384 * (4096 + uSq * (-768 + uSq * (320 - 175 * uSq)))B = uSq / 1024 * (256 + uSq * (-128 + uSq * (74 - 47 * uSq)))deltaSigma = B*sinSigma*(cos2SigmaM+B/4*(cosSigma*(-1+2*cos2SigmaM*cos2SigmaM)-B/6*cos2SigmaM*(-3+4*sinSigma*sinSigma)*(-3+4*cos2SigmaM*cos2SigmaM)))s = b * A * (sigma - deltaSigma)distance = sprint(distance)

240091.4562741062
240091.4562741062
240091.4562741062
240091.4562741062
240091.4562741062
240091.4562741062
…
Vincenty公式的计算结果大约是240.1km，这个数据和之前经投影变换之后的欧氏距离测量公式得到的结果的误差只有100m左右。

完全使用Python实现的geopy模块包含一个Vincenty公式的实现，并且它还支持地理位置编码，能够将地名转换为经纬度坐标。它的详情可以参考如下页
面：http://geopy.readthedocs.org/en/latest/。

总结

《Python地理空间分析指南第2版》学习笔记，仅供学习，如有侵权请联系删除。

《Python地理空间分析指南第2版》学习笔记-5.1 距离测量相关推荐

【笔记】《Python地理空间分析指南(第2版)》
转载地址:https://blog.csdn.net/jianbinzheng/article/details/80215228 概述部分地理空间数据地理空间技术概览 Python地理空间分析工具 ...
python地理空间分析指南pdf邓世超_Python地理空间分析指南（第2版）源代码.zip
[实例简介] Python地理空间分析指南(第2版)的随书源代码,需要的朋友可以下载一下~~ [实例截图] [核心代码] Python地理空间分析指南(第2版)源代码 └── Python地理空间分析 ...
python空间分析_读书笔记——《python地理空间分析指南》
本文为<Python地理空间分析指南(第2版)>的读书摘录,顺便挖个坑,进一步对python的几个包做学习整理. 本笔记的用途:了解python地理空间处理的技术框架和实现途径. 第三章 ...
《Python地理空间分析指南（第2版）》——1.9　地理信息系统基本概念
本节书摘来自异步社区<Python地理空间分析指南(第2版)>一书中的第1章,第1.9节,作者: [美]Joel Lawhead(莱哈德) 更多章节内容可以访问云栖社区"异步社区 ...
Python地理空间分析指南（第2版）学习笔记01
目录前言一.任务二.实现与解析 1.引入库 2.构造数据模型 3.渲染地图元素 4.执行查询操作以及完成绘图三.总结前言本书假定读者了解Pyhon.信息技术的基本知识,并且至少对地理空间分 ...
Python 地理空间分析
前文我们将快速浏览 Python 的(空间)数据科学生态系统,并了解如何使用一些基本的开源 Python 包,例如: pandas / geopandas shapely pysal pyproj ...
python网络爬虫权威指南百度云-分析《Python网络爬虫权威指南第2版》PDF及代码...
对那些没有学过编程的人来说,计算机编程看着就像变魔术.如果编程是魔术(magic),那么网页抓取(Web scraping)就是巫术(wizardry),也就是运用"魔术"来实现精 ...
一直在构建工作空间_国际资讯Python与地理空间分析
点击图片上方蓝色字体"慧天地"即可订阅英文原文来源:www.gislounge.com 英文原文链接:https://www.gislounge.com/python-and-g ...
Turf.js——用于地理空间分析的js库，处理各种地图算法
Turf.js--用于地理空间分析的js库,处理各种地图算法一.官网中文--https://turfjs.fenxianglu.cn/ 英文--https://turfjs.org/ npm地址- ...

《Python地理空间分析指南第2版》学习笔记-5.1 距离测量

第5章 Python与地理信息系统

5.1 距离测量

5.1.1勾股定理（也叫欧氏距离公式）——平面地球模型

5.1.2 半正矢量公式——球体地球模型

5.1.3Vincenty——椭球体地球模型

总结

《Python地理空间分析指南第2版》学习笔记-5.1 距离测量相关推荐

最新文章

热门文章

《Python地理空间分析指南 第2版》学习笔记-5.1 距离测量

第5章 Python与地理信息系统

5.1 距离测量

5.1.1勾股定理（也叫欧氏距离公式）——平面地球模型

5.1.2 半正矢量公式——球体地球模型

5.1.3Vincenty——椭球体地球模型

总结

《Python地理空间分析指南 第2版》学习笔记-5.1 距离测量相关推荐

最新文章

热门文章

《Python地理空间分析指南第2版》学习笔记-5.1 距离测量

《Python地理空间分析指南第2版》学习笔记-5.1 距离测量相关推荐