{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"import folium\n",
"import matplotlib\n",
"from matplotlib.pyplot import imread\n",
"from matplotlib.colors import Normalize\n",
"from folium import raster_layers"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"image = '/notebooks/resources/gpm/gpm_1d.20190531.tif'"
]
},
{
"cell_type": "code",
"execution_count": 22,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"10.0 -180.0\n",
"-60.0 -180.0\n",
"-60.0 -30.0\n",
"10.0 -30.0\n",
"ext = [[10.0, -180.0], [-60.0, -180.0], [-60.0, -30.0], [10.0, -30.0]]\n",
"geo_ext = [[0.17453292519943298, -3.141592653589794], [-1.0471975511965979, -3.141592653589794], [-6679169.447596415, -3503549.843504376], [1113194.9079327357, -3503549.843504376]]\n"
]
}
],
"source": [
"from osgeo import gdal,ogr,osr\n",
"\n",
"def GetExtent(gt,cols,rows):\n",
" ''' Return list of corner coordinates from a geotransform\n",
"\n",
" @type gt: C{tuple/list}\n",
" @param gt: geotransform\n",
" @type cols: C{int}\n",
" @param cols: number of columns in the dataset\n",
" @type rows: C{int}\n",
" @param rows: number of rows in the dataset\n",
" @rtype: C{[float,...,float]}\n",
" @return: coordinates of each corner\n",
" '''\n",
" ext=[]\n",
" xarr=[0,cols]\n",
" yarr=[0,rows]\n",
"\n",
" for px in xarr:\n",
" for py in yarr:\n",
" x=gt[0]+(px*gt[1])+(py*gt[2])\n",
" y=gt[3]+(px*gt[4])+(py*gt[5])\n",
" ext.append([y,x])\n",
" print (y,x)\n",
" yarr.reverse()\n",
" return ext\n",
"\n",
"def ReprojectCoords(coords,src_srs,tgt_srs):\n",
" ''' Reproject a list of x,y coordinates.\n",
"\n",
" @type geom: C{tuple/list}\n",
" @param geom: List of [[x,y],...[x,y]] coordinates\n",
" @type src_srs: C{osr.SpatialReference}\n",
" @param src_srs: OSR SpatialReference object\n",
" @type tgt_srs: C{osr.SpatialReference}\n",
" @param tgt_srs: OSR SpatialReference object\n",
" @rtype: C{tuple/list}\n",
" @return: List of transformed [[x,y],...[x,y]] coordinates\n",
" '''\n",
" trans_coords=[]\n",
" transform = osr.CoordinateTransformation( src_srs, tgt_srs)\n",
" for x,y in coords:\n",
" x,y,z = transform.TransformPoint(x,y)\n",
" trans_coords.append([x,y])\n",
" return trans_coords\n",
"\n",
"raster=image\n",
"ds=gdal.Open(raster)\n",
"\n",
"gt=ds.GetGeoTransform()\n",
"cols = ds.RasterXSize\n",
"rows = ds.RasterYSize\n",
"\n",
"ext=GetExtent(gt,cols,rows)\n",
"print(\"ext = \" + str(ext))\n",
"\n",
"src_srs=osr.SpatialReference()\n",
"src_srs.ImportFromWkt(ds.GetProjection())\n",
"# tgt_srs=osr.SpatialReference()\n",
"tgt_srs.ImportFromEPSG(3857)\n",
"# tgt_srs = src_srs.CloneGeogCS()\n",
"\n",
"geo_ext=ReprojectCoords(ext,src_srs,tgt_srs)\n",
"print(\"geo_ext = \" + str(geo_ext))"
]
},
{
"cell_type": "code",
"execution_count": 41,
"metadata": {
"scrolled": false
},
"outputs": [
{
"data": {
"text/html": [
""
],
"text/plain": [
""
]
},
"execution_count": 41,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"m = folium.Map(\n",
" location = [-22, -114]\n",
" , zoom_start = 1\n",
" , control_scale = True \n",
"# , tiles = 'Stamen Terrain'\n",
"# , crs = 'EPSG4326'\n",
")\n",
"\n",
"data = matplotlib.pyplot.imread(image)\n",
"\n",
"# Image bounds on the map in the form\n",
"# [[lat_min, lon_min], [lat_max, lon_max]]\n",
"m.add_child(raster_layers.ImageOverlay(\n",
" data\n",
" , opacity = 0.7\n",
"# , bounds = [[-60, -30], [10, -180]]\n",
"# , bounds = [[10, -180], [-60, -30]]\n",
"# , bounds = [[10, -30], [-60, -180]]\n",
" , bounds = [ext[2], ext[0]]\n",
"# , bounds = [[0.17453292519943298, -3.141592653589794], [-6679169.447596415, -3503549.843504376]]\n",
" , mercator_project = True\n",
")\n",
" )\n",
"\n",
"# folium.Marker(\n",
"# ext[2]\n",
"# , popup = '-60, -30'\n",
"# , tooltip = '-60, -30'\n",
"# ).add_to(m)\n",
"\n",
"# folium.Marker(\n",
"# ext[0]\n",
"# , popup = '10, -180'\n",
"# , tooltip = '10, -180'\n",
"# ).add_to(m)\n",
"\n",
"# folium.Marker(\n",
"# [10, -30]\n",
"# , popup = '10, -30'\n",
"# , tooltip = '10, -30'\n",
"# ).add_to(m)\n",
"\n",
"# folium.Marker(\n",
"# [-60, -180]\n",
"# , popup = '-60, -180'\n",
"# , tooltip = '-60, -180'\n",
"# ).add_to(m)\n",
"\n",
"folium.Marker(\n",
" [-8399737.89, -20037508.34]\n",
" , popup = 'uh...'\n",
" , tooltip = 'uh...'\n",
").add_to(m)\n",
"\n",
"\n",
"\n",
"m"
]
},
{
"cell_type": "code",
"execution_count": 42,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"{'_children': OrderedDict([('openstreetmap',\n",
" ),\n",
" ('image_overlay_029dcef611864ab69a5251c1e152a267',\n",
" ),\n",
" ('marker_1f4d013e5998403a9052985742ec7dd6',\n",
" )]),\n",
" '_env': ,\n",
" '_id': '264585f73b94452ca1f725f3a674fe99',\n",
" '_name': 'Map',\n",
" '_parent': ,\n",
" '_png_image': None,\n",
" 'control_scale': True,\n",
" 'crs': 'EPSG3857',\n",
" 'global_switches': ,\n",
" 'height': (100.0, '%'),\n",
" 'left': (0.0, '%'),\n",
" 'location': [-22, -114],\n",
" 'max_bounds': False,\n",
" 'max_lat': 90,\n",
" 'max_lon': 180,\n",
" 'min_lat': -90,\n",
" 'min_lon': -180,\n",
" 'no_wrap': False,\n",
" 'objects_to_stay_in_front': [],\n",
" 'png_enabled': False,\n",
" 'position': 'relative',\n",
" 'top': (0.0, '%'),\n",
" 'width': (100.0, '%'),\n",
" 'world_copy_jump': False,\n",
" 'zoom_control': True,\n",
" 'zoom_start': 1}"
]
},
"execution_count": 42,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"vars(m)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"print(data)\n",
"print(data.shape)\n",
"# data = data.transpose()\n",
"# print(data)\n",
"# print(data.shape)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"!gdalinfo '/notebooks/resources/gpm/gpm_1d.20190531.tif'"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# First: read the geotiff image with GDAL.\n",
"from osgeo import gdal, osr\n",
"\n",
"gdal.UseExceptions()\n",
"\n",
"ds = gdal.Open(image)\n",
"data = ds.ReadAsArray()\n",
"gt = ds.GetGeoTransform()\n",
"proj = ds.GetProjection()\n",
"\n",
"inproj = osr.SpatialReference()\n",
"inproj.ImportFromWkt(proj)\n",
"\n",
"print('\\n\\n## ds ##:\\n\\n' + str(ds))\n",
"print('\\n\\n## data ##:\\n\\n' + str(data))\n",
"print('\\n\\n## gt ##:\\n\\n' + str(gt))\n",
"print('\\n\\n## proj ##:\\n\\n' + str(proj))\n",
"print('\\n\\n## inproj ##:\\n\\n' + str(inproj))"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.5.2"
}
},
"nbformat": 4,
"nbformat_minor": 2
}