ICGEM官网下载 COST-G发布的40×40的球谐系数
close all;
clearvars -except;
% addpath('E:\Code\Tool\Function\GRACE_functions');
dir_degree_1 = 'E:\Code\GRACE_data\Degree_1\deg1_coef.txt';
dir_c20 = 'E:\Code\GRACE_data\Degree_2\C20_RL06.txt';
myDir_Swarm = 'E:\Data\Swarm\EGF\'; % COST-G
myFiles_gfc_Swarm = dir(fullfile(myDir_Swarm,'*.gfc'));
num_file = length(myFiles_gfc_Swarm);
lmax=40;name='Swarm';
cs= zeros(num_file,lmax+1,lmax+1);
cs_sgi= zeros(num_file,lmax+1,lmax+1);
cs_res= zeros(num_file,lmax+1,lmax+1);
cs_mss= zeros(num_file,lmax+1,lmax+1);
hwait=waitbar(0,'Waiting>>>>>>>>'); %加载等待对话框
for ii=1:num_file
str=['Processing...',num2str(ii),'/',num2str(num_file),' '];
hwait=waitbar(ii/num_file,hwait,str,'Name','SSM');
pathname=strcat(myDir_Swarm,myFiles_gfc_Swarm(ii,1).name);
[cs(ii,:,:),cs_sgi(ii,:,:),int_year(ii),int_month(ii),time(ii)] =gmt_readgfc(pathname);
end
cs_replace=cs;
[cs_replace] = gmt_replace_degree_1(dir_degree_1,cs_replace,int_year,int_month,num_file);
[cs_replace] = gmt_replace_C20(dir_c20,cs_replace,int_year,int_month,num_file);
cs_mean = mean(cs_replace(:,:,:));
for i=1:num_file
cs_res(i,:,:) = cs_replace(i,:,:)-cs_mean(1,:,:);%扣除平均值
end
CS_(:,:)=cs_res(1,:,:);
for i=1:num_file
cs_tmp(:,:) = cs_res(i,:,:);
cs_tmp(:,:)=gmt_gc2mc(cs_tmp);
cs_mss(i,:,:) = gmt_destriping(cs_tmp,'NONE');
end
eval(['grid_' name '=gmt_cs2grid(cs_mss,750,1);']);
ii=1;
eval(['grid_1(:,:)=grid_' name '(:,:,ii);']);
gmt_grid2map(grid_1.*100,-30,30,1,0,'cm',[name ' ' num2str(int_year(ii),'%02d') '-' num2str(int_month(ii),'%02d')],20)
未滤波处理结果:
750km高斯滤波结果:
从结果可以看出40阶的Swarm球谐系数反演的噪声太大,即使经过强约束的高斯滤波后效果仍然不理想,不太清楚是Swarm数据本身所致还是代码实现有误?还请有懂的友友帮忙解答下,不胜感激!
点赞支持,一键分享。欢迎评论交流!