Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
Helper file to create a geotiff from a given image or from an automatically downloaded satellite image from google maps which is centered around a given coordinate .
"""
BST License
# Description
File to create geotiffs depends on gdal installed on system. Can be used under two modes:
1) Download a satellite overview screenshot from google maps around a center
coordinate. Can give zoom (max is usually 20, but depends on area - will return
an image with error message if zoom is too high) and image size
(defaults to max 640x640)
2) Convert an existing image into a geotiff by giving UL and LR corner
coordinates of image.
Example:
Create geotiff from google around Barbados:
python create_geotiff.py --CLL "13.159717, -59.657345" -z 12
Create geotiff from existing map image around Barbados:
python create_geotiff.py -s ~/Desktop/map_image.png --ULLR "13.1601374805, -59.6577740797, 13.1592965187, -59.6569159218"
# Author: Johanna Hansen <jh1736@gmail.com>
# Edits:
created 19 Dec 2016
# Install (tested on Python 2.7)
pip install scipy LatLon urllib
"""
import os
from scipy.misc import imread, imsave
import math
def download_screenshot(img_path, LC, s_size_w=8000, s_size_h=8000, zoom=20):
import urllib
from io import BytesIO
import LatLon
latcenter = LC[0].strip()
loncenter = LC[1].strip()
# max zoom is 21, 1 is the whole earth
url = "http://maps.googleapis.com/maps/api/staticmap?maptype=satellite&center=%s,%s&size=%sx%s&zoom=%s" %(latcenter,loncenter,s_size_w, s_size_h, zoom)
buffer = BytesIO(urllib.urlopen(url).read())
img = imread(buffer)
# determine how many meters are represented in each pixel
km_per_pixel = (156543.03392 * math.cos((float(latcenter) * math.pi) / 180.0)) / (2.0**int(zoom)) *1e-3
h, w, _ = img.shape
hm = h/2
wm = w/2
imsave(img_path, img)
off_center_wkm = w/2.0 * km_per_pixel
off_center_hkm = h/2.0 * km_per_pixel
CLat = LatLon.Latitude(latcenter)
CLon = LatLon.Longitude(loncenter)
center = LatLon.LatLon(CLat, CLon)
# North is 0
# West is 270
CT = center.offset(0, off_center_hkm)
CB = center.offset(180, off_center_hkm)
UL = CT.offset(270, off_center_wkm).to_string('D')
LR = CB.offset(90, off_center_wkm).to_string('D')
print('Calculated UL, LR coordinates of image is:', UL, LR)
return UL, LR
def create_geotiff(img_path, UL, LR):
if not os.path.exists(img_path):
print("ERROR: image path to convert to geotiff does not exists: %s" %img_path)
geotiff_path = ''.join(img_path.split('.')[:-1]) + '_geo.tiff'
print("Creating geotiff named: %s" %geotiff_path)
# -a_ullr ulx uly lrx lry:
#Assign/override the georeferenced bounds of the output file.
#This assigns georeferenced bounds to the output file, ignoring what would have been derived from the source file.
#So this does not cause reprojection to the specified SRS.
cmd = "gdal_translate -a_nodata 0 -of GTiff -a_srs EPSG:4326 -a_ullr %s %s %s %s %s %s" %(UL[1], UL[0],
LR[1], LR[0],
img_path, geotiff_path)
print("Running gdal command:", cmd)
os.system(cmd)
if __name__ == '__main__':
import argparse
parser = argparse.ArgumentParser(
description="""
File to create geotiffs. Can be used under two modes:
1) Download a satellite overview screenshot from google maps around a center
coordinate. Can give zoom (max is usually 20, but depends on area - will return
an image with error message if zoom is too high) and image size
(defaults to max 640x640)
2) Convert an existing image into a geotiff by giving UL and LR corner
coordinates of image.
Example:
Create geotiff from google around Barbados:
python create_geotiff.py --CLL "13.159717, -59.657345" -z 12
Create geotiff from existing map image around Barbados:
python create_geotiff.py -s ~/Desktop/map_image.png --ULLR "13.1601374805, -59.6577740797, 13.1592965187, -59.6569159218"
""")
parser.add_argument('--ULLR',
type=str,
help='Coordinates of upper-left and lower right coordinates of provided screenshot: "ULlat,ULlon,LRlat,LRlon"')
parser.add_argument('--CLL',
type=str,
help='Center coordinates of google satellite imagery to be captured if we should create a new base image: "lat,lon"')
parser.add_argument('-s', type=str, default='map_image.png', help='filepath of image to be turned into geotiff or location to store google satellite image')
parser.add_argument('-z', type=str, default=20, help='zoom level of image to be downloaded from google, min zoom is 2 (whole earth), max zoom is 21')
parser.add_argument('--pixel_size', default="640,640", type=str, help='Pixel size of image to be downloaded: "pixel_width,pixel_height" maximum size is 640x640')
args = parser.parse_args()
map_image = args.s
zoom = args.z
pixel_size_w, pixel_size_h = args.pixel_size.split(',')
if args.CLL is not None:
LC = args.CLL.split(',')
print("Received Center Coordinates, creating geotiff from google at %s" %map_image)
UL, LR = download_screenshot(map_image, LC, pixel_size_w, pixel_size_h, zoom)
create_geotiff(map_image, UL, LR)
elif args.ULLR is not None:
print("Received Corner Coordinates, creating geotiff from provided image: %s" %map_image)
ULLR = [x.strip() for x in args.ULLR.split(',')]
UL = ULLR[:2]
LR = ULLR[:-2]
create_geotiff(map_image, UL, LR)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment