GEOHash

在GAE之中使用基于地理位置的查询

    当我决定将我的数据上传到GAE上的时候，我就预计到将来的数据调用将会是非常费劲的，果然，我现在仅仅是小小的使用了一下自己的数据就已经把自己折磨的够呛，因为GAE不能支持联合查询，在数据的排序和条件上也有诸多限制，因此，操作起来非常费劲，甚至有时候必须牺牲一些性能或效果。

    本文对那些问题不做详细的表述，仅仅谈谈我进行地理数据查询的实现。

    将地理的数据存储到Google的DataStore之中，可以直接用两个float的字段，代表经纬度，也可以直接使用db.GeoPtProperty,现在看来，直接使用float字段要直观容易的多，可是我当时上传的时候考虑到以后Google可能会在db.GeoPtProperty的基础是那个提供一些方便的功能，因此，就将数据上传为db.GeoPtProperty,可是上传之后遇到问题，原来db.GeoPtProperty不支持基于范围的查询，这就困难了，早知道我用两个float，建立好索引，应该可以实现周边查询等等功能的(后来补充：这样也不行，因为GQL在查询的时候仅仅支持一种"非等于"的条件，也就是说要实现lat>3 and lat<54 and lon>73 and lon<136这种查询也是不行的)，很怀恋以前直接使用SQL Server的时候，一个SQL语句，按距离排序就搞定了。

    那么GAE之中的经纬度是按照什么顺序存放的呢？我从网上查了一下，GAE的经纬度是严格按照(lat,lon)的顺序存放的，也就是说基本上是按照纬度排，纬度相同的再按照经度拍，这样的话，查找某一个纬度范围的点是很容易的，不过这个查询在实际应用之中用处不大。

    我从网上查了一下，基本上大家都是考虑用geohash来实现周边查询,这个东西我没有精力去做具体的研究（后来补充：准备研究，专文讲叙），只能大致的猜想：geohash将一个经纬度序列化成为一个字符串，最后能够大致做到“两个越相邻的经纬度，得到字符串开头部分相同的位数就越多”，也就是说将一个平面的2维坐标编码成为了一个线性的序列，在二维坐标之中离的近的点，在这个线性序列之中距离也近。这种算法让我觉得很奇怪，不知道怎么实现的，不过这个问题等有时间再去研究，我现在，且在网站上直接使用这个功能。

    1.首先要使用geohash库，从网上去下一个geohash.py文件，放到自己的应用程序根目录下面

    2.要将线上的数据包含经纬度的表增加一个geohash的字段，并且将geohash编码字符串存进去，这一点我是使用remote_data来实现的，这种方式来处理数据方便快捷的多，比以前一直刷新页面好多了，我的增加这个字段的代码如下：

 

 1# -*- coding: utf-8 -*- #   
 2import code
 3import getpass
 4import sys
 5#下面要改成自己的gae安装路径
 6sys.path.append("D:\Program Files\Google\google_appengine\lib\yaml\lib")
 7sys.path.append("D:\Program Files\Google\google_appengine")
 8#下面是我的应用程序的路径，刚才说过geohash是放在那里的
 9sys.path.append("D:\work\myapp")
10
11import geohash
12from google.appengine.ext.remote_api import remote_api_stub
13from google.appengine.ext import db
14#这是数据表的定义，注意，最下面一行被我加上了geohash字段
15class Train_stations(db.Model):
16  name = db.StringProperty()
17  latlng = db.GeoPtProperty()
18  superior = db.StringProperty()
19  address = db.StringProperty()
20  postcode = db.StringProperty()
21  regionCode = db.StringProperty()
22  level = db.StringProperty()
23  telephone = db.StringProperty()
24  oldName = db.StringProperty()
25  lineCount = db.IntegerProperty()
26  geohash = db.StringProperty()
27
28def auth_func():
29  return raw_input('Username:'), getpass.getpass('Password:')
30#这是应用程序ID
31app_id='myapp'
32host = '%s.appspot.com' % app_id
33
34remote_api_stub.ConfigureRemoteDatastore("myapp", '/remote_api',auth_func)
35
36key=''
37while True:
38    sql="select * from Train_stations"
39    if key:
40        sql+=" where __key__>KEY('%s') "%key
41    sql+=" order by __key__"
42    stations=db.GqlQuery(sql).fetch(100)
43    if len(stations)<=0:
44        break
45    for station in stations:
46        if station.latlng and station.latlng.lat and station.latlng.lon:
47            station.geohash=str(geohash.Geohash((float(station.latlng.lon), float(station.latlng.lat))))
48        else:
49            station.geohash=''
50        print "%s:%s" % (station.name,station.geohash)
51    db.put(stations)
52    key=str(stations[len(stations)-1].key())
53    print key
54print "OK"

 

    上面这段代码在本地运行即可，不需要上传到服务器上去（不过需要注意服务器必须打开了remote_api,并且路径对应），运行的时候，会看到将每行数据经纬度进行geohash之后的数据，等到运行显示"OK"，则说明所有的经纬度都hash完毕了，到GAE的后台可以看到数据，例如一个经纬度是"25.6138,109.484"，geohash之后是"he6nyfgqsbce4"。

    3.数据处理完毕之后，就可以进行查询了，先理清思路：要查询一个经纬度附近的点，应该先查询大于这个geohash的点之中最小的hash，再查小于这个geohash的点之中最大的hash，因此必须查两遍：

 

    tvs[ ' nearStations ' ] = db.GqlQuery( " select * from Train_stations where geohash < :1 order by geohash desc " ,(hash)).fetch( 10 )
    tvs[ ' nearStations ' ].extend(db.GqlQuery( " select * from Train_stations where geohash > :1 order by geohash " ,(hash)).fetch( 10 ))

 

    4.通过以上的方式，就可以查出临近的点，这只是geohash的一种使用而已，至少，你还可以通过查询和当前的geohash前几位字母相同的hash值，从而的得到这个经纬度周边一定范围内的所有点，这一点因为我的网站上没有用到，因此，没有相应的代码可以提供，不过从原理上是可以实现的。

    以上就是我对geohash的简单应用，现在，这些应用已经开始运行在http://www.dituren.cn上面，目前来看，效果不错，马上我会仔细研究一下geohash的实现，对于我这个非专业的人来讲，geohash真是太神奇了。

关于geohash的简单探讨

    在上次写了文章《在GAE之中使用基于地理位置的查询》，之后，我一直在奇怪geohash实现周边查询的原理是什么，毕竟地理数据可是二维的坐标，而geohash的结果只是一个简单的字符串，要说通过简单字符串的比较就能找到周边的点，无论如何我也不能相信，因此我研究了一下geohash的算法，果然发现我以前的做法其实是做不到精确的周边查询的，不过我也不得不承认geohash确实是一个很好的索引模式，下面简单的介绍一下。
    具体的geohash算法，应该去参考Geohash - Wiki，我仅仅做一下简单的原理介绍。
    geohash将一个二维的坐标转化为一个简单的可排序可比较的字符串的方法是这样的：先把经纬度范围(纬度范围-90到90，经度范围-180到180)当作一个纯平面的矩形，这样就得到一个二维的平面坐标系。怎么将这个二维坐标转化为一个一维的坐标呢，请看下面的示意图：
    如图所示，平面二维坐标(经纬度)和1维的坐标(geohash之后的字符串)转化方式如下(我绘图的本领不高，将就看看吧)：

 

    也就是说geohash先将平面坐标平分为4块，按照以下的顺序进行编号，之后再对每一个子块以同样的方式进行划分，得到另一个编号，一直按照这个规则划分下去，最终编号越来越多，而方格越来越小，直到能表达该坐标需要的精度为止。而在划分之中的得到的每个编号对应一维坐标之中的一段，这些编号最后可以编码成一个字符串，具体的编码过程请看前面提到的wiki的网址，我这里只是表述二维坐标向一维坐标的转化过程。
    知道了这个hash的原理之后，我们来进行下面的讨论：
    1.两个离的越近，geohash的结果相同的位数越多，对么？这一点是有些用户对geohash的误解，虽然geo确实尽可能的将位置相近的点hash到了一起，可是这并不是严格意义上的(实际上也并不可能，因为毕竟多一维坐标)，例如在上图大方格4的左下部分的点和大方格1的右下部分的点离的很近，可是它们的geohash值一定是相差的相当远，因为头一次的分块就相差太大了，很多时候我们对geohash的值进行简单的排序比较，结果貌似真的能够找出相近的点，并且似乎还是按照距离的远近排列的，可是实际上会有一些点被漏掉了。
    2.既然不能做到将相近的点hash值也相近，那么geohash的意义何在呢？我觉得geohash还是相当有用的一个算法，毕竟这个算法通过无穷的细分，能确保将每一个小块的geohash值确保在一定的范围之内，这样就为灵活的周边查找和范围查找提供了可能。
    我看到过一些其他实现周边查找的方法(这里不讨论直接通过组合sql语句进行查询的模式，因为这种模式实际上是很耗性能的)，基本上都是建立一个单一的块索引，例如将经纬度以0.1为单位将地图划分为若干个单元格，建立每个坐标所在单元格编号的索引，这样做虽然直观，可能只能支持一种特定的查询方式，如果需要更改为实现以0.01为单位的查询，就需要重新建立单元格索引，这显然会带来一些问题。
    如果使用geohash则不用担心这个问题，hash的结果可以支持多个级别范围的查询，每个级别之间查询的单元范围大小是2倍的关系，这实际已经可以实现灵活的查询了，例如要实现指定范围的查找，我们可以指定一个深度(深度越大，每个单元块就越小)，然后将此地理范围所覆盖的所有单元块的hash范围全部找出来，然后逐个在这些范围内查找。最后，假如有必要的话，再对结果进行筛选，去除不在此范围内的点。
    下面我贴出我通过geohash实现一个指定范围查找的例子，我的网站也是在这个例子的基础上实现了周边查找，如果不想研究具体的实现，可以考虑直接使用此代码，如果需要详细的明白实现的过程建议先去通过前面的wiki链接了解geohash的原理。

 1class GetNearStations(BaseService):
 2    def getIndexByTab(self,(t,n),bl):
 3        #根据单元格的二维序号计算1维序号
 4        num=0
 5        #先考虑经度多分一层的情况
 6        num+=n%2
 7        n=(n-num)/2
 8        #再逐次计算出num
 9        for i in range(0,bl):
10            num+=(int(t%2))<<(2*i+1)
11            t=(t-(t%2))/2
12            num+=(int(n%2))<<(2*i+2)
13            n=(n-(n%2))/2
14        return num
15    def getHashByIndex(self,index,bl):
16        #根据单元格的一维序号计算单元格的起始哈希值
17        if(index>=pow(2,bl*2+1)):
18            return "zzzzzzzzzzzzzz";
19        BASE_32 = "0123456789bcdefghjkmnpqrstuvwxyz"
20        index=index<<(64-bl*2)
21        return "".join([BASE_32[(long(index) >> ((13-i-1)*5)) & 31] for i in range(0,13)])
22    def searchBounds(self,bounds,bl):
23        #根据bounds(minLat,minLon,maxLat,maxLon)查找范围内的点
24        #此函数应确保在此范围内的点都被返回，但不确保返回的一定严格在此范围内
25        spans=[]
26        #单元格大小
27        precision=180*pow(2,-bl)
28        for t in range(int(floor((bounds[0]+90)/precision)),int(ceil((bounds[2]+90)/precision))):
29            if t<0 or t>=pow(2,bl): continue
30            for n in range(int(floor((bounds[1]+180)/precision)),int(ceil((bounds[3]+180)/precision))):
31                if n<0 or n>=pow(2,bl+1): continue
32                #调试时显示方格范围
33                #print "%.5f,%.5f,%.5f,%.5f" % (t*precision-90,n*precision-180,(t+1)*precision-90,(n+1)*precision-180)
34                index=self.getIndexByTab((t,n),bl)
35                spans.append([index,index+1])
36        #合并相邻的内容
37        spans.sort(lambda x,y: cmp(x[0],y[0]))
38        for i in range(len(spans)-1,0,-1):
39            if spans[i-1][1]==spans[i][0]:
40                spans[i-1][1]=spans[i][1]
41                spans.pop(i)
42        return [(self.getHashByIndex(span[0],bl),self.getHashByIndex(span[1],bl)) for span in spans]
43                
44    def get(self):
45        #此方法不完整，仅仅显示如何实现周边查找
46        #参数7表示深度，越大搜索单元范围越小，搜索的结果越精确，需要搜索的次数越多
47        spans=self.searchBounds((lat-0.5*span,lon-0.5*span,lat+0.5*span,lon+0.5*span),7)
48        stations=[]
49        for span in spans:
50            stations.extend(db.GqlQuery("select * from Train_stations where geohash >= '%s' and geohash<= '%s' order by geohash desc"%span).fetch(100))
51        self.outputJson(stations)

 

        实在抱歉，因为比较忙，以上文章和代码都写的比较仓促，希望没有浪费大家的时间。