hubble 发表于 2004-10-8 21:24

太阳位置计算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);



 


   
 
 

}

nngs 发表于 2004-10-9 01:45

请问hubble, 你是否打算发一个系列的计算程序?
如果那样,建议你把这些计算程序都跟贴在一个主贴里。方便大家检索查询:)

hubble 发表于 2004-10-9 09:37

好主意

我确实打算写一个系列,但是不一定什么时候完成,下次发到一个帖子里面。另外已经发的是否可以挪到一个帖子里面去?

nngs 发表于 2004-10-9 10:15

嗯,具体操作是要想想.

我建议你可以以这个贴为主贴,请大家不要跟贴.记下这主贴连接的地址.以后你再有新程序时,你可以新发一贴,在贴中注明这个主贴的连接.然后再在主贴中加入新贴的连接.

这样,网友不论读到这个系列的哪个贴子,都可以很容易找到主贴,进而找到其他分贴.

抱歉,这个方法麻烦了些,但即保证了论坛的特点(网友可以自由在这个系列的分贴里探讨),又容易检索.

主贴我给你加精华

hubble 发表于 2004-10-9 16:48

我已经把发的帖子都贴到了那个月球计算的帖子后面了。

福建ltycho 发表于 2004-11-1 19:48

这个程序其实不难,有计算的公式就差不多了。
不过开发成图形界面的回更好:)

macross11 发表于 2004-11-1 23:30

的确,强烈要求图形界面。

heihoo 发表于 2004-11-4 16:18

能不能,贴个具体计算,月亮的合朔,以及24节气的计算机的例子,也好叫我们这些门外汉,轻松入门。

霍曼 发表于 2010-7-31 01:14

太阳位置计算程序(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]
查看完整版本: 太阳位置计算c程序