MATLAB 2016计算NDVI

浩瀚地学

已于 2024-04-26 12:59:51 修改

阅读量751

点赞数 12

分类专栏： MATLAB 文章标签： matlab 人工智能经验分享笔记学习 arcgis

于 2024-04-23 08:00:00 首次发布

本文为博主原创文章，未经博主允许禁止转载。如有问题，欢迎指正。

本文链接：https://blog.youkuaiyun.com/zbh13859825167/article/details/138098904

版权

MATLAB 专栏收录该内容

1 篇文章

订阅专栏

MATLAB 2016计算NDVI

之前大二的一段MATLAB代码，突然找到，记录一下当时初学MATLAB的程序，用于读取TIFF计算归一化植被指数（NDVI）并将其保存为TIFF文件。

读取波段数据

[B4,R] = geotiffread('G:\data-2022\result_2_Red.tif');
[B5,R] = geotiffread('G:\data-2022\result_4_NIR.tif');

使用geotiffread函数分别读取两个GeoTIFF文件：result_2_Red.tif（红光波段，对应B4）和result_4_NIR.tif（近红外波段，对应B5）。每个波段数据被存储为一个二维数组（矩阵），有相同的地理参考信息R，包含图像的坐标系统、像素大小等空间属性。

初始化NDVI矩阵和转换数据类型

NDVI = zeros(10715,13980,'single');
B4 = single(B4);
B5 = single(B5);

根据波段数据的维度（此处为10715行×13980列，当时提前看属性的，这里改为size函数更通用），用zeros函数创建一个同样大小的零矩阵NDVI，用于存储计算出的NDVI值，并指定数据类型为'single'（单精度浮点数）。接着，将读取到的波段数据B4和B5转换为单精度浮点数类型。

计算NDVI

for i = 1:10715
    for j =  1:13980
        if B4(i,j) == 0 || B5(i,j) ==0
            continue;
        end
        NDVI(i,j) = (B5(i,j)-B4(i,j))/(B5(i,j)+B4(i,j));
    end
end

使用两层循环遍历波段数据的每一个像素位置（i为行索引，j为列索引）。在计算NDVI时，检查B4和B5是否为零，避免除以零。

获取参考图像信息并写入GeoTIFF文件

info = geotiffinfo('G:\data-2022\result_2_Red.tif');
geotiffwrite('G:\data-2022\result_NDVI.tif', NDVI, R, 'GeoKeyDirectoryTag', info.GeoTIFFTags.GeoKeyDirectoryTag);

使用geotiffinfo函数读取result_2_Red.tif的元数据信息，并将其存储在结构体变量info中。然后，使用geotiffwrite函数将计算得到的NDVI矩阵写入新的GeoTIFF文件。

完整代码如下

[B4,R] = geotiffread('G:\data-2022\result_2_Red.tif');
[B5,R] = geotiffread('G:\data-2022\result_4_NIR.tif');
NDVI = zeros(10715,13980,'single');
B4 = single(B4);
B5 = single(B5);
for i = 1:10715
    for j =  1:13980
        if B4(i,j) == 0 || B5(i,j) ==0
            continue;
        end
        NDVI(i,j) = (B5(i,j)-B4(i,j))/(B5(i,j)+B4(i,j));
    end
end
info = geotiffinfo('G:\data-2022\result_2_Red.tif');
geotiffwrite('G:\data-2022\result_NDVI.tif',NDVI,R,'GeoKeyDirectoryTag',info.GeoTIFFTags.GeoKeyDirectoryTag);