MODIS MOD13A3 ndvi数据sg批量滤波

[a,R]=geotiffread('D:\graduate\extracttotiff\extr2010152.tif');  %将数据内容放到A里(二维矩阵)，将地理参考放到R里面
info=geotiffinfo('D:\graduate\extracttotiff\extr2010152.tif');  %返回该tif文件的各种信息
[m,n]=size(a);%返回矩阵a的行列数，m是行，n是列
ndvilist= dir('D:\graduate\extracttotiff\*.tif');     % 获取所有tif文件名
qslist=dir('D:\graduate\qctotiff\*.tif');  %获取所有质量控制文件tif格式
qs = length(ndvilist);       % 总期数
ndvisum=zeros(m*n,qs)+NaN; %构建和NDVI同大小的NaN矩阵
qcsum=zeros(m*n,qs)+NaN;

%% 读取ndvi及qs并进行重排列
k = 1;                                  % 初始化为1
for i = 1:length(ndvilist)              %遍历所有文件
    
    filename=strcat(ndvilist(i).folder,'\',ndvilist(i).name); %ndvi影像文件名
    filename_pr=strcat(qslist(i).folder,'\',qslist(i).name); %ndvi质量控制文件文件名
    
    ndvi = double(importdata(filename));  % 获取ndvi数据
    ndvi = reshape(ndvi,m*n,1);           % 转为1列，方便处理
    ndvi(ndvi==-3000)&