0
点赞
收藏
分享

微信扫一扫

遥感技术入门与实践:从数据获取到图像分析

遥感技术作为我们提供了从空中俯瞰地球的能力,它通过卫星、飞机等平台获取地表信息,广泛应用于农业监测、环境评估、城市规划等领域。本文将带你了解遥感的基本概念,并通过 Python 代码实践遥感图像处理的核心技术。

什么是遥感?

遥感是指非接触式地获取物体或区域信息的技术。根据平台不同,可分为:

  • 卫星遥感:如 Landsat、Sentinel 系列卫星
  • 航空遥感:飞机、无人机搭载传感器
  • 地面遥感:高塔、三脚架等固定平台

遥感数据通常以数字图像形式存在,包含多个波段(如可见光、红外等),这使得我们能提取比普通照片更多的信息。

Python 遥感处理工具

处理遥感数据的常用 Python 库:

  • rasterio:读取和写入遥感图像
  • numpy:数值计算基础
  • matplotlib:图像可视化
  • scikit-image:图像处理算法
  • earthpy:简化遥感数据分析流程

实践:遥感图像处理基础

下面我们通过代码演示遥感图像的基本操作,包括读取数据、波段组合、计算植被指数等。

import rasterio
import numpy as np
import matplotlib.pyplot as plt
from rasterio.plot import show
import os

# 设置中文显示
plt.rcParams["font.family"] = ["SimHei", "WenQuanYi Micro Hei", "Heiti TC"]
plt.rcParams['axes.unicode_minus'] = False  # 解决负号显示问题

def read_remote_sensing_image(file_path):
    """读取遥感图像文件"""
    try:
        with rasterio.open(file_path) as src:
            # 打印图像信息
            print(f"图像文件: {file_path}")
            print(f"图像宽度: {src.width} 像素")
            print(f"图像高度: {src.height} 像素")
            print(f"波段数量: {src.count}")
            print(f"空间参考: {src.crs}")
            print(f"边界范围: {src.bounds}")
            
            # 读取所有波段
            img = src.read()
            # 转换形状为 (height, width, bands)
            img = np.transpose(img, (1, 2, 0))
            
            return img, src
    except Exception as e:
        print(f"读取图像出错: {e}")
        return None, None

def display_rgb_image(img, src, bands=[3,2,1], title="RGB合成图像"):
    """显示RGB合成图像"""
    # 确保波段索引从0开始
    bands = [b-1 for b in bands]
    
    # 提取RGB波段
    rgb_img = img[:, :, bands]
    
    # 归一化到0-1范围以便正确显示
    rgb_normalized = (rgb_img - np.min(rgb_img)) / (np.max(rgb_img) - np.min(rgb_img))
    
    plt.figure(figsize=(10, 8))
    plt.imshow(rgb_normalized)
    plt.title(title)
    plt.axis('off')
    plt.show()
    
    return rgb_normalized

def calculate_ndvi(img, nir_band=4, red_band=3):
    """计算归一化植被指数NDVI"""
    # 确保波段索引从0开始
    nir_idx = nir_band - 1
    red_idx = red_band - 1
    
    # 提取近红外和红波段
    nir = img[:, :, nir_idx].astype(float)
    red = img[:, :, red_idx].astype(float)
    
    # 计算NDVI,避免除零错误
    ndvi = np.where(
        (nir + red) == 0, 
        0, 
        (nir - red) / (nir + red)
    )
    
    # 显示NDVI图像
    plt.figure(figsize=(10, 8))
    plt.imshow(ndvi, cmap='Greens')
    plt.colorbar(label='NDVI值')
    plt.title('归一化植被指数(NDVI)')
    plt.axis('off')
    plt.show()
    
    return ndvi

def calculate_ndwi(img, green_band=2, nir_band=4):
    """计算归一化水体指数NDWI"""
    # 确保波段索引从0开始
    green_idx = green_band - 1
    nir_idx = nir_band - 1
    
    # 提取绿波段和近红外波段
    green = img[:, :, green_idx].astype(float)
    nir = img[:, :, nir_idx].astype(float)
    
    # 计算NDWI,避免除零错误
    ndwi = np.where(
        (green + nir) == 0, 
        0, 
        (green - nir) / (green + nir)
    )
    
    # 显示NDWI图像
    plt.figure(figsize=(10, 8))
    plt.imshow(ndwi, cmap='Blues')
    plt.colorbar(label='NDWI值')
    plt.title('归一化水体指数(NDWI)')
    plt.axis('off')
    plt.show()
    
    return ndwi

def classify_land_cover(ndvi, ndwi):
    """基于NDVI和NDWI进行简单土地覆盖分类"""
    # 创建分类图像
    classification = np.zeros_like(ndvi, dtype=np.uint8)
    
    # 水体:NDWI > 0.1
    classification[ndwi > 0.1] = 1
    
    # 植被:NDVI > 0.3 且 不是水体
    vegetation_mask = (ndvi > 0.3) & (classification != 1)
    classification[vegetation_mask] = 2
    
    # 裸地:NDVI <= 0.3 且 不是水体
    bare_ground_mask = (ndvi <= 0.3) & (classification != 1)
    classification[bare_ground_mask] = 3
    
    # 显示分类结果
    plt.figure(figsize=(10, 8))
    cmap = plt.cm.get_cmap('Set3', 4)
    im = plt.imshow(classification, cmap=cmap)
    plt.colorbar(im, ticks=[0, 1, 2, 3])
    plt.clim(-0.5, 3.5)
    plt.title('土地覆盖分类')
    plt.axis('off')
    
    # 添加图例
    from matplotlib.patches import Patch
    legend_elements = [
        Patch(facecolor=cmap(0), label='其他'),
        Patch(facecolor=cmap(1), label='水体'),
        Patch(facecolor=cmap(2), label='植被'),
        Patch(facecolor=cmap(3), label='裸地')
    ]
    plt.legend(handles=legend_elements, loc='lower right')
    
    plt.show()
    
    return classification

def main():
    # 这里使用Landsat 8样例数据,实际使用时替换为你的图像路径
    # 可以从USGS Earth Explorer下载免费的Landsat数据
    image_path = "LC08_L1TP_123032_20200515_20200527_01_T1_B4.TIF"
    
    # 检查文件是否存在,如果不存在则使用多个波段文件
    if not os.path.exists(image_path):
        print("找不到单个图像文件,尝试使用多波段文件模式")
        # 实际应用中,这里应该是你的多波段文件列表
        return
    
    # 读取遥感图像
    img, src = read_remote_sensing_image(image_path)
    if img is None:
        return
    
    # 显示真彩色合成图像 (红、绿、蓝波段)
    rgb_img = display_rgb_image(img, src, bands=[4,3,2], title="真彩色合成图像")
    
    # 计算并显示NDVI
    ndvi = calculate_ndvi(img)
    
    # 计算并显示NDWI
    ndwi = calculate_ndwi(img)
    
    # 进行简单的土地覆盖分类
    classification = classify_land_cover(ndvi, ndwi)
    
    print("遥感图像处理完成!")

if __name__ == "__main__":
    main()

代码解析

上面的代码实现了遥感图像处理的几个核心功能:

  1. 图像读取:使用rasterio库读取遥感图像文件,获取图像的基本信息(尺寸、波段数、空间参考等)。
  2. RGB 合成:遥感图像通常包含多个波段,通过选择合适的波段组合(如真彩色使用红、绿、蓝波段),可以生成人眼可见的图像。
  3. 植被指数计算:NDVI(归一化植被指数)是最常用的遥感指数之一,通过近红外波段和红波段计算得出,可用于评估植被生长状况。健康植被在近红外波段反射强烈,在红波段吸收强烈,因此 NDVI 值较高。
  4. 水体指数计算:NDWI(归一化水体指数)利用绿波段和近红外波段计算,可有效识别水体区域。
  5. 土地覆盖分类:基于 NDVI 和 NDWI 的阈值,进行简单的土地覆盖类型划分(水体、植被、裸地等)。

获取遥感数据

要运行上述代码,你需要实际的遥感数据。以下是一些获取免费遥感数据的网站:

  • USGS Earth Explorer:提供 Landsat 系列、Sentinel-2 等多种卫星数据
  • Sentinel Hub:专注于欧空局 Sentinel 系列卫星数据
  • NASA Earthdata:包含多种 NASA 卫星数据
  • 地理空间数据云:国内的遥感数据共享平台

遥感技术的应用场景

  1. 农业监测:通过 NDVI 监测作物生长状况,预测产量
  2. 环境监测:森林砍伐、湿地变化、荒漠化监测
  3. 灾害响应:洪水、火灾范围评估,灾后损失评估
  4. 城市规划:城市扩张监测,土地利用变化分析
  5. 水资源管理:水体监测,干旱评估

进阶学习方向

如果对遥感技术感兴趣,可以进一步学习:

  • 高级图像分类算法(如随机森林、深度学习)
  • 时间序列分析(监测长期变化趋势)
  • 三维遥感(激光雷达数据处理)
  • 定量遥感(反演地表参数)

遥感技术正不断发展,随着高分辨率卫星和无人机技术的普及,遥感数据的获取越来越容易,应用也越来越广泛。希望本文能为你打开遥感世界的大门!

举报

相关推荐

0 条评论