计算球面上任意两点间的球面距离（C++实现）

1 预备知识

在求解此问题之前首先要明确一下几点：

（1）两点间的球面距离： 球面上两点间的最短距离，即球心与球面上两点所确定的平面与球面相交，得到一截面圆，两点间的劣弧（圆心角<180°的圆弧）就是这两点间的球面距离。

（2）任意与球面相交的平面，球截面均为一个圆，不可能出现其他截面形状。

（3）余弦定理

（4） 弦长公式

角度制下：

已知圆弧的半径 $R$，圆心角 $n°$，则弧长计算公式为

$$l=\frac{n°\pi R}{180°}$$

弧度制下：

已知圆弧的半径 $R$，圆心角（弧度） $\theta$，则弧长计算公式为
$$l=R\cdot \theta$$

2 原理描述

已知球面上两点 $S_1(x_1,y_1,z_1)$、$S_2(x_2,y_2,z_2)$与球心 $O(x_0,y_0,z_0)$，可以唯一确定一个平面，且三点围成一个扇形 $⌔O{S}_1{S}_2$，半径为球半径 $R$。

根据余弦定理，可求得扇形的圆心角 $\alpha$

$$cos\alpha =\frac{2R^2-|S_1{S}_2{|}^2}{2R^2}=1-\frac{|S_1{S}_2{|}^2}{2R^2}$$
$$\alpha =arccos(1-\frac{|S_1{S}_2{|}^2}{2R^2})$$

其中，$|S_1{S}_2|$ 为两点构成的弦长，且
$$|S_1{S}_2|=\sqrt{(x_2-x_1{)}^2+(y_2-y_1{)}^2+(z_2-z_1{)}^2}$$
根据上面给出的弧长公式，即可计算出球面上任意两点间的球面距离。

3 代码实现

以球心位于坐标原点 $O(0,0,0)$，球半径 $R=1.0$ 的球面为例，已知球面上两点 $p_1(1,0,0)$，$p_2(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, p2);		//设置球面上两点
	vector<double> rec;					//存放弦长、圆心角、球面距离的向量
	rec = sdc.calculateSphereDistance();//执行计算，并将结果存放到rec中

	cout << "->球面两点的直线距离（圆弧弦长）：" << rec[0] << endl;
	cout << "->球面两点对应的圆心角（弧度）：" << rec[1] << endl
		 << "  球面两点对应的圆心角（角度）：" << rec[1] / PI * 180 << "°" << endl;
	cout << "->球面两点的球面距离：" << rec[2] << endl;

	return 0;
}

输出结果：

->球面两点的直线距离（圆弧弦长）：1.41421
->球面两点对应的圆心角（弧度）：1.5708
  球面两点对应的圆心角（角度）：90°
->球面两点的球面距离：1.5708