matlab计算岁差章动极移和地球自转矩阵

利用matlab计算岁差章动极移和地球自转矩阵，程序已验证

%%---------极移矩阵计算---%
function polarmotion=EP（MJD）
%这里的MJD是UTC的简约儒略日
%t为自J2000.0算起至t的儒略世纪数
% 206264.806247096 = 3600*DEG   用于秒化弧度

%MJD是约简的儒略世纪数
%MJD=JD-2400000.5
%t=（MJD-51544.0333）/36525.0;
%eop---地球定向文件 在LERS可以查到  注：尽量将文件补充完整
eop=importdata（‘地球定向参数_EOP.txt‘）;
%%
if MJD>eop（length（eop）1）||MJD% 时间超出表中时间范围
    xp = 0.0;
    yp = 0.0;
else
% 地极坐标 插值得到UTC时的XpYp
    xp = interp1（eop（:1）eop（:2）MJD）/3600*pi/180;
    yp = interp1（eop（:1）eop（:3）MJD）/3600*pi/180;
end
%% 极移矩阵
polarmotion=rotationy（-xp）*rotationx（-yp）;

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----

     文件        650  2015-03-18 20:17  EP.m

     文件       1407  2015-03-18 15:46  ER.m

     文件       1274  2015-03-18 15:40  NR.m

     文件        526  2015-03-17 09:56  PR.m

     文件      93131  2015-03-18 09:04  地球定向参数_EOP.txt

     文件       2317  2015-03-17 14:00  黄经和倾角章动序列表.txt

----------- ---------  ---------- -----  ----

                99305                    6