Python读取和处理TIFF数据教程 (由简入深)
阶段一:基础篇 - 读取和查看TIFF数据
1. 安装必要的库
首先安装处理TIFF文件常用的库:
pip install pillow numpy matplotlib rasterio
2. 使用Pillow读取TIFF文件
from PIL import Image
import matplotlib.pyplot as plt# 打开TIFF文件
img = Image.open('example.tif')# 显示基本信息
print(f"格式: {img.format}")
print(f"大小: {img.size}")
print(f"模式: {img.mode}")# 显示图像
plt.imshow(img)
plt.show()
3. 转换为NumPy数组进行处理
import numpy as np# 将TIFF图像转换为NumPy数组
img_array = np.array(img)# 打印数组信息
print(f"数组形状: {img_array.shape}")
print(f"数据类型: {img_array.dtype}")
print(f"最小值: {img_array.min()}, 最大值: {img_array.max()}")# 简单处理 - 反转图像
inverted_img = 255 - img_array # 假设是8位图像
plt.imshow(inverted_img, cmap='gray')
plt.show()
阶段二:中级篇 - 使用rasterio处理地理TIFF
1. 读取地理TIFF文件
import rasterio# 打开地理TIFF文件
with rasterio.open('geospatial.tif') as src:# 读取所有波段数据data = src.read()# 显示元数据print(f"波段数: {src.count}")print(f"宽度: {src.width}, 高度: {src.height}")print(f"坐标参考系统(CRS): {src.crs}")print(f"地理变换: {src.transform}")# 显示第一个波段plt.imshow(data[0], cmap='viridis')plt.colorbar()plt.show()
2. 多波段处理和统计
with rasterio.open('multiband.tif') as src:# 读取各个波段red = src.read(1) # 红波段green = src.read(2) # 绿波段blue = src.read(3) # 蓝波段# 计算NDVI (假设第4波段是近红外)nir = src.read(4)ndvi = (nir - red) / (nir + red + 1e-10) # 避免除以0# 显示NDVIplt.imshow(ndvi, cmap='RdYlGn', vmin=-1, vmax=1)plt.colorbar(label='NDVI')plt.title('NDVI Map')plt.show()
3. 裁剪和重投影
from rasterio.warp import calculate_default_transform, reproject, Resampling# 裁剪示例
with rasterio.open('large_image.tif') as src:# 定义裁剪窗口 (左上x, 左上y, 右下x, 右下y)window = rasterio.windows.Window(1000, 1000, 500, 500)clipped = src.read(window=window)# 显示裁剪结果plt.imshow(clipped[0])plt.show()# 重投影示例
dst_crs = 'EPSG:3857' # Web墨卡托投影
with rasterio.open('input.tif') as src:transform, width, height = calculate_default_transform(src.crs, dst_crs, src.width, src.height, *src.bounds)kwargs = src.meta.copy()kwargs.update({'crs': dst_crs,'transform': transform,'width': width,'height': height})with rasterio.open('reprojected.tif', 'w', **kwargs) as dst:for i in range(1, src.count + 1):reproject(source=rasterio.band(src, i),destination=rasterio.band(dst, i),src_transform=src.transform,src_crs=src.crs,dst_transform=transform,dst_crs=dst_crs,resampling=Resampling.nearest)
阶段三:高级篇 - 大型TIFF处理和优化
1. 分块处理大型TIFF
# 使用rasterio的分块读取功能
with rasterio.open('very_large.tif') as src:# 获取推荐的块大小block_shapes = src.block_shapesprint(f"块大小: {block_shapes}")# 迭代处理每个块for ji, window in src.block_windows(1):# 读取当前块block = src.read(window=window)# 处理块数据 (示例: 计算平均值)block_mean = np.mean(block)print(f"块 {ji} 的平均值: {block_mean}")# 这里可以添加你的处理逻辑# processed_block = some_processing(block)# 写入输出文件 (如果需要)# dst.write(processed_block, window=window)
2. 使用Dask处理超大型TIFF
import dask.array as da
import rasterio
from rasterio.windows import Window# 创建延迟加载的dask数组
with rasterio.open('huge.tif') as src:# 定义块大小 (例如 1024x1024)chunk_size = (1, 1024, 1024) # (波段, 高, 宽)# 创建dask数组dask_array = da.from_array(src, chunks=chunk_size)# 现在可以对dask数组执行各种操作 (延迟执行)mean_computation = dask_array.mean()# 实际计算结果 (触发计算)print(f"整个图像的平均值: {mean_computation.compute()}")
3. 并行处理多幅TIFF
from concurrent.futures import ThreadPoolExecutor
import globdef process_tif(file_path):"""处理单个TIFF文件的函数"""with rasterio.open(file_path) as src:data = src.read(1)# 这里添加你的处理逻辑result = np.mean(data)return file_path, result# 获取所有TIFF文件
tif_files = glob.glob('data/*.tif')# 使用线程池并行处理
with ThreadPoolExecutor(max_workers=4) as executor:results = list(executor.map(process_tif, tif_files))# 打印结果
for file_path, mean_val in results:print(f"{file_path}: {mean_val}")
实用技巧和小贴士
-
内存管理:
- 对于大型TIFF,使用
rasterio
的窗口读取或dask
进行分块处理 - 处理完成后及时关闭文件或使用
with
语句
- 对于大型TIFF,使用
-
元数据保存:
with rasterio.open('input.tif') as src:profile = src.profile # 保存元数据# 修改需要的参数 profile.update(dtype=rasterio.float32, count=1)# 写入新文件时使用保存的元数据 with rasterio.open('output.tif', 'w', **profile) as dst:dst.write(processed_data, 1)
-
TIFF格式选择:
- 地理空间数据通常使用GeoTIFF
- 科学数据可能使用BigTIFF (用于>4GB的文件)
- 考虑压缩选项 (DEFLATE, LZW等) 以减小文件大小
-
性能优化:
- 使用
rasterio
的window
参数只读取需要的区域 - 对于多波段操作,考虑使用
rasterio
的内存映射功能 - 对于重复读取相同文件,考虑转换为更适合分析的格式 (如Zarr)
- 使用