要在單位球面上隨機選取一個點,從 [-180,180) 中的均勻分布經度和 [-90,90] 中均勻分布的緯度中隨機選擇經緯度是不正確的,因為赤道和兩極的格點面積元並不相同,因此以這種方式選取的點將在極點附近“聚集”。(圖片引自 https://mathworld.wolfram.com/SpherePointPicking.html )

一個簡單的均勻采樣方法是:三維笛卡爾坐標系中生成隨機點,並投影到球面上,這樣生成的點才會均勻。即生成三個高斯隨機變量,x, y,和z,那么以下變量在三維球面上均勻分布的(Muller 1959, Marsaglia 1972):

在shell腳本中寫一個awk函數即可實現,以下代碼實現隨機選取兩個點,輸出它們的經緯度和大圓弧距離,輸出五列單位均為度。
1 seed=10000 2 sum=$(echo $seed*36/3.1416/3.1416 | bc) 3 awk 'BEGIN{srand(); 4 for ( i=0; i<='"$sum"'; i++ ){ 5 xa=2*rand()-1; 6 ya=2*rand()-1; 7 za=2*rand()-1; 8 xb=2*rand()-1; 9 yb=2*rand()-1; 10 zb=2*rand()-1; 11 norma=sqrt(xa*xa+ya*ya+za*za) 12 normb=sqrt(xb*xb+yb*yb+zb*zb) 13 if (norma<0.00001 || normb<0.00001) continue 14 if (xa*xa+ya*ya+za*za<1.0 && xb*xb+yb*yb+zb*zb<1.0){ 15 xa=xa/norma 16 ya=ya/norma 17 za=za/norma 18 xb=xb/normb 19 yb=yb/normb 20 zb=zb/normb 21 dis=sqrt((xb-xa)^2+(yb-ya)^2+(zb-za)^2) 22 lona=atan2(ya,xa)/0.017453 23 lonb=atan2(yb,xb)/0.017453 24 lata=90-atan2(sqrt(xa*xa+ya*ya),za)/0.017453 25 latb=90-atan2(sqrt(xb*xb+yb*yb),zb)/0.017453 26 arc=2*atan2(dis/2,sqrt(1-dis*dis/4))/0.017453 27 printf("%8.2f %7.2f %8.2f %7.2f %8.3f\n", 28 lona,lata,lonb,latb,arc) 29 } 30 } 31 }' > sphere_random_point.txt
備注:
seed代表最終生成文件大致行數
norma代表A點到球心距離
normb代表B點到球心距離
dis代表AB點之間弦長
arc代表AB點之間大圓弧長(度為單位)
References
Marsaglia, G. "Choosing a Point from the Surface of a Sphere." Ann. Math. Stat. 43, 645-646, 1972.
Muller, M. E. "A Note on a Method for Generating Points Uniformly on
-Dimensional Spheres." Comm. Assoc. Comput. Mach. 2, 19-20, Apr. 1959.
