查找地球坐标(经纬度),距离(米)和方位(角)
Finding Earth Coordinates (latitude,longitude), Distance (meters) and Bearing (angle)
我需要以各种方式处理地球坐标。在C/c++中没有一个函数可以直接做到这一点。
参考以下问题:
- Python:给定当前点,距离和轴承
- C:给定点(纬度,经度),距离和方位,如何得到新的纬度和经度
从第一个和活字脚本网站,我发现下面是公式:
查找两个坐标之间的方位(角度)
x = cos(lat1Rad)*sin(lat2Rad) - sin(lat1Rad)*cos(lat2Rad)*cos(lon2Rad-lon1Rad);
y = sin(lon2Rad-lon1Rad) * cos(lat2Rad);
bearing = atan2(y, x); // In radians;
// Convert to degrees and for -ve add 360
求两个坐标之间的距离(米)
PI = 3.14159265358979323846, earthDiameterMeters = 2*6371*1000;
x = sin((lat2Rad-lat1Rad) / 2);
y = sin((lon2Rad-lon1Rad) / 2);
meters = earthDiameterMeters * asin(sqrt(x*x + y*y*cos(lat1Rad)*cos(lat2Rad)));
从坐标+距离+角度查找坐标
meters *= 2 / earthDiameterMeters;
lat2Rad = asin(sin(lat1Rad)*cos(meters) + cos(lat1Rad)*sin(meters)*cos(bearing));
lon2Rad = lon1Rad + atan2(sin(bearing)*sin(meters)*cos(lat1Rad),
cos(meters) - sin(lat1Rad)*sin(lat2Rad));
下面的伪代码应相互验证以上3个方程:
struct Coordinate { double lat, lon; } c1, c2;
auto degree = FindBearing(c1, c2);
auto meters = FindDistance(c1, c2);
auto cX = FindCoordiante(c1, degree, meters);
现在实际上答案是几乎接近但不是正确的。即cX不等于c2!0.0005
的经度值总是存在差异。例如
c1 = (12.968460,77.641308)
c2 = (12.967862,77.653130)
angle = 92.97 ^^^
distance = 1282.74
cX = (12.967862,77.653613)
^^^
我没有太多的数学知识。但我知道的是,从fcc.gov网站上,答案总是正确的。
我做错了什么?
仅供参考的代码
虽然语法是c++的,但所有的数学函数都来自C,并且很容易移植到C中(因此标记为两者)
#include<iostream>
#include<iomanip>
#include<cmath>
// Source: http://www.movable-type.co.uk/scripts/latlong.html
static const double PI = 3.14159265358979323846, earthDiameterMeters = 6371.0 * 2 * 1000;
double degreeToRadian (const double degree) { return (degree * PI / 180); };
double radianToDegree (const double radian) { return (radian * 180 / PI); };
double CoordinatesToAngle (double latitude1,
const double longitude1,
double latitude2,
const double longitude2)
{
const auto longitudeDifference = degreeToRadian(longitude2 - longitude1);
latitude1 = degreeToRadian(latitude1);
latitude2 = degreeToRadian(latitude2);
using namespace std;
const auto x = (cos(latitude1) * sin(latitude2)) -
(sin(latitude1) * cos(latitude2) * cos(longitudeDifference));
const auto y = sin(longitudeDifference) * cos(latitude2);
const auto degree = radianToDegree(atan2(y, x));
return (degree >= 0)? degree : (degree + 360);
}
double CoordinatesToMeters (double latitude1,
double longitude1,
double latitude2,
double longitude2)
{
latitude1 = degreeToRadian(latitude1);
longitude1 = degreeToRadian(longitude1);
latitude2 = degreeToRadian(latitude2);
longitude2 = degreeToRadian(longitude2);
using namespace std;
auto x = sin((latitude2 - latitude1) / 2), y = sin((longitude2 - longitude1) / 2);
#if 1
return earthDiameterMeters * asin(sqrt((x * x) + (cos(latitude1) * cos(latitude2) * y * y)));
#else
auto value = (x * x) + (cos(latitude1) * cos(latitude2) * y * y);
return earthDiameterMeters * atan2(sqrt(value), sqrt(1 - value));
#endif
}
std::pair<double,double> CoordinateToCoordinate (double latitude,
double longitude,
double angle,
double meters)
{
latitude = degreeToRadian(latitude);
longitude = degreeToRadian(longitude);
angle = degreeToRadian(angle);
meters *= 2 / earthDiameterMeters;
using namespace std;
pair<double,double> coordinate;
coordinate.first = radianToDegree(asin((sin(latitude) * cos(meters))
+ (cos(latitude) * sin(meters) * cos(angle))));
coordinate.second = radianToDegree(longitude
+ atan2((sin(angle) * sin(meters) * cos(latitude)),
cos(meters) - (sin(latitude) * sin(coordinate.first))));
return coordinate;
}
int main ()
{
using namespace std;
const auto latitude1 = 12.968460, longitude1 = 77.641308,
latitude2 = 12.967862, longitude2 = 77.653130;
cout << std::setprecision(10);
cout << "(" << latitude1 << "," << longitude1 << ") --- "
"(" << latitude2 << "," << longitude2 << ")n";
auto angle = CoordinatesToAngle(latitude1, longitude1, latitude2, longitude2);
cout << "Angle = " << angle << endl;
auto meters = CoordinatesToMeters(latitude1, longitude1, latitude2, longitude2);
cout << "Meters = " << meters << endl;
auto coordinate = CoordinateToCoordinate(latitude1, longitude1, angle, meters);
cout << "Destination = (" << coordinate.first << "," << coordinate.second << ")n";
}
在CoordinateToCoordinate
中,您使用已经以度表示的sin(coordinate.first)
。使用sin(degreeToRadian(coordinate.first))
。
或者更干净:
... CoordinateToCoordinate (...)
{
...
coordinate.first = asin((sin(latitude) * cos(meters))
+ (cos(latitude) * sin(meters) * cos(angle)));
coordinate.second = longitude + atan2((sin(angle) * sin(meters) * cos(latitude)),
cos(meters) - (sin(latitude) * sin(coordinate.first)));
coordinate.first = radianToDegree(coordinate.first);
coordinate.second = radianToDegree(coordinate.second);
return coordinate;
}
这解决了这个问题。现场演示。
部分答案
将角度92.97°
转换为弧度后,调用sin/cos/tan
将有效地将角度更改为2.97°
。在度到弧度转换后以及在三角函数调用中发生周期缩减时,仅这一步就损失了6位精度。
具有较大角度的三角函数可以提高精度。使用幸运方向正好 360.0度在一个圆。执行"模数45°",可能使用remquo(angle, 45, &octant)
,并在使用radian参数调用三角函数之前将然后将转换为弧度。
示例sind()
77.641308和77.653130的答案相差约6600分之一(~13位精度)。这个答案可能不能完全解释这一点,但应该有所帮助。(如果某些地方出现了float
的用法,则应改为double
)
相关文章:
- 如何修复此错误:未定义对"距离(浮点数,浮点数,浮点数,浮点数,浮点数)"的引用
- 用C++程序计算圆锥体的体积、球体的体积、八边形的面积和两点之间的距离
- 计算所有对之间的曼哈顿距离
- 如何使用发送数据包所花费的时间计算两个节点之间的距离?
- Levenshtein 两个文件的距离花费了太多时间
- 两个有符号数字之间的距离
- 计算车辆之间的距离并设置速度,使距离保持不变,例如 5 米
- 计算两个迭代器之间的距离时"Vector Iterators Incompatible"
- 浮动距离为80位长双精度
- 查找彼此之间具有不同最小距离的随机数
- 相同字符之间的最小距离
- 如何根据距离矩阵确定节点?
- 如何排列二进制字符串以最小化它们之间的距离
- 曼哈顿距离比曼哈顿距离+线性冲突效果更好
- 使用C++中的链表对与原点的距离、x 和 y 坐标进行排序
- 应该如何编写用于计算最近点距离的C++函数?
- 有效地计算像素到其对应核线的距离
- 堆和堆栈之间的距离
- std::使用迭代器映射查找距离,程序不会终止
- 查找地球坐标(经纬度),距离(米)和方位(角)