function dood=doodsonArguments(mjd)

% compute gmst
Tu0   = (floor(mjd)-51544.5)/36525.0;
gmst0 = (6.0/24 + 41.0/(24*60) + 50.54841/(24*60*60));
gmst0 = gmst0 + (8640184.812866/(24*60*60))*Tu0;
gmst0 = gmst0 + (0.093104/(24*60*60))*Tu0*Tu0;
gmst0 = gmst0 + (-6.2e-6/(24*60*60))*Tu0*Tu0*Tu0;
r     = 1.002737909350795 + 5.9006e-11*Tu0 - 5.9e-15*Tu0*Tu0;
gmst  = mod(2*pi*(gmst0 + r * mod(mjd,1)), 2*pi);

t = (mjd-51544.5)/365250.0;

dood = zeros(6,1);
dood(2,1) = (218.31664562999 + (4812678.81195750 + (-0.14663889 + ( 0.00185140 + -0.00015355*t)*t)*t)*t) * pi/180;
dood(3,1) = (280.46645016002 + ( 360007.69748806 + ( 0.03032222 + ( 0.00002000 + -0.00006532*t)*t)*t)*t) * pi/180;
dood(4,1) = ( 83.35324311998 + (  40690.13635250 + (-1.03217222 + (-0.01249168 +  0.00052655*t)*t)*t)*t) * pi/180;
dood(5,1) = (234.95544499000 + (  19341.36261972 + (-0.20756111 + (-0.00213942 +  0.00016501*t)*t)*t)*t) * pi/180;
dood(6,1) = (282.93734098001 + (     17.19457667 + ( 0.04568889 + (-0.00001776 + -0.00003323*t)*t)*t)*t) * pi/180;
dood(1,1) = gmst + pi - dood(2,1);
