基于geohash的周边地理位置搜索

1.场景

随着智能手机和传感器技术的发展,LBS(Location based service)类的应用也逐渐多了起来。几乎每一个软件或多或少都要牵扯上点位置的事。这篇文章主要来讲解下怎么快速的搜寻某个位置点周边的数据。
将问题模型化就是,给一组数据,每个数据包括了这个数据的位置信息(经纬度)和其他信息,给一个搜寻点(经度和纬度)和搜寻范围(radius),在最少的计算和时间内将在搜寻点的搜寻范围内的数据找出来。

2.目前的解决方案

(1)圆形区域一一匹配搜索。
将list中的数据一个一个同搜索位置进行距离计算。几乎没有人这样用,每一个数据都需要进行开根号操作,数据量大的时候不可想象。


(2)先对数据库中的经度(longitude)和纬度(latitude)添加B+树索引,然后搜索以radius的两倍为边长的正方形,如下图:


(3)使用目前现成的拥有地理位置搜索的数据库,比如mongoDB,HBase等等。

(4)使用geohash数据结构,见下面介绍

3.geohash介绍

(1)感性认识GeoHash
首先来点感性认识,http://openlocation.org/geohash/geohash-js/ 提供了在地图上显示geohash编码的功能。
1)GeoHash将二维的经纬度转换成字符串,比如下图展示了北京9个区域的GeoHash字符串,分别是WX4ER,WX4G2、WX4G3等等,每一个字符串代表了某一矩形区域。也就是说,这个矩形区域内所有的点(经纬度坐标)都共享相同的GeoHash字符串,这样既可以保护隐私(只表示大概区域位置而不是具体的点),又比较容易做缓存,比如左上角这个区域内的用户不断发送位置信息请求餐馆数据,由于这些用户的GeoHash字符串都是WX4ER,所以可以把WX4ER当作key,把该区域的餐馆信息当作value来进行缓存,而如果不使用GeoHash的话,由于区域内的用户传来的经纬度是各不相同的,很难做缓存。

3.png

2)字符串越长,表示的范围越精确。如图所示,5位的编码能表示10平方千米范围的矩形区域,而6位编码能表示更精细的区域(约0.34平方千米)


4.png

3)字符串相似的表示距离相近(特殊情况后文阐述),这样可以利用字符串的前缀匹配来查询附近的POI信息。如下两个图所示,一个在城区,一个在郊区,城区的GeoHash字符串之间比较相似,郊区的字符串之间也比较相似,而城区和郊区的GeoHash字符串相似程度要低些。


城区

郊区

通过上面的介绍我们知道了GeoHash就是一种将经纬度转换成字符串的方法,并且使得在大部分情况下,字符串前缀匹配越多的距离越近,回到我们的案例,根据所在位置查询来查询附近餐馆时,只需要将所在位置经纬度转换成GeoHash字符串,并与各个餐馆的GeoHash字符串进行前缀匹配,匹配越多的距离越近。

(2)GeoHash算法的步骤
下面以北海公园为例介绍GeoHash算法的计算步骤


7.png

2.1. 根据经纬度计算GeoHash二进制编码地球纬度区间是[-90,90], 北海公园的纬度是39.928167,可以通过下面算法对纬度39.928167进行逼近编码:1)区间[-90,90]进行二分为[-90,0),[0,90],称为左右区间,可以确定39.928167属于右区间[0,90],给标记为1;2)接着将区间[0,90]进行二分为 [0,45),[45,90],可以确定39.928167属于左区间 [0,45),给标记为0;3)递归上述过程39.928167总是属于某个区间[a,b]。随着每次迭代区间[a,b]总在缩小,并越来越逼近39.928167;4)如果给定的纬度x(39.928167)属于左区间,则记录0,如果属于右区间则记录1,这样随着算法的进行会产生一个序列1011100,序列的长度跟给定的区间划分次数有关。

首先将纬度范围(-90, 90)平分成两个区间(-90, 0)、(0, 90), 如果目标纬度位于前一个区间,则编码为0,否则编码为1。
由于39.92324属于(0, 90),所以取编码为1。
然后再将(0, 90)分成 (0, 45), (45, 90)两个区间,而39.92324位于(0, 45),所以编码为0。
以此类推,直到精度符合要求为止,得到纬度编码为1011 1000 1100 0111 1001。

image.png

经度也用同样的算法,对(-180, 180)依次细分,得到116.3906的编码为1101 0010 1100 0100 0100

image.png

2.2. 组码通过上述计算,纬度产生的编码为10111 00011,经度产生的编码为11010 01011。偶数位放经度,奇数位放纬度,把2串编码组合生成新串:11100 11101 00100 01111。最后使用用0-9、b-z(去掉a, i, l, o)这32个字母进行base32编码,首先将11100 11101 00100 01111转成十进制,对应着28、29、4、15,十进制对应的编码就是wx4g。同理,将编码转换成经纬度的解码算法与之相反,具体不再赘述。


8.png

(3)GeoHash Base32编码长度与精度
下表摘自维基百科:http://en.wikipedia.org/wiki/Geohash可以看出,当geohash base32编码长度为8时,精度在19米左右,而当编码长度为9时,精度在2米左右,编码长度需要根据数据情况进行选择。

9.png

(4)GeoHash算法
上文讲了GeoHash的计算步骤,仅仅说明是什么而没有说明为什么?为什么分别给经度和维度编码?为什么需要将经纬度两串编码交叉组合成一串编码?本节试图回答这一问题。如图所示,我们将二进制编码的结果填写到空间中,当将空间划分为四块时候,编码的顺序分别是左下角00,左上角01,右下脚10,右上角11,也就是类似于Z的曲线,当我们递归的将各个块分解成更小的子块时,编码的顺序是自相似的(分形),每一个子快也形成Z曲线,这种类型的曲线被称为Peano空间填充曲线。这种类型的空间填充曲线的优点是将二维空间转换成一维曲线(事实上是分形维),对大部分而言,编码相似的距离也相近, 但Peano空间填充曲线最大的缺点就是突变性,有些编码相邻但距离却相差很远,比如0111与1000,编码是相邻的,但距离相差很大。


10.png

除Peano空间填充曲线外,还有很多空间填充曲线,如图所示,其中效果公认较好是Hilbert空间填充曲线,相较于Peano曲线而言,Hilbert曲线没有较大的突变。为什么GeoHash不选择Hilbert空间填充曲线呢?可能是Peano曲线思路以及计算上比较简单吧,事实上,Peano曲线就是一种四叉树线性编码方式。


(5)使用注意点
1)由于GeoHash是将区域划分为一个个规则矩形,并对每个矩形进行编码,这样在查询附近POI信息时会导致以下问题,比如红色的点是我们的位置,绿色的两个点分别是附近的两个餐馆,但是在查询的时候会发现距离较远餐馆的GeoHash编码与我们一样(因为在同一个GeoHash区域块上),而较近餐馆的GeoHash编码与我们不一致。这个问题往往产生在边界处。


12.png

解决的思路很简单,我们查询时,除了使用定位点的GeoHash编码进行匹配外,还使用周围8个区域的GeoHash编码,这样可以避免这个问题。2)我们已经知道现有的GeoHash算法使用的是Peano空间填充曲线,这种曲线会产生突变,造成了编码虽然相似但距离可能相差很大的问题,因此在查询附近餐馆时候,首先筛选GeoHash编码相似的POI点,然后进行实际距离计算。

4.实践

geohash

<?php
  class Geohash{
private $coding="0123456789bcdefghjkmnpqrstuvwxyz";
private $codingMap=array();

public function Geohash()
{
    for($i=0; $i<32; $i++)
    {
        $this->codingMap[substr($this->coding,$i,1)]=str_pad(decbin($i), 5, "0", STR_PAD_LEFT);
    }

}

public function decode($hash)
{
    $binary="";
    $hl=strlen($hash);
    for($i=0; $i<$hl; $i++)
    {
        $binary.=$this->codingMap[substr($hash,$i,1)];
    }

    $bl=strlen($binary);
    $blat="";
    $blong="";
    for ($i=0; $i<$bl; $i++)
    {
        if ($i%2)
            $blat=$blat.substr($binary,$i,1);
        else
            $blong=$blong.substr($binary,$i,1);

    }

     $lat=$this->binDecode($blat,-90,90);
     $long=$this->binDecode($blong,-180,180);
    // $lat=$this->binDecode($blat,2,54);
    // $long=$this->binDecode($blong,72,136);

    $latErr=$this->calcError(strlen($blat),-90,90);
    $longErr=$this->calcError(strlen($blong),-180,180);

    $latPlaces=max(1, -round(log10($latErr))) - 1;
    $longPlaces=max(1, -round(log10($longErr))) - 1;

    $lat=round($lat, $latPlaces);
    $long=round($long, $longPlaces);

    return array($lat,$long);
}

public function encode($lat,$long)
{
    $plat=$this->precision($lat);
    $latbits=1;
    $err=45;
    while($err>$plat)
    {
        $latbits++;
        $err/=2;
    }

    $plong=$this->precision($long);
    $longbits=1;
    $err=90;
    while($err>$plong)
    {
        $longbits++;
        $err/=2;
    }

    $bits=max($latbits,$longbits);

    $longbits=$bits;
    $latbits=$bits;
    $addlong=1;
    while (($longbits+$latbits)%5 != 0)
    {
        $longbits+=$addlong;
        $latbits+=!$addlong;
        $addlong=!$addlong;
    }

    $blat=$this->binEncode($lat,-90,90, $latbits);

    $blong=$this->binEncode($long,-180,180,$longbits);

    $binary="";
    $uselong=1;
    while (strlen($blat)+strlen($blong))
    {
        if ($uselong)
        {
            $binary=$binary.substr($blong,0,1);
            $blong=substr($blong,1);
        }
        else
        {
            $binary=$binary.substr($blat,0,1);
            $blat=substr($blat,1);
        }
        $uselong=!$uselong;
    }

    $hash="";
    for ($i=0; $i<strlen($binary); $i+=5)
    {
        $n=bindec(substr($binary,$i,5));
        $hash=$hash.$this->coding[$n];
    }

    return $hash;
}

private function calcError($bits,$min,$max)
{
    $err=($max-$min)/2;
    while ($bits--)
        $err/=2;
    return $err;
}

private function precision($number)
{
    $precision=0;
    $pt=strpos($number,'.');
    if ($pt!==false)
    {
        $precision=-(strlen($number)-$pt-1);
    }

    return pow(10,$precision)/2;
}

private function binEncode($number, $min, $max, $bitcount)
{
    if ($bitcount==0)
        return "";
    $mid=($min+$max)/2;
    if ($number>$mid)
        return "1".$this->binEncode($number, $mid, $max,$bitcount-1);
    else
        return "0".$this->binEncode($number, $min, $mid,$bitcount-1);
}

private function binDecode($binary, $min, $max)
{
    $mid=($min+$max)/2;

    if (strlen($binary)==0)
        return $mid;

    $bit=substr($binary,0,1);
    $binary=substr($binary,1);

    if ($bit==1)
        return $this->binDecode($binary, $mid, $max);
    else
        return $this->binDecode($binary, $min, $mid);
}
}

测试

    <?php

require_once("db_config.php");
require_once('geohash.class.php');

$geohash=new Geohash;
//经纬度转换成Geohash
//获取附近的信息

$n_latitude  =  34.236080797698;
$n_longitude = 109.0145193757;

echo "当前位置为:经度108.7455,纬度34.3608<br/><br/>
以下网点离我最近:";

//开始
$b_time = microtime(true);

//方案A,直接利用数据库存储函数,遍历排序

//方案B geohash求出附近,然后排序

//当前 geohash值
 $n_geohash = $geohash->encode($n_latitude,$n_longitude);

//附近
$n = 3;
$like_geohash = substr($n_geohash, 0, $n);

$sql = 'select * from retailersinfotable where geohash like         
 "'.$like_geohash.'%"';



$query = mysql_query($sql);
if(mysql_num_rows($query))
{
    while($row=mysql_fetch_array($query))
    {
        $data[] = array (
                        "latitude" =>$row["Latitude"],
                        "longitude"=>$row["Longitude"],
                        "name"           
  =>$row["RetailersName"],
                );
    }
}

//算出实际距离
 foreach($data as $key=>$val)
{
    $distance =    getDistance($n_latitude,$n_longitude,$val['latitude'],$val['longitude']
 );
  $data[$key]['distance'] = $distance;

   //排序列
 $sortdistance[$key] = $distance;
 } 

  //距离排序
  array_multisort($sortdistance,SORT_ASC,$data);

 //结束
 $e_time = microtime(true);

 echo "(计算耗时:" ;
 echo $e_time - $b_time; 
 echo "秒)<br/>";

//var_dump($data);

foreach($data as $key=>$val)
{
echo "</br>";
echo $val['distance']. " 米-------".$val['name'];
} 

/**
*  @desc 根据两点间的经纬度计算距离
*  @param float $latitude 纬度值
*  @param float $longitude 经度值
*/
function getDistance($latitude1, $longitude1, $latitude2, $longitude2) 
{
$earth_radius = 6371000;   //approximate radius of earth in meters

$dLat = deg2rad($latitude2 - $latitude1);
$dLon = deg2rad($longitude2 - $longitude1);
 /*
   Using the
   Haversine formula

   http://en.wikipedia.org/wiki/Haversine_formula
   http://www.codecodex.com/wiki/Calculate_Distance_Between_Two_Points_on_a_Globe
   验证:百度地图  http://developer.baidu.com/map/jsdemo.htm#a6_1
   calculate the distance
 */ 
$a = sin($dLat/2) * sin($dLat/2) + cos(deg2rad($latitude1)) * cos(deg2rad($latitude2)) * sin($dLon/2) * sin($dLon/2);
$c = 2 * asin(sqrt($a));
$d = $earth_radius * $c;

return round($d);   //四舍五入
}

 ?>

5.备注

  在纬度相等的情况下:
    经度每隔0.00001度,距离相差约1米;
    每隔0.0001度,距离相差约10米;
    每隔0.001度,距离相差约100米;
    每隔0.01度,距离相差约1000米;
    每隔0.1度,距离相差约10000米。

  在经度相等的情况下:
    纬度每隔0.00001度,距离相差约1.1米;
    每隔0.0001度,距离相差约11米;
    每隔0.001度,距离相差约111米;
    每隔0.01度,距离相差约1113米;
    每隔0.1度,距离相差约11132米。

6.参考

查找附近点--Geohash方案讨论
Haversine formula球面距离公式
球面距离公式代码实现
球面距离公式验证

?>

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 199,393评论 5 467
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 83,790评论 2 376
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 146,391评论 0 330
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 53,703评论 1 270
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 62,613评论 5 359
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,003评论 1 275
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 37,507评论 3 390
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,158评论 0 254
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 40,300评论 1 294
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,256评论 2 317
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,274评论 1 328
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 32,984评论 3 316
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 38,569评论 3 303
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,662评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 30,899评论 1 255
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 42,268评论 2 345
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 41,840评论 2 339

推荐阅读更多精彩内容