基于三次多項式擬合三角函數(shù)的地理空間距離計算算法

李雪佳，封紅旗，梅　宇，趙玉寧

（1.常州大學(xué)信息科學(xué)與工程學(xué)院，江蘇常州　213000；2.江蘇省南京市鼓樓醫(yī)院信息科，南京　210000）

基于三次多項式擬合三角函數(shù)的地理空間距離計算算法

李雪佳1，封紅旗1，梅宇1，趙玉寧2

（1.常州大學(xué)信息科學(xué)與工程學(xué)院，江蘇常州213000；2.江蘇省南京市鼓樓醫(yī)院信息科，南京210000）

移動互聯(lián)網(wǎng)技術(shù)的快速發(fā)展，現(xiàn)有的APP提供了非常人性化的操作，用戶可以進(jìn)行對商家的篩選，APP默認(rèn)的選項“智能排序”，“離我最近”；這兩個選擇項，系統(tǒng)都會去計算當(dāng)前用戶的位置和各個商家之間的距離；這種大量計算距離的場景十分消耗資源，按照現(xiàn)在的統(tǒng)計數(shù)據(jù)100 W條數(shù)據(jù)時候，延遲將會達(dá)到140 ms；且隨著數(shù)據(jù)的增長，服務(wù)器的性能堪憂；當(dāng)前計算空間距離的方法是使用Lucene算法進(jìn)行計算，但是開銷比較大，不利于用戶一眼，文章提出了一種基于三次多項式擬合三角函數(shù)的優(yōu)化算法對空間距離進(jìn)行計算，文章從距離計算的準(zhǔn)確性和快速性兩個方面與原有方法進(jìn)行比較；最后總結(jié)出優(yōu)化算法的優(yōu)勢之處。

余弦函數(shù)；空間距離；三次多項式；準(zhǔn)確性

0　引言

互聯(lián)網(wǎng)技術(shù)發(fā)展越來越快，以團購為代表OTO業(yè)務(wù)的移動終端受到用戶的喜愛，OTO平臺最大的特點是基于本地化的服務(wù)，在各大OTO平臺中，都存在一個非常重要功能，商家的篩選，按距離、智能排序（智能排序中距離因素也是比較重要的因素）。其中對商家的篩選，主要依靠的是距離的因素。如果獲取商家與用戶之間的距離，是解決這一問題的關(guān)鍵所在。這里涉及到了地理空間距離計算的問題。在該算法中，其中的原理是比較簡單的，最關(guān)鍵之處在于對算法的優(yōu)化，在用戶進(jìn)行篩選時候系統(tǒng)可以在最快的時間給出最準(zhǔn)確的答案。這樣用戶的體驗性好。

本文提出的基于三次多項式擬合三角函數(shù)的地理空間距離計算算法，通過使用三次多項式來抵消原公式中三角函數(shù)的計算，使得算法在效率上提高很多。但是本文提出的算法是在同城范圍內(nèi)，主要是由于論文中約定了在同城范圍內(nèi)，經(jīng)緯度可以看作是相互垂直的直線。最后文章給出了初始方法、基于坐標(biāo)轉(zhuǎn)化方法、以及文中所提出的方法之間相互比較。

1.1原理

地理空間距離計算方法比較多，現(xiàn)在主要使用的是以下兩種模型：一個是球面模型，該模型將地球看做一個標(biāo)準(zhǔn)的球體，球面上任意兩點之間的距離表示了圓最大的弧長。這個方法使用的比較廣泛。另一個是橢球模型，這個模型更加貼近地球的真實環(huán)境，但是橢圓的計算方法非常復(fù)雜，但是在實際的工程中對精度的要求沒有那么高。

現(xiàn)在主要介紹基于球面模型的地理空間距離計算公式，

假設(shè)地球上有兩個點A（ja，wa），B（jb，wb），其中ja，jb表示A、B兩點的精度，wa，wb表示A、B兩點的維度。所以A、B兩點的距離就是AB之間的弧長，所以AB弧長=R*角AOB（注：角AOB是A跟B的夾角，O是地球的球心，R是地球半徑，約為6367000米）。為了求解角AOB的角度，可以先對AOB角的最大邊AB的長度進(jìn)行求解，然后根據(jù)余弦定理算出其夾角。

1）由經(jīng)緯度，地球半徑R，將A、B兩點的經(jīng)緯度轉(zhuǎn)化

圖1　球面模式示意圖

成球體的三維坐標(biāo)：

2）求出AB之間的長度：

AB2=（Xa-Xb）2+（Ya-Yb）2+（Za-Zb）2=...=

2R2（1-cos（wa）cos（wb）cos（jb-ja）-

sin（wa）sin（wb））3）求出角AOB：

4）AB的弧長：

AB=Rarc cos（cos（wa）cos（wb）-cos（ja-jb）+ sin（wa）sin（wb））

以上便是空間距離計算的原理。在這里只是簡單對如果進(jìn)行論述。

1.2使用Lucene進(jìn)行距離的計算

現(xiàn)有的APP后臺基本使用Lucene來對商家進(jìn)行篩選，Lucene使用了spatial4j工具包進(jìn)行對距離的計算，在spatial4j工具包中同時提供了基于多種球面的地理空間距離公式，上面描述的就是其中的一種。在這個工具包中有一中常用的工具Haversine公式，這也是spatial4j工具包中默認(rèn)的公式。引入這個公式的好處：前面的余弦公式中cos（jb-ja），當(dāng)系統(tǒng)的浮點運算精度比較低，在計算AB兩點距離比較近的時候產(chǎn)生的誤差比較大，Haversine方法中通過某種變化的方式對cos（jb-ja）進(jìn)行消除。這樣解決了近距離誤差較大的問題。

（1）Haversine公式調(diào)用：

public static double dist HaversineRAD（double lat1，double lon1，double lat2，double lon2）｛

double hsin X=Math.sin（（lon1-lon2）*0.5）；

double hsin Y=Math.sin（（lat1-lat2）*0.5）

double h=hsin Y*hsin Y+

（Math.cos（lat1）*Math.cos（lat2）*hsin X* hsinX）；

return 2*Math.atan2（Math.sqrt（h），Math.sqrt（1 -h））*6367000；

｝

（2）Haversine公式性能：

本文的測試環(huán)境處理器：2.9 GHz Intel Core i7，內(nèi)存為8 GB 1 600 MHz DDR3，操作系統(tǒng)為OSX10.8.3，實驗在單線程環(huán)境下運行。測試的數(shù)據(jù)等級5 w，10 w，100 w。測試結(jié)果如下圖，在下圖中可以發(fā)現(xiàn)當(dāng)數(shù)據(jù)只有5 w個時候，計算一遍距離7 ms，但當(dāng)數(shù)據(jù)增加到100 W的時候，系統(tǒng)只是進(jìn)行計算距離一項就需要144 ms，還不包括篩選過程中其他的運行項目，性能十分堪憂。

表1　Lucene性能分析

2　優(yōu)化方案

為了提高計算的效率，提供更好的用戶體驗，在實驗過程中進(jìn)行了抓棧實驗，在實現(xiàn)過程中發(fā)現(xiàn)了消耗cpu較多線程的地方大多數(shù)存在計算執(zhí)行距離公式中的三角函數(shù)。正是這個給我們提出了一些優(yōu)化方案，去消除或者消減三角函數(shù)。

2.1三次多項式擬合三角函數(shù)的地理空間距離計算算法

1）思路：

業(yè)務(wù)場景大多是僅限于在同一個城市進(jìn)行距離的計算，也就是說兩點之間的直線距離一般不會超過300 KM，針對于范圍比較小，則可以近似的認(rèn)為經(jīng)緯度是垂直的，要求A（116.8，39.78）和B（116.9，39.68）兩點的距離，我們可以先求出南北方向距離AM，然后求出東西方向距離BM，最后求矩形對角線距離，即sqrt（AMAM+BMBM）。

圖2　近距離下經(jīng)緯度示意圖

南北方向AM=R緯度差 Math.PI/180.0；

東西方向BM=R經(jīng)度差Cos＜當(dāng)?shù)鼐暥葦?shù)*Math.PI/ 180.0＞

這種方式僅僅需要計算一次cos函數(shù)。

public static double distanceSimplify（double lat1，double lng1，double lat2，double lng2，double［］a）｛

double dx=lng1-lng2；//經(jīng)度差值

double dy=lat1-lat2；//緯度差值

double b=（lat1+lat2）/2.0；//平均緯度

double Lx=toRadians（dx）*6367000.0*Math.cos（toRadians（b））；//東西距離

double Ly=6367000.0*toRadians（dy）；//南北距離

return Math.sqrt（Lx*Lx+Ly*Ly）；//用平面的矩形對角距離公式計算總距離

｝

｝

2）有效性驗證：

我們首先檢驗這種簡化是否能滿足我們應(yīng)用的精度，如果精度較差將不能用于實際生產(chǎn)環(huán)境。

我們的方法叫distanceSimplify，lucene的方法叫dist HaversineRAD。下表是在不同尺度下兩個方法的相差情況。

可以看到兩者在百米、千米尺度上幾乎沒有差別，在萬米尺度上也僅有分米的差別，此外由于我們的業(yè)務(wù)是在一個城市范圍內(nèi)進(jìn)行篩選排序，所以我們選擇了北京左下角和右上角兩點進(jìn)行比較，兩點相距有260多千米，兩個方法差別17 m。從精度上看該優(yōu)化方法能滿足我們應(yīng)用需求。

表2　改進(jìn)算法的有效性

3）性能測試：

表3　改進(jìn)性能分析

2.2進(jìn)一步優(yōu)化方案

在上面的過程中，可以看出進(jìn)行了一次的cos三角函數(shù)的運算，如何將三角函數(shù)完全的進(jìn)行消除，這樣便于進(jìn)一步提高優(yōu)化的效率。

在這里主要使用利用多項式來擬合三角函數(shù)，這是高等代數(shù)里面高斯不等式的相關(guān)知識。在一定的范圍內(nèi)，當(dāng)然等式的次數(shù)越高，獲得到的準(zhǔn)確率將越來越精確。在這里選擇了三次多項式進(jìn)行擬合。

這里對原有算法進(jìn)行改進(jìn)，中國的緯度范圍在10～60之間，即我們將此區(qū)間離散成Length份作為我們的訓(xùn)練集。

public static double［］trainPoly Fit（int degree，int Length）｛

PolynomialCurveFitter polynomialCurveFitter=PolynomialCurve-Fitter.create（degree）；

double minLat=10.0；//中國最低緯度

double max Lat=60.0；//中國最高緯度

double interv=（max Lat-min Lat）/（double）Length；

List＜WeightedObserved Point＞weightedObservedPoints=new Array List＜WeightedObservedPoint＞（）；

for（int i=0；i＜Length；i++）｛

WeightedObservedPoint weightedObservedPoint=newWeightedObservedPoint（1，minLat+（double）i*interv，Math.cos（toRadians（x［i］）））；

weighted ObservedPoints.add（weightedObservedPoint）；

｝

return polynomialCurveFitter.fit（weightedObservedPoints）；

｝

public static double distanceSimplify More（double lat1，double lng1，double lat2，double lng2，double［］a）｛

//1）計算3個參數(shù)

double dx=lng1-lng2；//經(jīng)度差值

double dy=lat1-lat2；//緯度差值

double b=（lat1+lat2）/2.0；//平均緯度

//2）計算東西方向距離和南北方向距離（單位：米），東西距離采用三階多項式

double Lx=（a［3］*b*b*b+a［2］*b*b+a［1］*b+a［0］）*toRadians（dx）*6367000.0；//東西距離

double Ly=6367000.0*toRadians（dy）；//南北距離

//3）用平面的矩形對角距離公式計算總距離

return Math.sqrt（Lx*Lx+Ly*Ly）；

｝

1）有效性驗證：

我們的優(yōu)化方法叫distanceSimplify More，lucene的方法叫dist HaversineRAD，下表是在不同尺度下兩個方法的相差情況。

表4　最終優(yōu)化的有效性

可以看到在百米尺度上兩者幾乎未有差別，在千米尺度上僅有分米的區(qū)別，在更高尺度上如72千米僅有5.6 M米別，在264千米也僅有8.1米區(qū)別，因此該優(yōu)化方法的精度能滿足我們的應(yīng)用需求。

2）性能驗證：

表5　最終改進(jìn)算法性能分析

3　三種方案的對比

在總結(jié)前，還將簡單論述了一種比較簡單的轉(zhuǎn)化方法，最后將文章中所提出的3種方法進(jìn)行比較詳細(xì)的總結(jié)：

我們的計算場景是計算用戶位置與所有篩選出來的商戶的距離，這里會涉及到大量三角函數(shù)計算。一種優(yōu)化思路是商戶數(shù)據(jù)不保存經(jīng)緯度而保存球面模型下的三維坐標(biāo)（x，y，z），映射方法如下：

x=Math.cos（lat）Math.cos（lon）；

y=Math.cos（lat）Math.sin（lon）；

z=Math.sin（lat）；

那么當(dāng)我們求夾角AOB時，只需要做一次點乘操作。比如求（lon1，lat1）和（lon2，lat2）的夾角，只需要計算x1x 2 +y 1y2+z1*z2，這樣避免了大量三角函數(shù)的計算。

在得到夾角之后，還需要執(zhí)行arccos函數(shù)，將其轉(zhuǎn)換成角度，AB弧長=角AOB*R（R是地球半徑）。上面的這種方法成為坐標(biāo)轉(zhuǎn)換方法。

坐標(biāo)轉(zhuǎn)換方法和基于三次多項式擬合三角函數(shù)方法性能都非常高，相比lucene使用的Haversine算法大大提高了計算效率，然而坐標(biāo)轉(zhuǎn)換方法存在一些缺點：

1）坐標(biāo)轉(zhuǎn)換后的數(shù)據(jù)不能被直接用于空間索引。lucene可以直接對經(jīng)緯度進(jìn)行g(shù)eohash空間索引，而通過空間轉(zhuǎn)換變成三維數(shù)據(jù)后不能直接使用。我們的應(yīng)用有附近范圍篩選功能（例如附近5 km的團購單子），通過geohash空間索引可以提高范圍篩選的效率；

2）坐標(biāo)轉(zhuǎn)換方法增大內(nèi)存開銷。我們會將坐標(biāo)寫入倒排索引中，之前坐標(biāo)是2列（經(jīng)度和緯度），現(xiàn)在變成3列（x，y，z），在使用中我們往往會將這數(shù)據(jù)放入到cache中，因此會增大內(nèi)存開銷；

3）坐標(biāo)轉(zhuǎn)換方法增大建索引開銷。此方法本質(zhì)上是將計算從查詢階段放至到索引階段，因此提高了建索引的開銷。

所以，文章中所提出的基于三次多項式擬合三角函數(shù)方法計算地球空間距離，在同城范圍內(nèi)的優(yōu)勢非常明顯，在效率和有效性方面想比較與以前的方法都有比較大的提升。這是今后研究的一個方向。

4　結(jié)束語

文中通過利用三次多項式擬合三角函數(shù)的方案改進(jìn)了原有的地理空間距離計算的算法，從實驗結(jié)果證明了，改進(jìn)算法所具備了搜索目標(biāo)點的快速性，同時也保證了相應(yīng)的準(zhǔn)確率。當(dāng)然這個實驗的結(jié)果是基于同城范圍的基礎(chǔ)上，完全符合OTO業(yè)務(wù)對于搜索場景的要求。同時方案中也存在有待改進(jìn)的地方，在方案中是用三次多項式去擬合三角函數(shù)，根據(jù)數(shù)學(xué)知識可知，多項式的次數(shù)越高擬合的三角函數(shù)越精確，但是隨之的計算的復(fù)雜度也在增加。所以如何在計算的復(fù)雜度和多項式擬合三角函數(shù)的精確度之間做出比較好的選擇，可以進(jìn)一步提高搜索的準(zhǔn)確率，這個方向是今后研究的重點方向。

［1］Hua M C，Lou D C，Chang M C.Dual-Wrapped digital watermarking scheme for image copyright protection［J］.Computers&Security，2007，16（1）：1-12.

［2］司少林，關(guān)永.三角函數(shù)曲線擬合最佳次數(shù)的確定［J］.計算機工程與設(shè)計，2006.24（1）：30-35.

［3］Xiao Y Q，Ji Q.A robust content-based digital image watermarking scheme［J］.Signal Processing，2007，7：1264-1280.

［4］陳娛，徐君.考慮地理距離的復(fù)雜網(wǎng)絡(luò)社區(qū)挖據(jù)算法［J］.地球信息科學(xué)，2013（3）：58-65.

［5］李耀彬，曾祥斌，沈鋮武.基于拋物線插值的正弦波擬合算法［J］.計算機工程與設(shè)計，2009（11）：100-105.

［6］薛國新，孫玉強.正弦曲線三點擬合問題的一種新方法［J］.計算機仿真，2006（2）：99-105.

［7］羅成漢，劉小山.曲線擬合法的Matlab實現(xiàn)［J］.現(xiàn)代電子技術(shù)，2003（20）：54-59.

［8］田崢，徐成，米超，等.基于消失點和主方向估計的道路分割算法［J］.計算機研究與發(fā)展，2014（4）：99-104.

［9］譚方勇，于復(fù)生，吳建平.基于消失點的坐標(biāo)校準(zhǔn)算法［J］.計算機應(yīng)用，2011（1）：101-105.

［10］羅小松，房斌，楊維斌.采用韋伯局部特征的道路消失點檢測［J］.計算機應(yīng)用，2014（S1）：219-222.

Geographical Space Distance Calculation Algorithm Based on Three Time Polynomial Fitting Trigonometric Function

Li Xuejia1，F(xiàn)eng Hongqi1，Mei Yu1，Zhao Yuning2
（1.School of information science&Engineering，Changzhou213000，China；2.Nanjing Gulou Hospital，Nanjing210000，China）

With the rapid development of mobile Internet technology，the existing APP provides a very user-friendly operation，the user can make a selection of businesses，APP default option“smart sort”，“from me recently.”These two options，the system will calculate the current position of the user and the distance between the various businesses.This is a large number of computing distance of the scene is very consumption of resources，in accordance with the current statistical data 100 W data when the delay will reach 140 ms.And with the growth of data，the performance of the server is worrying.At present，the method of calculating the space distance is to use Lucene algorithm to calculate，but the cost is relatively large，and is not conducive to the user one eye，the paper presents a method based on the three polynomial fitting trigonometric function optimization algorithm to calculate the space distance，the accuracy and speed of the distance from the two aspects of the original method.Finally，the advantages of the optimization algorithm are summarized

cosine function；space distance；three degree polynomial；accuracy

1671-4598（2016）05-0199-03

10.16526/j.cnki.11-4762/tp.2016.05.056

TP3

A

2015-11-09；

2015-12-15。

國家自然科學(xué)基金資助項目（61103172）。

李雪佳（1990-），湖北荊州人，碩士研究生，主要從事數(shù)據(jù)挖掘、地理信息系統(tǒng)。

封紅旗（1976-），江蘇常州人，博士，副教授，主要從事數(shù)據(jù)挖掘、地理信息系統(tǒng)方向的研究。