1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
|
// Encoding: UTF-8
#include "calendar.h"
#include <cstdio>
#include <ctime>
#include <iomanip>
/*以下为"calendar.h规定之量*/
double day = 86400;
double delta = 1e-11;
Julian julian;
parameter p;
pair<double, double> getSunPos(time_t time) {
double t = julian.getJulianKiloYear(time);
pair<double, double> sun = p.sun_longitude(t);
radToDMS("Sun longitude:\t\t", 0, 360, sun.first);
radToDMS("Sun latitude:\t\t", -90, 90, sun.second);
double alpha, delta; // 太阳赤经赤纬
double lambda = sun.first, beta = sun.second; // 太阳黄经黄纬
// 黄赤交角epsilon
double epsilon = p.get_epsilon(t);
radToDMS("Ecliptic Obliquity:\t\t", 0, 90, epsilon);
delta =
asin(sin(epsilon) * sin(lambda) * cos(beta) + cos(epsilon) * sin(beta));
alpha = atan2(
(cos(epsilon) * sin(lambda) * cos(beta) - sin(epsilon) * sin(beta)),
(cos(lambda) * cos(beta)));
radToDMS("Sun right ascension:\t\t", 0, 360, alpha);
radToDMS("Sun declination:\t\t", -90, 90, delta);
// 时角
double H = p.getHourAngle(t, LONGITUDE, LATITUDE, alpha);
double A, h; //方位角和高度角
h = asin(sin(LATITUDE) * sin(delta) + cos(LATITUDE) * cos(delta) * cos(H));
A = atan2(-cos(delta) * sin(H),
cos(LATITUDE) * sin(delta) - sin(LATITUDE) * cos(delta) * cos(H));
A += M_PI;
// 转换为角度制输出
radToDMS("Hour angle:\t\t", 0, 360, H);
radToDMS("Azimuth:\t\t", 0, 360, A);
radToDMS("Elevation:\t\t", -90, 90, h);
}
int main(int argc, char *argv[]) {
Date date;
if (argc != 2) {
printf("Input the time you want to calculate in <YYYY-MM-DD,HH:MM:SS> "
"format:\t\t");
scanf("%d-%d-%d,%d:%d:%d", &date.tm_year, &date.tm_mon, &date.tm_mday,
&date.tm_hour, &date.tm_min, &date.tm_sec);
} else {
sscanf(argv[1], "%d-%d-%d,%d:%d:%d", &date.tm_year, &date.tm_mon,
&date.tm_mday, &date.tm_hour, &date.tm_min, &date.tm_sec);
}
date.tm_year -= 1900;
date.tm_mon -= 1;
date.tm_isdst = -1;
time_t time = mktime(&date);
getSunPos(time);
return 0;
}
|
