下载代码 仅供参考.瓦片服务是我们自己搭的服务器，参考geoserver。
import io
import math
import multiprocessing
import time
import urllib.request as ur
from threading import Thread
import traceback
import PIL.Image as pil
from pathlib import Path
from numpy.lib.type_check import imag
import sys
# from PIL import Image
# import Image
import cv2
import numpy as np
from osgeo import gdal, osr
from tools import *
# from tqdm import tqdm, trange
# -----------------------------------------------------------
# ---------------------------------------------------------
class Downloader(Thread):
# multiple threads downloader
def __init__(self, index, count, urls, datas):
# index represents the number of threads
# count represents the total number of threads
# urls represents the list of URLs nedd to be downloaded
# datas represents the list of data need to be returned.
super().__init__()
self.urls = urls
self.datas = datas
self.index = index
self.count = count
def download(self, url):
print("下载图片地址",url,)
HEADERS = {
'User-Agent': 'Mozilla/5.0 (Windows NT 10.0; Win64; x64) AppleWebKit/537.36 (KHTML, like Gecko) Chrome/88.0.4324.150 Safari/537.36 Edg/88.0.705.68'}
header = ur.Request(url, headers=HEADERS)
err = 1
while err < 11:
try:
data = ur.urlopen(header, timeout=10).read()
except:
traceback.print_exc()
print("下载错误 重新下载...")
err += 1
else:
img = cv2.imdecode(np.frombuffer(data, np.uint8), cv2.IMREAD_COLOR)
print("最小像素点", np.min(img))
if np.min(img) == 255:
print(f"图片是空白的,开始重新下载第{err}次")
err += 1
else:
return data
raise Exception("网络连接异常.")
def run(self):
for i, url in enumerate(self.urls):
print("下载任务", i, "/", len(self.urls))
tile_x, tile_y, z = url.split("&xyz=")[1].split(",")
# lon, lat = tile2lonlat(tile_x, tile_y, z)
# 如果i除线程总数1取余 不等于 第几个线程 总是0
if i % self.count != self.index:
print("warning!!!跳过该url，记录一下", url)
continue
self.datas[str(tile_x)+"_"+str(tile_y)] = self.download(url)
# ---------------------------------------------------------
# ---------------------------------------------------------
def getExtent(x1, y1, x2, y2, z, source="Google China"):
pos1x, pos1y = wgs_to_tile(x1, y1, z)
pos2x, pos2y = wgs_to_tile(x2, y2, z)
Xframe = pixls_to_mercator(
{"LT": (pos1x, pos1y), "RT": (pos2x, pos1y), "LB": (pos1x, pos2y), "RB": (pos2x, pos2y), "z": z})
for i in ["LT", "LB", "RT", "RB"]:
Xframe[i] = mercator_to_wgs(*Xframe[i])
if source == "Google":
pass
elif source == "ZKLF":
pass
elif source == "Google China":
for i in ["LT", "LB", "RT", "RB"]:
Xframe[i] = gcj_to_wgs(*Xframe[i])
else:
raise Exception("Invalid argument: source.")
return Xframe
def saveTiff(r, g, b, gt, filePath):
fname_out = filePath
driver = gdal.GetDriverByName('GTiff')
# Create a 3-band dataset
dset_output = driver.Create(fname_out, r.shape[1], r.shape[0], 3, gdal.GDT_Byte)
dset_output.SetGeoTransform(gt)
try:
proj = osr.SpatialReference()
proj.ImportFromEPSG(4326)
dset_output.SetSpatialRef(proj)
except:
print("Error: Coordinate system setting failed")
dset_output.GetRasterBand(1).WriteArray(r)
dset_output.GetRasterBand(2).WriteArray(g)
dset_output.GetRasterBand(3).WriteArray(b)
dset_output.FlushCache()
dset_output = None
print("Image Saved tif图片生成成功")
# ---------------------------------------------------------
# ---------------------------------------------------------
MAP_URLS = {
"onwer_server": "http://ip:port/geoserver/ows?service=WMS&version=1.3.0&transparent=true&request=GetMap&layers={style}&width=256&height=256&srs=EPSG%3A3857&format=image/png&bbox={bbox}&xyz={xyz}",
"Google": "http://mts0.googleapis.com/vt?lyrs={style}&x={x}&y={y}&z={z}",
"Google China": "http://mt2.google.cn/vt/lyrs={style}&hl=zh-CN&gl=CN&src=app&x={x}&y={y}&z={z}&s=Galile"}
def tile2lonlat(x, y, z):
"""
@description :切图坐标转lonlat
---------
@param x :瓦片的x轴
@param y :瓦片的y轴
@param z :瓦片的层级
-------
@Returns :[经度, 纬度]
-------
"""
x = float(x)
y = float(y)
n = math.pow(2, int(z))
lon = x / n * 360.0 - 180.0
lat = math.atan(math.sinh(math.pi * (1 - 2 * y / n)))
lat = lat * 180.0 / math.pi
return [lon, lat]
def get_url(source, x, y, z, style): #
if source == 'Google China':
url = MAP_URLS["Google China"].format(x=x, y=y, z=z, style=style)
elif source == 'Google':
url = MAP_URLS["Google"].format(x=x, y=y, z=z, style=style)
elif source == "onwer_server":
min_xy_list = tile2lonlat(x, y + 1, z)
max_xy_list = tile2lonlat(x + 1, y, z)
print("min_xy_list:", min_xy_list)
print("max_xy_list:", max_xy_list)
lng_min, lat_min = min_xy_list[0], min_xy_list[1]
lng_max, lat_max = max_xy_list[0], max_xy_list[1]
# gcj02转wgs84
# lng_min, lat_min = GCJ02_to_WGS84(min_xy_list[0], min_xy_list[1])
lng_min, lat_min = WGS84_to_WebMercator(lng_min, lat_min)
# lng_max, lat_max = GCJ02_to_WGS84(max_xy_list[0], max_xy_list[1])
lng_max, lat_max = WGS84_to_WebMercator(lng_max, lat_max)
bbox = ",".join([str(lng_min), str(lat_min), str(lng_max), str(lat_max)])
xyz = ",".join([str(x), str(y), str(z)])
url = MAP_URLS["ZKLF"].format(bbox=bbox, style=style, xyz=xyz)
else:
raise Exception("Unknown Map Source ! ")
return url
def get_urls(x1, y1, x2, y2, z, source, style):
"""
@description :左上角x1,y1右下角x2,y2
---------
@param :
-------
@Returns :
-------
"""
pos1x, pos1y = wgs_to_tile(x1, y1, z) # 左上角的瓦片图坐标
pos2x, pos2y = wgs_to_tile(x2, y2, z)
print("pos1x, pos1y:", pos1x, pos1y)
print("pos2x, pos2y:", pos2x, pos2y)
lenx = abs(pos2x - pos1x) + 1
leny = abs(pos2y - pos1y) + 1
print("Total tiles number：{x} X {y}".format(x=lenx, y=leny))
print("pos1y, pos1y + leny:", pos1y, pos1y + leny)
print("pos1x, pos1x + lenx:", pos1x, pos1x + lenx)
urls = [get_url(source, i, j, z, style) for j in range(pos1y, pos1y + leny) for i in range(pos1x, pos1x + lenx)]
return urls
# ---------------------------------------------------------
# ---------------------------------------------------------
def merge_tiles(datas, x1, y1, x2, y2, z):
pos1x, pos1y = wgs_to_tile(x1, y1, z)
pos2x, pos2y = wgs_to_tile(x2, y2, z)
lenx = pos2x - pos1x + 1
leny = pos2y - pos1y + 1
outpic = pil.new('RGBA', (lenx * 256, leny * 256))
for i, data in enumerate(datas):
picio = io.BytesIO(data)
small_pic = pil.open(picio)
y, x = i // lenx, i % lenx
outpic.paste(small_pic, (x * 256, y * 256))
print('Tiles merge completed')
return outpic
def download_tiles(urls, multi=1):
url_len = len(urls)
# datas = [None] * url_len
datas = dict()
if multi < 1 or multi > 20 or not isinstance(multi, int):
raise Exception("multi of Downloader shuold be int and between 1 to 20.")
tasks = [Downloader(i, multi, urls, datas) for i in range(multi)]
for i in tasks:
i.start()
for i in tasks:
i.join()
return datas
# ---------------------------------------------------------
# ---------------------------------------------------------
def main(left, top, right, bottom, zoom, filePath, style='s', server="Google China", zone_id=""):
"""
Download images based on spatial extent.
East longitude is positive and west longitude is negative.
North latitude is positive, south latitude is negative.
Parameters
----------
left, top : left-top coordinate, for example (100.361,38.866)
right, bottom : right-bottom coordinate
z : zoom
filePath : File path for storing results, TIFF format
style :
m for map;
s for satellite;
y for satellite with label;
t for terrain;
p for terrain with label;
h for label;
source : Google China (default) or Google
"""
# ---------------------------------------------------------
# Get the urls of all tiles in the extent
urls = get_urls(left, top, right, bottom, zoom, server, style)
print("瓦片图总数:", len(urls))
# Group URLs based on the number of CPU cores to achieve roughly equal amounts of tasks
urls_group = [urls[i:i + math.ceil(len(urls) / 2)] for i in
range(0, len(urls), math.ceil(len(urls) / 2))]
# urls_group = [urls[i:i + math.ceil(len(urls) / multiprocessing.cpu_count())] for i in
# range(0, len(urls), math.ceil(len(urls) / multiprocessing.cpu_count()))]
print(urls_group)
# return False
# Each set of URLs corresponds to a process for downloading tile maps
print('Tiles downloading......瓦片图下载中')
# results = []
pool = multiprocessing.Pool(2)
# pool = multiprocessing.Pool(multiprocessing.cpu_count())
print("cpu_count:", multiprocessing.cpu_count())
print("pool", pool)
results = pool.map(download_tiles, urls_group)
pool.close()
pool.join()
# print("results:", type(results[0]))
image_number = 1
for res in results:
for key in res.keys():
print(f"开始保存第{image_number}张图片")
print("image_name:", key)
x = str(key).split("_")[0]
y = str(key).split("_")[1]
parent_dir = Path(f"D:\\tiles\\{zone_id}\\tiles\\{zoom}\\{x}") # 父级目录
if not parent_dir.exists():
parent_dir.mkdir(exist_ok=True, parents=True)
with open(f"D:\\tiles\\{zone_id}\\tiles\\{zoom}\\{x}\\{y}.png", "wb") as f:
# print(_)
print("sys.getsizeof(res[key]):", sys.getsizeof(res[key]))
f.write(res[key])
image_number += 1
# print(res.get())
# result = [x for j in results for x in j]
print('Tiles download complete 瓦片图 下载成功')
# # break
# Combine downloaded tile maps into one map
# outpic = merge_tiles(result, left, top, right, bottom, zoom)
# outpic = outpic.convert('RGB')
# r, g, b = cv2.split(np.array(outpic))
# # Get the spatial information of the four corners of the merged map and use it for outputting
# extent = getExtent(left, top, right, bottom, zoom, server)
# gt = (extent['LT'][0], (extent['RB'][0] - extent['LT'][0]) / r.shape[1], 0, extent['LT'][1], 0,
# (extent['RB'][1] - extent['LT'][1]) / r.shape[0])
# saveTiff(r, g, b, gt, filePath)
# ---------------------------------------------------------
if __name__ == '__main__':
from sqls import select_zone_info
import json
start_time = time.time()
zone_id = "1547212382202232832"
# main(118.901, 32.066,118.934, 32.109, 18, r'.\Temp\test.tif', server="Google China")
# main(100.361, 38.866, 100.386, 38.839, 18, r'.\Temp\test.tif', server="Google China")
# main(97.376955,37.133309,97.4074156,37.118448, 5, r'.\Temp\test-one.tif', server="onwer_server", style="onwer_server_AREA:onwer_server_super_hd")
# path = [{"lng": 118.909207, "lat": 32.081909}, {"lng": 118.912005, "lat": 32.081703}, {"lng": 118.911776, "lat": 32.081512}, {"lng": 118.909180, "lat": 32.080421}]
left_top_x, left_top_y, right_buttom_x, right_buttom_y = get_zone_gps_max(zone_id)
# print(float(left_top[0]), float(left_top[1]), float(right_buttom[0]), float(right_buttom[1]))
for zoom in range(18, 19):
parent_dir = Path(f"D:\\tiles\\{zone_id}\\tiles\\{zoom}") # 父级目录
if not parent_dir.exists():
parent_dir.mkdir(exist_ok=True, parents=True)
if 1:
main(float(left_top_x),float(left_top_y),float(right_buttom_x),float(right_buttom_y), zoom, r'.\\Temp\\test-one.tif', server="onwer_server", style="onwer_server_AREA:onwer_server_super_hd", zone_id=zone_id)
#main(97.376955,37.133309,97.376955,37.133309, 19, r'.\Temp\test-one.tif', server="Google China")
end_time = time.time()
print('lasted a total of {:.2f} seconds'.format(end_time - start_time))2.合并瓦片图 谷歌的瓦片图长这样.瓦片图是金字塔类型的，这里就不多做解释了.
合并代码 就是创建一个大的画布，然后把下载好的瓦片图一张张的贴上去，没有难度
import glob
import re
from PIL import Image
files = glob.glob('D:\\tiles\\1547212382202232832\\tiles\\17\\*\\*.png')
# print(re.findall(r"\d+", files[0])[-2:])
# print(tuple(int(i) for i in re.findall(r"\d+", files[0])[-2:]))
files.sort(key=lambda x: tuple(int(i) for i in re.findall(r"\d+", x)[-2:]))
# print(files)
imagefiles = {}
for item in files:
match = re.findall(r'\d+', item)
# print(match[-2])
pre = int(match[-2])
if not imagefiles.get(pre):
imagefiles[pre] = []
imagefiles[pre].append(item)
imagefiles = sorted(zip(imagefiles.keys(), imagefiles.values()))
print(imagefiles)
total_width = len(imagefiles) * 256
total_height = len(imagefiles[0][1]) * 256
new_image = Image.new('RGB', (total_width, total_height))
x_offset = 0
for item in imagefiles:
y_offset = 0
# print(item[1])
images = map(Image.open, item[1])
# print(list(images))
for subitem in images:
new_image.paste(subitem, (x_offset, y_offset))
y_offset += subitem.size[0]
_x = subitem.size[0]
# print(list(images))
x_offset += _x
new_image.save('merge.jpg', quality=90)3.在合并好的瓦片图上绘制我们的mark点和多边形。
思路：首先我们合并好的瓦片图上只有像素一个计量单位，如果要化gps点上去的话，就要找到一个全局的参考坐标。最好的参考坐标就是左上角第一块瓦片。找到左上角的点坐标。因为切出来的瓦片像素是我们自定义的，我用的是256*256,同时可以获取到瓦片的实际长度和宽度(就是bbox参数/墨卡托投影),由此我们可以算出单位像素对应的实际长度(单位是米)。 接下来我们只需要找到左上角第一块瓦片的左上角的点坐标即可。 通过左上角瓦片图的gps可以算出对应的瓦片图坐标，根据瓦片图坐标既可以算出 瓦片的左下角坐标和右上角坐标，既得左上角坐标。再转成墨卡托投影即可作为全局的参考坐标了。 ps:瓦片的范围用最小外接矩形来算。
def get_zone_gps_max(zone_id):
path_info = select_zone_info(zone_id)
path = json.loads(path_info[0]["path"])
path = [list(WGS84_to_WebMercator(_["lng"], _["lat"])) for _ in path]
print(path)
# left_top, right_buttom = get_maxarea_box(path)
min_x, max_x, min_y, max_y = get_maxarea_box(path)
print("右下→左下→左上→右上:", min_x, max_x, min_y, max_y)
left_top_x, left_top_y = WebMercator_to_WGS84(min_x, max_y)
right_buttom_x, right_buttom_y = WebMercator_to_WGS84(max_x, min_y)
return (float(left_top_x), float(left_top_y), float(right_buttom_x), float(right_buttom_y))
def get_first_tile_poi(zone_id, z=17):
"""
@description :获取第一块瓦片图的坐标
---------
@param :
-------
@Returns :左下角和右上角
-------
"""
left_top_x, left_top_y, right_buttom_x, right_buttom_y = get_zone_gps_max(zone_id)
pos1x, pos1y = wgs_to_tile(left_top_x, left_top_y, z) # 左上角的瓦片图坐标
min_xy_list = tile2lonlat(pos1x, pos1y + 1, z)
max_xy_list = tile2lonlat(pos1x + 1, pos1y, z)
lng_min, lat_min = min_xy_list[0], min_xy_list[1]
lng_max, lat_max = max_xy_list[0], max_xy_list[1]
# gcj02转wgs84
# lng_min, lat_min = GCJ02_to_WGS84(min_xy_list[0], min_xy_list[1])
lng_min, lat_min = WGS84_to_WebMercator(lng_min, lat_min)
# lng_max, lat_max = GCJ02_to_WGS84(max_xy_list[0], max_xy_list[1])
lng_max, lat_max = WGS84_to_WebMercator(lng_max, lat_max)
bbox = [lng_min, lat_min, lng_max, lat_max]
return bbox有疑问或者错误的地方 可以留言讨论 互相学习
本文章为转载内容，我们尊重原作者对文章享有的著作权。如有内容错误或侵权问题，欢迎原作者联系我们进行内容更正或删除文章。
