关于AltAz转换成RaDec的问题
下面是MatLab的AltAz转换成RaDec代码（转自MatLab官方用户交流），我用skymap验证了一下，Dec正确，Ra有点问题，可能是我对赤经、赤纬以及坐标系概念不清楚，小弟M币有限，只能奖励1个，求求各位高手能指点下，谢谢！function = AltAz2RaDec(Altitude,Azimuth,Latitude,Longitude,time)
% By: Christopher Wilcox and Ty Martinez
% Jan 22 2010
% Naval Research Laboratory
%
% Description:Convert Altitude/Azimuth angles in degrees to Right
% Ascension/Declination in degrees for an Alt/Az telescope
% mount.
%
% Input: Altitude - Telescope Altitude in degrees
% Azimuth - Telescope Azimuth in degrees
% Latidute - Observer's Latitude (Negative for South) in degrees
% Longitude - Observer's Longiture (Negative for West) in degrees
% (optional) time - Date vector, as returned from 'clock.m', if
% not supplied, the current date/time is used.
% Output: RA - Right Ascension in degrees
% Dec - Declination in degrees
if nargin == 4
time = clock;
year = time(1);
month = time(2);
day = time(3);
hour = time(4);
min = time(5);
sec = time(6);
else
year = time(1);
month = time(2);
day = time(3);
hour = time(4);
min = time(5);
sec = time(6);
end
JD = floor(365.25*(year + 4716.0)) + floor(30.6001*( month + 1.0)) + 2.0 - ...
floor(year/100.0) + floor(floor(year/100.0 )/4.0) + day - 1524.5 + ...
(hour + min/60 + sec/3600)/24;
T = (JD - 2451545)/36525;
ThetaGMST = 67310.54841 + (876600*3600 + 8640184.812866)*T + .093104*(T^2) - (6.2*10^-6)*(T^3);
ThetaGMST = mod((mod(ThetaGMST,86400*(ThetaGMST/abs(ThetaGMST)))/240),360);
ThetaLST = ThetaGMST + Longitude;
ThetaLST = mod(ThetaLST,360);
Dec = asind(sind(Altitude)*sind(Latitude) + cosd(Altitude)*cosd(Latitude)*cosd(Azimuth));
cos_RA = (sind(Altitude) - sind(Dec)*sind(Latitude))/(cosd(Dec)*cosd(Latitude));
RA = acosd(cos_RA);
if sind(Azimuth) > 0
RA = 360 - RA;
end
RA = ThetaLST - RA;
if RA < 0
RA = RA + 360;
end
if Dec >= 0
Dec = abs(Dec);
else
Dec = -abs(Dec);
end
本帖最后由 goldenham 于 2012-5-4 18:06 编辑
可能是计算机时间的问题，你的电脑上的时间是本地时间吧，那么由于clock函数使用本地时间计算ThetaLST，之后就得把多算的减去。将其中这一行代码
ThetaLST = ThetaGMST + Longitude;
改成
ThetaLST = ThetaGMST + Longitude - 120;%GMT+8
就行了。
goldenham 发表于 2012-5-4 18:02 static/image/common/back.gif
可能是计算机时间的问题，你的电脑上的时间是本地时间吧，那么由于clock函数使用本地时间计算ThetaLST，之 ...
问题已经解决了！
谢谢你，我查看matlab里有计算儒略时间的函数，同时还得把输入时间改成UTC时间。已经过skymap验证啦！
下面是matlab里计算儒略时间函数的部分代码，我下载的代码主要就是缺少前面部分，输入时间应为UTC时间。
juliandate：
if ( month <= 2 ) % january & february
year= year - 1.0;
month = month + 12.0;
endjd = floor( 365.25*(year + 4716.0)) + floor( 30.6001*( month + 1.0)) + 2.0 - ...
floor( year/100.0 ) + floor( floor( year/100.0 )/4.0 ) + day - 1524.5 + ...
(hour + min/60 + sec/3600)/24;
页:
[1]