summaryrefslogtreecommitdiffstats
path: root/src/parameters.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'src/parameters.cpp')
-rw-r--r--src/parameters.cpp44
1 files changed, 20 insertions, 24 deletions
diff --git a/src/parameters.cpp b/src/parameters.cpp
index 57d98c2..2e2758a 100644
--- a/src/parameters.cpp
+++ b/src/parameters.cpp
@@ -1,5 +1,6 @@
1#include "parameters.h"
2#include "calendar.h" 1#include "calendar.h"
2#include "parameters.h"
3#include <utility>
3using namespace std; 4using namespace std;
4 5
5//凡无特殊说明的函数,全部返回弧度制 6//凡无特殊说明的函数,全部返回弧度制
@@ -11,7 +12,7 @@ double parameter::get_longitude(vector<double> l, double t) {
11 L += l[i] * pow(t, i); 12 L += l[i] * pow(t, i);
12 } 13 }
13 //返回弧度制,太阳地心黄经=地球日心黄经+180° 14 //返回弧度制,太阳地心黄经=地球日心黄经+180°
14 return L + pi; 15 return L + M_PI;
15} 16}
16 17
17// 计算地球日心黄纬 18// 计算地球日心黄纬
@@ -36,19 +37,19 @@ double parameter::get_R(vector<double> r, double t) {
36// 转换FK5误差 37// 转换FK5误差
37double parameter::delta_FK5(double L, double B, double T) { 38double parameter::delta_FK5(double L, double B, double T) {
38 // L转角度制 39 // L转角度制
39 L *= 180 / pi; 40 L *= 180 / M_PI;
40 //计算L',角度制 41 //计算L',角度制
41 double L_p = L - 1.397 * T - 0.00031 * T * T; 42 double L_p = L - 1.397 * T - 0.00031 * T * T;
42 //计算L',弧度制 43 //计算L',弧度制
43 L_p *= pi / 180; 44 L_p *= M_PI / 180;
44 while (L_p < 0) 45 while (L_p < 0)
45 L_p += 2 * pi; 46 L_p += 2 * M_PI;
46 while (L_p >= 2 * pi) 47 while (L_p >= 2 * M_PI)
47 L_p -= 2 * pi; 48 L_p -= 2 * M_PI;
48 //计算delta_L,单位为角秒 49 //计算delta_L,单位为角秒
49 double delta_L = -0.09033 + 0.03916 * (cos(L_p) + sin(L_p)) * tan(B); 50 double delta_L = -0.09033 + 0.03916 * (cos(L_p) + sin(L_p)) * tan(B);
50 //转换为弧度制 51 //转换为弧度制
51 delta_L *= pi / (180 * 3600); 52 delta_L *= M_PI / (180 * 3600);
52 return delta_L; 53 return delta_L;
53} 54}
54 55
@@ -88,14 +89,14 @@ double parameter::nutation(double T) {
88 nutation_parameters[i][2] * ang[2] + 89 nutation_parameters[i][2] * ang[2] +
89 nutation_parameters[i][3] * ang[3] + 90 nutation_parameters[i][3] * ang[3] +
90 nutation_parameters[i][4] * ang[4]; 91 nutation_parameters[i][4] * ang[4];
91 sita *= pi / 180; 92 sita *= M_PI / 180;
92 s += (nutation_parameters[i][5] + nutation_parameters[i][6] * T) * 93 s += (nutation_parameters[i][5] + nutation_parameters[i][6] * T) *
93 sin(sita); 94 sin(sita);
94 // c+=(nutation_parameters[i][7]+nutation_parameters[i][8]*T)*cos(sita); 95 // c+=(nutation_parameters[i][7]+nutation_parameters[i][8]*T)*cos(sita);
95 } 96 }
96 // 计算得到的s和c单位0.0001秒,转换为弧度制 97 // 计算得到的s和c单位0.0001秒,转换为弧度制
97 s *= 0.0001; 98 s *= 0.0001;
98 s *= pi / (180 * 3600); 99 s *= M_PI / (180 * 3600);
99 return s; 100 return s;
100} 101}
101 102
@@ -104,7 +105,7 @@ double parameter::aberration(double R) {
104 //公式,单位0.0001角秒 105 //公式,单位0.0001角秒
105 double delta = -20.4898 / R; 106 double delta = -20.4898 / R;
106 //转换为弧度制 107 //转换为弧度制
107 delta *= pi / (180 * 3600); 108 delta *= M_PI / (180 * 3600);
108 return delta; 109 return delta;
109} 110}
110 111
@@ -157,7 +158,7 @@ vector<vector<double>> parameter::get_parameters(double t) {
157} 158}
158 159
159// 计算太阳地心黄经,外部调用,返回角度制 160// 计算太阳地心黄经,外部调用,返回角度制
160double parameter::sun_longitude(double t) { 161pair<double, double> parameter::sun_longitude(double t) {
161 //从文件中计算地日运行参数 162 //从文件中计算地日运行参数
162 vector<vector<double>> p = get_parameters(t); 163 vector<vector<double>> p = get_parameters(t);
163 164
@@ -177,13 +178,8 @@ double parameter::sun_longitude(double t) {
177 //光行差修正,返回弧度制 178 //光行差修正,返回弧度制
178 L += aberration(R); 179 L += aberration(R);
179 180
180 // L转角度制 181 // 返回太阳地心黄经和黄纬
181 L *= 180 / pi; 182 return make_pair(L, -B);
182 while (L < 0)
183 L += 360;
184 while (L >= 360)
185 L -= 360;
186 return L;
187} 183}
188 184
189// 计算月球地心黄经,外部调用,返回角度制 185// 计算月球地心黄经,外部调用,返回角度制
@@ -202,7 +198,7 @@ double parameter::moon_longitude(double t) {
202 moon_longitude_parameters[i][1] * angles[1] + 198 moon_longitude_parameters[i][1] * angles[1] +
203 moon_longitude_parameters[i][2] * angles[2] + 199 moon_longitude_parameters[i][2] * angles[2] +
204 moon_longitude_parameters[i][3] * angles[3]; 200 moon_longitude_parameters[i][3] * angles[3];
205 sita *= pi / 180; 201 sita *= M_PI / 180;
206 EI += moon_longitude_parameters[i][4] * sin(sita) * 202 EI += moon_longitude_parameters[i][4] * sin(sita) *
207 pow(E, fabs(moon_longitude_parameters[i][1])); 203 pow(E, fabs(moon_longitude_parameters[i][1]));
208 } 204 }
@@ -212,12 +208,12 @@ double parameter::moon_longitude(double t) {
212 a1 = 119.75 + 131.849 * T; 208 a1 = 119.75 + 131.849 * T;
213 a2 = 53.09 + 479264.290 * T; 209 a2 = 53.09 + 479264.290 * T;
214 210
215 EI += 3958.0 * sin(a1 * pi / 180); 211 EI += 3958.0 * sin(a1 * M_PI / 180);
216 EI += 1962.0 * sin((L - angles[3]) * pi / 180); 212 EI += 1962.0 * sin((L - angles[3]) * M_PI / 180);
217 EI += 318.0 * sin(a2 * pi / 180); 213 EI += 318.0 * sin(a2 * M_PI / 180);
218 214
219 L += EI / 1000000; 215 L += EI / 1000000;
220 L += 180 / pi * nutation(T); //地球章动修正 216 L += 180 / M_PI * nutation(T); //地球章动修正
221 217
222 while (L < 0) 218 while (L < 0)
223 L += 360; 219 L += 360;