计算日落日出时间算法与代码

admin

下面是一种随经纬度变化的日出日落时间计算方法,我成功运用在一智能路灯控制器中,希望对需要的朋友有帮助。   
已知：日出日落时太阳的位置h＝-0.833°,要计算地的地理位置,经度Long,纬度G1at,时区zone,UTo为上次计算的日出日落时间,第一次计算时UTo＝180°。
(1)先计算出从格林威治时间公元2000年1月1日到计算日天数days;
(2)计算从格林威治时间公元2000年1月1日到计算日的世纪数t, 则t＝(days+UTo／360)／36525;
(3)计算太阳的平黄径 ： L=280.460+36000.770×t;
(4)计算太阳的平近点角 ：G＝357.528+35999.050×t
(5)计算太阳的黄道经度 ：λ＝L+1.915×sinG+0.020xsin(2G);
(6)计算地球的倾角     ε＝23.4393-0.0130×t;
(7)计算太阳的偏差     δ＝arcsin(sinε×sinλ);
(8)计算格林威治时间的太阳时间角GHA： GHA=UTo-180-1.915×sinG-0.020×sin(2G)   +2.466×sin(2λ)-0.053×sin(4λ)
(9)计算修正值e： e=arcos{[   sinh-sin(Glat)sin(δ)]/cos(Glat)cos(δ)}
(10)计算新的日出日落时间 ：UT＝UTo-(GHA+Long±e); 其中“+”表示计算日出时间,“-”表示计算日落时间;
(11)比较UTo和UT之差的绝对值,如果大于0.1°即0.007小时,把UT作为新的日出日落时间值,重新从第(2)步开始进行迭代计算,如果UTo和UT之差的绝对值小于0.007小时,则UT即为所求的格林威治日出日落时间;
(12)上面的计算以度为单位,即180°=12小时,因此需要转化为以小时表示的时间,再加上所在的时区数Zone,即要计算地的日出日落时间为 ：T=UT/15+Zone
上面的计算日出日落时间方法适用于小于北纬60°和南纬60°之间的区域,如果计算位置为西半球时,经度Long为负数。
Java代码：

1. public class SunRiseSet {
   
2. private static int[] days_of_month_1={31,28,31,30,31,30,31,31,30,31,30,31};
   
3. private static int[] days_of_month_2={31,29,31,30,31,30,31,31,30,31,30,31};
   
4. private final static double h= -0.833;//日出日落时太阳的位置
   
5. private final static double UTo=180.0;//上次计算的日落日出时间，初始迭代值180.0
   
6.     //输入日期
   
7. //输入经纬度
   
8. //判断是否为闰年：若为闰年，返回1；若不是闰年,返回0
   
9. public static boolean leap_year(int year){
   
10.     if(((year%400==0) || (year%100!=0) && (year%4==0))) return true;
    
11.     else return false;
    
12. }
    
13. //求从格林威治时间公元2000年1月1日到计算日天数days
    
14. public static int days(int year, int month, int date){  
    
15.     int i,a=0;  
    
16.     for(i=2000;i<year;i++){  
    
17.         if(leap_year(i)) a=a+366;  
    
18.         else a=a+365;  
    
19.     }  
    
20.     if(leap_year(year)){  
    
21.         for(i=0;i<month-1;i++){
    
22.             a=a+days_of_month_2[i];  
    
23.         }
    
24.     }else {  
    
25.         for(i=0;i<month-1;i++){
    
26.             a=a+days_of_month_1[i];
    
27.         }  
    
28.     }  
    
29.     a=a+date;  
    
30.     return a;  
    
31. }
    
32. //求格林威治时间公元2000年1月1日到计算日的世纪数t
    
33. public static double t_century(int days,  double UTo){  
    
34.     return ((double)days+UTo/360)/36525;  
    
35. }
    
36. //求太阳的平黄径
    
37. public static double L_sun(double t_century){  
    
38.     return (280.460+36000.770*t_century);  
    
39. }
    
40. //求太阳的平近点角
    
41. public static double G_sun(double t_century){  
    
42.     return (357.528+35999.050*t_century);  
    
43. }
    
44. //求黄道经度
    
45. public static double ecliptic_longitude(double L_sun,double G_sun){  
    
46.     return (L_sun+1.915*Math.sin(G_sun*Math.PI/180)+0.02*Math.sin(2*G_sun*Math.PI/180));
    
47. }
    
48. //求地球倾角
    
49. public static double earth_tilt(double t_century){  
    
50.     return (23.4393-0.0130*t_century);  
    
51. }
    
52. //求太阳偏差
    
53. public static double sun_deviation(double earth_tilt,double ecliptic_longitude){  
    
54.     return (180/Math.PI*Math.asin(Math.sin(Math.PI/180*earth_tilt)*Math.sin(Math.PI/180*ecliptic_longitude)));  
    
55. }
    
56. //求格林威治时间的太阳时间角GHA
    
57. public static double GHA(double UTo,double G_sun,double ecliptic_longitude){  
    
58.     return (UTo-180-1.915*Math.sin(G_sun*Math.PI/180)-0.02*Math.sin(2*G_sun*Math.PI/180)+2.466*Math.sin(2*ecliptic_longitude*Math.PI/180)-0.053*Math.sin(4*ecliptic_longitude*Math.PI/180));
    
59. }
    
60. //求修正值e
    
61. public static double e(double h,double glat,double sun_deviation){  
    
62.     return 180/Math.PI*Math.acos((Math.sin(h*Math.PI/180)-Math.sin(glat*Math.PI/180)*Math.sin(sun_deviation*Math.PI/180))/(Math.cos(glat*Math.PI/180)*Math.cos(sun_deviation*Math.PI/180)));  
    
63. }
    
64. //求日出时间
    
65. public static double UT_rise(double UTo,double GHA,double glong,double e){  
    
66.     return (UTo-(GHA+glong+e));  
    
67. }
    
68. //求日落时间
    
69. public static double UT_set(double UTo,double GHA,double glong,double e){  
    
70.     return (UTo-(GHA+glong-e));
    
71. }
    
72. //判断并返回结果（日出）
    
73. public static double result_rise(double UT,double UTo,double glong, double glat, int year, int month, int date){  
    
74.     double d;  
    
75.     if(UT>=UTo) d=UT-UTo;  
    
76.     else d=UTo-UT;  
    
77.     if(d>=0.1) {  
    
78.         UTo=UT;  
    
79.         UT=UT_rise(UTo,
    
80.             GHA(UTo,G_sun(t_century(days(year,month,date),UTo)),
    
81.                     ecliptic_longitude(L_sun(t_century(days(year,month,date),UTo)),
    
82.                     G_sun(t_century(days(year,month,date),UTo)))),
    
83.         glong,
    
84.         e(h,glat,sun_deviation(earth_tilt(t_century(days(year,month,date),UTo)),
    
85.             ecliptic_longitude(L_sun(t_century(days(year,month,date),UTo)),
    
86.             G_sun(t_century(days(year,month,date),UTo))))));  
    
87.         result_rise(UT,UTo,glong,glat,year,month,date);
    
88. 
    
89.     }  
    
90.     return UT;  
    
91. }
    
92. //判断并返回结果（日落）
    
93. public static double result_set(double UT,double UTo,double glong,double glat, int year, int month, int date){  
    
94.     double d;  
    
95.     if(UT>=UTo) d=UT-UTo;  
    
96.     else d=UTo-UT;  
    
97.     if(d>=0.1){  
    
98.         UTo=UT;  
    
99.         UT=UT_set(UTo,
    
100.           GHA(UTo,G_sun(t_century(days(year,month,date),UTo)),
     
101.           ecliptic_longitude(L_sun(t_century(days(year,month,date),UTo)),
     
102.           G_sun(t_century(days(year,month,date),UTo)))),
     
103.           glong,
     
104.           e(h,glat,sun_deviation(earth_tilt(t_century(days(year,month,date),UTo)),
     
105.           ecliptic_longitude(L_sun(t_century(days(year,month,date),UTo)),
     
106.           G_sun(t_century(days(year,month,date),UTo))))));  
     
107.         result_set(UT,UTo,glong,glat,year,month,date);  
     
108.     }  
     
109.     return UT;  
     
110. }
     
111. //求时区
     
112. public static int Zone(double glong){  
     
113.     if(glong>=0) return (int)((int)(glong/15.0)+1);  
     
114.     else return (int)((int)(glong/15.0)-1);
     
115. }
     
116. //打印结果
     
117. // public static void output(double rise, double set, double glong){  
     
118. //     if((int)(60*(rise/15+Zone(glong)-(int)(rise/15+Zone(glong))))<10)  
     
119. //         System.out.println("The time at which the sunrise is: "+(int)(rise/15+Zone(glong))+":"+(int)(60*(rise/15+Zone(glong)-(int)(rise/15+Zone(glong))))+" .\n");  
     
120. //     else System.out.println("The time at which the sunrise is: "+(int)(rise/15+Zone(glong))+":"+(int)(60*(rise/15+Zone(glong)-(int)(rise/15+Zone(glong))))+" .\n");
     
121. //      
     
122. //     if((int)(60*(set/15+Zone(glong)-(int)(set/15+Zone(glong))))<10)
     
123. //         System.out.println("The time at which the sunset is: "+(int)(set/15+Zone(glong))+": "+(int)(60*(set/15+Zone(glong)-(int)(set/15+Zone(glong))))+" .\n");  
     
124. //     else System.out.println("The time at which the sunset is: "+(int)(set/15+Zone(glong))+":"+(int)(60*(set/15+Zone(glong)-(int)(set/15+Zone(glong))))+" .\n");
     
125. // }
     
126. public static String getSunrise(GeoPoint geoPoint,Time sunTime){
     
127. double sunrise,glong,glat;
     
128. int year,month,date;
     
129. year=sunTime.year;
     
130. month=sunTime.month;
     
131. date=sunTime.monthDay;
     
132. glong=geoPoint.getLongitude();
     
133. glat=geoPoint.getLatitude();
     
134. 
     
135. sunrise=result_rise(UT_rise(UTo,
     
136.                  GHA(UTo,G_sun(t_century(days(year,month,date),UTo)),
     
137.                         ecliptic_longitude(L_sun(t_century(days(year,month,date),UTo)),
     
138.                         G_sun(t_century(days(year,month,date),UTo)))),
     
139.             glong,
     
140.             e(h,glat,sun_deviation(earth_tilt(t_century(days(year,month,date),UTo)),
     
141.             ecliptic_longitude(L_sun(t_century(days(year,month,date),UTo)),
     
142.         G_sun(t_century(days(year,month,date),UTo)))))),UTo,glong,glat,year,month,date);
     
143. //System.out.println("Sunrise is: "+(int)(sunrise/15+Zone(glong))+":"+(int)(60*(sunrise/15+Zone(glong)-(int)(sunrise/15+Zone(glong))))+" .\n");
     
144. Log.d("Sunrise", "Sunrise is: "+(int)(sunrise/15+8)+":"+(int)(60*(sunrise/15+8-(int)(sunrise/15+8)))+" .\n");
     
145.         //return "Sunrise is: "+(int)(sunrise/15+Zone(glong))+":"+(int)(60*(sunrise/15+Zone(glong)-(int)(sunrise/15+Zone(glong))))+" .\n";
     
146.         return "Sunrise is: "+(int)(sunrise/15+8)+":"+(int)(60*(sunrise/15+8-(int)(sunrise/15+8)))+" .\n";
     
147. }
     
148. public static String getSunset(GeoPoint geoPoint,Time sunTime){
     
149. double sunset,glong,glat;
     
150. int year,month,date;
     
151. year=sunTime.year;
     
152. month=sunTime.month;
     
153. date=sunTime.monthDay;
     
154. glong=geoPoint.getLongitude();
     
155. glat=geoPoint.getLatitude();
     
156. sunset=result_set(UT_set(UTo,
     
157.                 GHA(UTo,G_sun(t_century(days(year,month,date),UTo)),
     
158.                   ecliptic_longitude(L_sun(t_century(days(year,month,date),UTo)),
     
159.                G_sun(t_century(days(year,month,date),UTo)))),
     
160.                 glong,
     
161.                 e(h,glat,sun_deviation(earth_tilt(t_century(days(year,month,date),UTo)),
     
162.                   ecliptic_longitude(L_sun(t_century(days(year,month,date),UTo)),
     
163.                G_sun(t_century(days(year,month,date),UTo)))))),UTo,glong,glat,year,month,date);
     
164. //System.out.println("The time at which the sunset is: "+(int)(sunset/15+Zone(glong))+":"+(int)(60*(sunset/15+Zone(glong)-(int)(sunset/15+Zone(glong))))+" .\n");
     
165. Log.d("Sunset", "Sunset is: "+(int)(sunset/15+8)+":"+(int)(60*(sunset/15+8-(int)(sunset/15+8)))+" .\n");
     
166.        //return "Sunset is: "+(int)(sunset/15+Zone(glong))+":"+(int)(60*(sunset/15+Zone(glong)-(int)(sunset/15+Zone(glong))))+" .\n";
     
167.        return "Sunset is: "+(int)(sunset/15+8)+":"+(int)(60*(sunset/15+8-(int)(sunset/15+8)))+" .\n";
     
168. }
     
169. }

复制代码

//SunRiseSet.getSunRise(GeoPoint geoPoint,Time sunTime),SunRiseSet.getSunSet(GeoPoint geoPoint,Time sunTime)算出的是北京时间下各地区日出日落时间表，计算当地时间的话，把时区8改为Zone(glong)即可。GeoPoint类是地理位置坐标点类，构造器输入是经纬度。

admin

文章来自：
http://blog.csdn.net/gebitan505/article/category/1727369