太阳位置计算c程序
根据Paul Schlyter, Stockholm, Sweden中等精度的行星位置计算方法,实现了其中太阳位置的计算方法。/////////////////////////////////////////////////////////////////////////////////////////////////////
// 名称:太阳位置计算程序
// 作者:胡铂(http://hubble.lamost.org)
// 日期 2004-09-30
// 说明:根据Paul Schlyter, Stockholm, Sweden的中等精度计算方法实现
////////////////////////////////////////////////////////////////////////////////////////////////////
#include "stdio.h"
#include "math.h"
#define Pi 3.14159265
#define DE Pi/180
/////////////////////////////////////////////////////////////////
float jde(int Y,int M,int D,int hour,int min,int sec)
{
int f,g;
double mid1,mid2,J,JDE,A;
if(M>=3)
{
f=Y;
g=M;
}
if(M==1||M==2)
{
f=Y-1;
g=M+12;
};
mid1=floor(365.25*f);
mid2=floor(30.6001*(g+1));
A=2-floor(f/100)+floor(f/400);
J=mid1+mid2+D+A+1720994.5;
JDE=J+hour/24+min/1440+sec/86400;
return JDE;
}
/////////////////////////////////////变量声明///////////////////////////////
void main()
{ int i,year,month,day,hour,min,sec;
double d,w,a,e,M,oblecl,L,E,xe,ye,r,v,lon,x,y,z,xequt,yequt,zequt,dist,RA,Decl;
// scanf("%d,%d,%d,%d,%d,%d",&year,&month,&day,&hour,&min,&sec);
// printf("\n");
///////////////////////////////测试数据////////////////////////////////////////////////////
year=1990;
month=4;
day=19;
hour=0;
min=0;
sec=0;
////////////////////////////////////////////轨道根数//////////////////////////////////////
d=jde(year,month,day,hour,min,sec)-2451543.5;//相对儒略日;
w=282.9404+4.70935*0.00001*d;//升交点经度
a=1;
e=0.016709-1.151*0.000000001*d;//偏心率
M = 356.0470 + 0.9856002585 * d;//平近点角
oblecl = 23.4393-3.563*0.0000001 * d;//黄赤交角
////////////////////////////////////////////////////////////////////////////////////////////
L=w+M;//太阳的平均精度
L=fmod(L,360);
E=M + (180/Pi) * e * sin(M*DE) * (1 + e * cos(M*DE));//(开普勒方程的近似解)
E=fmod(E,360);
xe=cos(E*DE) - e;//椭圆轨道上的直角坐标x
ye=sin(E*DE) * sqrt(1 - e*e);//椭圆轨道上的直角坐标y
r=sqrt(xe*xe+ye*ye);//距离
v=atan2(ye,xe)*180/Pi;//真近点角
//lon=fmod((v+w),360);//太阳的精度
//太阳的黄道直角坐标
lon=v+w;
x=r*cos(lon*Pi/180);
y=r*sin(lon*Pi/180);
z=0;
//太阳的赤道直角坐标
xequt=x;
yequt=y*cos(oblecl*Pi/180);
zequt = y* sin(oblecl*Pi/180);
//日地距离、赤经赤纬
dist=sqrt(xequt*xequt+yequt*yequt);
RA=atan2(yequt,xequt)*180/Pi;
RA=fmod(RA,360);
Decl=asin(zequt/r)*180/Pi;
Decl=fmod(Decl,360);
////////////////////////////////////////////////////////////////////
printf("d=%f\n",d);
printf("w=%f\n",w);
printf("e=%f\n",e);
printf("M=%f\n",M);
printf("L=%f\n",L);
printf("E=%f\n",E);
printf("oblecl=%f",oblecl);
printf("xe=%f,ye=%f\n",xe,ye);
printf("r=%f,v=%f\n",r,v);
printf("lon=%f\n",lon);
printf("x=%f,y=%f\n",x,y);
printf("xequt=%f,yequt=%f,zequt=%f\n",xequt,yequt,zequt);
printf("dist=%f\n",dist);
printf("RA=%f,Decl=%f\n",RA,Decl);
} 请问hubble, 你是否打算发一个系列的计算程序?
如果那样,建议你把这些计算程序都跟贴在一个主贴里。方便大家检索查询:)
好主意
我确实打算写一个系列,但是不一定什么时候完成,下次发到一个帖子里面。另外已经发的是否可以挪到一个帖子里面去? 嗯,具体操作是要想想.我建议你可以以这个贴为主贴,请大家不要跟贴.记下这主贴连接的地址.以后你再有新程序时,你可以新发一贴,在贴中注明这个主贴的连接.然后再在主贴中加入新贴的连接.
这样,网友不论读到这个系列的哪个贴子,都可以很容易找到主贴,进而找到其他分贴.
抱歉,这个方法麻烦了些,但即保证了论坛的特点(网友可以自由在这个系列的分贴里探讨),又容易检索.
主贴我给你加精华 我已经把发的帖子都贴到了那个月球计算的帖子后面了。 这个程序其实不难,有计算的公式就差不多了。
不过开发成图形界面的回更好:) 的确,强烈要求图形界面。 能不能,贴个具体计算,月亮的合朔,以及24节气的计算机的例子,也好叫我们这些门外汉,轻松入门。 太阳位置计算程序(2006-12-24 15:07:26)
转载 分类:格物偶得
/////////////////////////////////////////////////////////////////////////////////////////////////////
// 名称:太阳位置计算程序
// 作者:胡铂(http://leohubble.go.3322.org)
// 日期 2004-09-30
// 说明:根据Paul Schlyter, Stockholm, Sweden的中等精度计算方法实现
////////////////////////////////////////////////////////////////////////////////////////////////////
#include "stdio.h"
#include "math.h"
#define Pi 3.14159265
#define DE Pi/180
/////////////////////////////////////////////////////////////////
float jde(int Y,int M,int D,int hour,int min,int sec)
{
int f,g;
double mid1,mid2,J,JDE,A;
if(M>=3)
{
f=Y;
g=M;
}
if(M==1||M==2)
{
f=Y-1;
g=M+12;
};
mid1=floor(365.25*f);
mid2=floor(30.6001*(g+1));
A=2-floor(f/100)+floor(f/400);
J=mid1+mid2+D+A+1720994.5;
JDE=J+hour/24+min/1440+sec/86400;
return JDE;
}
/////////////////////////////////////变量声明///////////////////////////////
void main()
{ int i,year,month,day,hour,min,sec;
double d,w,a,e,M,oblecl,L,E,xe,ye,r,v,lon,x,y,z,xequt,yequt,zequt,dist,RA,Decl;
// scanf("%d,%d,%d,%d,%d,%d",&year,&month,&day,&hour,&min,&sec);
// printf("\n");
///////////////////////////////测试数据////////////////////////////////////////////////////
year=1990;
month=4;
day=19;
hour=0;
min=0;
sec=0;
////////////////////////////////////////////轨道根数//////////////////////////////////////
d=jde(year,month,day,hour,min,sec)-2451543.5;//相对儒略日;
w=282.9404+4.70935*0.00001*d;//升交点经度
a=1;
e=0.016709-1.151*0.000000001*d;//偏心率
M = 356.0470 + 0.9856002585 * d;//平近点角
oblecl = 23.4393-3.563*0.0000001 * d;//黄赤交角
////////////////////////////////////////////////////////////////////////////////////////////
L=w+M;//太阳的平均精度
L=fmod(L,360);
E=M + (180/Pi) * e * sin(M*DE) * (1 + e * cos(M*DE));//(开普勒方程的近似解)
E=fmod(E,360);
xe=cos(E*DE) - e;//椭圆轨道上的直角坐标x
ye=sin(E*DE) * sqrt(1 - e*e);//椭圆轨道上的直角坐标y
r=sqrt(xe*xe+ye*ye);//距离
v=atan2(ye,xe)*180/Pi;//真近点角
//lon=fmod((v+w),360);//太阳的精度
//太阳的黄道直角坐标
lon=v+w;
x=r*cos(lon*Pi/180);
y=r*sin(lon*Pi/180);
z=0;
//太阳的赤道直角坐标
xequt=x;
yequt=y*cos(oblecl*Pi/180);
zequt = y* sin(oblecl*Pi/180);
//日地距离、赤经赤纬
dist=sqrt(xequt*xequt+yequt*yequt);
RA=atan2(yequt,xequt)*180/Pi;
RA=fmod(RA,360);
Decl=asin(zequt/r)*180/Pi;
Decl=fmod(Decl,360);
////////////////////////////////////////////////////////////////////
printf("d=%f\n",d);
printf("w=%f\n",w);
printf("e=%f\n",e);
printf("M=%f\n",M);
printf("L=%f\n",L);
printf("E=%f\n",E);
printf("oblecl=%f",oblecl);
printf("xe=%f,ye=%f\n",xe,ye);
printf("r=%f,v=%f\n",r,v);
printf("lon=%f\n",lon);
printf("x=%f,y=%f\n",x,y);
printf("xequt=%f,yequt=%f,zequt=%f\n",xequt,yequt,zequt);
printf("dist=%f\n",dist);
printf("RA=%f,Decl=%f\n",RA,Decl);
}
请教一下 我想得到wgs-84的坐标 怎么转换》
页:
[1]