python-visualization / folium

Python Data. Leaflet.js Maps.
https://python-visualization.github.io/folium/
MIT License
6.9k stars 2.22k forks source link

ImageOverlay alignment issue #1631

Closed harshmistry96 closed 2 years ago

harshmistry96 commented 2 years ago

I am trying to overlay a GeoTiff image on folium map but I think there is an alignment issue. Already tried with mercator_project=True, but still there is an alignment issue.

Link to Geotiff: https://sedac.ciesin.columbia.edu/downloads/data/gpw-v4/gpw-v4-population-density-rev11/gpw-v4-population-density-rev11_2020_2pt5_min_tif.zip

Code is as follow:

import xarray as xr
import numpy as np
import numpy.ma as ma
import os
path = os.getcwd()
world_path =os.path.join(path,'World_data')
file_path = world_path+'/gpw_v4_population_density_rev11_2020_2pt5_min_3857.tif''
#%%
xr_img = xr.open_rasterio(file_path)
nodata  = xr_img.nodatavals[0]

xr_img  = xr_img.where(xr_img>nodata, np.nan)
arr_dem = xr_img.values.astype(np.float64)
#Calucate center and bounds 
clat = xr_img.y.values.mean()
clon = xr_img.x.values.mean()

mlat = xr_img.y.values.min()
mlon = xr_img.x.values.min()

xlat = xr_img.y.values.max()
xlon = xr_img.x.values.max()

print(clat, clon, mlat, mlon, xlat, xlon)
#%%
def colorize(array, cmap='viridis'):
    normed_data = (array - array.min()) / (array.max() - array.min())    
    # normed_data = np.log10(array)
    cm = plt.cm.get_cmap(cmap,20)    

    return cm(normed_data) 

abb=arr_dem[0]
data         = ma.masked_invalid(abb)
colored_data = colorize(np.log10(data), cmap='jet')

#%%
mapa = folium.Map(location=[clat, clon], 
                   tiles="cartodbpositron",
                   zoom_start=2)

folium.raster_layers.ImageOverlay(colored_data,
                                    [[mlat,mlon], [xlat, xlon]],
                                    mercator_project=True,
                                  opacity=0.5).add_to(mapa)

Result is as follow:

image

ocefpaf commented 2 years ago

Link to Geotiff: https://sedac.ciesin.columbia.edu/downloads/data/gpw-v4/gpw-v4-population-density-rev11/gpw-v4-population-density-rev11_2020_2pt5_min_tif.zip

Cannot download that. It requires a sign in. I'd like to check the image CRS, if there is one, and try warp it. Here is an example:

import cartopy.crs as ccrs
from cartopy.img_transform import warp_array

source_extent = [lon.min(), lon.max(), lat.min(), lat.max()]

new_data = warp_array(
    colored_data,
    target_proj=ccrs.GOOGLE_MERCATOR,
    source_proj=ccrs.PlateCarree(),  # change here for the source projection you have
    target_res=data.shape,
    source_extent=source_extent,
    target_extent=None,
    mask_extrapolated=False,
)

m = folium.Map(location=[lat.mean(), lon.mean()], zoom_start=1)

folium.raster_layers.ImageOverlay(
    image=new_data[0],
    bounds=[[lat.min(), lon.min()], [lat.max(), lon.max()]],
    opacity=0.25,
).add_to(m)
harshmistry96 commented 2 years ago

Here is the file: The crs for geotiff is epsg 3857. Please let me know your thoughts.

gpw_v4_population_density_rev11_2020_2pt5_min_3857.zip

ocefpaf commented 2 years ago

Can't promise but I'll play with it and see if I can make it work. I'll report back if successful.

ocefpaf commented 2 years ago

You can, and you should probably do it, with warp_array like in my example above. However, this image was quite big for me to handle on my laptop so I used the sub-sampled cartopy image instead (it is virtually the same thing b/c cartopy is using warp_array under the hood). I went that route just to be sure I had the proper metadata.

https://nbviewer.org/gist/ocefpaf/5e570a1b11d51a9f65eae8592290f4c2

BTW, closing this b/c it is now a folium problem but data handling.