原文地址: https://blog.fanscore.cn/a/51/
背景
最近接到一個需求，開發中使用了redis georadius命令取附近給定距離內的點。完工後對服務進行壓測後發現georadius的效能比預期要差，因此我分析了georadius的原始碼，並對原始的實現方案進行了最佳化，總結成了本文。
我們生產環境使用的redis版本為4.0.13，因此本文redis原始碼皆為4.0.13版本的原始碼
redis geo原理
往redis中新增座標的命令是GEOADD key longitude latitude member [longitude latitude member ...]，實際上redis會將經緯度轉成一個52bit的整數作為zset的score，然後新增到zset中，所以實際上redis geo底層就是個zset，你甚至可以直接使用zset的命令來操作一個geo型別的key。
那麼經緯度是如何轉成52bit整數的呢？業內廣泛使用的方法是首先對經緯度分別按照二分法編碼，然後將各自的編碼交叉組合成最後的編碼。我們以116.505021, 39.950898這個座標為例看下如何編碼：
- 第一次二分操作，把經度分為兩個區間：
[-180,0)和[0,180]，116.505021落在右區間，因此用1表示第一次編碼後的值 - 第二次二分操作，把
[0,180]分為兩個區間[0,90)和[90,180]，116.505021落在右區間，因此用1表示第二次編碼後的值 - 第三次二分操作，把
[90,180]分為兩個區間[90,135)和[135,180]，116.505021落在左區間，因此用0表示第二次編碼後的值 - 按照這種方法依次處理，做完5次後，得到經度值的5位編碼值：
11010
| 分割槽次數 | 左區間 | 右區間 | 經度116.505021在區間 | 編碼值 |
|---|---|---|---|---|
| 1 | [-180, 0) | [0, 180] | [0, 180] | 1 |
| 2 | [0, 90) | [90, 180] | [90, 180] | 1 |
| 3 | [90, 135) | [135, 180] | [90, 135]) | 0 |
| 4 | [90, 112.5) | [112.5, 135] | [112.5, 135] | 1 |
| 5 | [112.5, 123.75) | [123.75, 180] | [112.5, 123.75] | 0 |
- 按照同樣的方法對緯度值進行編碼，得到緯度值的5位編碼值：
10111
| 分割槽次數 | 左區間 | 右區間 | 緯度39.950898在區間 | 編碼值 |
|---|---|---|---|---|
| 1 | [-90, 0) | [0, 90] | [0, 90] | 1 |
| 2 | [0, 45) | [45, 90] | [0, 45] | 0 |
| 3 | [0, 22.5) | [22.5, 45] | [22.5, 45]) | 1 |
| 4 | [22.5, 33.75) | [33.75, 45] | [33.75, 45] | 1 |
| 5 | [33.75, 39.375) | [39.375, 45] | [39.375, 45] | 1 |
然後將經度編碼11010和緯度編碼值10111交叉得到最終geohash值1110011101
通常會使用base32將編碼值轉成字串表示的hash值，與本文無關這裡不多做介紹
根據如上的演算法通常可以直觀的寫出如下的程式碼：
// 該程式碼來源於https://github.com/HDT3213/godis/blob/master/lib/geohash/geohash.go
func encode0(latitude, longitude float64, bitSize uint) ([]byte, [2][2]float64) {
box := [2][2]float64{
{-180, 180}, // lng
{-90, 90}, // lat
}
pos := [2]float64{longitude, latitude}
hash := &bytes.Buffer{}
bit := 0
var precision uint = 0
code := uint8(0)
for precision < bitSize {
for direction, val := range pos {
mid := (box[direction][0] + box[direction][1]) / 2
if val < mid {
box[direction][1] = mid
} else {
box[direction][0] = mid
code |= bits[bit]
}
bit++
if bit == 8 {
hash.WriteByte(code)
bit = 0
code = 0
}
precision++
if precision == bitSize {
break
}
}
}
if code > 0 {
hash.WriteByte(code)
}
return hash.Bytes(), box
}
可以看到基本就是上述演算法的實際描述，但是redis原始碼中卻是另外一種演算法：
int geohashEncode(const GeoHashRange *long_range, const GeoHashRange *lat_range,
double longitude, double latitude, uint8_t step,
GeoHashBits *hash) {
// 引數檢查此處程式碼省略
...
double lat_offset =
(latitude - lat_range->min) / (lat_range->max - lat_range->min);
double long_offset =
(longitude - long_range->min) / (long_range->max - long_range->min);
lat_offset *= (1 << step);
long_offset *= (1 << step);
// lat_offset與long_offset交叉
hash->bits = interleave64(lat_offset, long_offset);
return 1;
}
那麼該如何理解redis的這種演算法呢？我們假設經度用3位來編碼
可以看到編碼值從左到右實際就是從000到111依次加1遞進的，給定的經度值在這條線的位置(偏移量)就是其編碼值。假設給定經度值為50，那麼它在這條線的偏移量就是(50 - -180) / (180 - -180) * 8 = 5即101
georadius原理
georadius命令格式為GEORADIUS key longitude latitude radius m|km|ft|mi [WITHCOORD] [WITHDIST] [WITHHASH] [COUNT count] [ASC|DESC] [STORE key] [STOREDIST key]，以給定的經緯度為中心， 返回鍵包含的位置元素當中， 與中心的距離不超過給定最大距離的所有位置元素。
首先需要明確一點的是並非兩個座標點編碼相近其距離越近，以上圖為例，雖然A所在區塊的編碼與C所在區塊編碼較之B更相近，但實際B點距離A點更近。為了避免這種問題redis中會先計算出給定點東南西北以及東北、東南、西北、西南八個區塊以及自己身所在的區塊即九宮格區域內所有座標點，然後計算與當前點的距離，再進一步篩選出符合距離條件的點。
假設要查附近100km的點，那麼要保證矩形的邊長要大於100km，才能保證能獲取到所有符合條件的點，地球半徑約6372.797km，第一次分割後可以得到四個東西長6372.797*π，南北長3186.319*π，繼續切割：
| 分割次數 | 東西長(km) | 南北長(km) |
|---|---|---|
| 1 | 6372.797*π | 3186.319*π |
| 2 | 3186.319*π | 1593.160*π |
| 3 | 1593.160*π | 796.58*π |
| 4 | 796.58*π | 398.29*π |
| 5 | 398.29*π | 199.145*π |
| 6 | 199.145*π | 99.573*π |
| 7 | 99.573*π | 49.787*π |
分割到第七次時南北長49.787*π，如果再切分長度為24.894*π，長度小於100km，因此停止分割，所以如果要查附近100km的點，我們需要的精度為7
redis中根據給定的距離估算出需要的精度的程式碼如下
const double MERCATOR_MAX = 20037726.37;
uint8_t geohashEstimateStepsByRadius(double range_meters, double lat) {
if (range_meters == 0) return 26;
int step = 1;
while (range_meters < MERCATOR_MAX) {
range_meters *= 2;
step++;
}
step -= 2;
// 高緯度地區地球半徑小因此適當降低精度
if (lat > 66 || lat < -66) {
step--;
if (lat > 80 || lat < -80) step--;
}
if (step < 1) step = 1;
if (step > 26) step = 26;
return step;
}
呼叫encode0函式就能計算出給定點在step = geohashEstimateStepsByRadius()精度級別所在矩形區域的geohash值。接下來計算該矩形區域附近的八個區域。
...
// 呼叫encode0函式計算geohash
geohashEncode(&long_range,&lat_range,longitude,latitude,steps,&hash);
// 計算出附近八個區域
geohashNeighbors(&hash,&neighbors);
...
一個區域的東側區域只要將經度的編碼值+1即可，反之西側區域只要將經度編碼值-1即可，北側區域只要將緯度的編碼值+1即可，南側區域只要將緯度的編碼值-1即可。對應redis原始碼如下:
void geohashNeighbors(const GeoHashBits *hash, GeoHashNeighbors *neighbors) {
neighbors->east = *hash;
neighbors->west = *hash;
neighbors->north = *hash;
neighbors->south = *hash;
neighbors->south_east = *hash;
neighbors->south_west = *hash;
neighbors->north_east = *hash;
neighbors->north_west = *hash;
// 緯度加1就是東側區域
geohash_move_x(&neighbors->east, 1);
geohash_move_y(&neighbors->east, 0);
// 緯度減1就是西側區域
geohash_move_x(&neighbors->west, -1);
geohash_move_y(&neighbors->west, 0);
// 精度減1就是南側區域
geohash_move_x(&neighbors->south, 0);
geohash_move_y(&neighbors->south, -1);
geohash_move_x(&neighbors->north, 0);
geohash_move_y(&neighbors->north, 1);
geohash_move_x(&neighbors->north_west, -1);
geohash_move_y(&neighbors->north_west, 1);
geohash_move_x(&neighbors->north_east, 1);
geohash_move_y(&neighbors->north_east, 1);
geohash_move_x(&neighbors->south_east, 1);
geohash_move_y(&neighbors->south_east, -1);
geohash_move_x(&neighbors->south_west, -1);
geohash_move_y(&neighbors->south_west, -1);
}
如上圖所示，當給定點在中心區域的東北側時，西北、西、西南、南、東南五個方向的區域中的所有點距離給定點肯定超過了給定距離，所以可以過濾掉，redis程式碼如下所示：
if (steps >= 2) {
if (area.latitude.min < min_lat) {
GZERO(neighbors.south); // 南側區域置零，過濾南側區域
GZERO(neighbors.south_west);
GZERO(neighbors.south_east);
}
if (area.latitude.max > max_lat) {
GZERO(neighbors.north);
GZERO(neighbors.north_east);
GZERO(neighbors.north_west);
}
if (area.longitude.min < min_lon) {
GZERO(neighbors.west);
GZERO(neighbors.south_west);
GZERO(neighbors.north_west);
}
if (area.longitude.max > max_lon) {
GZERO(neighbors.east);
GZERO(neighbors.south_east);
GZERO(neighbors.north_east);
}
}
計算出區塊後下一步就需要將九宮格區域中的所有座標點拿出來，依次計算與給定點的距離，然後過濾出符合給定距離的點
// 遍歷九宮格內所有點，依次計算與給定點的距離，然後過濾出符合給定距離的點新增到ga中
int membersOfAllNeighbors(robj *zobj, GeoHashRadius n, double lon, double lat, double radius, geoArray *ga) {
GeoHashBits neighbors[9];
unsigned int i, count = 0, last_processed = 0;
int debugmsg = 1;
neighbors[0] = n.hash;
neighbors[1] = n.neighbors.north;
neighbors[2] = n.neighbors.south;
neighbors[3] = n.neighbors.east;
neighbors[4] = n.neighbors.west;
neighbors[5] = n.neighbors.north_east;
neighbors[6] = n.neighbors.north_west;
neighbors[7] = n.neighbors.south_east;
neighbors[8] = n.neighbors.south_west;
// 遍歷九宮格
for (i = 0; i < sizeof(neighbors) / sizeof(*neighbors); i++) {
...
// 當給定距離過大時，區塊可能會重複
if (last_processed &&
neighbors[i].bits == neighbors[last_processed].bits &&
neighbors[i].step == neighbors[last_processed].step)
{
continue;
}
// 取出宮格內所有點，依次計算距離，符合條件後新增到ga中
count += membersOfGeoHashBox(zobj, neighbors[i], ga, lon, lat, radius);
last_processed = i;
}
return count;
}
int membersOfGeoHashBox(robj *zobj, GeoHashBits hash, geoArray *ga, double lon, double lat, double radius) {
GeoHashFix52Bits min, max;
// 根據區塊的geohash值計算出對應的zset的score的上下限[min,max]
scoresOfGeoHashBox(hash,&min,&max);
// 取出底層的zset中的[min,max]範圍內的元素，依次計算距離，符合條件後新增到ga中
return geoGetPointsInRange(zobj, min, max, lon, lat, radius, ga);
}
georadius最佳化
從上一節中可以看到，給定距離範圍越大，則九宮格區域越大，九宮格區域內的點就越多，而每個點都需要計算與中間點的距離，距離計算又涉及到大量的三角函式計算，所以這部分計算是十分消耗CPU的。又因為redis工作執行緒是單執行緒的，因此無法充分利用多核，無法透過增加redis server的CPU核數來提升效能，只能新增從庫。
距離計算演算法及最佳化可以看下美團的這篇文章: https://tech.meituan.com/2014/09/05/lucene-distance.html
對於這個問題，我們可以將九宮格以及距離計算部分提升到我們的應用程式即redis客戶端來進行，步驟如下：
- 在客戶端計算出九宮格區域，然後轉為zset score的範圍
- 使用zrangebyscore命令從redis取出score範圍內的所有點
- 遍歷所有點依次計算與給定點的距離，篩選出符合距離條件的點
陌陌好像也是使用了這種方案：https://mp.weixin.qq.com/s/DL2P49y4R1AE2MIdkxkZtQ
由於我們使用golang進行開發，因此我將redis中的georadius部分程式碼轉為了golang程式碼，並整理成一個庫開源在了github：https://github.com/Orlion/go-georadius
原本的寫法是:
client.GeoRadius(key, longitude, latitude, &redis.GeoRadiusQuery{
Radius: 1000,
Unit: "m", // 距離單位
Count: 1, // 返回1條
WithCoord: true, // 將位置元素的經緯度一併返回
WithDist: true, // 一併返回距離
})
改造後：
ga := make([]redis.Z, 0)
ranges := geo.NeighborRanges(longitude, latitude, 1000)
for _, v := range ranges {
zs, _ := client.ZRangeByScoreWithScores(key, redis.ZRangeBy{
Min: strconv.Itoa(int(v[0])),
Max: strconv.Itoa(int(v[1])),
}).Result()
for _, z := range zs {
dist := geox.GetDistanceByScore(longitude, latitude, uint64(z.Score))
if dist < 1000 {
ga = append(ga, z)
}
}
}
壓測結果對比
43w座標點，取附近50km（九宮格內有14774點，符合條件的點約6000個）
50km最佳化前
Concurrency Level: 5
Time taken for tests: 89.770 seconds
Complete requests: 5000
Failed requests: 0
Write errors: 0
Total transferred: 720000 bytes
HTML transferred: 0 bytes
Requests per second: 55.70 [#/sec] (mean)
Time per request: 89.770 [ms] (mean)
Time per request: 17.954 [ms] (mean, across all concurrent requests)
Transfer rate: 7.83 [Kbytes/sec] received
Connection Times (ms)
min mean[+/-sd] median max
Connect: 0 0 0.0 0 0
Processing: 23 90 10.7 90 159
Waiting: 23 89 10.7 89 159
Total: 23 90 10.7 90 159
Percentage of the requests served within a certain time (ms)
50% 90
66% 93
75% 96
80% 97
90% 102
95% 107
98% 111
99% 116
100% 159 (longest request)
50km最佳化後
Concurrency Level: 5
Time taken for tests: 75.447 seconds
Complete requests: 5000
Failed requests: 0
Write errors: 0
Total transferred: 720000 bytes
HTML transferred: 0 bytes
Requests per second: 66.27 [#/sec] (mean)
Time per request: 75.447 [ms] (mean)
Time per request: 15.089 [ms] (mean, across all concurrent requests)
Transfer rate: 9.32 [Kbytes/sec] received
Connection Times (ms)
min mean[+/-sd] median max
Connect: 0 0 0.0 0 0
Processing: 21 75 14.2 75 159
Waiting: 21 75 14.1 75 159
Total: 21 75 14.2 75 159
Percentage of the requests served within a certain time (ms)
50% 75
66% 80
75% 84
80% 86
90% 92
95% 98
98% 104
99% 111
100% 159 (longest request)
可以看到效能並沒有巨大的提升，我們減小距離範圍到5km（符合條件的點有130個）再看下壓測結果
5km最佳化前
Concurrency Level: 5
Time taken for tests: 14.006 seconds
Complete requests: 5000
Failed requests: 0
Write errors: 0
Total transferred: 720000 bytes
HTML transferred: 0 bytes
Requests per second: 356.99 [#/sec] (mean)
Time per request: 14.006 [ms] (mean)
Time per request: 2.801 [ms] (mean, across all concurrent requests)
Transfer rate: 50.20 [Kbytes/sec] received
Connection Times (ms)
min mean[+/-sd] median max
Connect: 0 0 0.0 0 0
Processing: 2 14 5.5 12 33
Waiting: 2 14 5.5 12 33
Total: 2 14 5.5 12 34
Percentage of the requests served within a certain time (ms)
50% 12
66% 16
75% 19
80% 20
90% 22
95% 23
98% 27
99% 28
100% 34 (longest request)
5km最佳化後
Concurrency Level: 5
Time taken for tests: 16.661 seconds
Complete requests: 5000
Failed requests: 0
Write errors: 0
Total transferred: 720000 bytes
HTML transferred: 0 bytes
Requests per second: 300.11 [#/sec] (mean)
Time per request: 16.661 [ms] (mean)
Time per request: 3.332 [ms] (mean, across all concurrent requests)
Transfer rate: 42.20 [Kbytes/sec] received
Connection Times (ms)
min mean[+/-sd] median max
Connect: 0 0 0.0 0 0
Processing: 3 17 5.8 16 66
Waiting: 3 16 5.8 16 66
Total: 3 17 5.8 16 66
Percentage of the requests served within a certain time (ms)
50% 16
66% 20
75% 21
80% 22
90% 24
95% 26
98% 28
99% 30
100% 66 (longest request)
可以看到當最佳化後效能更差了
猜測造成這個結果的原因應該是附近5km九宮格內的點比較少，所以最佳化後實際沒減少多少距離計算，但多了n(n<=9)倍的請求數，多了額外的命令解析與響應內容的消耗，因此這種最佳化方案僅僅適用於附近點特別多的情況。