Geotiff影像与Folium地图叠加时内部像素偏移问题求助

问题根源 #

你遇到的整体位置正确但内部像素偏移的问题，核心原因是GDAL读取的栅格行顺序与Leaflet（Folium基于Leaflet）的渲染顺序相反：GDAL返回的数组第一行对应影像的顶部，而Leaflet的ImageOverlay是从底部开始向上绘制，导致像素行错位，从而出现轻微偏移。

解决方案 #

1. 反转栅格数组行顺序

直接调整GDAL读取的数组，匹配Leaflet的渲染逻辑：
将img = ds.GetRasterBand(1).ReadAsArray()修改为：

img = np.flipud(ds.GetRasterBand(1).ReadAsArray())

2. 严谨处理坐标转换（可选但推荐）

如果你的GeoTIFF不是EPSG:4326投影，手动计算边界可能因投影不匹配产生误差，建议用GDAL的坐标转换工具获取准确的WGS84边界：

from osgeo import osr

# 获取影像原始投影
srs = osr.SpatialReference()
srs.ImportFromWkt(ds.GetProjection())
# 定义目标投影（WGS84，Folium的bounds需要此坐标）
target_srs = osr.SpatialReference()
target_srs.ImportFromEPSG(4326)
transform = osr.CoordinateTransformation(srs, target_srs)

# 计算影像四至点并转换为WGS84
ul_x, ul_y = GT[0], GT[3]
lr_x, lr_y = GT[0] + width*xRes, GT[3] + height*yRes
ul_lon, ul_lat, _ = transform.TransformPoint(ul_x, ul_y)
lr_lon, lr_lat, _ = transform.TransformPoint(lr_x, lr_y)

# 设置正确的bounds格式：[[西南纬度, 西南经度], [东北纬度, 东北经度]]
bounds = [[lr_lat, ul_lon], [ul_lat, lr_lon]]

3. 正确启用Mercator投影转换

如果影像原始投影是EPSG:4326，给ImageOverlay添加Mercator_project=True参数，确保自动转换为地图使用的EPSG:3857投影。

修改后的完整代码 #

#!/usr/bin/env python3
import numpy as np
from osgeo import gdal, osr
import folium
import branca.colormap as cm
from matplotlib import colors as colors

def mapvalue2color(value, cmap): 
    if (value == 0):
        return (1, 0, 0, 0)
    else:
        return colors.to_rgba(cmap(value), 0.7)   # alpha channel=0.7 tuple (r, g, b, a)

colormap = cm.LinearColormap(["red", "magenta"], vmin=6, vmax=8) 

# 初始化地图，EPSG:3857是Leaflet默认CRS，无需修改
folium_map = folium.Map(location=[46.3, 3.5], tiles='openstreetmap', zoom_start=10, max_zoom=18)

ds = gdal.Open("./data.tif")
# 反转栅格行顺序，匹配Leaflet渲染逻辑
img = np.flipud(ds.GetRasterBand(1).ReadAsArray())
GT = ds.GetGeoTransform()
xRes = GT[1]
yRes = GT[5]
width = ds.RasterXSize
height = ds.RasterYSize 

# 处理坐标转换，确保bounds使用WGS84坐标
srs = osr.SpatialReference()
srs.ImportFromWkt(ds.GetProjection())
target_srs = osr.SpatialReference()
target_srs.ImportFromEPSG(4326)
transform = osr.CoordinateTransformation(srs, target_srs)

# 计算影像四至点并转换为WGS84
ul_x, ul_y = GT[0], GT[3]
lr_x, lr_y = GT[0] + width*xRes, GT[3] + height*yRes
ul_lon, ul_lat, _ = transform.TransformPoint(ul_x, ul_y)
lr_lon, lr_lat, _ = transform.TransformPoint(lr_x, lr_y)

# 正确设置bounds：[[西南角(lat, lon)], [东北角(lat, lon)]]
bounds = [[lr_lat, ul_lon], [ul_lat, lr_lon]]

folium.raster_layers.ImageOverlay(
    image=img,
    opacity=0.5,
    bounds=bounds,
    colormap=lambda value: mapvalue2color(value, colormap),
    Mercator_project=True  # 影像为EPSG:4326时开启
).add_to(folium_map)

# 转换标记点坐标（如果原始投影非WGS84）
marker_x = GT[0] + xRes*9
marker_y = GT[3] + yRes*4
marker_lon, marker_lat, _ = transform.TransformPoint(marker_x, marker_y)
folium.Marker(
    location=[marker_lat, marker_lon],
    popup=folium.Popup(f"lat : {marker_lat:.6f}<br>lon : {marker_lon:.6f}", max_width=500),
).add_to(folium_map)

map_path = "./map.html"
folium_map.save(map_path)

说明 #

修改后，栅格行顺序与Leaflet渲染逻辑一致，坐标转换也更准确，能彻底解决像素偏移问题。如果你的GeoTIFF本身就是EPSG:4326投影，坐标转换部分可以简化，但保留该代码能兼容更多投影场景。

内容的提问来源于stack exchange，提问作者vincentL