GAMP: An open-source software of multi-GNSS precise point positioning using undifferenced and uncombined observations PPPH: A MATLAB-based software for multi-GNSS precise point positioning analysis goGPS: open-source MATLAB software -- by Antonio M. Herrera, Hendy F. Suhandri, Eugenio Realini, Mirko Reguzzoni, and M. Clara de Lacy
Доброе утро. Что то не сходится у меня расчет в это коде. Запускаю с теми исходными данными что заданы. результат сходится. считываю навигационное сообщение на этот же спутник и на это же время но на следующий день. результат сходится. Ввожу данные и время другого спутника, - результат не совпадает. вот 15 спутник - эфемериды: * 2019 2 12 8 0 0.00000000 PG15 323.842837 -24393.481094 -10268.435720 -328.217094 8 9 5 84 навигационное сообщение 15 19 2 12 8 0 0.0-0.328210648149D-03 0.170530256582D-11 0.000000000000D+00 0.160000000000D+02-0.331250000000D+01 0.523093213189D-08-0.124802397881D+01 -0.182539224625D-06 0.112132804934D-01 0.928156077862D-05 0.515363610840D+04 0.201600000000D+06 0.296160578728D-06 0.897675791800D+00 0.726431608200D-07 0.928038015777D+00 0.181906250000D+03 0.763414028836D+00-0.831177437988D-08 -0.243938730327D-09 0.100000000000D+01 0.204000000000D+04 0.000000000000D+00 0.200000000000D+01 0.000000000000D+00-0.107102096081D-07 0.160000000000D+02 0.197370000000D+06 0.000000000000D+00 0.000000000000D+00 0.000000000000D+00 расшифровка: double Crs = -0.331250000000E+01; double dn = 0.523093213189E-08; double M0 = -0.124802397881E+01; double Cuc = -0.182539224625E-06; double ec = 0.112132804934E-01; double Cus = 0.928156077862E-05; double a = 0.515363610840E+04 * 0.515363610840E+04; double Toe = 0.201600000000E+06; double Cic = 0.296160578728E-06; double W0 = 0.897675791800E+00; double Cis = 0.726431608200E-07; double i0 = 0.928038015777E+00; double Crc = 0.181906250000E+03; double w = 0.763414028836E+00; double Wdot = -0.831177437988E-08; double idot = -0.243938730327E-09; время T - double T = 28800; результат X = 23761731.699 Y = 4267023.204 Z = -11757519.146
Сравниваем с данными в SP3 на момент * 2019 2 12 8 0 0.00000000 … PG15 323.842837 -24393.481094 -10268.435720 Берём из навигационного файла последовательно на моменты (выделено жирным) 15 19 2 12 5 59 44.0-0.328223221004D-03 0.170530256582D-11 0.000000000000D+00 … 15 19 2 12 8 0 0.0-0.328210648149D-03 0.170530256582D-11 0.000000000000D+00 … 15 19 2 12 10 0 0.0-0.328198540956D-03 0.170530256582D-11 0.000000000000D+00 Получим: T = 7216.0000000 X = 323855.110 Y = -24393486.501 Z = -10268423.228 T = 0.0000000 X = 323854.697 Y = -24393486.153 Z = -10268423.669 T = -7200.0000000 X = 323854.624 Y = -24393486.029 Z = -10268423.571
Добрый день. Вот время Toe, которое в эфемеридах, - оно там представлено на начало или конец интервала действия навигационного сообщения? вот я смотрю toe на 19 число с 2 до 4 UTC, там 1800000 у меня в это время в 9-45 (новосибирск, разница 7 часов), 187200. (это данные от модуля). Отчего вопрос? У меня время (Врем¤ недели измерени¤ по локальному времени приемника. rcvTow \ в мс ) меньше чем toe. цифры такие, - если toe180000, то Tow 179135, в другом случае, - 194400 - 190784. Мне кажется, так быть не должно.
Все таки не понимаю....... везде пишут: t_oe - опорная эпоха (время навигационного сообщения) в секундах от начала недели смотрю часы: вот время на этот момент - 265444s. беру данные от навигатора (в сообщении псевдодальности) примерно в это же время - сообщение Tow 265ххх секунд, Toe из навигатора с сообщении эфемерид - 266400. Но ведь не должна же быть опорная эпоха больше, чем настоящее время?
Почему? --- Сообщения объединены, 20 фев 2019, Оригинальное время сообщения: 20 фев 2019 --- Воспользуйтесь советом Александр Яковченко, Обратите внимание на раздел 20.3.4.5 Reference Times. и фразу "Note that a user applying current navigation data will normally be working with negative values of (t-toc) and (t-toe) in evaluating the expansions."
1. Эпоха в астрономии (от греч. έποχή — «остановка») — момент времени, для которого определены астрономические координаты или элементы орбиты. 2. t_oe - опорная эпоха (время навигационного сообщения) в секундах от начала недели То есть в моем понимании, - от начала, - до начала следующей эпохи. другого объяснения не знаю. --- Сообщения объединены, 20 фев 2019, Оригинальное время сообщения: 20 фев 2019 --- спасибо. перечитаю. до этого я не дочитал...... значит все правильно пока.........
Отличие астрономического термина эпоха от "исторического" (и бытового, без всякого отрицательного смысла) в том, что в астрономии эпоха – это просто конкретный момент времени (любой), а в истории эпоха – интервал времени между моментами. И как момент времени его совсем не обязательно связывать с координатами или любым фазовым вектором состояния любого тела. Поэтому совсем нет смысла говорить о "начале следующей эпохи". Это начало следующей эпохи совсем не нужно, да и его может просто не быть. P.S. Русская вики по своему обыкновению врёт, чаще в деталях.
Добрый вечер. Подскажите пожалуйста, - какой точности вычисления координат спутников достаточно для расчета координат точки на местности, с точностью в пределах плюс/минус 100 метров.
Доброе утро. применил поправки. Код: void calk_SV1(int x) // расчет положения спутника по его эфемеридам { double E; a_bp = pow(GP[x].square_big_axis,2); // вычмсляем большую полуось n_o = sqrt(gp_earth/(pow(a_bp,3))); // n = n_o + GP[x].delta_n; // находим среднее движение // t = t - ( PR_sv[x] / c );// t_k = t - GP[x].t_oe ; // вычисляем промежуток времени от опорной эпохи эфемерид if ( t_k > 302400) // если t больше 302400 секунд, вычтите 604800 секунд из T. t_k = t_k - 604800; else if ( t_k < -302400) // Если T меньше -302400 секунд, добавьте 604800 секунд к T. t_k = t_k + 604800; M_sr = GP[x].Mo + (n * t_k); // определяем среднюю аномалию E1 = M_sr; // начальное условие для вычислений уравнения Кеплера for (;;) // находим эксцентрическую аномалию через уравнение Кеплера { E1 = M_sr + (GP[x].E_a * sin(E)); // Kepler_equation(x); if(fabs(E1 - E) > 1e-19) E = E1; // E это по документу Ek else break; } ////////////////////////коррекция синхроимпульсов double F,d_tr, d_tsv; F = - 4.442807633E-10; // (-2 * sqrt( gp_earth)) / (c * c); d_tr = F * GP[x].E_a * GP[x].square_big_axis * sin (E); // E_a - это бывшая e d_tsv = GP[x].af0 + (GP[x].af1*(t - GP[x].t_oc)) + (GP[x].af2 * pow (t -GP[x].t_oc,2)) + d_tr - GP[x].t_gd; ; // t = t - ( PR_sv[x] / c );// t = t - d_tsv; // вычисляем t с учетом поправки t_k = t - GP[x].t_oe ; // вычисляем промежуток времени от опорной эпохи эфемерид if ( t_k > 302400) // если t больше 302400 секунд, вычтите 604800 секунд из T. t_k = t_k - 604800; else if ( t_k < -302400) // Если T меньше -302400 секунд, добавьте 604800 секунд к T. t_k = t_k + 604800; M_sr = GP[x].Mo + (n * t_k); // определяем среднюю аномалию E1 = M_sr; // начальное условие для вычислений уравнения Кеплера for (;;) // находим эксцентрическую аномалию через уравнение Кеплера { E1 = M_sr + (GP[x].E_a * sin(E)); // Kepler_equation(x); if(fabs(E1 - E) > 1e-19) E = E1; // E_a - это бывшая e else break; } v = (atan((sqrt((1 + GP[x].E_a)/(1 - GP[x].E_a))) * tan(E/2)))*2; //v = atan((sqrt(1-pow(GP[x].E_a,2))*sin(E)/(1-GP[x].E_a*cos(E)))/((cos(E)-GP[x].E_a)/(1-GP[x].E_a*cos(E)))); printf("v : %.12f \n",v); nash = v + GP[x].omega; // находим невозмущенный аргумент широты sigma_u = GP[x].Cus * sin(2 * nash) + GP[x].Cuc * cos(2 * nash); sigma_r = GP[x].Crs * sin(2 * nash) + GP[x].Crc * cos(2 * nash); sigma_i = GP[x].Cis * sin(2 * nash) + GP[x].Cic * cos(2 * nash); u_eph = nash + sigma_u; r_eph = a_bp * (1 - GP[x].E_a * cos(E)) + sigma_r; i_eph = GP[x].io + sigma_i + GP[x].IDOT * t_k; X1 = r_eph * cos(u_eph); Y1 = r_eph * sin(u_eph); OMEGA_r = GP[x].omega_0 + (GP[x].OMEGADOT - OMEGA_zemli)*t_k - OMEGA_zemli*GP[x].t_oe; X_wgs84 = X1 * cos(OMEGA_r)- Y1 * cos(i_eph) * sin(OMEGA_r); Y_wgs84 = X1 * sin(OMEGA_r)+ Y1 * cos(i_eph) * cos(OMEGA_r); Z_wgs84 = Y1 * sin(i_eph); } Вот такие результаты: X23 : 6582196.1679999996 из файла sp3 X : 6582196.4182939604 результат без поправок X : 6582196.1691133399 результат с поправками Y23 : 17468760.1230000000 Y : 17468759.9883824740 Y : 17468760.3942986840 Z23 : 19054289.1160000000 Z : 19054289.1614533770 Z : 19054288.8624326470 Мне хотелось бы понять, - правильно ли я применяю поправки, и до какой точности нужно добиваться совпадения результата с файлом SP3/ желательно не отправлять меня в англоязычный сегмент интернета, а просто ответить.......
Добрый день. кажется написал расчет координат по методу наименьших квадратов, - для 4 спутников. получились вот такие результаты: координаты спутников: //24/02/19 10:30 SV: 1 X : -13511943.4625868720 Y : 22537777.8291581910 Z : -2579689.0497690514 - * ---- end ---- * - SV: 7 X : -2361138.8937401767 Y : 21700294.4831329210 Z : 15045482.4396653310 - * ---- end ---- * - SV: 8 X : -10765179.0955109730 Y : 11387996.6558811390 Z : 21462376.4363597040 - * ---- end ---- * - SV: 30 X : 7438780.3448691871 Y : 15023634.4812066070 Z : 20578692.6607901340 - * ---- end ---- * - //Time gps: 37800 псевдодальности: //SV: 1 Q: 7 Psr: 25109839.64387868 //SV: 7 Q: 4 Psr: 21056670.78789259 //SV: 8 Q: 6 Psr: 21555412.29389584 //SV: 30 Q: 4 Psr: 20686354.37529096 результат: L : 83.229887854074036 B : 54.286477535936783 H : 107614.049749834040000 Координаты: E : 83,13,47 N : 54,17,11 нет ли возможности проверить правильность моих вычислений?
Вот я выше немного как раз об этом и спрашивал, - и код приложил. Все ли поправки я учел, если нет, - что упустил?
Добрый день. В рекомендованном исходнике программы на Паскале, есть процедура, вычисляющая азимут и угол возвышения. Не мог кто нибудь написать формулу, по которой здесь угол возвышения считается. Формул в интернете много, - но хотелось бы опознать эту. Вопрос возник оттого, что как мне кажется, я неверно перевожу этот, на первый взгляд простенький алгоритм на Си. Код: procedure calcAzEl(Xs, {satellite ECEF XYZ} Xu : vec3; {user ECEF XYZ} var Az, {azimuth [rad]} El : real; {elevation [rad]} var stat: boolean); {calculation succeeded: stat = true} var R, p, x, y, z, s: real; e: array[1..3,1..3] of real; i, k: integer; d: vec3; begin x := Xu[1]; y := Xu[2]; z := Xu[3]; p := sqrt(x*x + y*y); if p = 0.0 then stat := false else stat := true; if stat then begin R := sqrt(x*x + y*y + z*z); e[1,1] := - y / p; e[1,2] := + x / p; e[1,3] := 0.0; e[2,1] := - x*z / (p*R); e[2,2] := - y*z / (p*R); e[2,3] := p / R; e[3,1] := x / R; e[3,2] := y / R; e[3,3] := z / R; for k := 1 to 3 do begin d[k] := 0.0; for i := 1 to 3 do d[k] := d[k] + (Xs - Xu) * e[k,i] end; s := d[3] / sqrt(d[1]*d[1] + d[2]*d[2] + d[3]*d[3]); if s = 1.0 then El := 0.5 * pi else El := arctan(s / sqrt(1.0 - s*s)); if (d[2] = 0.0) and (d[1] > 0.0) then Az := 0.5 * pi else if (d[2] = 0.0) and (d[1] < 0.0) then Az := 1.5 * pi else begin Az := arctan(d[1] / d[2]); if d[2] < 0.0 then Az := Az + pi else if (d[2] > 0.0) and (d[1] < 0.0) then Az := Az + 2.0 * pi end; end; end; {procedure calcAzEl}
Elevation and azimuth computation З.Ы. То, что малограмотные переводчики с англицкого называют углом возвышения в русскоязычной астрономической литературе называется просто высотой. Высота и азимут – это две координаты в горизонтной (горизонтальной) сферической системе координат.
Доброе утро. Докопался до коррекции псевдодальностей, - {correct pseudorange}Pcor[prn] := Praw[prn] + dTclck * c - dTiono * c - dRtrop + Cr Praw[prn] - псевдодальность полученная от модуля. dTiono - ионосферная поправка, - которую сейчас считаю, с - скорость света. вот где мне взять еще эти три поправки. dTclck , dRtrop, Cr В основном документе я их не увидел.
Тропосферную модель для расчета поправки мне подсказали, dTclck вроде бы высчитывается при расчетах методом наименьших квадратов, одновременно с координатами, А вот Cr - так и не понятно, что это и для чего.
