% ’n‹…Ž©“]‘¬“x•Ï“®‹ßŽ—�¸“xŠm”F 2004/06/04
function testgoeddist

dirs='e:\gps\erp\iers';

% station position lat,lon,h
p=[
  0,  0, 0
  0, 90, 0
 60, 90, 0
-60,180, 0
 80,270, 0
];
tu=caltomjd([2004,10,1,6,0,0]); % time
dt=0.1; % delta t

% ecef to eci transformation matrix
U0=ecsftoecef(tu-dt/86400,readerp(tu-dt/86400,dirs),-32)';   % strict model
U1=ecsftoecef(tu,readerp(tu,dirs),-32)'*Rz(dt*0.7292115E-4); % approx model

disp(sprintf('%04d/%02d/%02d %02d:%02d:%02d dt=%.2fsec',mjdtocal(tu),dt));
for n=1:size(p,1)
    pos=geodtoecef(p(n,:));
    dp=U0*pos-U1*pos;
    disp(sprintf('lat=%5.1fdeg,lon=%5.1fdeg,h=%.0fm : dp=%8.5fm,%8.5fm,%8.5fm (%7.5fm)',...
                 p(n,:),dp,norm(dp)));
end

function R=Rz(t), R=[cos(t),sin(t),0;-sin(t),cos(t),0;0,0,1];
