MATLAB对TIFF遥感影像批量辐射定标

将HDF格式的MODIS数据批量拼接、投影及格式转换、重采样、裁剪到研究区后，想批量对数据进行辐射定标（把DN值计算为具有物理意义的值），这里以MOD13A2产品为例

inpath='E:\1SMAPDATA\NDVIREF2016_2019\05CATEGORY\02CLASS\NDVI';%NDVI原始数据 NDVI文件夹下都是TIFF格式的NDVI遥感影像
outpath='E:\1SMAPDATA\NDVIREF2016_2019\06SCALE';%输出路径
min=-2000; %最小值
max=10000;%最大值
fill=-3000;%填充值
back=-32768;%背景值 这个需要注意 需要先在MATLAB中先读取一幅影像后确定 有时为65535(LST)
scale=0.0001;%缩放因子
off=0;%偏移量
EqualFillBackSca(inpath,outpath,fill,back,min,max,scale,off);%统一辐射定标函数

% fill back min max scale off  MOD13A2
% REF: -1000 -32678 0 10000 0.0001
% NDVI EVI :-3000 -32678 -2000 10000 0.0001

% min=0; 
% max=10000;
% fill=-1000;
% back=-32768;
% scale=0.0001;
% off=0;

%读取TIFF数据
[A,R]=geotiffread('E:\02CLASS\blue\A2015353.1_km_16_days_blue_reflectance.tif');

统一辐射定标函数：

function [] = EqualFillBackSca(imagepath,outpath,fill,background,min,max,scale,offset)
%EqualFillBackSca  
%imagepath表示含有某一种类数据的文件夹 
%outpath 影像输出路径
%有些影像没有offset值,相应位置就为0
%统一至：背景值 填充值均为nan 这样在ArcGIS中显示时就是正常范围
%背景值可能是最大的LST：65535 也可能是0 建议先读取数据确认
%Fill值一般是产品说明中有的，还是再确定一些比较好
% 要在有效值范围内
fileFolder=fullfile(imagepath);
dirOutput=dir(fullfile(fileFolder,'*.tif')); 
fileNames={dirOutput.name}'; l
for i=1:length(fileNames)
filepath=[imagepath,'\',fileNames{i}];
[A,R]=geotiffread(filepath);%注意这里读取的数据矩阵A是不是double格式,这个格式不能为负数或者是nan
A=double(A);
A(A==fill)=nan;
A(A==background)=nan;
A((A<min))=nan;
A((A>max))=nan;
A=A*scale+offset;
info = geotiffinfo(filepath);
geotiffwrite([outpath,'\',fileNames{i}],A,R,'GeoKeyDirectoryTag', info.GeoTIFFTags.GeoKeyDirectoryTag);
end
end

最后结果：