计算球面上任意两点间的球面距离（C++实现）
Posted 没事就要敲代码
tags:
篇首语：本文由小常识网(cha138.com)小编为大家整理，主要介绍了计算球面上任意两点间的球面距离（C++实现）相关的知识，希望对你有一定的参考价值。
1 预备知识
在求解此问题之前首先要明确一下几点：
（1）两点间的球面距离： 球面上两点间的最短距离，即球心与球面上两点所确定的平面与球面相交，得到一截面圆，两点间的劣弧（圆心角<180°的圆弧）就是这两点间的球面距离。
（2）任意与球面相交的平面，球截面均为一个圆，不可能出现其他截面形状。
（4） 弦长公式
角度制下：
已知圆弧的半径 R R R，圆心角 n ° n° n°，则弧长计算公式为
l = n ° π R 180 ° l=\\cfrac{n°\\pi R}{180°} l=180°n°πR
弧度制下：
已知圆弧的半径
R
R
R，圆心角（弧度）
θ
\\theta
θ，则弧长计算公式为
l
=
R
⋅
θ
l=R·\\theta
l=R⋅θ
2 原理描述
已知球面上两点
S
1
(
x
1
,
y
1
,
z
1
)
S_1(x_1,y_1,z_1)
S1(x1,y1,z1)、
S
2
(
x
2
,
y
2
,
z
2
)
S_2(x_2,y_2,z_2)
S2(x2,y2,z2)与球心
O
(
x
0
,
y
0
,
z
0
)
O(x_0,y_0,z_0)
O(x0,y0,z0)，可以唯一确定一个平面，且三点围成一个扇形
⌔
O
S
1
S
2
⌔OS_1S_2
⌔OS1S2，半径为球半径
R
R
R。
根据余弦定理，可求得扇形的圆心角 α \\alpha α
c
o
s
α
=
2
R
2
−
∣
S
1
S
2
∣
2
2
R
2
=
1
−
∣
S
1
S
2
∣
2
2
R
2
cos\\alpha=\\cfrac{2R^2-|S_1S_2|^2}{2R^2}=1-\\cfrac{|S_1S_2|^2}{2R^2}
cosα=2R22R2−∣S1S2∣2=1−2R2∣S1S2∣2
α
=
a
r
c
c
o
s
(
1
−
∣
S
1
S
2
∣
2
2
R
2
)
\\alpha=arccos(1-\\cfrac{|S_1S_2|^2}{2R^2})
α=arccos(1−2R2∣S1S2∣2)
其中，
∣
S
1
S
2
∣
|S_1S_2|
∣S1S2∣ 为两点构成的弦长，且
∣
S
1
S
2
∣
=
(
x
2
−
x
1
)
2
+
(
y
2
−
y
1
)
2
+
(
z
2
−
z
1
)
2
|S_1S_2|=\\sqrt{(x_2-x_1)^2+(y_2-y_1)^2+(z_2-z_1)^2}
∣S1S2∣=(x2−x1)2+(y2−y1)2+(z2−z1)2
根据上面给出的弧长公式，即可计算出球面上任意两点间的球面距离。
3 代码实现
以球心位于坐标原点 O ( 0 , 0 , 0 ) O(0,0,0) O(0,0,0)，球半径 R = 1.0 R=1.0 R=1.0 的球面为例，已知球面上两点 p 1 ( 1 , 0 , 0 ) p_1(1,0,0) p1(1,0,0)， p 2 ( 0 , 0 , 1 ) p_2(0,0,1) p2(0,0,1)，求这两点的球面距离。
代码：
#include <iostream>
#include <vector>
#include <cmath>
using namespace std;
//XYZ点类型结构体
struct PointXYZ
{
public:
double x;
double y;
double z;
};
//球面距离求解类
class SphereDistanceCalculator
{
public:
/**
* @brief ：设置球半径
* @param[I]：R（球半径）
* @param[O]：none
* @return ：none
* @note ：
**/
void setSphereRadius(double R);
/**
* @brief ：设置球面两点
* @param[I]：p1（球面第一点）
* @param[I]：p2（球面第二点）
* @param[O]：none
* @return ：none
* @note ：
**/
void setSpherePoints(PointXYZ p1, PointXYZ p2);
/**
* @brief ：计算球面距离
* @param[I]：none
* @param[O]：none
* @return ：vector<double>，依次为弦长、圆心角、球面距离
* @note ：
**/
vector<double> calculateSphereDistance();
private:
double m_R; //球半径
bool is_setSphereRadius = false;
PointXYZ m_p1; //球面第一点
PointXYZ m_p2; //球面第二点
bool is_setSpherePoints = false;
};
/**
* @brief ：设置球半径
* @param[I]：R（球半径）
* @param[O]：none
* @return ：none
* @note ：
**/
void SphereDistanceCalculator::setSphereRadius(double R)
{
if (R > 0)
{
m_R = R;
is_setSphereRadius = true;
}
else
{
cerr << "\\a->球体半径应为一个正数！\\n";
system("pause");
abort();
}
}
/**
* @brief ：设置球面两点
* @param[I]：p1（球面第一点）
* @param[I]：p2（球面第二点）
* @param[O]：none
* @return ：none
* @note ：
**/
void SphereDistanceCalculator::setSpherePoints(PointXYZ p1, PointXYZ p2)
{
m_p1 = p1;
m_p2 = p2;
is_setSpherePoints = true;
}
/**
* @brief ：计算球面距离
* @param[I]：none
* @param[O]：none
* @return ：vector<double>，依次为弦长、圆心角、球面距离
* @note ：
**/
vector<double> SphereDistanceCalculator::calculateSphereDistance()
{
if (!is_setSphereRadius)
{
cerr << "\\a->请输入球半径！\\n";
system("pause");
abort();
}
if (!is_setSpherePoints)
{
cerr << "\\a->请输入球面上两点三维坐标！\\n";
system("pause");
abort();
}
vector<double> res; //存放返回值的向量
double S1S2; //弦长
S1S2 = sqrt(pow((m_p2.x - m_p1.x), 2) + pow((m_p2.y - m_p1.y), 2) + pow((m_p2.z - m_p1.z), 2));
res.push_back(S1S2);
double alpha; //圆心角
alpha = acos(1 - pow(S1S2 / m_R, 2) / 2);
res.push_back(alpha);
double sphereDistance; //球面距离
sphereDistance = m_R * alpha;
res.push_back(sphereDistance);
return res;
}
int main()
{
const double PI = 3.14159265359; //常量PI
PointXYZ p1, p2; //球面两点坐标
p1.x = 1.0;
p1.y = 0.0;
p1.z = 0.0;
p2.x = 0.0;
p2.y = 0.0;
p2.z = 1.0;
SphereDistanceCalculator sdc; //创建球面距离求解对象
sdc.setSphereRadius(1.0); //设置球半径
sdc.setSpherePoints(p1关于已知两点经纬度求球面最短距离的公式推导