一种Python爬取百度地图瓦片的方式

2019-09-29 11:39:07 浏览数 (2)

瓦片是互联网地图webGIS组织地图数据的一种方式,最近的一个项目需要获取一定区域内的百度地图瓦片;ArcGIS的一个插件ArcBruTile支持很多地图源(如OSM、Bing)的瓦片获取,但是没有百度地图的,Github上的一些项目(如pyMap)也不支持百度瓦片的下载,

因此打算自己去获取。开始是想用百度API的静态图服务,但是我只有一个区域的边界矩形(MBR),而它接受的参数是中心点坐标和zoom的尺度,其请求URL的核心参数为center=116.403,39.914&width=1024&height=1024&zoom=11,要分块下载很麻烦,想了很多办法去算坐标和进行坐标距离换算,走了一些弯路,后来发现百度地图有一个隐晦的支持参数为x,y,z的调用,而x,y的变化是比较有规律的,因此只需要有边界的x和y再叠加就行。我是下载区域内的小图片再拼接为大图使用的,最后有上万张小图片,合并为一个200多MB的大图。

下载下来的瓦片

具体实现过程如下:

1,获取图片

百度坐标拾取系统(可以用关键词搜索得到网址)网页,先按F12调出控制台,用坐标反查定位到左下坐标,并且调好层级,我要用的是17级的地图,然后找到一张图片手动确定x1和y1;(虽然有逻辑可以根据坐标和层级算x和y,但百度地图版的实现起来还是费些时间的,我没找到现成的代码,如果读者有发现或写过实现这个功能的代码欢迎在评论里告知,非常感谢),再定位到边界的右上角,同样确定好x2,y2,然后用requests库写获取图片的代码,比较建议先定y,改x,我用先循环y的方式保存的图片合并起来更复杂些,循环下去;可以得到整个区域的图片;

百度坐标拾取系统采用瓦片底图

(先验知识:百度地图的瓦片是从左下角算的,而不是Google Map的左上角开始;)

代码语言:javascript复制
import requests

def getTileByXYZ(): #根据x,y,z参数获取瓦片
   z=17
   xidx=[22568,22676]
   yidx=[6898,7008]   
    for y in range(yidx[0],yidx[1] 1):        
        for x in range(xidx[0],xidx[1] 1):
           url="http://online3.map.bdimg.com/tile/?qt=tile&x={x}&y={y}&z=17&styles=pl"  
                "&scaler=1&udt=20180601".format(x=x,y=y)
           savePngByXYZ(url, x, y, z)
       print(y)    
                    
def savePngByXYZ(url,x,y,z=17): #保存图片
   r = requests.get(url)
   sname="./cdZoomImg/cd_{x}_{y}.png".format(x=x,y=y) #这里建议保存编码是y_x 这样下面合并图片也要适当改代码
   with open(sname, 'ab') as pngf:       
        for chunk in r.iter_content(chunk_size=1024):           
            if chunk:
               pngf.write(chunk)
               pngf.flush()

2,合并瓦片

一共爬了1万多张瓦片,花了40多分钟(这个要看计算机性能),每张瓦片是256*256像素的方形区域;接下来就是合并的过程了,本来打算一次性合并完毕,但是弄得电脑内存不足了,于是先合并x相同的图片到另一个文件夹,形成一个个长条形的图片,每张尺寸变成了256*28672像素,再跑一遍把这些图片合并到一起,就形成了一张27648*28672的地图图片,可以用来作为ArcGIS一些空间分析的底图。

逐步合并瓦片

代码语言:javascript复制
import os
import glob
from PIL import Image
def complieImg():
   #命名规则:cd_x_y.png 左下坐标系
   #同一个x 同1列,y增加,图片在上面
   #假设输入排好序了
   p = "./chengduImg"
   plst = glob.glob(os.path.join(p, '*.png'))   xmin=((plst[0].split("\")[1]).split(".")[0]).split('_')[1]
   alst=[] #3维
   qlst=[]    
    for f in plst:
       w=((f.split("\")[1]).split(".")[0]).split('_') #['cd', '22568', '6898']
       w[0] = f        
        if w[1]==xmin:
           qlst.append(w.copy())        
        else:
           alst.append(qlst.copy())
           xmin=w[1]
           qlst=[]
   m2 = [256*len(alst[0]), 256 * len(alst)]  #im2=Image.new('RGBA', (m2[0], m2[1]))
   print(m2)
   psave = "./complexLevel"
   iw=0
   for k in alst:#k里面装的是x相同的值,y应该递增
       plen=len(k)
       msize = [256, 256 * (plen 1)]
       print(msize)
       toImage = Image.new('RGBA', (msize[0], msize[1]))       
        for i in range(plen):
           fromImage = Image.open(k[plen - i - 1][0])
           toImage.paste(fromImage, (0 * msize[0], i * msize[0]))       sname="/m_{x}.png".format(x=k[0][1])
       iw =1
       
       toImage.save(psave sname)    def complieImgInY():
   #合并长条形图片,x变化,y不变 长图是complieImg()里生成的
   p = "./complexLevel" #chengduImg
   plst = glob.glob(os.path.join(p, '*.png'))   xmin=((plst[0].split("\")[1]).split(".")[0]).split('_')[1]
   ima21=Image.open(plst[0])
   w=np.array(ima21).shape
   print(w)
   psave = "D:/wexcel"
   plen=len(plst)
   msize = [w[1]*plen/2, w[0]/2]
   print(msize)
   toImage = Image.new('RGBA', (int(msize[0]), int(msize[1])))    
    for i in range(plen):
       fromImage = Image.open(plst[i])
       fromImage=fromImage.resize((int(256/2),int(msize[1])), Image.ANTIALIAS)
       toImage.paste(fromImage, (int(i * 256/2), 0))   sname="/chengduMap.png"   toImage.save(psave sname) #保存图片

最后得到一张200多MB的整合图片。最后代码更新于https://github.com/QLWeilcf/LcfGeoProject/tree/master/WebGISLyn。

OutputImg

0 人点赞