pip install rasterio

conda install -c conda-forge rasterio

import rasterio

with rasterio.open('sentinel2.tif') as src:
    # 獲取數據矩陣（NumPy數組）
    data = src.read()
    
    # 查看元數據
    print(f"數據集信息：\n{src.meta}")
    print(f"坐標系：{src.crs}")
    print(f"影像尺寸：{src.width}x{src.height}")
    print(f"空間分辨率：{src.res}")
    print(f"地理范圍：{src.bounds}")

# 創(chuàng)建新柵格文件
profile = src.profile
profile.update(
    dtype=rasterio.float32,
    count=3,
    compress='lzw'
)

with rasterio.open('output.tif', 'w', **profile) as dst:
    dst.write(data.astype(rasterio.float32))

from rasterio.mask import mask
import geojson

# 按像素范圍裁剪
window = Window(100, 200, 800, 600)
clip_data = src.read(window=window)

# 按地理范圍裁剪，GeoJSON格式
geometry = {
    "type": "Polygon",
    "coordinates": [
        [
            [116.0, 39.5],
            [116.5, 39.5],
            [116.5, 40.0],
            [116.0, 40.0],
            [116.0, 39.5]
        ]
    ]
}

clip_data, clip_transform = rasterio.mask.mask(src, [geometry], crop=True)

from rasterio.enums import Resampling

# 選擇重采樣方法：最近鄰、雙線性、立方卷積等
resampled_data = src.read(
    out_shape=(src.height//2, src.width//2),
    resampling=Resampling.bilinear  # 雙線性插值法
)

# 地理坐標 ? 像素坐標
row, col = src.index(116.3974, 39.9093)  # 天安門坐標轉像素位置
lon, lat = src.xy(500, 500)  # 中心點坐標

# 坐標系統(tǒng)轉換
from rasterio.warp import transform
dst_crs = 'EPSG:3857'  # Web墨卡托
transformed = transform(src.crs, dst_crs, [lon], [lat])

# 波段合并
with rasterio.open('RGB.tif', 'w', **profile) as dst:
    for i in range(3):
        dst.write(data[i], i+1)

# 計算NDVI
red = src.read(3).astype(float)
nir = src.read(4).astype(float)
ndvi = (nir - red) / (nir + red + 1e-8)

# 可視化NDVI
import matplotlib.pyplot as plt

plt.figure(figsize=(12, 8))
plt.imshow(ndvi, cmap='RdYlGn', vmin=-1, vmax=1)
plt.colorbar(label='NDVI')
plt.title('植被指數分布', fontsize=16)
plt.axis('off')
plt.savefig('ndvi_map.png', dpi=300, bbox_inches='tight')

# 示例：批量處理Landsat8數據提取水體信息

# 1. 讀取熱紅外波段
with rasterio.open('landsat8_b10.tif') as src:
    thermal_band = src.read(1)

# 2. 計算地表溫度（簡單的示例，實際需要進行輻射校正）
temperature = (thermal_band - 273.15)  # 假設熱紅外波段已經是開爾文溫度

# 3. 應用閾值分割提取水體（基于溫度或植被指數）
water_mask = temperature < 20  # 假設溫度低于20°C為水體

# 4. 輸出二值化結果
with rasterio.open('water_mask.tif', 'w', **src.profile) as dst:
    dst.write(water_mask.astype(rasterio.uint8), 1)

# 5. 生成統(tǒng)計報告
import numpy as np
water_area = np.sum(water_mask) * src.res[0] * src.res[1]  # 計算水體面積
print(f'水體面積: {water_area} 平方米')

from rasterio import Env

with Env(GDAL_NUM_THREADS=4):
    with rasterio.open('large_image.tif') as src:
        data = src.read()