文章目录
背景
图像重采样是指改变图像的分辨率,通常指分辨率变小,即像素的宽高变长,是科研常用的操作。比如我们手头只有300m分辨率的图像,但是为了统一不同数据源的数据精度,需要统一为1000m分辨率,这时候就需要重采样。
实战
import numpy as np
import rasterio
from rasterio.warp import calculate_default_transform, reproject
from rasterio import crs
dataset = rasterio.open('D:/A_2021寒假/城市群相关/Data/ESACCI土地覆盖数据1992-2015/final/mask_proj_con_mask_ESACCI-LC-L4-LCCS-Map-300m-P1Y-1992-v2.0.7.tif')
src_img = 'input.tif' # 输入图像
dst_img = 'output.tif' # 输出图像
#dst_crs = crs.CRS.from_proj4('+proj=lcc +lat_1=15 +lat_2=65 +lat_0=30 +lon_0=95 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs') # 输出图像坐标系
with rasterio.open(src_img) as src_ds:
# 计算在新空间参考系下的仿射变换参数,图像尺寸
dst_transform, dst_width, dst_height = calculate_default_transform(
src_ds.crs, # 输入坐标系
src_ds.crs, # 输出坐标系
src_ds.width, # 输入图像宽
src_ds.height, # 输入图像高
resolution=(1000,1000), # 输出图像分辨率,单位m
*src_ds.bounds) # 输入数据源的图像范围
# 更新数据集的元数据信息
profile = src_ds.meta.copy()
profile.update({
'crs': src_ds.crs,
'transform': dst_transform,
'width': dst_width,
'height': dst_height
})
# 重投影并写入数据
with rasterio.open(dst_img, 'w', **profile) as dst_ds:
for i in range(1, src_ds.count + 1): # 遍历每个图层,通常只需要第一层即可
src_array = src_ds.read(i)
dst_array = np.empty((dst_height, dst_width), dtype=profile['dtype']) # 初始化输出图像数据
# 重投影
reproject(
# 源文件参数
source=src_array,
src_crs=src_ds.crs,
src_transform=src_ds.transform,
# 目标文件参数
destination=dst_array,
dst_transform=dst_transform,
dst_crs=src_ds.crs,
# 其它配置
resampling=Resampling.average,
num_threads=2)
# 写入图像
dst_ds.write(dst_array, i)
这里需要注意的是,我提供的示例是输入和输出坐标系一致,仅分辨率不同的情况,如果输入和输出坐标系不一致,是需要改变dst_crs的。投影变换我会在下一篇博文具体展开,这里我将两个功能拆开。
另外,重采样的图像必须有投影坐标系,因为只有投影后的平面坐标才有米这种长度单位,地理坐标系是用度分秒这样的角度单位来度量的。