將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
最后結果：
