Файл: Учебное пособие ( Лабораторный практикум на компьютере ) Київ 2008 1.pdf

Скачать файл
Заказать решение

ВУЗ: Не указан

Категория: Не указан

Дисциплина: Не указана

Добавлен: 11.01.2024

Просмотров: 338

Скачиваний: 5

ВНИМАНИЕ! Если данный файл нарушает Ваши авторские права, то обязательно сообщите нам.

СОДЕРЖАНИЕ

Функция Read_GL_Alm function [alm,max_kol] = Read_GL_Alm(Dat); %Имя функции: Read_GL_Alm %Функция читает данные альманаха спутников ГЛОНАСС. %Входные данные: %Dat - имя файла, содержащего альманах спутников ГЛОНАСС , например, Dat = 'GLN_all.alm'; %Выходные данные: %alm - данные альманаха спутников ГЛОНАСС (структура) , %max_kol- максимальное количество спутников в альманахе fid =fopen(Dat,'rt');%открытие файла max_kol = 0; %чтение данных из файла while not(feof(fid)) s1=fscanf(fid,'%s',6); if not(feof(fid)) lenstr = length(s1); max_kol = max_kol + 1; % while (fscanf(fid,'%s',1) == '-') end str1 = fscanf(fid,'%s',1); % lenstr = length(str1); str2 = fscanf(fid,'%s',1); lenstr = length(str2); n_sv = sscanf(str2,'%d'); strID = str2(1:lenstr); ID = sscanf(strID,'%d'); alm(ID).ID = ID; %номер спутника (1...24) %%%%%%%%%%% t_1=fscanf(fid,'%s',1); alm(ID).Hn=fscanf(fid,'%d'); % номер частоты спутника t_2=fscanf(fid,'%s',1); alm(ID).Health=fscanf(fid,'%d',1);%здоровье (0- спутник здоров) t_3=fscanf(fid,'%s',1); alm(ID).ecc = fscanf(fid,'%g',1);% эксцентриситет орбиты спутника t_4=fscanf(fid,'%s',1); 108 alm(ID).Na =fscanf(fid,'%g',1 );% номер суток к которым относятся данные альманаха, отсчиты- ваемых %от ближайшего високосного года t_5=fscanf(fid,'%s',3); alm(ID).deltai = fscanf(fid,'%g',1); % наклонение орбиты спутника (радианы) t_6=fscanf(fid,'%s',2); alm(ID).LambdaN = fscanf(fid,'%g',1); % долгота восходящего узла орбиты спутника (радианы) while not(fscanf(fid,'%c',1) == ':') end % t_7=fscanf(fid,'%s',1) alm(ID).TLambdaN = fscanf(fid,'%g',1) - 10800;% время прохождения восходящего узла орбиты спутника, %к которому относятся данные альманаха, приведенное к времени UTC (секунда) t_8=fscanf(fid,'%s',2); alm(ID).omegan = fscanf(fid,'%g',1); %аргумент перигея (радианы) t_11=fscanf(fid,'%s',1); alm(ID).Tdr=fscanf(fid,'%g',1); % драконический период (секунда/виток) t_12=fscanf(fid,'%s',1); alm(ID).dTdr=fscanf(fid,'%g',1); % скорость изменения драконического периода (секунда/виток в квадрате) t_13=fscanf(fid,'%s',2); alm(ID).tau_n=fscanf(fid,'%g',1) / 1000.0;%коэффициет коррекции шкалы времени (сдвиг времени спутника %относительно системного времени ГЛОНАСС в секундах) end; end fclose(fid); Функция Gln_data_from_NA function [time_UTC] = Gln_data_from_NA(leap_year, day_from_leap); %Имя функции:Gln_data_from_NA %Функция предназначена для преобразования номера дня NA (день привязки альманаха % от ближайшего високосного года) в текущую дату %function [time_UTC] = Gln_data_from_NA(leap_year, day_from_leap); %Входные данные: % leap_year - ближайший високосный год %day_from_leap (NA) - день привязки альманаха %Выходные данные: %структура timeUTC (год, месяц, день) - текущая дата DnMon = [31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31]; %количество дней в месяцах 109 n4 = mod(leap_year, 4); n100 = mod(leap_year, 100); n400 = mod(leap_year, 400); if (n4 == 0) n_leap = 1; else n_leap = 0; end if ((n100 == 0) & (n400 > 0)) n_leap = 0; end if (day_from_leap > (365 + n_leap)) day = day_from_leap - (365 + n_leap); k =fix (day / 365); day = mod(day, 365); else day = day_from_leap; k = -1; end; god = leap_year + k + 1; if (god > leap_year) n_leap = 0; end; mon = 1; mon_day = 31; while (day > mon_day) day = day - mon_day; mon = mon + 1; % LK mon_day = DnMon(mon); % mon = mon+1; if (mon == 2) mon_day = mon_day + n_leap; end end time_UTC.year = god; time_UTC.mon = mon; time_UTC.day = day; Функция JD_data function [JD, day_year] = JD_data(timeUTC) %Имя:JD_data % Функция JD_data(timeUTC) вычисляет : 110 %JD - номер юлианского дня, day_year - номер дня года. %Входные данные: %Структура timeUTC %timeUTC.year - год, % timeUTC.mon - месяц, % timeUTC.day - день. %Выходные данные: %JD - юлианский день; %day_year- день от начала года. %количество дней в месяцах DnMon = [31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31]; %Вычисление номера юлианского дня опорной эпохи jd0 = JD_epohi(timeUTC.year); %Учет високосного года nfebr = 0; if mod(timeUTC.year,4) == 0 nfebr = 1; end; %Расчетномера дня года k = 0; for i = 2 : timeUTC.mon k = k + DnMon(i - 1); if (i == 2) k = k + nfebr; end; end; day_year = k + timeUTC.day; %Расчет номера юлианского дня JD = jd0 + day_year; Функция LLH_to_ECEF function [XYZ] = LLH_to_ECEF( a, b, llh) %Имя функции:LLH_to_ECEF.m %Функция преобразования географических координат в прямоугольную геоцентрическую систему координат (ECEF) %Входные данные: %Структура llh %llh.lon-долгота (радиан), %llh.lat-широта (радиан), %llh.h-высота (метр); %a, b-большая и малая полуоси земного эллипсоида в WGS-84 (метр); %Выходные данные: 111 %Структура XYZ %XYZ.x - координата x в ECEF; %XYZ.y - координата y в ECEF; %XYZ.z - координата z в ECEF; a2=a * a; b2=b * b; r = a2 / sqrt(a2 * cos(llh.lat) * cos(llh.lat) + b2 * sin(llh.lat) * sin(llh.lat)); XYZ.x = (r + llh.h) * cos(llh.lat) * cos(llh.lon); XYZ.y = (r + llh.h) * cos(llh.lat) * sin(llh.lon); XYZ.z = (b2 / a2 * r + llh.h) * sin(llh.lat); Функция WGS84_to_PZ90 function [pos_pz90] = WGS84_to_PZ90( pos_wgs84) %Имя функции: WGS84_to_PZ90.m %Функция преобразует координаты из системы WGS 84 в систему ПЗ 90 p = (1.0 - 0.12e-6); x = pos_wgs84.x; y = pos_wgs84.y; z = pos_wgs84.z; pos_pz90.x = (p * x + p * 0.82e-6 * y + 1.4 ); pos_pz90.y = (-p * 0.82e-6 * x + p * y + 1.4); pos_pz90.z = (p * z + 0.9); Функция ECEF_to_LLH function [llh] = ECEF_to_LLH(a, b, XYZ) %Имя функции:ECEF_to_LLH %Функция преобразует координаты из системы ECEF в географическую систему %Входные данные: %a, b-большая и малая полуоси земного эллипсоида (метр); ( a=6378137.0; %b=6356752.314; метры) %Структура XYZ %XYZ.x - координата x в ECEF; %XYZ.y - координата y в ECEF; %XYZ.z- координата z в ECEF; %Выходные данные: %Структура llh %llh.lon - долгота (радиан); %llh.lat - широта (радиан); %llh.h-высота (метр); a2 = a * a; b2 = b * b; xy = sqrt(XYZ.x * XYZ.x + XYZ.y * XYZ.y); thet = atan2(XYZ.z * a , (xy * b)); 112 %thet = atan(XYZ.z * a / (xy * b)); esq = 1.0 - b2 / a2; epsq = a2 / b2 - 1.0; llh.lat = atan((XYZ.z + epsq * b * (sin(thet) ^3)) / (xy - esq * a * (cos(thet) ^3))); llh.lon = atan2(XYZ.y, XYZ.x);%! while (llh.lon < 0) llh.lon = 2*pi + llh.lon; end; %if llh.lon < 0 % llh.lon = 2*pi + llh.lon; %end ; r = a2 / sqrt(a2 * cos(llh.lat) * cos(llh.lat) + b2 * sin(llh.lat) * sin(llh.lat)); llh.h = xy / cos(llh.lat) - r; end Функция GLN_satfind function [kol_gln_a, satvis_gln, satpos_gln_ns, satpos_gln_ecef] = GLN_satfind(a, b, timeUTC, cur- rent_loc_pz90, alm_gln); %Имя функции:GLN_satfind %Функция вычисляет координаты спутников ГЛОНАСС, углы видимости и количество видимых спутников %Входные данные: %a, b- большая и малая полуоси земного эллипсоида; %timeUTC- параметры времени; %current_loc_pz90- широта, долгота и высота точки, из которой наблюдаются спутники (структура); % alm_gln-альманах спутников ГЛОНАСС (структура) %Выходные данные: %kol_gln_a- количество видимых спутников ГЛОНАСС, %satvis_gln- углы видимости спутников (структура), %satpos_gln_ns-координаты спутников (структура) в абсолютной системе координат (ECI) %satpos_gln_ecef- координаты спутников (структура) в подвижной системе координат (ECEF) RAD_TO_DEG = 180.00 / pi ; KOL_GLN = 24; nom_ns = 1; Min_Elev = 0; % минимальный угол видимости kol_gln_a = 0; [nom_sat_gln_a, trac_sat_gln_a, sat_pos, satvis_gln] = init_data(KOL_GLN); % число дней от фундаментальной эпохи 2000г от 12h, 0, January: % jd_up_epoh = JD_from_epohi( 2000, timeUTC) - 0.5; % from 12h UTC 2000 nut = 0; 113 S0 = s0_Nut( timeUTC, nut); time_s0 = S0.s0_m_mod ; %time_s0 - истинное звездное время в день обсервации, 0h UTC year = timeUTC.year; leap_year = fix( year / 4) * 4; % ближайший к текущему (предыдущий) високосный год ti = timeUTC.ti; % текущее время обсервации от начала дня n00 = fix(ti / 86400); % n0 - номер текущих суток внутри 4-х летнего периода (от ближайшего високосного года) n0 = JD_from_epohi(leap_year, timeUTC) + n00 + 1; [eci_current_loc, eci_rec_pos_xyz] = llh_to_eci(a, b, ti, time_s0, current_loc_pz90); satpos_eci = init_satpos_gln(); for ( i = 1 : KOL_GLN) prn = alm_gln (i).ID; health = alm_gln (i).Health; if ( (prn > 0) & (health == 0)) satpos_eci = gln_a_1(i, n0, ti, time_s0, alm_gln); [satpos_eci, satpos_gln_a] = satpos_eci_in_metr(satpos_eci); %Transformation coordinates from "Km" in "m" satpos_gln_ns(prn) = satpos_gln_a; %----- [satpos_ecef] =eci_to_ecef(time_s0, ti, satpos_eci); satpos_gln_ecef(prn) = satpos_ecef; % calculates angles of visibility by almanac: [top] = top_coord(prn, eci_current_loc, eci_rec_pos_xyz,satpos_gln_a); satvis_gln(prn) = top; if (satvis_gln(prn).r > 0.0) el = RAD_TO_DEG * satvis_gln(prn).el; if ( el >= Min_Elev) % fprintf('prn=%2i el=%f az=%f \n ', prn, el, RAD_TO_DEG * satvis_gln(prn).az); satpos_gln( prn) = satpos_eci; [satpos_ecef] =eci_to_ecef(time_s0, ti, satpos_eci); [trac_sat_gln_a] = rewrite_satpos(nom_ns, satpos_ecef); trac_sat_gln_a( nom_ns).r = satvis_gln( prn).r; trac_sat_gln_a (nom_ns).prn = prn; nom_sat_gln_a(nom_ns) = prn; nom_ns = nom_ns + 1; % number of visible satellites end; % if ( el >= Min_Elev) end; % if ( satvis_gln[prn-1].r > 0.0) 114 end; % f ( (prn > 0) & (health == 1)) end; % for ( i = 1 : KOL_GLN) kol_gln_a = nom_ns - 1; Функция ris_vis_sat function ris_vis_sat(max_n, nj, ti_start, step_t, ris_vid) %Имя функции:ris_vis_sat %Функция строит графики спутников, видимых из точки наблюдения mm = step_t/3600.0; j_color = 0; color7(1:7) = ['b' 'g' 'r' 'c' 'm' 'y' 'k']; for i = 1 : max_n str1 = ris_vid(i, 1: nj ); k = find(str1); if k > 0 len = length(k); if len == 1 hPlot = plot(k(1):k(1),ris_vid(i,k(1):k(1)), s); set(hPlot, 'LineWidth', 3); hold on; else n0 = k(1); j_color = j_color + 1; if (j_color > 7 ) j_color = 1; end s = color7(j_color); for j = 2 : len if (k(j) > (k(j-1) + 1)) | (j == len) k0 = k(j - 1); hPlot = plot(n0*mm : mm : k0*mm,ris_vid(i,n0:k0), s); set(hPlot, 'LineWidth', 3); hold on; n0 = k(j); end; % if end; % for j = 2 : len end; % else end; % if k > 0 end; % for i = 1 : max_n 115 grid on hAxes=gca; My = [0:2:25]; kol_x = ti_start + nj * step_t; Mx = [0 : 2 : (kol_x - ti_start) / 3600]; set(hAxes, 'ytick', My,'xtick',Mx) ; hText = gca; set(hText,'FontSize',12,'FontName','TimesNewRoman') xlabel('Время наблюдения спутников ГЛОНАСС, час') ylabel('Номер навигационного спутника') Функция JD_epohi function jden = JD_epohi(epoha) %Имя: JD_epohi %Фукция JD_epohi(epoha) рассчитывает номер юлианского дня %до начала года (epoha) на 12h, 0 день, январь. %Входные данные: epoha, размерность-год %Выходные данные: % jden- номер юлианского дня на 12h, 0 день, январь ( размерность -дни) rk = mod(epoha,4); if ( rk == 0 ) rk = 1.0; else rk = 2.0 - rk * 0.25; end; n100 = floor(epoha / 100); n400 = floor(epoha / 400); jden = (4712 + epoha) * 365.25 + n400 - n100 + rk; % fprintf('epoha=%d rk=%f jden=%6.2f \n', epoha, rk, jden); Функция init_data function[nom_sat_gln_a, trac_sat, sat_pos, satvis_gln] = init_data(kol) %Имя функции: init_data %Функция предназначена для инициализации массивов sat_pos.prn = 0; sat_pos.x = -1.0; sat_pos.y = -1.0; sat_pos.z = -1.0; sat_pos.vx = -1.0; sat_pos.vy = -1.0; sat_pos.vz = -1.0; 116 sat_pos.dvx = -1.0; sat_pos.dvy = -1.0; sat_pos.dvz = -1.0; sat_pos.range_m = -1.0; for( i = 1 : kol) nom_sat_gln_a(i) = 0; satvis_gln(i).el = -1; satvis_gln(i).az = -1; satvis_gln(i).r = 0; trac_sat(i).prn = 0; trac_sat(i).x = -1.0; trac_sat(i).y = -1.0; trac_sat(i).z = -1.0; trac_sat(i).vx = -1.0; trac_sat(i).vy = -1.0; trac_sat(i).vz = -1.0; trac_sat(i).range_m = -1.0; end; Функция JD_from_epohi function [jd] =JD_from_epohi( epoha, timeUTC); %Имя функции:JD_from_epohi.m %Функция вычисляет jd - количество дней от указанного года (epoha) % до текущей даты, указанной в структуре timeUTC, представленной в виде % (timeUTC.year, timeUTC.mon, timeUTC.day) jd0 = JD_epohi(epoha) + 1; % 12h, 1 den January [day, day_year] = JD_data(timeUTC); jd = day - jd0; Функция s0_Nut function [S0] = s0_Nut( timeUTC, nut) %Имя функции: s0_Nut %Функция рассчитывает истинное или среднее звездное время на 0ч UTC %Входные данные: %timeUTC - дата, на которую требуется рассчитать истинное или среднее звездное время %nut- признак (если nut= 0, то вычисляется среднее звездное время без учета нутации, иначе вычис- ляется %истинное звездное время) %Выходные данные: %S0 - истинное или среднее звездное время на 0ч UTC jd2000 = 2451545; % 12h UTC 1 января % Применить функцию JD_data [jd, day_year] = JD_data( timeUTC); if (jd == NaN) 117 s0_mod = NaN; h = NaN; min = NaN; sec = NaN; fprintf('function s0_m - end0 \n'); return; end; jd = jd - 0.5; d = jd - jd2000; t = d / 36525.0; % 36525 - юлианский период 100 лет t2 = t * t; h1 = 24110.54841; %h1=6.0*3600.0+41.0*60.0+50.54841; % h2 = 236.555367908 * d; h2 = 8640184.812866 * t ; h3 = 0.093104 * t2; h4 = t2 * t * 6.2E-6; if ( nut == 0) na = 0; else na = utc_nut(t); end; s0_m = h1 + h2 + h3 - h4; S0.s0_nut = s0_m + na; S0.s0_m_mod = mod(s0_m, 86400); s0_day = floor(s0_m / 86400); S0.s0_m_hour = S0.s0_m_mod / 3600.0; S0.s0_m_hour = floor(S0.s0_m_mod / 3600); sec_min = S0.s0_m_mod - S0.s0_m_hour * 3600; S0.s0_m_min = floor(sec_min / 60); S0.s0_m_sec = sec_min - S0.s0_m_min * 60; S0.s0_nut_mod = mod(S0.s0_nut, 86400); s0_day = floor(S0.s0_nut / 86400); S0.s0_nut_hour = S0.s0_nut_mod / 3600.0; S0.s0_nut_hour = floor(S0.s0_nut_mod / 3600); sec_min = S0.s0_nut_mod - S0.s0_nut_hour * 3600; S0.s0_nut_min = floor(sec_min / 60); S0.s0_nut_sec = sec_min - S0.s0_nut_min * 60; Функция llh_to_eci function [eci_llh, eci_xyz] = llh_to_eci(a, b, ti, time_s0, llh_loc) ; %Имя функции:llh_to_eci %Функция вычисляет позицию приемника в абсолютной геоцентрической системе координат (ECI) %Входные данные: %a, b-большая и малая полуоси земного эллипсоида (метр); 118 %ti- текущее время (секунды) , %time_s0- истиное звездное время, %Структура llh_loc - координаты приемника; { %llh_loc.lon-долгота (радиан); %llh_loc.lat-широта (радиан); %llh_loc.h-высота (метр); %Выходные данные: %Структура eci_llh - географические координаты приемника в абсолютной геоцентрической системе координат (ECI) %eci_llh.lon - долгота (радиан); %eci_llh.lat - широта (радиан); %eci_llh.h - высота (метр); %Структура eci_xyz- координаты приемника в абсолютной прямоугольной геоцентрической системе координат (ECI) % eci_xyz.x - координата x; % eci_xyz.y - координата y; % eci_xyz.z - координата z; OMEGA_Z = 0.7292115e-4; % угловая скорость вращения Земли, (рад/cek) SEC_IN_RAD = 7.2722052166430e-005; % (PI / 43200.0 ) // Number radian eci_llh.lon = llh_loc.lon + ti * OMEGA_Z + SEC_IN_RAD * time_s0; eci_llh.lat = llh_loc.lat; eci_llh.h = llh_loc.h; eci_xyz = llh_to_ecef( a, b, eci_llh); Функция init_satpos_gln function [sat_pos] = init_satpos_gln() ; %Имя функции:init_satpos_gln %Функция предназначена для инициализации структуры sat_pos sat_pos.prn = 0; sat_pos.x = -1.0; sat_pos.y = -1.0; sat_pos.z = -1.0; sat_pos.vx = -1.0; sat_pos.vy = -1.0; sat_pos.vz = -1.0; sat_pos.dvx = -1.0; sat_pos.dvy = -1.0; sat_pos.dvz = -1.0; sat_pos.range_m = -1.0; Функция gln_a_1 function [satpos_xyz_gln] = gln_a_1(ns, n0, ti_current, time_s0, alm_gln); 119 %Имя функции:gln_a_1 %Функция рассчитывает координаты и скорости спутников ГЛОНАСС в соответствии с интерфейс- ным %контрольным документом ГЛОНАСС %Входные данные: %ns- номер спутника , %n0 - номер текущих суток внутри 4-х летнего периода (от ближайшего високосного года), %ti_current - текущее время обсервации от начала дня, %time_s0 - истинное звездное время в текущий момент обсервации, %alm_gln - альманах спутников ГЛОНАСС %Выходные данные: %Структура satpos_xyz_gln %satpos_xyz_gln.x - координата по оси x; %satpos_xyz_gln.y- координата по оси y; %satpos_xyz_gln.z- координата по оси z; %satpos_xyz_gln.vx- скорость по оси x; %satpos_xyz_gln.vy - скорость по оси y; %satpos_xyz_gln.vz- скорость по оси z ; % I_MID - Mean value of an inclination of a plane of orbit of a satellite I_MID = 1.0995574287564; %(double)(PI * 63.0 / 180.0) T_DR_MID = 43200.0; %Mean value a dragon of cycle time of a satellite A_PZ90_KM = 6378.136; %(Km) Equatorial radius of the Earth MU = 398600.44; % (Km^3/cek^2) constant of a gravitational OMEGA_Z = 0.7292115e-4; %скорость вращения Земли (Angular speed of rotation of the Earth, рад/cek) SEC_IN_RAD = 7.2722052166430e-005; %(PI / 43200.0 ) Number radian in second of time 2*PI/(24*3600)=PI/43200 HALF_PI = pi / 2; D7_3 = 7.0 / 3.0; D7_4 = 7.0 / 4.0; D7_6 = 7.0 / 6.0; D7_24 = 7.0 / 24.0; D49_72 = 49.0 / 72.0; nn = fix(ti_current / 86400); ti = ti_current - nn * 86400; % i_incl = I_MID + alm_gln(ns).deltai; i_incl = alm_gln(ns).deltai; ecc2_1 = 1.0 - alm_gln(ns).ecc * alm_gln(ns).ecc; % t_dr = T_DR_MID + alm_gln(ns).Tdr; t_dr = alm_gln(ns).Tdr; 120 n_dr = pi * 2 / t_dr; %1. Вычисление полуоси a_n: a_n = semi_axis_1( t_dr, i_incl, alm_gln(ns).ecc, alm_gln(ns).omegan); alm_gln(ns).a = a_n; %2. Вычисление t_lambda_k - время пересечения восходящего узла орбиты спутника sin_i = sin(i_incl); cos_i = cos(i_incl); sin_i2 = sin_i * sin_i; cos_i2 = cos_i * cos_i; v = - alm_gln(ns).omegan; % omega_n - angle of a perigee J = -0.00162393855 ; % J = 3/2 *C20; C20=1082.6257 * 10-6; ae_a = A_PZ90_KM / a_n; ae_a2 = ae_a * ae_a; j_ae_a2 = J * ae_a2; omega1 = j_ae_a2 * n_dr * cos_i / (ecc2_1 * ecc2_1); n0_na = (n0 - alm_gln(ns).Na); tz = ti - alm_gln(ns).TLambdaN + 86400.0 * n0_na; wk = tz / t_dr; wi = fix(wk); wi2 = wi * wi; t_lambda_kk = alm_gln(ns).TLambdaN + t_dr * wi + alm_gln(ns).dTdr * wi2; t_lambda_k = t_lambda_kk - n0_na * 86400.0 ; lambda_k = alm_gln(ns).LambdaN + (omega1 - OMEGA_Z) * (wi * t_dr + alm_gln(ns).dTdr * wi2) ; %time_s0 - a true sidereal time to Greenwich midnight of date n0 to which time ti concerns time_s = time_s0 * SEC_IN_RAD + OMEGA_Z * t_lambda_k ; omega_0 = lambda_k + time_s; % Auxiliary values: d1 = 1.0 - 1.5 * sin_i2; j_ae_a2_d1 = j_ae_a2 * d1; j_ae_a2_d = j_ae_a2 * (1.0 - 1.5 * sin_i); j_ae_a2_sin_i = j_ae_a2 * sin_i; j_ae_a2_sin_i2 = j_ae_a2 * sin_i2; j_ae_a2_cos_i2 = j_ae_a2 * cos_i2; % 3. Calculation of constants: tau = 0; l = cos(alm_gln(ns).omegan) * alm_gln(ns).ecc; h = alm_gln(ns).ecc * sin(alm_gln(ns).omegan); dop_x = (sqrt(1.0 - alm_gln(ns).ecc) * tan( v / 2.0)); dop_y = sqrt(1.0 + alm_gln(ns).ecc); 121 ee = 2.0 * atan2(dop_x, dop_y ); ma = ee - alm_gln(ns).ecc * sin(ee); ll = ma + alm_gln(ns).omegan; % for (j = 0; j <= 1; j++) for ( j= 1 : 2) sin_l = sin(ll); sin_2l = sin(2 * ll); sin_3l = sin(3 * ll); sin_4l = sin(4 * ll); cos_l = cos(ll); cos_2l = cos(2 * ll); cos_3l = cos(3 * ll); cos_4l = cos(4 * ll); l_cosl = l * cos_l; delta_a(j) = 2.0 * a_n * j_ae_a2_d1 * l_cosl +... j_ae_a2 * sin_i * (0.5 * h * sin_l - 0.5 * l_cosl +... cos(2.0 * alm_gln(ns).LambdaN) + 3.5 * l * cos_3l + 3.5 * h * sin_3l); delta_h(j) = j_ae_a2_d1 * (l * n_dr * tau + sin_3l +... 1.5 * l * sin_2l - 1.5 * h * cos_2l) -... 0.25 * j_ae_a2_sin_i2 * (sin_l - D7_3 * sin_3l +... 5.0 * l * sin_2l - 8.5 * l * sin_4l +... 8.5 * h * cos_4l + h * cos_2l) + j_ae_a2_cos_i2 *... (tau * l * n_dr - 0.5 * l * sin_2l); delta_l(j) = j_ae_a2_d1 * ( - tau * h * n_dr + cos_l +... 1.5 * l * cos_2l + 1.5 * h * sin_2l) -... 0.25 * j_ae_a2_sin_i2 * ( - cos_l - D7_3 * cos_3l -... 5.0 * h * sin_2l - 8.5 * l * cos_4l -... 8.5 * h * sin_4l + l * cos_2l) + j_ae_a2_cos_i2 *... ( - tau * h * n_dr + 0.5 * h * sin_2l); d_omega(j) = j_ae_a2 * cos_i * (tau * n_dr + 3.5 * l * sin_l -... 3.5 * h * cos_l - 0.5 * sin_2l - D7_6 * sin_3l + ... D7_6 * h * cos_3l); delta_i(j) = 0.5 * j_ae_a2_sin_i * cos_i * ( - l * cos_l +... h * sin_l + cos_2l + D7_3 * l * cos_3l + D7_3 * h * sin_3l); delta_ll(j) = 2.0 * j_ae_a2_d *... % 2.0 * j_ae_a2_d1 * (tau * n_dr + D7_4 * l * sin_l -... D7_4 * h * cos_l) + 3 * j_ae_a2_sin_i *... ( - D7_24 * h * cos_l - D7_24 * l * sin_l -... D49_72 * h * cos_3l + D49_72 * l * sin_3l +... 0.25 * sin_2l) + j_ae_a2 * cos_i *... 122 (tau * n_dr + 3.5 * l * sin_l - 2.5 * h * cos_l -... 0.5 * sin_2l - D7_6 * l * sin_3l + D7_6 * h * cos_3l); tau = ti - t_lambda_k; ll = ma + alm_gln(ns).omegan + n_dr * tau; end; % 4. Corrections to orbit elements of a satellite by second zone harmonic J20 % influence at the moments of time ti dda = delta_a(2) - delta_a(1); ddh = delta_h(2) - delta_h(1); ddl = delta_l(2) - delta_l(1); dd_omega = d_omega(2) - d_omega(1); ddi = delta_i(2) - delta_i(1); dd_ll = delta_ll(2) - delta_ll(1); % 5. calculation of influenced elements of orbits at the moment of ti time ai = a_n + dda; hi = h + ddh; li = l + ddl; % ecc_i = alm_gln[ns].ecc; ecc_i = sqrt(hi * hi + li * li); if (ecc_i == 0) w_i = 0; else if (li == ecc_i) w_i = HALF_PI; else if (li == - ecc_i) w_i = - HALF_PI; else if (li = 0) w_i = atan2(hi,li); else w_i = HALF_PI; end; end; end; end; omega_i = omega_0 + dd_omega; ii_incl = i_incl + ddi; ll_z = ma + alm_gln(ns).omegan + n_dr * tau + dd_ll; 123 % ll_z = ma + alm_gln(ns).omegan + n_dr * (ti - t_lambda_k) + dd_ll; mm_i = ll_z - w_i; % 6. Coordinates and components of satellite velocity vector % in absolute coordinate system OXaYaZa at the moment of ti time : % Ei(n) = Mi + ei * sin Ei(n-1), Ei(0) = Mi, | Ei(n) - Ei(n-1) | <= 10-8, eps_kepler = 1E-8; ee_i = kepler( mm_i, ecc_i, eps_kepler);% ee_i - эксцентрическая аномалия % alm_mca(ns).ek = ee_i; v_i = 2.0 * atan2(sqrt(1.0 + ecc_i) * tan(ee_i / 2.0) , sqrt(1.0 - ecc_i) ); u_i = v_i + w_i; ri = ai * (1.0 - ecc_i * cos(ee_i)) ; ecc_i2 = 1.0 - ecc_i * ecc_i; sqrt_mu_ai = sqrt( MU / ai); vr_i = sqrt_mu_ai * (ecc_i - sin(v_i) ) * ecc_i2; vu_i = sqrt_mu_ai * (1.0 + ecc_i * cos(v_i) ) * ecc_i2; cos_ui = cos(u_i); sin_ui = sin(u_i); cos_ii = cos(ii_incl); sin_ii = sin(ii_incl); cos_omega_i = cos(omega_i); sin_omega_i = sin(omega_i); x_dop = (cos_ui * cos_omega_i - sin_ui * sin_omega_i * cos_ii); y_dop = (cos_ui * sin_omega_i + sin_ui * cos_omega_i * cos_ii); satpos_xyz_gln.x = ri * x_dop; satpos_xyz_gln.y = ri * y_dop; satpos_xyz_gln.z = ri * sin_ui * sin_ii; satpos_xyz_gln.vx = vr_i * x_dop - vu_i * (sin_ui * cos_omega_i +... cos_ui * sin_omega_i * cos_ii); satpos_xyz_gln.vy = vr_i * y_dop - vu_i * (sin_ui * sin_omega_i -... cos_ui * cos_omega_i * cos_ii); satpos_xyz_gln.vz = vr_i * sin_ui * sin_ii + vu_i * cos_ui * sin_ii; Функция satpos_eci_in_metr function [satpos_eci, satpos_gln] = satpos_eci_in_metr(satpos_eci); %Имя функции: satpos_eci_in_metr %Функция преобразует координаты satpos_eci (структура), заданные в километрах в координаты %satpos_eci, satpos_gln в метрах satpos_eci.x = satpos_eci.x * 1000.0; satpos_eci.y = satpos_eci.y * 1000.0; 124 satpos_eci.z = satpos_eci.z * 1000.0; satpos_eci.vx = satpos_eci.vx * 1000.0; satpos_eci.vy = satpos_eci.vy * 1000.0; satpos_eci.vz = satpos_eci.vz * 1000.0; satpos_gln.x = satpos_eci.x; satpos_gln.y = satpos_eci.y; satpos_gln.z = satpos_eci.z; Функция init_satpos_gln function [sat_pos] = init_satpos_gln() ; %Имя функции:init_satpos_gln %Функция предназначена для инициализации структуры sat_pos sat_pos.prn = 0; sat_pos.x = -1.0; sat_pos.y = -1.0; sat_pos.z = -1.0; sat_pos.vx = -1.0; sat_pos.vy = -1.0; sat_pos.vz = -1.0; sat_pos.dvx = -1.0; sat_pos.dvy = -1.0; sat_pos.dvz = -1.0; sat_pos.range_m = -1.0; Функция eci_to_ecef function [satpos_ecef] =eci_to_ecef(s0, ti, satpos_eci) %Имя функции:eci_to_ecef %Функция преобразования координат %Входные данные: s0 - истинное звездное время в текущий момент обсервации , %ti - текущее время; satpos_eci %Структура satpos_eci %satpos_eci.x - координата x в абсолютной неподвижной системе координат (ECI); %satpos_eci.y - координата y в абсолютной неподвижной системе координат (ECI); %satpos_eci.z - координата z в абсолютной неподвижной системе координат (ECI); %satpos_eci.vx - скорость vx в абсолютной неподвижной системе координат (ECI); %satpos_eci.vy - скорость vy в абсолютной неподвижной системе координат (ECI); % satpos_eci.vz - скорость vz в абсолютной неподвижной системе координат (ECI); %Выходные данные: % Структура satpos_ecef %satpos_ecef.x - координата x в подвижной системе координат (ECEF); %satpos_ecef.y - координата y в подвижной системе координат (ECEF); %satpos_ecef.z - координата z в подвижной системе координат (ECEF); 125 %satpos_ecef.vx - скорость по оси x в подвижной системе координат (ECEF); %satpos_ecef.vy - скорость по оси z в подвижной системе координат (ECEF); %satpos_ecef.vz- скорость по оси z в подвижной системе координат (ECEF); %Коэффициенты % SEC_IN_RAD - коэффициент преобразования секунд в радианы % s0(radian) = s0 (sek) * SEC_IN_RAD, where % SEC_IN_RAD = 2 * pi / (24 * 3600) = pi / 43200 SEC_IN_RAD = pi / 43200; OMEGA_Z = 0.7292115e-4; %( скорость вращения Земли (angular speed of rotation of the Earth, рад/cek) s_zv = s0 * SEC_IN_RAD + OMEGA_Z * ti; cos_s = cos(s_zv); sin_s = sin(s_zv); satpos_ecef.x = satpos_eci.x * cos_s + satpos_eci.y * sin_s; satpos_ecef.y = -satpos_eci.x * sin_s + satpos_eci.y * cos_s; satpos_ecef.z = satpos_eci.z; satpos_ecef.vx = satpos_eci.vx * cos_s + satpos_eci.vy * sin_s + OMEGA_Z * satpos_ecef.y; satpos_ecef.vy = -satpos_eci.vx * sin_s + satpos_eci.vy * cos_s - OMEGA_Z * satpos_ecef.x; satpos_ecef.vz = satpos_eci.vz; Функция top_coord function [satvis] = top_coord(prn, rec_llh, rec_xyz, nlo_xyz) % Имя функции: top_coord % Функция выполняет расчет топоцентрических координат объекта по заданным %географическим (долгота, широта, высота) и геоцентрическим (x, y, z) %координатам приемника, а также геоцентрическим координатам объекта(x, y, z) %Входные данные: %prn - номер спутника %структура % rec_llh.lat - широта (рад) приемника %rec_llh.lon - долгота (рад) приемника %rec_llh.h - высота (м) приемника %структура %прямоугольные геоцентрические координаты приемника (м) % rec_xyz.x %rec_xyz.y %rec_xyz.z %структура 126 % прямоугольные геоцентрические координаты объекта (м) %nlo_xyz.x %nlo_xyz.y %nlo_xyz.z %Выходные данные: %структура satvis % top.s - проекция вектора дальности на ось, направленную на Юг (South) % top.e - проекция вектора дальности на ось, направленную на Восток (East) % top.z - проекция вектора дальности на ось, направленную в Зенит % top.daln - дальность до объекта % top.az - угол азимута объекта % top.el - угол видимости объекта rx = nlo_xyz.x - rec_xyz.x; ry = nlo_xyz.y - rec_xyz.y; rz = nlo_xyz.z - rec_xyz.z; r_sat = sqrt( rx * rx + ry * ry + rz * rz); r_rec = sqrt((rec_xyz.x)^2 + (rec_xyz.y)^2+ (rec_xyz.z)^2); top.r = r_sat; rx1 = rx; ry1 = ry; rz1 = rz; sin_lat = sin(rec_llh.lat); cos_lat = cos(rec_llh.lat); sin_lon = sin(rec_llh.lon); cos_lon = cos(rec_llh.lon); % Projections of vector of range in topocentric coordinate system: top.e = -sin_lon * rx1 + cos_lon * ry1; top.s = cos_lon * sin_lat * rx1 + sin_lon * sin_lat * ry1 - cos_lat * rz1; top.z = cos_lat * cos_lon * rx1 + cos_lat * sin_lon * ry1 + sin_lat * rz1; % azimut: отсчет по часовой стрелке от оси направленной на Север (N or -S) (-top.s) eps = 10e-10; if ( (abs(top.e) < eps) || (abs(top.s) < eps)) top.az = 0.0; else top.az = atan2(top.e,-top.s); end; if (top.az < 0.0) top.az = top.az + pi * 2; end; % elevation: cos_el_top = (rec_xyz.x * rx + rec_xyz.y * ry + rec_xyz.z * rz) / (r_sat * r_rec); if ( cos_el_top >= 1.00 ) 127 el = 0.0; else if ( cos_el_top <= -1.00 ) el = pi; else el = acos(cos_el_top); end; end; top.el = pi / 2.0 - el; satvis.el = top.el; satvis.az = top.az; satvis.r = r_sat; Функция rewrite_satpos function [result] = rewrite_satpos(nom, satpos) %Имя функции:rewrite_satpos %Функция предназначена для перезаписи структуры satpos в массив структур result result(nom).x = satpos.x; result(nom).y = satpos.y; result(nom).z = satpos.z; result(nom).vx = satpos.vx; result(nom).vy = satpos.vy; result(nom).vz = satpos.vz; Функция utc_nut function nut = utc_nut(t) %Имя: utc_nut %Функция предназначена для расчета нутации %Входные данные: %t =6.023472005475702e+002; %Выходные данные: %nut - нутация R = 1296000; % ( 1r=360grad=1 296 000 cek) RAD_SEK_ANGL = pi/(3600*180); t2 = t * t; t3 = t2 * t; l = 485866.733 + (1325.0 * R + 715922.633) * t + 31.310 * t2 + 0.064 * t3;%1.034807679476340e+012 l1 = 1287099.804 + (99 * R + 1292581.224) * t - 0.577 * t2 - 0.012 * t3; f = 335778.877 + (1342 * R + 295263.137) * t - 13.257 * t2 + 0.011 * t3; dd = 1072261.307 + (1236 * R + 1105601.328) * t - 6.891 * t2 + 0.019 * t3; omega = 450160.280 - (5 * R + 482890.539) * t + 7.455 * t2 + 0.008 * t3; eps0 = 84381.448 - 46.8150 * t - 0.00059 * t2 + 0.001813 * t3; % eps0 = 84381.447996 - 46.8150 * t - 0.00059 * t2 + 0.001813 * t3; 128 eps_d = utc_nut_fi_eps(t, l, l1, f, dd, omega, 'd','e'); eps_k = utc_nut_fi_eps(t, l, l1, f, dd, omega, 'k','e'); eps = eps0 + eps_d + eps_k; cos_eps = cos(RAD_SEK_ANGL * eps ) / 15.0; d_fi = utc_nut_fi_eps(t, l, l1, f, dd, omega, 'd', 'f'); k_fi = utc_nut_fi_eps(t, l, l1, f, dd, omega, 'k', 'f'); nut1 = d_fi * cos_eps; nut2 = k_fi * cos_eps; % nut3 = 0.00264 * sin(omega) + 0.000063 * sin(2.0 * omega) nut3 = 0; nut1_dop = nut1; nut2_dop = nut2; nut3_dop = nut3; nut = nut1 + nut2 + nut3; Функция utc_nut_fi_ep function nut_fi_eps = utc_nut_fi_eps(t, l, l1, f, dd, omega, typ_nut, fi_eps) % Имя функции:utc_nut_fi_eps %Функция предназначена для расчета нутации по епсилон и фи %Входные данные: %t - значение параметра приведено в функции utc_nut; %l - значение параметра приведено в функции utc_nut; %l1 - значение параметра приведено в функции utc_nut; %f -значение параметра приведено в функции utc_nut; %dd - значение параметра приведено в функции utc_nut; % оmega - значение параметра приведено в функции utc_nut; %typ_nut - признак параметра приведен в функции utc_nut; %fi_eps - признак параметра приведен в функции utc_nut; %Выходные данные: %nut_fi_eps - значение нутации по по епсилон и фи % применить функцию koef [koef_id, koef_abd, koef_ik, koef_abk] = koef; RAD_SEK_ANGL = pi/(3600*180); if (typ_nut == 'd') n = 30; else n = 76; end; sum_a = 0; sum_b = 0; for i = 1 : n if (typ_nut == 'd') 129 s1 = koef_id(i,1) * l + koef_id(i,2) * l1 + koef_id(i,3) * f + koef_id(i,4) * dd + koef_id(i,5) * omega; if (fi_eps == 'f') a = koef_abd(i,1) * 1e-4; bt = koef_abd(i,2) * 1e-4; else a = koef_abd(i,3) * 1e-4; bt = koef_abd(i,4) * 1e-4; end; else s1 = koef_ik(i,1) * l + koef_ik(i,2) * l1 + koef_ik(i,3) * f + koef_ik(i,4) * dd + koef_ik(i,5) * omega; if (fi_eps == 'f') a = koef_abk(i,1) * 1e-4; bt = koef_abk(i,2) * 1e-4; else a = koef_abk(i,3) * 1e-4; bt = koef_abk(i,4) * 1e-4; end; end; if (fi_eps == 'f') sin_s1 = sin(RAD_SEK_ANGL * s1); sa = a * sin_s1; sb = bt * sin_s1; else cos_s1 = cos(RAD_SEK_ANGL * s1); sa = a * cos_s1; sb = bt * cos_s1; end; arg = RAD_SEK_ANGL * s1; sum_a = sum_a + sa; sum_b = sum_b + sb; end; nut_fi_eps = sum_a + sum_b * t; Функция koef function[koef_id, koef_abd, koef_ik, koef_abk] = koef % Имя функции:koef % Функция предназначена для ввода коэффициентов при расчете нутации koef_id = [ 0, 0, 0, 0, 1 ; % 1 0, 0, 0, 0, 2 ; % 2 -2, 0, 2, 0, 1 ; % 3 2, 0,-2, 0, 0 ; % 4 -2, 0, 2, 0, 2 ; % 5 130 1,-1, 0,-1, 0 ; % 6 0,-2, 2,-2, 1 ; % 7 2, 0,-2, 0, 1 ; % 8 0, 0, 2,-2, 2 ; % 9 0, -1, 0, 0, 0 ; % 10 % исправлено 0, 1, 2,-2, 2 ; % 11 0,-1, 2,-2, 2 ; % 12 0, 0, 2,-2, 1 ; % 13 -2, 0, 0,2, 0 ; % 14 % исправлено 0, 0, 2,-2, 0 ; % 15 0, 2, 0, 0, 0 ; % 16 0, 1, 0, 0, 1 ; % 17 0, 2, 2,-2, 2 ; % 18 0,-1, 0, 0, 1 ; % 19 -2, 0, 0, 2, 1 ; % 20 0,-1, 2,-2, 1 ; % 21 2, 0, 0,-2, 1 ; % 22 0, 1, 2,-2, 1 ; % 23 1, 0, 0,-1, 0 ; % 24 2, 1, 0,-2, 0 ; % 25 0, 0,-2, 2, 1 ; % 26 0, 1,-2, 2, 0 ; % 27 0, 1, 0, 0, 2 ; % 28 -1, 0, 0, 1, 1 ; % 29 0, 1, 2,-2, 0 ];% 30 koef_abd = [ -171996.0,-174.2, 92025.0, 8.9; % 1 2062.0, 0.2, -895.0, 0.5; % 2 46.0, 0.0, -24.0, 0.0; % 3 11.0, 0.0, 0.0, 0.0; % 4 -3.0, 0.0, 1.0, 0.0; % 5 -3.0, 0.0, 0.0, 0.0; % 6 -2.0, 0.0, 1.0, 0.0; % 7 1.0, 0.0, 0.0, 0.0; % 8 -13187.0, -1.6, 5736.0,-3.1; % 9 -1426.0, 3.4, 54.0,-0.1; % 10 % исправлено -517.0, 1.2, 224.0,-0.6; % 11 217.0, -0.5, -95.0, 0.3; % 12 129.0, 0.1, -70.0, 0.0; % 13 -48.0, 0.0, 1.0, 0.0; % 14 % исправлено -22.0, 0.0, 0.0, 0.0; % 15 17.0, -0.1, 0.0, 0.0; % 16 131 -15.0, 0.0, 9.0, 0.0; % 17 -16.0, 0.1, 7.0, 0.0; % 18 -12.0, 0.0, 6.0, 0.0; % 19 -6.0, 0.0, 3.0, 0.0; % 20 -5.0, 0.0, 3.0, 0.0; % 21 4.0, 0.0, -2.0, 0.0; % 22 4.0, 0.0, -2.0, 0.0; % 23 -4.0, 0.0, 0.0, 0.0; % 24 1.0, 0.0, 0.0, 0.0; % 25 1.0, 0.0, 0.0, 0.0; % 26 -1.0, 0.0, 0.0, 0.0; % 27 1.0, 0.0, 0.0, 0.0; % 28 1.0, 0.0, 0.0, 0.0; % 29 -1.0, 0.0, 0.0, 0.0]; % 30 koef_ik = [ 0, 0, 2, 0, 2; % 31 1, 0, 0, 0, 0; % 32 0, 0, 2, 0, 1; % 33 1, 0, 2, 0, 2; % 34 1, 0, 0,-2, 0; % 35 -1, 0, 2, 0, 2; % 36 0, 0, 0, 2, 0; % 37 1, 0, 0, 0, 1; % 38 -1, 0, 0, 0, 1; % 39 -1, 0, 2, 2, 2; % 40 1, 0, 2, 0, 1; % 41 0, 0, 2, 2, 2; % 42 2, 0, 0, 0, 0; % 43 1, 0, 2,-2, 2; % 44 2, 0, 2, 0, 2; % 45 0, 0, 2, 0, 0; % 46 -1, 0, 2, 0, 1; % 47 -1, 0, 0, 2, 1; % 48 1, 0, 0,-2, 1; % 49 -1, 0, 2, 2, 1; % 50 1, 1, 0,-2, 0; % 51 0, 1, 2, 0, 2; % 52 0,-1, 2, 0, 2; % 53 1, 0, 2, 2, 2; % 54 1, 0, 0, 2, 0; % 55 2, 0, 2,-2, 2; % 56 0, 0, 0, 2, 1; % 57 0, 0, 2, 2, 1; % 58 132 1, 0, 2,-2, 1; % 59 0, 0, 0,-2, 1; % 60 1,-1, 0, 0, 0; % 61 2, 0, 2, 0, 1; % 62 0, 1, 0,-2, 0; % 63 1, 0,-2, 0, 0; % 64 0, 0, 0, 1, 0; % 65 1, 1, 0, 0, 0; % 66 1, 0, 2, 0, 0; % 67 1,-1, 2, 0, 2; % 68 -1,-1, 2, 2, 2; % 69 -2, 0, 0, 0, 1; % 70 3, 0, 2, 0, 2; % 71 0,-1, 2, 2, 2; % 72 1, 1, 2, 0, 2; % 73 -1, 0, 2,-2, 1; % 74 2, 0, 0, 0, 1; % 75 1, 0, 0, 0, 2; % 76 3, 0, 0, 0, 0; % 77 0, 0, 2, 1, 2; % 78 -1, 0, 0, 0, 2; % 79 1, 0, 0,-4, 0; % 80 -2, 0, 2, 2, 2; % 81 -1, 0, 2, 4, 2; % 82 2, 0, 0,-4, 0; % 83 1, 1, 2,-2, 2; % 84 1, 0, 2, 2, 1; % 85 -2, 0, 2, 4, 2; % 86 -1, 0, 4, 0, 2; % 87 1,-1, 0,-2, 0; % 88 2, 0, 2,-2, 1; % 89 2, 0, 2, 2, 2; % 90 1, 0, 0, 2, 1; % 91 0, 0, 4,-2, 2; % 92 3, 0, 2,-2, 2; % 93 1, 0, 2,-2, 0; % 94 0, 1, 2, 0, 1; % 95 -1,-1, 0, 2, 1; % 96 0, 0,-2, 0, 1; % 97 0, 0, 2,-1, 2; % 98 0, 1, 0, 2, 0; % 99 1, 0,-2,-2, 0; % 100 133 0,-1, 2, 0, 1; % 101 1, 1, 0,-2, 1; % 102 1, 0,-2, 2, 0; % 103 2, 0, 0, 2, 0; % 104 0, 0, 2, 4, 2; % 105 0, 1, 0, 1, 0]; % 106 koef_abk = [-2274.0, -0.2, 977.0, -0.5; % 31 712.0, 0.1, -7.0, 0.0; % 32 -386.0, -0.4, 200.0, 0.0; % 33 -301.0, 0.0, 129.0, -0.1; % 34 -158.0, 0.0, -1.0, 0.0; % 35 123.0, 0.0, -53.0, 0.0; % 36 63.0, 0.0, -2.0, 0.0; % 37 63.0, 0.1, -33.0, 0.0; % 38 -58.0, -0.1, 32.0, 0.0; % 39 -59.0, 0.0, 26.0, 0.0; % 40 -51.0, 0.0, 27.0, 0.0; % 41 -38.0, 0.0, 16.0, 0.0; % 42 29.0, 0.0, -1.0, 0.0; % 43 29.0, 0.0, -12.0, 0.0; % 44 -31.0, 0.0, 13.0, 0.0; % 45 26.0, 0.0, -1.0, 0.0; % 46 21.0, 0.0, -10.0, 0.0; % 47 16.0, 0.0, -8.0, 0.0; % 48 -13.0, 0.0, 7.0, 0.0; % 49 -10.0, 0.0, 5.0, 0.0; % 50 -7.0, 0.0, 0.0, 0.0; % 51 7.0, 0.0, -3.0, 0.0; % 52 -7.0, 0.0, 3.0, 0.0; % 53 -8.0, 0.0, 3.0, 0.0; % 54 6.0, 0.0, 0.0, 0.0; % 55 6.0, 0.0, -3.0, 0.0; % 56 -6.0, 0.0, 3.0, 0.0; % 57 -7.0, 0.0, 3.0, 0.0; % 58 6.0, 0.0, -3.0, 0.0; % 59 -5.0, 0.0, 3.0, 0.0; % 60 5.0, 0.0, 0.0, 0.0; % 61 -5.0, 0.0, 3.0, 0.0; % 62 -4.0, 0.0, 0.0, 0.0; % 63 4.0, 0.0, 0.0, 0.0; % 64 -4.0, 0.0, 0.0, 0.0; % 65 -3.0, 0.0, 0.0, 0.0; % 66 134 3.0, 0.0, 0.0, 0.0; % 67 -3.0, 0.0, 1.0, 0.0; % 68 -3.0, 0.0, 1.0, 0.0; % 69 -2.0, 0.0, 1.0, 0.0; % 70 -3.0, 0.0, 1.0, 0.0; % 71 -3.0, 0.0, 1.0, 0.0; % 72 2.0, 0.0, -1.0, 0.0; % 73 -2.0, 0.0, 1.0, 0.0; % 74 2.0, 0.0, -1.0, 0.0; % 75 -2.0, 0.0, 1.0, 0.0; % 76 2.0, 0.0, 0.0, 0.0; % 77 2.0, 0.0, -1.0, 0.0; % 78 1.0, 0.0, -1.0, 0.0; % 79 -1.0, 0.0, 0.0, 0.0; % 80 1.0, 0.0, -1.0, 0.0; % 81 -2.0, 0.0, 1.0, 0.0; % 82 -1.0, 0.0, 0.0, 0.0; % 83 1.0, 0.0, -1.0, 0.0; % 84 -1.0, 0.0, 1.0, 0.0; % 85 -1.0, 0.0, 1.0, 0.0; % 86 1.0, 0.0, 0.0, 0.0; % 87 1.0, 0.0, 0.0, 0.0; % 88 1.0, 0.0, -1.0, 0.0; % 89 -1.0, 0.0, 0.0, 0.0; % 90 -1.0, 0.0, 0.0, 0.0; % 91 1.0, 0.0, 0.0, 0.0; % 92 1.0, 0.0, 0.0, 0.0; % 93 -1.0, 0.0, 0.0, 0.0; % 94 1.0, 0.0, 0.0, 0.0; % 95 1.0, 0.0, 0.0, 0.0; % 96 -1.0, 0.0, 0.0, 0.0; % 97 -1.0, 0.0, 0.0, 0.0; % 98 -1.0, 0.0, 0.0, 0.0; % 99 -1.0, 0.0, 0.0, 0.0; % 100 -1.0, 0.0, 0.0, 0.0; % 101 -1.0, 0.0, 0.0, 0.0; % 102 -1.0, 0.0, 0.0, 0.0; % 103 1.0, 0.0, 0.0, 0.0; % 104 -1.0, 0.0, 0.0, 0.0; % 105 1.0, 0.0, 0.0, 0.0]; % 106 Функция semi_axis_1 135 % double semi_axis(double t_dr, double i_incl, double ecc, double omega_n) function [a_npp] = semi_axis_1(t_dr, i_incl, ecc, omega_n); %Имя функции:semi_axis_1 %Функция вычисляет радиус орбиты спутника ГЛОНАСС в соответствии с интерфейсным контроль- ным %документом ГЛОНАСС %Входные данные: %t_dr - драконический период обращения спутника ГЛОНАСС (секунды) %i_incl - наклонение орбиты спутника ГЛОНАСС (радиан) %ecc - эксцентриситет %omega_n - аргумент перигея орбиты спутника ГЛОНАСС (радиан) %Выходные данные: %a_npp - большая полуось орбиты спутника ГЛОНАСС (километр) A_PZ90_KM = 6378.136; %(Km) Equatorial radius of the Earth J20 = -1082625.7e-9; % Factor at the second zone harmonicm % B_PZ90_KM = 6356.75136174571344; % AP_LAND (Km) Polar radius of the Earth */ %a_npp = 1; MU = 398600.44; % (Km^3/cek^2) constant of a gravitational epsilon = 1.0e-3; sin_i = sin(i_incl); sin_i2 = sin_i * sin_i; v = -omega_n;% omega_n - angle of a perigee ecc2_1 = 1.0 - ecc * ecc; b1 = 2.0 - 2.5 * sin_i2; b2 = sqrt(ecc2_1 * ecc2_1 * ecc2_1 ); b3 = 1.0 + ecc * cos(omega_n); b3 = b3 * b3; b4 = 1.0 + ecc * cos(v); b5 = b4 * b4 * b4 / ecc2_1; b = b1 * b2 / b3 + b5; t_ock = t_dr; tock_2pi = t_ock /(pi * 2); p1_3 = 1.0 / 3.0; a_dop = MU * tock_2pi * tock_2pi; a_n = a_dop^(1/3); % a_n = pow(a_dop, p1_3); nn = 0; dda = epsilon + 1; while ( (dda > epsilon) & (nn < 50) ) p = a_n * ecc2_1; % Focal parameter 136 ae_p = A_PZ90_KM / p; b0 = 1.0 + (1.5 * J20 * ae_p * ae_p) * b; t_ock = t_dr / b0; tock_2pi = t_ock / (pi * 2); a_dop = MU * tock_2pi * tock_2pi; a_npp = a_dop^(1/3); dda = abs(a_n - a_npp); a_n = a_npp; nn = nn + 1; end; Функция kepler function [Ek] = kepler(Mk, ecc, eps); %Имя функции:kepler %Функция предназначена для решения уравнения Кеплера % eps = 1.0E-15; y = ecc * sin(Mk); x1 = Mk; x = y; for k = 0 : 15 x2 = x1; x1 = x; y1 = y; y = Mk - (x - ecc * sin(x)); if (abs(y - y1) < eps) break end x = (x2 * y - x * y1) / (y - y1); end %k Ek = x; Функция map function map(N) %Имя функции:map %Применена функция MatLab для внесения в графики орбитального движения изображения Земли load('topo.mat','topo','topomap1'); [x,y,z] = sphere(50); cla reset %axis square off props.AmbientStrength = 0.1; 137 props.DiffuseStrength = 1; props.SpecularColorReflectance = .5; props.SpecularExponent = 20; props.SpecularStrength = 1; props.FaceColor= 'texture'; props.EdgeColor = 'none'; props.FaceLighting = 'phong'; props.Cdata = topo; surface(x*N,y*N,z*N,props); light('position',[-1 0 1]); light('position',[-1.5 0.5 -0.5], 'color', [.6 .2 .2]); view(3) ; grid on 4.3.5 Примеры выполнения комплекса программ ORBITA_GL_NAVIOR Некоторые результаты, иллюстрирующие работу программ, изображены на рис. 4.8 - рис. 4.11. Рис. 4.8. Орбиты спутников ГЛОНАСС в системе координат ECEF 138 Рис. 4.9. Орбиты спутников ГЛОНАСС в системе координат ECI Рис. 4.10. Орбита спутника ГЛОНАСС № 24 за 7 суток 139 0 24 68 10 12 14 16 18 20 22 24 02 46 810 12 14 16 18 20 22 24Время наблюдения спутников ГЛОНАСС, часНо ме р нав игац ионн ого сп ут ник аРис. 4.11. Видимость спутников ГЛОНАСС Программы имеют и другие возможности в зависимости от задач, которые может поставить специалист, исследующий орбитальное движение спутников. 140 1411   ...   4   5   6   7   8   9   10   11   ...   14

Функция gln_a function [ ] = gln_a(A_PZ90_KM, ns, n0, ti_current, time_s0, alm_gln, timeUTC_leap, fid); %Имя функции:gln_a %Функция формирует альманах ГЛОНАСС в формате YUMA. Данные альманаха ГЛОНАСС соот- ветствуют %данным альманаха ГЛОНАСС, записанным приемником СН 4701 и преобразованным в строгом со- ответствии %интерфейсному контрольному документу ГЛОНАСС %Постоянные: I_MID = 1.0995574287564; %(double)(PI * 63.0 / 180.0). I_MID - Mean value of an inclination of a plane of orbit of a satellite T_DR_MID = 43200.0; %Mean value a dragon of cycle time of a satellite MU = 398600.44; % (Km^3/cek^2) гравитационная постоянная (constant of a gravitational) OMEGA_Z = 0.7292115e-4; %(рад/cek)- Угловая скорость вращения Земли (Angular speed of rotation of the Earth ) SEC_IN_RAD = 7.2722052166430e-005; %(PI / 43200.0 ) Number radian in second of time 2*PI/(24*3600)=PI/43200 HALF_PI = pi / 2; D7_3 = 7.0 / 3.0; D7_4 = 7.0 / 4.0; D7_6 = 7.0 / 6.0; D7_24 = 7.0 / 24.0; D49_72 = 49.0 / 72.0; nn = fix(ti_current / 86400); ti = ti_current - nn * 86400;%время в секундах i_incl = I_MID + alm_gln(ns).deltai; %наклонение орбиты спутника ГЛОНАСС с номером ns 156 ecc2_1 = 1.0 - alm_gln(ns).ecc * alm_gln(ns).ecc; % t_dr = T_DR_MID + alm_gln(ns).Tdr; % GG24 t_dr = alm_gln(ns).Tdr; %драконический период спутника ГЛОНАСС (данные приемника СН 4701) n_dr = pi * 2 / t_dr; %1. Вычисление полуоси a_n: a_n = semi_axis( A_PZ90_KM, t_dr, i_incl, alm_gln(ns).ecc, alm_gln(ns).omegan); alm_gln(ns).a = a_n; %2. Вычисление t_lambda_k - время пересечения восходящего узла орбиты спутника ГЛОНАСС sin_i = sin(i_incl); cos_i = cos(i_incl); sin_i2 = sin_i * sin_i; cos_i2 = cos_i * cos_i; v = - alm_gln(ns).omegan; % omega_n - angle of a perigee J = -0.00162393855 ; % J = 3/2 *C20; C20=1082.6257 * 10-6; ae_a = A_PZ90_KM / a_n; ae_a2 = ae_a * ae_a; j_ae_a2 = J * ae_a2; omega1 = j_ae_a2 * n_dr * cos_i / (ecc2_1 * ecc2_1); n0_na = (n0 - alm_gln(ns).Na) ; tz = ti - alm_gln(ns).TLambdaN + 86400.0 * n0_na; wk = tz / t_dr; wi = fix(wk); wi2 = wi * wi; t_lambda_kk = alm_gln(ns).TLambdaN + t_dr * wi + alm_gln(ns).dTdr * wi2; t_lambda_k = t_lambda_kk - n0_na * 86400.0 ; lambda_k = alm_gln(ns).LambdaN + (omega1 - OMEGA_Z) * (wi * t_dr + alm_gln(ns).dTdr * wi2) ; %time_s0 - a true sidereal time to Greenwich midnight of date n0 to which time ti concerns %time_s0 - истинное звездное время в текущий момент обсервации time_s0 = 0; time_s = time_s0 * SEC_IN_RAD + OMEGA_Z * t_lambda_k ; omega_0 = lambda_k + time_s; % Auxiliary values: d1 = 1.0 - 1.5 * sin_i2; j_ae_a2_d1 = j_ae_a2 * d1; j_ae_a2_d = j_ae_a2 * (1.0 - 1.5 * sin_i); j_ae_a2_sin_i = j_ae_a2 * sin_i; j_ae_a2_sin_i2 = j_ae_a2 * sin_i2; j_ae_a2_cos_i2 = j_ae_a2 * cos_i2; % 3. Calculation of constants: tau = 0; 157 l = cos(alm_gln(ns).omegan) * alm_gln(ns).ecc; h = alm_gln(ns).ecc * sin(alm_gln(ns).omegan); dop_x = (sqrt(1.0 - alm_gln(ns).ecc) * tan( v / 2.0)); dop_y = sqrt(1.0 + alm_gln(ns).ecc); ee = 2.0 * atan2(dop_x, dop_y ); ma = ee - alm_gln(ns).ecc * sin(ee); ll = ma + alm_gln(ns).omegan; % for (j = 0; j <= 1; j++) for ( j= 1 : 2) sin_l = sin(ll); sin_2l = sin(2 * ll); sin_3l = sin(3 * ll); sin_4l = sin(4 * ll); cos_l = cos(ll); cos_2l = cos(2 * ll); cos_3l = cos(3 * ll); cos_4l = cos(4 * ll); l_cosl = l * cos_l; delta_a(j) = 2.0 * a_n * j_ae_a2_d1 * l_cosl +... j_ae_a2 * sin_i * (0.5 * h * sin_l - 0.5 * l_cosl +... cos(2.0 * alm_gln(ns).LambdaN) + 3.5 * l * cos_3l + 3.5 * h * sin_3l); delta_h(j) = j_ae_a2_d1 * (l * n_dr * tau + sin_3l +... 1.5 * l * sin_2l - 1.5 * h * cos_2l) -... 0.25 * j_ae_a2_sin_i2 * (sin_l - D7_3 * sin_3l +... 5.0 * l * sin_2l - 8.5 * l * sin_4l +... 8.5 * h * cos_4l + h * cos_2l) + j_ae_a2_cos_i2 *... (tau * l * n_dr - 0.5 * l * sin_2l); delta_l(j) = j_ae_a2_d1 * ( - tau * h * n_dr + cos_l +... 1.5 * l * cos_2l + 1.5 * h * sin_2l) -... 0.25 * j_ae_a2_sin_i2 * ( - cos_l - D7_3 * cos_3l -... 5.0 * h * sin_2l - 8.5 * l * cos_4l -... 8.5 * h * sin_4l + l * cos_2l) + j_ae_a2_cos_i2 *... ( - tau * h * n_dr + 0.5 * h * sin_2l); d_omega(j) = j_ae_a2 * cos_i * (tau * n_dr + 3.5 * l * sin_l -... 3.5 * h * cos_l - 0.5 * sin_2l - D7_6 * sin_3l + ... D7_6 * h * cos_3l); delta_i(j) = 0.5 * j_ae_a2_sin_i * cos_i * ( - l * cos_l +... h * sin_l + cos_2l + D7_3 * l * cos_3l + D7_3 * h * sin_3l); delta_ll(j) = 2.0 * j_ae_a2_d *... % 2.0 * j_ae_a2_d1 * (tau * n_dr + D7_4 * l * sin_l -... 158 D7_4 * h * cos_l) + 3 * j_ae_a2_sin_i *... ( - D7_24 * h * cos_l - D7_24 * l * sin_l -... D49_72 * h * cos_3l + D49_72 * l * sin_3l +... 0.25 * sin_2l) + j_ae_a2 * cos_i *... (tau * n_dr + 3.5 * l * sin_l - 2.5 * h * cos_l -... 0.5 * sin_2l - D7_6 * l * sin_3l + D7_6 * h * cos_3l); tau = ti - t_lambda_k; ll = ma + alm_gln(ns).omegan + n_dr * tau; end; % 4. Corrections to orbit elements of a satellite by second zone harmonic J20 % influence at the moments of time ti dda = delta_a(2) - delta_a(1); ddh = delta_h(2) - delta_h(1); ddl = delta_l(2) - delta_l(1); dd_omega = d_omega(2) - d_omega(1); ddi = delta_i(2) - delta_i(1); dd_ll = delta_ll(2) - delta_ll(1); % 5. calculation of influenced elements of orbits at the moment of ti time ai = a_n + dda; hi = h + ddh; li = l + ddl; % ecc_i = alm_gln[ns].ecc; ecc_i = sqrt(hi * hi + li * li); if (ecc_i == 0) w_i = 0; else if (li == ecc_i) w_i = HALF_PI; else if (li == - ecc_i) w_i = - HALF_PI; else if (li = 0) w_i = atan2(hi,li); else w_i = HALF_PI; end; end; end; end; 159 omega_i = omega_0 + dd_omega; ii_incl = i_incl + ddi; ll_z = ma + alm_gln(ns).omegan + n_dr * tau + dd_ll; % ll_z = ma + alm_gln(ns).omegan + n_dr * (ti - t_lambda_k) + dd_ll; mm_i = ll_z - w_i; %-------------------------------------------------- %timeUTC_leap Na = alm_gln(ns).Na; [wGPS,weekday] = week_GLONAS_gps(Na, timeUTC_leap); Yuma(ns).week = wGPS;% - 1024 ;%варианты записи недели Yuma(ns).id = alm_gln(ns).ID + 37; %Запись альманаха ГЛОНАСС в формате YUMA %Заголовок альманаха ГЛОНАСС. В заголовок вынесен номер (литера) частоты- Hn fprintf(fid,'GLO*Week %i almanac for prn-%i NAUHn-%i***** \n',Yuma(ns).week, Yuma(ns).id,alm_gln(ns).Hn); %Номер спутника ГЛОНАСС 1...24. В заголовке номера спутников ГЛОНАСС 38...61. fprintf(fid,'ID: %i\n',Yuma(ns).id-37);% для Planning %fprintf(fid,'ID: %i\n',Yuma(ns).id);% вариант записи номера спутника %Здоровье спутников ГЛОНАСС: 000- здоров Yuma(ns).Health = alm_gln (ns).Health ; if (Yuma(ns).Health > 100) fprintf(fid,'Health: %i\n', Yuma(ns).Health); else if (Yuma(ns).Health > 10) fprintf(fid,'Health: 0%i\n', Yuma(ns).Health); else fprintf(fid,'Health: 00%i\n', Yuma(ns).Health); end; end; %Ексцентриситет орбиты спутника strdop = e_norm( alm_gln(ns).ecc, 9); fprintf(fid,'Eccentricity: %s \n', strdop); %Время, к которому относятся данные орбита спутника ГЛОНАСС, приведенное к началу недели GPS. %При необходимости нужно учесть поправку смещения времени GPS от UTC (в 2006 году- 14 се- кунд) strdop = e_norm((alm_gln(ns).TLambdaN + 86400 * (weekday)),9); fprintf(fid,'Time of Applicability(s): %s \n',strdop); %Наклонение орбиты спутника ГЛОНАСС Yuma(ns).inc = alm_gln(ns).deltai; 160 fprintf(fid,'Orbital Incluation(rad): %11.9e \n',Yuma(ns).inc); %Скорость изменения восходящего узла орбиты спутника ГЛОНАСС strdop = e_norm(omega1,9); fprintf(fid,'Rate of Right Ascen(r/s): %s\n', strdop); %Корень квадратный из большой полуоси орбиты спутника ГЛОНАСС Yuma(ns).sqrt_a = sqrt(a_n * 1000); fprintf(fid,'SQRT(A) (m^1/2): %11.6f \n',Yuma(ns).sqrt_a); %Долгота восходящего узла орбиты спутника ГЛОНАСС strdop = e_norm(omega_i,9); fprintf(fid,'Right Ascen at Week(rad): %s\n', strdop); %Аргумент перигея орбиты спутника ГЛОНАСС strdop = e_norm(w_i,9); fprintf(fid,'Argument of Perigee(rad): %s\n',strdop); %Средняя аномалия спутника ГЛОНАСС strdop = e_norm(mm_i,9); fprintf(fid,'Mean Anom(rad): %s\n', strdop); %Коэффициенты полинома для учета поправок времени strdop = e_norm(alm_gln(ns).tau_n,9); fprintf(fid,'Af0(s): %s\n', strdop); Yuma(ns).Af1 = 0.0; fprintf(fid,'Af1(s/s): %11.9e\n', Yuma(ns).Af1); %Номер недели fprintf(fid,'week: %i \n\n',Yuma(ns).week); %Номер дня от начала недели, к которому относятся данные альманаха Yuma(ns).day = weekday; Функция week_GLONAS_gps function [wGPS,weekday] = week_GLONAS_gps(Na, timeUTC) % Имя функции: week_GLONAS_gps % Функция вычисляет : % Выходные данные: wGPS, weekday. %wGPS -номер GPS- недели, % weekday - день недели (0 - воскресенье, 1 - понедельник, 2 - вторник, 3 - среда, %4 - четверг, 5 - пятница, 6 - суббота. % Входные данные: %Na- номер суток, к которым относятся данные альманаха спутника ГЛОНАСС, отсчитываемые от %ближайшего високосного года % timesUTC -время UTC c начала текущей недели (секунда), %Используемые константы: % сдвиг времени между UTC и системным временем GPS на начало 2004 года %diff_UTC_GPS = 13;(до конца 2005 года) diff_UTC_GPS = 14;%(с 2006 года) 161 %количество дней в месяцах DnMon = [31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31]; % юлианский день начала отсчета недели GPS на ноль часов UTC c 5 на 6 % января 1980 года: wGPS_0 = 2444244.5; % Расчет номеров юлианского дня и дня года, от которого отсчитываются дни ГЛОНАСС [JD, day_year] = JD_data(timeUTC); %Сдвиг на ноль часов UTC на гривчском меридиане JD = JD - 0.5; %Расчет номера GPS неделе стандартной функцией MatLab "floor" wGPS = floor((JD + Na - wGPS_0) / 7); %Расчет дня недели стандартной функцией MatLab "mod" weekday = mod((JD+Na - wGPS_0) , 7); Функция semi_axis function [a_npp] = semi_axis(A_PZ90_KM, t_dr, i_incl, ecc, omega_n); %Имя функции:semi_axis %Функция вычисляет радиус орбиты спутника ГЛОНАСС в соответствии с интерфейсным контроль- ным %документом ГЛОНАСС % A_WGS84 = 6378.1370 ; %km %A_PZ90_KM = A_WGS84; % A_PZ90_KM = 6378.136; %(Km) Equatorial radius of the Earth J20 = -1082625.7e-9; % Factor at the second zone harmonicm % B_PZ90_KM = 6356.75136174571344; % AP_LAND (Km) Polar radius of the Earth */ %a_npp = 1; MU = 398600.44; % (Km^3/cek^2) constant of a gravitational epsilon = 1.0e-3; sin_i = sin(i_incl); sin_i2 = sin_i * sin_i; v = -omega_n;% omega_n - angle of a perigee ecc2_1 = 1.0 - ecc * ecc; b1 = 2.0 - 2.5 * sin_i2; b2 = sqrt(ecc2_1 * ecc2_1 * ecc2_1 ); b3 = 1.0 + ecc * cos(omega_n); b3 = b3 * b3; b4 = 1.0 + ecc * cos(v); b5 = b4 * b4 * b4 / ecc2_1; b = b1 * b2 / b3 + b5; t_ock = t_dr; 162 tock_2pi = t_ock /(pi * 2); p1_3 = 1.0 / 3.0; a_dop = MU * tock_2pi * tock_2pi; a_n = a_dop^(1/3); % a_n = pow(a_dop, p1_3); nn = 0; dda = epsilon + 1; while ( (dda > epsilon) & (nn < 50) ) p = a_n * ecc2_1; % Focal parameter ae_p = A_PZ90_KM / p; b0 = 1.0 + (1.5 * J20 * ae_p * ae_p) * b; t_ock = t_dr / b0; tock_2pi = t_ock / (pi * 2); a_dop = MU * tock_2pi * tock_2pi; a_npp = a_dop^(1/3); dda = abs(a_n - a_npp); a_n = a_npp; nn = nn + 1; end; Функция JD_epohi function jden = JD_epohi(epoha) %Имя: JD_epohi %Фукция JD_epohi(epoha) рассчитывает номер юлианского дня %до начала года (epoha) на 12h, 0 день, январь. %Входные данные: epoha, размерность-год %Выходные данные: % jden- номер юлианского дня на 12h, 0 день, январь ( размерность -дни) rk = mod(epoha,4); if ( rk == 0 ) rk = 1.0; else rk = 2.0 - rk * 0.25; end; n100 = floor(epoha / 100); n400 = floor(epoha / 400); jden = (4712 + epoha) * 365.25 + n400 - n100 + rk; % fprintf('epoha=%d rk=%f jden=%6.2f \n', epoha, rk, jden); 163 Функция JD_data function [JD, day_year] = JD_data(timeUTC) %Имя:JD_data % Функция JD_data(timeUTC) вычисляет : %JD - номер юлианского дня, day_year - номер дня года. %Входные данные: %timeUTC.year - год, % timeUTC.mon - месяц, % timeUTC.day - день. %Выходные данные: %JD - юлианский день, day_year- день от начала года. %Оригинальные функции: function jd0 = JD_epohi(year), %(смотри комментарий). % число дней от начала периода до 12ч. всемирного времени (полдень) % указанной даты ( по Гринвичу) %Входные данные для контрольного примера: %year = 2000; mon = 1; day = 1; %количество дней в месяцах DnMon = [31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31]; %Вычисление номера юлианского дня опорной эпохи jd0 = JD_epohi(timeUTC.year); %Учет високосного года nfebr = 0; if mod(timeUTC.year,4) == 0 nfebr = 1; end; %Расчетномера дня года k = 0; for i = 2 : timeUTC.mon k = k + DnMon(i - 1); if (i == 2) k = k + nfebr; end; end; day_year = k + timeUTC.day; %Расчет номера юлианского дня JD = jd0 + day_year; %2451545 - номер юлианского дня; 1- первый день года. 164 5.2 Декодирование данных альманаха спутников GPS 5.2.1 Краткие сведения из теории Спутниковые навигационные приемники обрабатывают данные навигационных спутников и вырабатывают информации в определенных форматах, как правило, в форма- тах приемников характерных для фирм-производителей. Рассмотрим программу преобра- зования данных альманаха, получаемых с одного из широко распространенных приемни- ков фирмы Novatel, ProPak-G2, выполненного на базе платы OEM 4 [10, 11].. Информация, вырабатываемая приемником, записывается в виде файла в бинарной форме и в форме кода ASCII. В программе декодирования альманаха используется бинар- ная форма. Декодируется сырой альманах (RAWALM, Message ID: 74) [10, 11]. Формат данных из оригинала [10, 11].имеет вид Таблица Формат данных альманаха RAWALM Field Field type Data Descrsption Format Binary Bytes Binary Offset 1 header Log header H 0 2 ref week Almanac reference week number Ulong 4 H 3 ref secs Almanac reference time (seconds.) Ulong 4 H+4 4 subframes Number of subframes to follow Ulong 4 H+8 5 svid SV ID (satellite vehicle ID)a UShort 2 H+12 6 data Subframe page data. Hex 30 H+14 7... Next subframe offset = H + 12 + (subframe x 32) variable xxxx 32-bit CRC (ASCII and Binary only) Hex 4 H + 12 + (32 x subframes) variable [CR][LF] Sentence terminator (ASCII only) - - - На CD-диске программа записана в папке ALM_ProPakG2 Цель лабораторной работы изучение процедур и правил декодирования данных навигационного приемника. 5.2.2 Лабораторная работа 5. 2 «Декодирование данных альманаха навигационных приемников на базе плат OEM- 4» Пакет программ для выполнения лабораторной работы расположен в папке ALM_ProPakG2 Для выполнения работы в качестве входных данных потребуется альма- нах спутников GPS в виде бинарного файла, который можно получить с навигационных приемников на базе плат OEM 4. Для тестирования программы такой альманах приводит- ся в приложении. 165 Рекомендуется следующий порядок (не обязательный) выполнения лабораторной работы. 1. Создайте папки ALM_ProPakG2_My, RAW_ALM_BIN_DAT, RAW_ALM_BIN_ALM. 2. Скопируйте в папку ALM_ProPakG2_My из ALM_ProPakG2 папку RAW_ALM_PRG и RAW_ALM_BIN_DAT, RAW_ALM_BIN_ALM. 3. Запишите в папку RAW_ALM_BIN_DAT бинарный файл альманаха, который требу- ется декодировать. 4. Откроете папку RAW_ALM_PRG и изучите файлы и функции по комментария и про- граммным процедурам. 5. В папке RAW_ALM_PRG откройтефайл Raw_Alm.m и в соответствующие строки запишите имя файла альманаха, который требуется декодировать и имя файла альма- наха, под которым запишется декодированный альманах в формате YUMA. 6. Задание. Сравните данные декодированного альманаха с содержанием альманаха YUMA на сайте7. Выполните файл Raw_Alm.m и в папке RAW_ALM_BIN_ALM прочитайте декодиро- ванный альманах. 5.2.3 Задание для самоподготовки В руководстве пользователя на приемники, выполненные на базе плат OEM 4 [10, 11]. Изучите соответствующие разделы. 5.2.4 Функции и файлы из папки RAW_ALM_PRG Файл Raw_Alm.m %Имя файла:Raw_Alm.m %Программа преобразует данные бинарного файла альманаха спутников GPS в %формат YUMA %Входные данные: %Бинарный файл "сырого" альманаха навигационных приемников ProPak G2, ProPak G2 plus, DL 4 %для приведенного примера имя файла (2006_05_09_16.31_G2_RAW).bin %Выходные данные: %Файл альманаха GPS в формате YUMA, %для приведенного примера имя файла (2006_05_09_16.31_G2_RAW_OK).alm %начало цикла для 63 спутников (в том числе и спутников-нагрузок) for i=1:63 %структура согласно формата альманаха YUMA и ICD GPS-200 C alm(i)= struct('SV_ID',0, 'week',0, 'health_0',0, 'health',0,'ecc',0,... 't0a',0, 'incl_angle',0, 'omega_dot',0,... ЧЧЧЧЧЧЧЧЧЧЧЧЧЧЧЧЧЧЧЧЧ.166 'sqrtA',0, 'omega0',0, 'omega', 0, 'M0',0,'af0',0, 'af1',0); end; %диск и последовательность папок, из которых считывается файл с записью альманаха в бинарном %формате. При необходимости можно указать любой другой путь к файлу name_in = 'E:\MATLAB7\work\RAW_ALM_BIN_DAT\(2006_05_09_16.31_G2_RAW).bin'; %открытие данных для чтения и записи fw = fopen('E:\MATLAB7\work\RAW_ALM_BIN_ALM\(2006_05_09_16.31_G2_RAW_OK).alm','Wt'); %формирование CRC-кода для проверки считываемых данных fid_crc = fopen(name_in,'rb'); fseek(fid_crc, 0, 'eof'); size1 = ftell(fid_crc)- 4; fseek(fid_crc, size1, 'bof'); CRC = fread(fid_crc,1,'uint32'); fprintf('CRC_hex = %X \n',CRC); fseek(fid_crc, 0, 'bof'); CRC_own = CalculateBlockCRC32(fid_crc, size1); fprintf('CRC_own_hex= %X \n',CRC_own); fclose(fid_crc); %Чтение "сырых" данных альманаха fid = fopen(name_in,'rb'); % Заголовок данных- Header Structure (Binary Message) Syn_c = fread(fid,3,'char'); %В OEM-4 'Char' Hear_d= fread(fid,1,'uchar');%В OEM-4 'Uchar' Message_ID = fread(fid,1,'uint16');%В OEM-4 'Long' Message_Type=fread(fid,1,'char'); %В OEM-4 'Char' Port_Address = fread(fid,1,'char'); %В OEM-4 'Char' Message_Length = fread(fid,1,'uint16'); %В OEM-4 'Ushort' Sequenc_e = fread(fid,1,'uint16'); %В OEM-4 'Ushort' Idle_Time = fread(fid,1,'char'); %В OEM-4 'Char' Time_Status = fread(fid,1,'char'); %В OEM-4 'Enum' week_1 = fread(fid,1,'uint16'); %В OEM-4 'Ushort' Time = fread(fid,1,'ulong'); %В OEM-4 'Double' N_ul = fread(fid,1,'uint'); d_fcc=fread(fid,1,'uint16'); d_fcc = dec2hex(d_fcc); e1 = fread(fid,1,'uint16'); week =fread(fid,1,'uint');% неделя GPS T0a= fread(fid,1,'ulong');% время привязки альманаха kol_str= fread(fid,1,'ulong'); %количество строк %Чтение данных по каждому спутнику for i=1:kol_str 167 Nsv= fread(fid,1,'ushort'); % номер спутника согласно табл.20-V ICD GPS-200 C TLM= fread(fid,3,'char'); TL_M = dec2hex(TLM,2); HOW= fread(fid,3,'char'); HO_W = dec2hex(HOW,2); e6= fread(fid,1,'char'); maska = bin2dec('11000000'); ID_dop = bitand(e6, maska); data_ID = bitshift(ID_dop, -6); maska = bin2dec('00111111'); SV_ID = bitand(e6, maska); ecc1 = fread(fid,1,'char'); ecc2 = fread(fid,1,'char'); ecc = (bitshift(ecc1, 8) + ecc2)*2^(-21);%ексцентриситет t0a= fread(fid,1,'char')*2^12 ;% время привязки альманаха d1= fread(fid,1,'bit8'); d2=fread(fid,1,'char'); incl_angle= ((bitshift(d1, 8) + d2)*2^(-19) + 0.3) * pi;%наклонение орбиты спутника d1= fread(fid,1,'bit8'); d2=fread(fid,1,'char'); omega_dot= (bitshift(d1, 8) + d2)*2^(-38) * pi;%скорость изменения восходящего узла maska = bin2dec('00111111'); %здоровье спутника GPS health_1 = fread(fid,1,'char'); if health_1 > 255 health_1 = health_1; end; health = bitand(health_1,maska); if health > 63 health = health; end; maska = bin2dec('11000000'); health_0 = bitand(health_1,maska); health_0 = bitshift(health_0, -6); if health_0 > 3 health_0 = health_0; end; d1 = fread(fid,1,'char'); d2 = fread(fid,1,'char'); 168 d3 = fread(fid,1,'char'); sqrtA=(bitshift(d1, 16) + bitshift(d2, 8) + d3)*2^(-11);%корень квадратный из большей полуоси орбиты d1 = fread(fid,1,'bit8'); sign2 = sign(d1); d2 = fread(fid,1,'char'); d3 = fread(fid,1,'char'); omega0 =((bitshift(d1, 16) + bitshift(d2, 8) + d3)*2^(-23)) * pi; %долгота восходящего узла орбиты спут- ника d1 = fread(fid,1,'bit8'); sign1 = sign(d1); d2 = fread(fid,1,'char'); d3 = fread(fid,1,'char'); omega =((bitshift(d1, 16) + bitshift(d2, 8) + d3)*2^(-23)) * pi; %аргумент перигея d1 = fread(fid,1,'bit8'); d2 = fread(fid,1,'char'); d3 = fread(fid,1,'char'); M0 =((bitshift(d1, 16) + bitshift(d2, 8) + d3)*2^(-23)) * pi;%средняя аномалия на время привязки t0a af0_1 = fread(fid,1,'char'); af0_1_bin = dec2bin(af0_1); d2 = fread(fid,1,'char'); d2bin = dec2bin(d2); d3 = fread(fid,1,'char'); d3bin = dec2bin(d3); maska = bin2dec('00011100'); d4 = bitand(d3,maska); af0_2 = bitshift(d4, -2); af0_dop = bitshift(af0_1, 3) + af0_2; koef = 2^(-20); af0 = var_sign(af0_dop,11,koef);%коэффициент коррекции времени af0 % af1: maska = bin2dec('11100000'); d4 = bitand(d3,maska); af1_2 = bitshift(d2, 3); af1_1 = bitshift(d4, -5); 169 af1_dop = bitor(af1_2, af1_1); koef = 2^(-38); af1 = var_sign(af1_dop,11,koef);%коэффициент коррекции времени af1 %-------------------------------------------------- if SV_ID > 0 alm(SV_ID)= struct('SV_ID',SV_ID, 'week',week, 'health_0', health_0, 'health',health, 'ecc',ecc,... 't0a',t0a, 'incl_angle',incl_angle, 'omega_dot',omega_dot,... 'sqrtA',sqrtA, 'omega0',omega0, 'omega', omega, 'M0',M0,... 'af0',af0, 'af1',af1); end; end; %Формирование вывода данных альманаха в формате YUMA for i=1:30 %63 if alm(i).SV_ID > 0 if alm(i).SV_ID < 10 fprintf(fw,'**** Week %i almaNAU for PRN-0%i **********\n',alm(i).week, alm(i).SV_ID); fprintf(fw,'ID: 0%i\n',alm(i).SV_ID); else fprintf(fw,'**** Week %i almaNAU for PRN-%i **********\n',alm(i).week, alm(i).SV_ID); fprintf(fw,'ID: %i\n',alm(i).SV_ID); end; if alm(i).health < 10 fprintf(fw,'Health: %i0%i\n', alm(i).health_0, alm(i).health); else fprintf(fw,'Health: %i%i\n', alm(i).health_0, alm(i).health); end; strdop = e_norm(alm(i).ecc, 10); fprintf(fw,'Eccentricity: %s\n', strdop); fprintf(fw,'Time of Applicability(s): %6.4f\n',alm(i).t0a); fprintf(fw,'Orbital Incluation(rad): %0.10f \n',alm(i).incl_angle); strdop = e_norm(alm(i).omega_dot, 10); fprintf(fw,'Rate of Right Ascen(r/s): %s\n', strdop); 170 fprintf(fw,'SQRT(A) (m^1/2): %4.7f \n',alm(i).sqrtA); strdop = e_norm(alm(i).omega0, 10); fprintf(fw,'Right Ascen at Week(rad): %s\n', strdop); if alm(i).omega < 0 fprintf(fw,'Argument of Perigee(rad): %1.10f \n',alm(i).omega); else fprintf(fw,'Argument of Perigee(rad): %1.10f \n',alm(i).omega); end; strdop = e_norm(alm(i).M0, 10); fprintf(fw,'Mean Anom(rad): %s\n', strdop); strdop = e_norm(alm(i).af0, 10); fprintf(fw,'Af0(s): %s\n', strdop); strdop = e_norm(alm(i).af1, 10); fprintf(fw,'Af1(s/s): %s\n', strdop); fprintf(fw,'week: %i \n',alm(i).week); fprintf(fw,' \n'); end; %if alm(i).SV_ID > 0 end; %i %CR_C = fread(fid,1,'uint32'); %CR_C = dec2hex(CR_C); fclose(fid);%закрытия файла для чтения fclose(fw);%закрытия файла для записи 1   ...   6   7   8   9   10   11   12   13   14

Цель лабораторной работы изучение орбитального движения навигационных спутников GPS, ГЛОНАСС, GALILEO с последующей визуализацией спутников на любой момент времени. 5.4.2 Лабораторная работа 5. 4 «Модель движения и визуализация спутников GPS, ГЛОНАСС, GALILEO » Комплекс программ для выполнения настоящей работы основан на пакете Vsion_GLONASS_GPS, вкотором файл Vision_GLONASS_GPS.m заменяется файлом Vision_GL_GPS_GAL.m, а входными данными является альманах спутников трех сис- тем- GPS, ГЛОНАСС, GALILEO приведенный в приложении. Рекомендуется следующий порядок выполнения лабораторной работы. 1. Создайте папку Vision_GLONASS_GPS_GALILEO _My. 2. Скопируйте в папку Vision_GLONASS_GPS_GALILEO _My изVsion_GLONASS_GPS функции и файлы, из папки Vi-sion_GLONASS_GPS_GALILEO файл Vision_GL_GPS_GAL.m, из приложения файл альманаха AlmGGG.yum. Комплекс готов к работе. Выполните задания 1, 2. 3. Задание 1 .Откройте файл Vision_GL_GPS_GAL.m, убедитесь в правильности введен- ных входных данных ( альманах , дата и время) исполните файл и результаты визуали- зации по каждой спутниковой системе внесите в отчет. 4. Задание 2. Измените в Vision_GL_GPS_GAL.m дату на 1 сутки, выполните файл, вне- сите результаты в отчет. Объясните причины, изменения местоположения спутников по каждой из систем. 5.4.3 Задание для самоподготовки Подготовьте сравнительную характеристику спутниковых радионавигационных систем GPS, ГЛОНАСС, GALILEO 186 5.4.4 Файл из папки Vision_GLONASS_GPS_GALILEO Файла Vision_GL_GPS_GAL.m clear all %Имя m-файла:Vision_GL_GPS_GAL.m %Программа рассчитывает углы видимости навигационных спутников GPS, ГЛОНАСС, GALLILEO %(азимута и места) спутников на заданный момент времени %Входные данные: %файл альманаха в формате Yuma,имя файла альманаха присваивается %переменной "Dat",например,Dat = 'имя файла альманаха'; %данные о начале отсчета "d2",d2='месяц/день/год';h=час;min=минута;s=секунда; %координаты позиции приемника -lat(широта в радианах),lon (долгота в радианах, %hr (высота в метрах); %шаг с каким будут рассчитываться параметры (step,секунды); %количество точек (L), в которых будут рассчитываться параметры орбит %L=12*3600/step,L читается так: количество часов (например,12) %число секунд в часе (3600) деленное на шаг (step) %Постоянные: %скорость вращения Земли OMEGAeDOT=7.2921151467e-005; %радиусы земного эллипсоида A_WGS84=6378137.0; B_WGS84=6356752.314; %константы mu=398600500000000; F_CONST = 4.442807633E-10; kt=1;%установка времени на титульной надписи графика, определяется параметрами d2'; h; min; s и j или L; %Задание цветов для графики color6(1:14) = ['k' 'b' 'g' 'r' 'c' 'm' 'r' ':' 'g' ':' 'b' ':' 'k' 'h']; %Входные данные Dat = 'AlmGGG.yum';%альманах спутников GPS, ГЛОНАСС, GALLILEO d2='12/21/2006'; h=14; min=52.0; s=58.0; %координаты точки, из которой наблюдаются спутники lat = 0.88032730015257;%50 градусов 26минут 20.54 секунд lon = 0.53109641675259;%30 градусов 25 минут 46.4995секунд hr=187.488;% метров X_label=['Широта' ':' num2str(lat) ';' 'долгота' ':' num2str(lon) ';' 'высота' ':' num2str(hr)]; 187 step=300;%шаг отсчета времени в секундах (300=5 минутам);шаг можно изменять step=0; %L=(24*3600) / step;% убрать % перед L для вывода таблицы улов видимости и азимута L=1;% установить % перед L для вывода таблицы улов видимости и азимута %Чтение альманаха [alm,max_kol] = Yuma_GPS_GLONAS_Alm(Dat); %max_kol = 32; nom = 1; i = 0; k = 0; i = 0; % for i = 1 : max_kol while ( i <91) id = alm(nom).ID; Health = alm(nom).Health; %fprintf('1: i=%i k=%i nom=%i id=%i Health=%i \n', i, k, nom, id, Health); if ( id > 0) Health = alm(nom).Health; if ( Health == 0) k = k + 1; nom_ns(k) = id; % fprintf('2: i=%i k=%i nom=%i id=%i Health=%i \n', i, k, nom, id, Health); nom = nom + 1; else nom = nom + 1; end; else nom = nom + 1; end; i = i + 1; end; % i kol = k; fprintf('kol=%i \n', kol); nom_ns % - номер навигационного спутника [Rx,Ry,Rz] = ECEFLLH(lon, lat,hr); %Rx=0;Ry=0;,Rz=0;%центр масс Земли %Начало отсчета текущего времени [week,modeweek,d,dweek,weeks]=Tim(d2,h,min,s); for j = 1% 0:L 188 % for j = 1:L % 0:L t( j )=weeks+step*(j); %-step; %t1(j) = t(j)/60; %изменение дискретности текущего времени %d_wn = (week - alm(i).Week); %d_wn = 0; for k = 1 : kol i = nom_ns(k) ; % input "i" !!! % fprintf('i=%i alm(i).Week= %i \n',i, alm(i).Week) d_wn=(modeweek - alm(i).Week);%если в альманахе не учтено 1024 tk = t(j) + d_wn * 604800 - alm(i).TOA; d_wn = abs(modeweek - alm(i).Week); dd = 302400.0 + d_wn * 604800; if ( ( alm(i).A05 > 0 ) & ( alm(i).Health == 0 ) ) while (abs(tk) > dd) if tk > dd tk = tk - 604800; else if tk < -dd tk = tk + 604800; end end % if end % while %Справочник по альманаху- цифра в скобках обозначает порядковый номер %параметра альманаха в формате YUMA %alm(ID).ID(1); alm(ID).Health(2); alm(ID).e(3); alm(ID).TOA(4); alm(ID).deltai(5); %alm(ID).OMEGADOT(6); alm(ID).A05(7); alm(ID).omega0(8); alm(ID).omega(9); %alm(ID).M0(10); alm(ID).Af0(11); alm(ID).Af1(12); alm(ID).Week(13); n0=sqrt((mu) / (alm(i).A05^6)); j2 = 1082.68E-6; re = (A_WGS84 + B_WGS84) / 2.; sin55 = sin(55.0 * pi / 180.0); dn = 1.5 * j2 * re * re / (alm(i).A05^4 ) * (1. - 1.5 * sin55 * sin55); %dn = 0; n=n0 * (1 + dn); 189 Mk = alm(i).M0 + n*tk; e=alm(i).e; %решение уравнения Кеплера eps = 1.0E-15; y = e * sin(Mk); x1 = Mk; x = y; for k = 0 : 15 % количество итераций x2 = x1; x1 = x; y1 = y; y = Mk - (x - e * sin(x)); if (abs(y - y1) < eps) break end x = (x2 * y - x * y1) / (y - y1); end % kepler Ek = x; deltr = F_CONST * alm(i).e * alm(i).A05 * sin(Ek); dt1 = alm(i).Af0 + alm(i).Af1 * tk + deltr; tk = tk - dt1; vd = 1. - alm(i).e * cos(Ek); nuk =atan2(sqrt(1-alm(i).e^2)*sin(Ek) / vd,(cos(Ek)-alm(i).e) / vd); Ek = acos((alm(i).e+cos(nuk))/(1+alm(i).e*cos(nuk))); Fk =nuk + alm(i).omega; uk =Fk; ik=alm(i).deltai; % alm(i) rk =(alm(i).A05^2)*(1.0-alm(i).e*cos(Ek)); xkk =rk*cos(uk); ykk =rk*sin(uk); OMEGAk =alm(i).omega0+(alm(i).OMEGADOT-OMEGAeDOT)*tk-OMEGAeDOT*alm(i).TOA; % fprintf('i=%i j=%i \n ',i, j); Xk(j,i) = xkk*cos(OMEGAk)-ykk*cos(ik)*sin(OMEGAk); Yk(j,i) = xkk*sin(OMEGAk)+ykk*cos(ik)*cos(OMEGAk); Zk(j,i) = ykk*sin(ik); %Дальности до спутников PR(j,i) = sqrt((Xk(j,i) - Rx)^2 + (Yk(j,i) - Ry)^2 + (Zk(j,i) - Rz)^2); 190 %Перевод в географическую систему если требуется %[lons,lats,hrs] = LLHECEF(Xk,Yk,Zk); %(Llon(j),Llat(j),Hhr(j)) = [lons,lats,hrs]; %расчет угла видимости спутника xls = Xk(j,i) - Rx; yls = Yk(j,i) - Ry; zls = Zk(j,i) - Rz; range1 = sqrt(xls*xls+yls*yls+zls*zls); if j>1 doppler(j-1) = (range1 - range2) * 5.2514 / step; end range2 = range1; P = sqrt(Rx * Rx + Ry * Ry + Rz * Rz); tdot = ( Rx*xls+Ry*yls+Rz*zls)/range1/P; xll = xls/range1; yll = yls/range1; zll = zls/range1; if tdot >= 1.00 b = 0.0; elseif tdot <= -1.00 b = pi; else b = acos( tdot); end satang = pi/2.0 - b; TT =satang; TT(j,i) =TT;%угол видимости спутников %расчет угла азимута спутников xn =-cos(lon)*sin(lat); yn =-sin(lon)*sin(lat); zn = cos(lat); xe =-sin(lon); ye = cos(lon); xaz = xe*xll + ye*yll; yaz = xn*xll + yn*yll + zn*zll; if (xaz == 0) or (yaz == 0) az(j)= 0; else az(j,i) = atan2(xaz,yaz); end 191 if az(j,i) < 0 az(j,i) = az(j,i) + pi*2; end AZ(j,i) =az(j,i) *180/pi;%угла азимута спутников в градусах EL(j,i) = TT(j,i) *180/pi;%угла видимости спутников в градусах % ПЕРЕСЧЕТ ВРЕМЕНИ A(j)=mod(t(j),86400); her(j)=floor(A(j)/3600); m(j)=floor((A(j)-her(j)*3600)/60); sek(j)=A(j)-her(j)*3600-m(j)*60; %Построение полярной системы координат if EL(j,i) < 0 elp = 180; else elp = (EL(j,i) - 90); end; azp = (AZ(j,i) + 90.0); rad = pi / 180; x0 = 0; y0 = 0; xt(j,i) = (elp * cos(azp * rad)); yt(j,i) = -(elp * sin(azp * rad)); % fprintf('i=%i j=%i \n' , i, j); end % i = ns end; % if ( alm(i).A05 > 0 ) end % j = time %ВНИМАНИЕ. Для вывода времени визуализации спутников на график установите kt t_itle=[d2 ' ' num2str(her(kt)) ':' num2str(m(kt)) ':' num2str(sek(kt))]; %X_label=['Широта' ':' num2str(lat) ';' 'долгота' ':' num2str(lon) ';' 'высота' ':' num2str(hr)]; %num2ctr(lat) %num2str(her(kt)) %X_label=['66' ':']; n = 6; max_n = max(nom_ns); n_end = mod(max(nom_ns),n); n_end = mod(kol, n); n2 = fix(kol / n) * n - n +1; 192 %Формирование таблицы вывода времени UTC (Time), GPS (Tgps в секундах), номера спутника (nns), % углов видимости и азимута от времени и номера спутника for i=1:n:kol fprintf(' Time Tgps '); for k=1: n nns = nom_ns(i+k-1); fprintf(' %2i ', nns); end; fprintf(' \n'); for j=1:L fprintf('%2i:%2i:%2i %i ',her(j),m(j),sek(j), t(j)); for k=1: n nns = nom_ns(i+k-1); fprintf('%6.1f *%6.1f ', EL(j,nns), AZ(j,nns) ); end; fprintf(' \n'); end ; % J=1:L if (i) == (n2) n = n_end; end; end% i hold on %Окружности уровней на круговой диаграмме видимости спутников k1 = 10; k2 = 30; k3 = 50; k4 = 70; k5 = 85; k6=90; n=0; for k=1:5:365 n=n+1; m1 = pi / 180; x(n)=cos(k*m1); y(n)=sin(k*m1); end; %График круговой диаграммы plot(k1*x(:),k1*y(:),'k:', k2*x(:),k2*y(:),'k:', k3*x(:),k3*y(:),'k:',k4*x(:),k4*y(:),'k:', k5*x(:),k5*y(:),'r', k6*x(:),k6*y(:),'r:'); text(5, 10,'80','FontSize',12,'FontName','TimesNewRoman'); 193 text(18, 23,'60','FontSize',12,'FontName','TimesNewRoman'); text(32, 37,'40','FontSize',12,'FontName','TimesNewRoman'); text(45, 50,'20','FontSize',12,'FontName','TimesNewRoman'); text(55, 60,'5','FontSize',12,'FontName','TimesNewRoman'); text(62, 67,'0','FontSize',12,'FontName','TimesNewRoman'); grid on; %Построение изображений видимых спутников на круговой диаграмме i=1; for k=1:kol i = nom_ns(k) ; if (i < 31) plot(xt(kt,i),yt(kt,i), 'Marker' ,'o','MarkerSize',20,'MarkerFaceColor','g' ) ; end; if ((i > 31) & ( i < 62)) plot(xt(kt,i),yt(kt,i), 'Marker' ,'o','MarkerSize',20,'MarkerFaceColor','r' ) ; end; if (i > 61) plot(xt(kt,i),yt(kt,i), 'Marker' ,'o','MarkerSize',20,'MarkerFaceColor','y' ) ; end; title(t_itle); xlabel(X_label,'FontSize',12,'FontName','TimesNewRoman') set(get(gcf,'CurrentAxes'),'FontSize',14,'FontName','TimesNewRoman') hold on str1 = num2str( i, 2); %{ if i<32 str1 = num2str( i, 2); else str1 = num2str( (i), 2); end %} text(xt(kt,i), yt(kt,i),str1,'FontSize',14,'FontName','TimesNewRoman','HorizontalAlignment','center' ); axis( [-100 100 -100 100]); %axis( [-90 90 -90 90]); end clear Результат выполнения программы изображен на рис. 5.5. 194 -100-50 050 100-100-80-60-40-20 020 40 60 80 100 80 60 40 20 50 12/21/2006 14:52:58Широта:0.88033;долгота:0.5311;высота:187.488 37 14 16 18 19 21 22 26 29 38 44 45 59 61 64 65 66 72 73 74 75 80 81 82 83 91Рис. 5.5. Видимость спутников GPS, ГЛОНАСС, GALILEO 5.4.5 Лабораторная работа 5. 5 «Орбиты спутников GPS, ГЛОНАСС, GALILEO » Цель лабораторной работы моделирование орбит спутников GPS, ГЛОНАСС, GALILEOДля одновременного расчета орбит спутников трех систем GPS, ГЛОНАСС, GALI-LEO пакет программ, сформированный в лабораторной работе 5. 4 может быть дополнен файлами ORBITA_GGG.m и map. m из папки ORBIT_GGG. Рекомендуется следующий порядок выполнения лабораторной работы. 1. Создайте папку ORBIT_GGG _My. 2. Из папки Vision_GLONASS_GPS_GALILEO скопируйте в папку ORBIT_GGG _My функции ECEFLLH, LLHECEF, Tim; файлы Yuma_GPS_GLONAS_Alm.m, AlmGGG.yum; из папки ORBIT_GGG файл ORBITA_GGG.m. 3. Запустите MatLab [7, 8]. 4. Обратитесь к папке ORBIT_GGG _My , откройте ее, изучите функции, файлы и про- граммные процедуры по комментариям и выполните задания 1, 2. 5. Задание 1. Откройте файл ORBITA_GGG.m. Сформируйте желаемую конфигурацию спутников, орбиты которых требуется рассчитать (см. комментарий к файлу). Устано- вите скорость вращения Земли, равную 0 (OMEGAeDOT=0). Выполните файл. Дайте описание, полученного графического изображения орбит. Результаты внесите в отчет. 6. Задание 2. В файле ORBITA_GGG.m сформируйте конфигурацию спутников: 1 спут- ник GPS, 1 спутник ГЛОНАСС, 1 спутник GALILEO. Установите 195 OMEGAeDOT=7.2921151467e-005. Выполните файл. Дайте описание, полученного графического изображения орбит. Результаты внесите в отчет. 5.4.6 Задания и вопросы для самоподготовки 1. Чем различаются орбиты спутников GPS, ГЛОНАСС, GALILEO? 2. По каким параметрам можно определить принадлежность спутника той или иной ор- бите? 3. Сформулируйте алгоритм определения дальности до спутника из точки наблюдения. 4. В какой системе координат рассчитываются углы видимости спутника? 5.4.7 Листинг файла ORBITA_GGG.m Файл ORBITA_GGG.m clear all %Имя m-файла:ORBITA_GGG.m %Программа рассчитывает орбиты спутников GPS, ГЛОНАСС, GALLILEO %Входные данные: %файл альманаха в формате Yuma,имя файла альманаха присваивается %переменной "Dat",например,Dat = 'имя файла альманаха'; %данные о начале отсчета "d2",d2='месяц/день/год';h=час;min=минута;s=секунда; %координаты позиции приемника -lat(широта в радианах),lon (долгота в радианах, %hr (высота в метрах); %шаг с каким будут рассчитываться параметры (step,секунды); %количество точек (L), в которых будут рассчитываться параметры орбит %L=12*3600/step,L читается так: количество часов (например,12) %число секунд в часе (3600) деленное на шаг (step) %Выходные данные: %Xk(j,i) - координата X спутника с номером i на момент времени j; %Yk(j,i) - координата Y спутника с номером i на момент времени j; %Zk(j,i) - координата Z спутника с номером i на момент времени j; %Примечание: входные и выходные данные могут быть дополнены / изменены закомментированны- ми %параметрами, помеченными строкой %%%%%%%%%%%%%%%%%%%%%%%; %Постоянные: %скорость вращения Земли %OMEGAeDOT=7.2921151467e-005; OMEGAeDOT=0;%%%%%%%%%%%%%%для орбитальных плоскостей %радиусы земного эллипсоида A_WGS84=6378137.0; B_WGS84=6356752.314; %константы mu=398600500000000; 196 F_CONST = 4.442807633E-10; %Задание цветов для графики j_color = 0; color6(1:11) = ['B' 'B' 'B' 'B' 'B' 'B' 'o' 'o' '+' '+' '+'];%%%%%%%%%%%%%%%%%%%%%%% %Входные данные Dat = 'AlmGGG.yum';%альманах спутников GPS, ГЛОНАСС, GALLILEO N=6378136;% радиус Земли (используется, как нормирующий коэффициент map(N);%функция выводит на графики Землю d2='12/21/2006'; %%%%%%%%%%%%%%%%%%%%%%% h=14; min=52.0; s=58.0;%%%%%%%%%%%%%%%%%%%%%%% %координаты точки, из которой наблюдаются спутники lat = 0.88032730015257;%50 градусов 26минут 20.54 секунд %%%%%%%%%%%%%%%%%%%%%% lon = 0.53109641675259;%30 градусов 25 минут 46.4995секунд %%%%%%%%%%%%%%%%%%%%% hr=187.488;% метров %%%%%%%%%%%%%%%%%%%%%%% step=600;%шаг отсчета времени в секундах (300=5 минутам);шаг можно изменять %step=0; L=(12*3600) / step;% убрать перед L для вывода таблицы улов видимости и азимута %%%%%%%%%% %L=1;% установить % перед L для вывода таблицы улов видимости и азимута %Чтение альманаха [alm,max_kol] = Yuma_GPS_GLONAS_Alm(Dat);%%%%%%%%%%%%%%%%%%%%%%% nom = 1; i = 0; k = 0; i = 0; while ( i <91) id = alm(nom).ID; Health = alm(nom).Health; if ( id > 0) Health = alm(nom).Health; if ( Health == 0) k = k + 1; nom_ns(k) = id; nom = nom + 1; else nom = nom + 1; end; else 197 nom = nom + 1; end; i = i + 1; end; % i kol = k; nom_ns ; % - номер навигационного спутника [Rx,Ry,Rz] = ECEFLLH(lon, lat,hr);%%%%%%%%%%%%%%%%%%%%%%% [week,modeweek,d,dweek,weeks]=Tim(d2,h,min,s); for j = 1:L % for j = 1:L % 0:L t( j )=weeks+step*(j); %-step;%%%%%%%%%%%%%%%%%%%%%%% %t1(j) = t(j)/60; %изменение дискретности текущего времени %d_wn = (week - alm(i).Week); %d_wn = 0; kol=11;%%%%%%%%%%%%%%%%%%%%%%% nom_ns(1:kol) =[1 3 4 5 9 10 38 58 62 72 82 ];%%%%%%%%%%%%%%%%%%%%%%% for k = 1 : kol i = nom_ns(k) ; % input "i" !!! d_wn=(modeweek - alm(i).Week);%если в альманахе не учтено 1024 tk = t(j) + d_wn * 604800 - alm(i).TOA; d_wn = abs(modeweek - alm(i).Week); dd = 302400.0 + d_wn * 604800; if ( ( alm(i).A05 > 0 ) & ( alm(i).Health == 0 ) ) while (abs(tk) > dd) if tk > dd tk = tk - 604800; else if tk < -dd tk = tk + 604800; end end % if end % while %Справочник по альманаху- цифра в скобках обозначает порядковый номер %параметра альманаха в формате YUMA %alm(ID).ID(1); alm(ID).Health(2); alm(ID).e(3); alm(ID).TOA(4); alm(ID).deltai(5); %alm(ID).OMEGADOT(6); alm(ID).A05(7); alm(ID).omega0(8); alm(ID).omega(9); %alm(ID).M0(10); alm(ID).Af0(11); alm(ID).Af1(12); alm(ID).Week(13); n0=sqrt((mu) / (alm(i).A05^6)); j2 = 1082.68E-6; re = (A_WGS84 + B_WGS84) / 2.; sin55 = sin(55.0 * pi / 180.0); dn = 1.5 * j2 * re * re / (alm(i).A05^4 ) * (1. - 1.5 * sin55 * sin55); 198 %dn = 0; n=n0 * (1 + dn); Mk = alm(i).M0 + n*tk; e=alm(i).e; %решение уравнения Кеплера eps = 1.0E-15; y = e * sin(Mk); x1 = Mk; x = y; for k = 0 : 15 % количество итераций x2 = x1; x1 = x; y1 = y; y = Mk - (x - e * sin(x)); if (abs(y - y1) < eps) break end x = (x2 * y - x * y1) / (y - y1); end % kepler Ek = x; deltr = F_CONST * alm(i).e * alm(i).A05 * sin(Ek); dt1 = alm(i).Af0 + alm(i).Af1 * tk + deltr; tk = tk - dt1; vd = 1. - alm(i).e * cos(Ek); nuk =atan2(sqrt(1-alm(i).e^2)*sin(Ek) / vd,(cos(Ek)-alm(i).e) / vd); Ek = acos((alm(i).e+cos(nuk))/(1+alm(i).e*cos(nuk))); Fk =nuk + alm(i).omega; uk =Fk; ik=alm(i).deltai; rk =(alm(i).A05^2)*(1.0-alm(i).e*cos(Ek)); xkk =rk*cos(uk); ykk =rk*sin(uk); OMEGAk =alm(i).omega0+(alm(i).OMEGADOT-OMEGAeDOT)*tk-OMEGAeDOT*alm(i).TOA; Xk(j,i) = xkk*cos(OMEGAk)-ykk*cos(ik)*sin(OMEGAk); Yk(j,i) = xkk*sin(OMEGAk)+ykk*cos(ik)*cos(OMEGAk); Zk(j,i) = ykk*sin(ik); %%%%%%%%%%%%%%%%%%%%%%% %{ %Дальности до спутников PR(j,i) = sqrt((Xk(j,i) - Rx)^2 + (Yk(j,i) - Ry)^2 + (Zk(j,i) - Rz)^2); %Перевод в географическую систему если требуется %[lons,lats,hrs] = LLHECEF(Xk,Yk,Zk); 199 %(Llon(j),Llat(j),Hhr(j)) = [lons,lats,hrs]; %расчет угла видимости спутника xls = Xk(j,i) - Rx; yls = Yk(j,i) - Ry; zls = Zk(j,i) - Rz; range1 = sqrt(xls*xls+yls*yls+zls*zls); if j>1 doppler(j-1) = (range1 - range2) * 5.2514 / step; end range2 = range1; P = sqrt(Rx * Rx + Ry * Ry + Rz * Rz); tdot = ( Rx*xls+Ry*yls+Rz*zls)/range1/P; xll = xls/range1; yll = yls/range1; zll = zls/range1; if tdot >= 1.00 b = 0.0; elseif tdot <= -1.00 b = pi; else b = acos( tdot); end satang = pi/2.0 - b; TT =satang; TT(j,i) =TT;%угол видимости спутников %расчет угла азимута спутников xn =-cos(lon)*sin(lat); yn =-sin(lon)*sin(lat); zn = cos(lat); xe =-sin(lon); ye = cos(lon); xaz = xe*xll + ye*yll; yaz = xn*xll + yn*yll + zn*zll; if (xaz == 0) or (yaz == 0) az(j)= 0; else az(j,i) = atan2(xaz,yaz); end if az(j,i) < 0 az(j,i) = az(j,i) + pi*2; 200 end AZ(j,i) =az(j,i) *180/pi;%угла азимута спутников в градусах EL(j,i) = TT(j,i) *180/pi;%угла видимости спутников в градусах % ПЕРЕСЧЕТ ВРЕМЕНИ A(j)=mod(t(j),86400); her(j)=floor(A(j)/3600); m(j)=floor((A(j)-her(j)*3600)/60); sek(j)=A(j)-her(j)*3600-m(j)*60; %Построение полярной системы координат if EL(j,i) < 0 elp = 180; else elp = (EL(j,i) - 90); end; azp = (AZ(j,i) + 90.0); rad = pi / 180; x0 = 0; y0 = 0; xt(j,i) = (elp * cos(azp * rad)); yt(j,i) = -(elp * sin(azp * rad)); %} %%%%%%%%%%%%%%%%%%%%%%%% end % i = ns end; % if ( alm(i).A05 > 0 ) end % j = time for (i=1:kol) j_color = j_color + 1; if (j_color > 11 ) j_color = 1; end S = color6(j_color); prn = nom_ns(i); %%%%%%%%%%%%%%%%%%%%%%% hold on h_F1 = gca; plot3(Xk(:,prn),Yk(:,prn),Zk(:,prn),S,'LineWidth',0.5); axis([ -2.552*10^(7) 2.552*10^(7) -2.552*10^(7) 2.552*10^(7) -2.552*10^(7) 2.552*10^(7)]); set(get(gcf,'CurrentAxes'),'FontSize',14,'FontName','TimesNewRoman'); set(h_F1,'Position',[0.1 0.1 0.85 0.9]) ; xlabel('Координата X') 201 ylabel('Координата Y'), zlabel('Координата Z'),grid on str1 = num2str( prn); text(Xk(j,prn), Yk(j,prn),Zk(j,prn),str1,'FontSize',14,'FontName','TimesNewRoman','HorizontalAlignment','center' ); hold on %%%%%%%%%%%%%%%%%%%%%%% end clear 5.4.8 Пример выполнения файла:ORBITA_GGG.m Результат выполнения программы изображен на рис. 5.6. Рис. 5.6. Орбиты спутников GPS, ГЛОНАСС, GALILEO На рис. 5.6 сплошными линиями показаны 6 орбит спутников GPS, маркером «о» 2 орбиты ГЛОНАСС, маркером «+» 3 орбиты GALILEO; цифрами обозначены номера спутников соответствующих систем. 5.5 Декодирование и расшифровка данных спутников ГЛОНАСС 5.5.1 Краткие сведения из теории В данном подразделе приводится программа декодирования данных спутников ГЛОНАСС. Данные, которые декодируются получены экспериментально с реального на-202 вигационного приемника после корреляционной обработки в виде последовательности символов «0» и «1». Поскольку объем данных относительной большой около 60000 сим- волов, то приводятся два вспомогательных m-файла. Эти вспомогательные программы по- зволяют сформировать данные на бумажном носителе в виде последовательности шестна- дцатеричных чисел, из которых можно восстановить m-файл для декодирования. Экспе- риментальные данные рассматриваемые далее представляют собой символы длительно- стью 10 миллисекунд и являются составляющими информационных символов сложенных по модулю 2 с меандром имеющим частоту 100 Гц. Примененный алгоритм декодирования основан на строгом соответствии данных, передаваемых с навигационного спутника ГЛОНАСС, которые соответствуют формату сообщений. Более подробно формат сообщений и данные, передаваемые в кадрах сообще- ний приводятся в интерфейсном контрольном документе ГЛОНАССи в книге [1] под- раздел 2. 2, стр. 77- 87). Комплекс программ декодирования и расшифровки находится в папке Decod_GLONASS. Цель лабораторной работы освоение и применение методики декодирования и расшифровки данных спутниковой системы ГЛОНАСС. 5.5.2 Лабораторная работа 5. 6 «Декодирование данных спутников ГЛОНАСС в на- вигационном приемнике» Рекомендуется следующий порядок выполнения лабораторной работы. 1. Создайте папку Decod_GLONASS _My и скопируете в ее файлы и функции из папки Decod_GLONASS. 2. Запустите MatLab. 3. Обратитесь к папке Decod_GLONASS _My , откройте ее, изучите функции, файлы и программные процедуры по комментариям и выполните задания 1, 2. 4. Задание 1. Откройте файл Hex_GLONASS.m. Изучите программные процедуры фай- ла, используя его текст и расширенный комментарий.Файл преобразования экспериментальных данных в шестнадцатеричные симво-лы Hex_GLONASS.m %Имя m- файла: Hex_GLONASS.m %Файл Hex_GLONASS.m предназначен для преобразования экспериментальных данных, записанных %в фйле с именем raw_bits_SV18J.dat в файл SV1.dat, записанный в шестнадцатеричных символах name_in = 'raw_bits_SV18J.dat'; % имя экспериментального файла fw = fopen('SV1.dat','Wt');%имя файла в шестнадцатеричных символах M=dlmread(name_in);%чтение экспериментального файла в матрицу M %чтение двоичных символов в 1875 строк по 16 символов в строке [ ](203 for i=1:1875 str_danD(i,:)=sprintf('%d',M(i+15*(i-1):(15+i+15*(i-1)))); end bin_dec=bin2dec(str_danD);%преобразование строк символов в десятичные числа dec_hex=dec2base(bin_dec,16)';%преобразование десятичных чисел в шестнадцатеричные fprintf(fw,'%c',dec_hex);%запись файла в шестнадцатеричном представлении fclose(fw);%закрытие записанного файла Расширенный комментарий к файлу Hex_GLONASS.m.При условии, что имеется файл данных спутника с именем raw_bits_SV18J.dat ре- зультат выполнения m- файла: Hex_GLONASS.m, записанный под именем SV1.dat и про- читанный клавишей F3 файла в Total Commander в шестнадцатеричном представлении имеет вид AAD34B52CAAB3534B52CAB5332ACD4AD3F1BA84B2D4AB4D4B54ACD54AAD555555555555555555552B2BF1BA84B2D34ACB4D5555555554B54AAAAAAAAAAAAAAAAACB33F1BA84B2CB52D3555334ACAAD4B55552B32CACAAB4AB2B34D3F1BA84B2CCB4B54D52AAAB2B2CAACB535554AB4ACD54B52CB3F1BA84B354AD3554ACCB2D2AAAACAAAD4CB2D35554ACAD4D2BF1BA84B352CCCCD4CD52D534D534CB535554D53532AD2AB34BF1BA84B34B52CD52B4CD52AD54B52AAD4CB34CAAD4B352D2B3F1BA84B34D32AAACCD4B5532CB5334ACAAAAD2D532AD34AACBF1BA84B32B5354AAB4CD2B2CAB552AAD534B35555532CD2D4BF1BA84B32CAAB2B2AAB32B55554CCB535554D535352ACAB4ABF1BA84B33555555555 5555555555555555555555555554D2ABF1BA84B332AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABF1BA84B2AD2B534B54B4B34CB5334CAB54D2CB2AD32D2D2B33F1BA84B2B554ACB554AB3 534B4D4D4AAACAB4ACAD4B54D32D3F1BA84B2B32AAAB55552CB4AD3554CACB4AD354ACCD532CB4BF1BA84B2D4AB4D4B54ACD54AAD555555555555555555552B2BF1BA84B2D34ACB4D5555555554B54AAAAAAAAAAAAAAAAACB33F1BA84B2CAAAAD55555555555555555555555555555554D4B3F1BA84B2CD555555555555555555555555555555555554D2B3F1BA84B355554AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAB534BF1BA84B352AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAB3F1BA84B3 4B554CACAAB4D2D4B4AD3552D4D32D54AD4B54B4ABF1BA84B34CCB34CD33552B2AAACD34AD3333 552D4CACAACCCBF1BA84B32AAAD555555555555555555555555555555554D4B3F1BA84B32D5555555555 555555555555555555555555552D53F1BA84B334AAD34B2AAB52D4CCB32AAD2B2D55554D2D2D54D3F1BA84B332CACB534AAB2CCACAAD34AD3333552D4D2AD4ACD3F1BA84B2AD2B534B534B4CB34ACCB35 4AB2D34D52CD2D332ABF1BA84B2B554ACB554AB3534B4D4D4AAACAB4ACAD4B54D32D3F1BA84B2B32AAAB55552CB4AD3554CACB4AD354ACCD532CB4BF1BA84B2D4AB4D4B54ACD54AAD555555555555 555555552B2BF1BA84B2D34ACB4D5555555554B54AAAAAAAAAAAAAAAAACB33F1BA84B2CAAACAAAD32AD3554D2CCAAD2B2ACAAAD4D2B33333F1BA84B2CD4CCAD32D4AACCAB4B2CB52CCCD334D4D2D4AB4CBF1BA84B355553354B4D552CCB54ACAACAAD333554B4CB4CAD3F1BA84B352D4B4B32D4D34B532D4CB52CCCD4B52B2D4CAB54BF1BA84B34AAB554B55354ACAD3553553554CACAACAACACD2ABF1BA84B34D4B2B32B52B4D55354D34AD3332AAD2B2D2AACD4BF1BA84B32AAB52AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAB554BF1BA84B32D5555555555555555555555555555555555552D53F1BA84B33554B5555555555555555555555555555555552D4BF1BA84B332AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABF1BA84B2AD2B534B4B4B4CB34ACCB354AB2D34D52CD2D354CBF1BA84B2B554ACB554AB3534B4D4D4AAACAB4ACAD4B54D32D3F1BA84B2B32AAAB55552CB4AD3554CACB4AD354ACCD532CB4BF1BA84B2D4AB4D4B54ACD54AAD555555555555555555552B2BF1BA84B2D34ACB4D5555555554B54AAAAAAAAAAAAAAAAACB33F1BA84B2CAAB4D55555555555555555555555555555555 2B53F1BA84B2CD555555555555555555555555555555555554D2B3F1BA84B35554D555555555555555555555 5555555555552D4BF1BA84B352AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAB3F1BA84B34AAB2D55555555555555555555555555555554D4B3F1BA84B34D555555555555555555555555555555555555 2D4BF1BA84B32AAB35555555555555555555555555555555552B53F1BA84B32D555555555555555555555555 5555555555552D53F1BA84B33554CD55555555555555555555555555555554D2ABF1BA84B332AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABF1BA84B2AD2B534B4CB4B34CB5334CAB54D2CB2AD32D2D4D53F1BA84B2B554ACB554AB3534B4D4D4AAACAB4ACAD4B54D32D3F1BA84B2B32AAAB555 52CB4AD3554CACB4AD354ACCD532CB4BF1BA84B2D4AB4D4B54ACD54AAD555555555555555555552B204 2BF1BA84B2D34ACB4D5555555554B54AAAAAAAAAAAAAAAAACB33F1BA84B2CAAD55555555555555 555555555555555555552B33F1BA84B2CD555555555555555555555555555555555554D2B3F1BA84B354AD52AAAD4CAB54CD54AAAD534AB4AAACB2B2B533F1BA84B352ACD2B52D2D5554AAACCB535554B2D53 52B4B4B53F1BA84B34B52B4AB4D332D54CAD34AAD52CD4D52AAD34D4CB3F1BA84B34CCAAD554ACACCACCD3334AD3334ACCAAB2D2AAD4BF1BA84B32AAD4D555555555555555555555555555555552B2BF1BA84B32D5555555555555555555555555555555555552D53F1BA84B33552D55555555555555555555555555 5555554D2CBF1BA84B332AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABF1BA84B2AD2B534B2B4B4CB34ACCB354AB2D34D52CD2D354D3F1BA84B2B554ACB554AB3534B4D4D4AAACAB4ACAD4B54D32D3F1BA84B2B2D5554AAAAD34B52CAAB3534B52CAB5332ACD4AD3F1BA84B2D4AB4D4B54ACD54AAD555555555555555555552B2BF1BA84B2D34ACB4D5555555554B54AAAAAAAAAAAAAAAAACB33F1BA84B2CB52D3555334ACAAD4B55552B32CACAAB4AB2B34D3F1BA84B2CCB4B54D52AAAB2B2CAACB535554AB4ACD54B52CB3F1BA84B354AD3554ACCB2D2AAAACAAAD4CB2D35554ACAD4D2BF1BA84B352CCCCD4CD52D534D534CB535554D53532AD2AB34BF1BA84B34B52CD52B4CD52AD54B52AAD4CB34CAAD4B352D2B3F1BA84B34D32AAACCD4B5532CB5334ACAAAAD2D532AD34AACBF1BA8 4B32B5354AAB4CD2B2CAB552AAD534B35555532CD2D4BF1BA84B32CAAB2B2AAB32B55554CCB53555 4D535352ACAB4ABF1BA84B335555555555555555555555555555555555554D2ABF1BA84B332AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABF1BA84B2AD2B534B2CB4B34CB5334CAB54D2CB2AD32D2D4D4BF1BA84B2B554ACB554AB3534B4D4D4AAACAB4ACAD4B54D32D3F1BA84B2B32AAAB555 52CB4AD3554CACB4AD354ACCD532CB4BF1BA84B2D4AB4D4B54ACD54AAD555555555555555555552B2BF1BA84B2D34ACB4D5555555554B54AAAAAAAAAAAAAAAAACB33F1BA84B2CAAAAD55555555555 555555555555555555554D4B3F1BA84B2CD555555555555555555555555555555555554D2B3F1BA84B355554AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAB534BF1BA84B352AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAB3F1BA84B34B554CACAAB4D2D4B4AD3552D4D32D54AD4B54B4ABF1BA84B34CCB34CD33552B2AAACD34AD3333552D4CACAACCCBF1BA84B32AAAD555555555555555555555555 555555554D4B3F1BA84B32D5555555555555555555555555555555555552D53F1BA84B334AAD34B2AAB52D4CCB32AAD2B2D55554D2D2D54D3F1BA84B332CACB534AAB2CCACAAD34AD3333552D4D2AD4ACD3F1BA84B2AD2B534B34B4B34CB5334CAB54D2CB2AD32D2D2B2BF1BA84B2B554ACB554AB3534B4D4D4AAACAB4ACAD4B54D32D3F1BA84B2B32AAAB55552CB4AD3554CACB4AD354ACCD532CB4BF1BA84B2D4AB4D4B54ACD54AAD555555555555555555552B2BF1BA84B2D34ACB4D5555555554B54AAAAAAAAAAAAAAAAACB33F1BA84B2CAAACAAAD32AD3554D2CCAAD2B2ACAAAD4D2B33333F1BA84B2CD4CCAD32D4AACCAB4B2CB52CCCD334D4D2D4AB4CBF1BA84B355553354B4D552CCB54ACAACAAD333 554B4CB4CAD3F1BA84B352D4B4B32D4D34B532D4CB52CCCD4B52B2D4CAB54BF1BA84B34AAB554B5 5354ACAD3553553554CACAACAACACD2ABF1BA84B34D4B2B32B52B4D55354D34AD3332AAD2B2D2AACD4BF1BA84B32AAB52AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAB554BF1BA84B32D55555555 55555555555555555555555555552D53F1BA84B33554B5555555555555555555555555555555552D4BF1BA84B332AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABF1BA84B2AD2B534B334B4CB34ACCB354AB2D34D52CD2D332B3F1BA84B2B554ACB554AB3534B4D4D4AAACAB4ACAD4B54D32D3F1BA8 4B2B32AAAB55552CB4AD3554CACB4AD354ACCD532CB4BF1BA84B2D4AB4D4B54ACD54AAD5555555 55555555555552B2BF1BA84B2D34ACB4D5555555554B54AAAAAAAAAAAAAAAAACB33F1BA84B2CAAB4D555555555555555555555555555555552B53F1BA84B2CD555555555555555555555555555555555554D2B3F1BA84B35554D5555555555555555555555555555555552D4BF1BA84B352AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAB3F1BA84B34AAB2D55555555555555555555555555555554D4B3F1BA84B34D5555555555555555555555555555555555552D4BF1BA84B32AAB35555555555555555555555555555555552B53F1BA84B32D5555555555555555555555555555555555552D53F1BA84B33554CD55555555555555555555555 555555554D2ABF1BA84B332AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABF1BA84B2AD2B534D54B4B34CB5334CAB54D2CB2AD32D2CD4ABF1BA84B2B554ACB554AB3534B4D4D4AAACAB4ACAD4B54D32D3F1BA84B2B32AAAB55552CB4AD3554CACB4AD354ACCD532CB4BF1BA84B2D4AB4D4B54ACD54AAD555555555555555555552B2BF1BA84B2D34ACB4D5555555554B54AAAAAAAAAAAAAAAAACB33F1BA84B2CAAD55555555555555555555555555555555552B33F1BA84B2CD5555555555555555555 55555555555555554D2B3F1BA84B354AD52AAAD4CAB54CD54AAAD534AB4AAACB2B2B533F1BA84B35 2ACD2B52D2D5554AAACCB535554B2D5352B4B4B53F1BA84B34B52B4AB4D332D54CAD34AAD52CD4D52AAD34D4CB3F1BA84B34CCAAD554ACACCACCD3334AD3334ACCAAB2D2AAD4BF1BA84B32AAD4D555555555555555555555555555555552B2BF1BA84B32D5555555555555555555555555555555555552D53F1BA84B33552D555555555555555555555555555555554D2CBF1BA84B332AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABF1BA84B2AD2B534D534B4CB34ACCB354AB2D34D52CD2D2CD33F1BA84B2B554ACB4000000000000000000000000000000000000000000000000000 Этот результат можно было получить при условии, что имеется эксперименталь- ный файл raw_bits_SV18J.dat. Если экспериментального файла нет, то из его шестнадцате- ричного представления, приведенного выше можно сформировать файл для последующе-205 го декодирования. Для этого скопируем приведенное шестнадцатеричное представление файла и запишем его в виде файла с именем FilGL1.txt с помощью текстового редактора. В форму для декодирования (эта форма в точности соответствует файлу raw_bits_SV18J.dat файл FilGL1.txt преобразуется с помощью программы HexBin_GLONASS.m. Файл HexBin_GLONASS.m %Имя m- файла: HexBin_GLONASS.m %Файл HexBin_GLONASS предназначен для преобразования файла FilGL1.txt %в шестнадцатеричных символах в файл FilGL1_bin.txt в символах "0" и "1" name_in1='FilGL1.txt';%имя файла в шестнадцатеричном представлении fid = fopen(name_in1,'rt'); %открытие файла fw_bin=fopen('FilGL1_bin.txt','wb');%запись преобразованного файла %чтение и переформатирование файла MAS = fscanf(fid,'%s'); for i=1:1875 str_MAS(1:4,i)=sprintf('%s',MAS(i+3*(i-1):(3+i+3*(i-1)))); end MATR=reshape(MAS,4,1875)'; hex_dec= base2dec(MATR,16);%преобразование шестнадцатеричных чисел в десятичные dec_bin= dec2bin(hex_dec)';%преобразование десятичных чисел в двоичные MATR_Char = reshape(dec_bin,1,1875*16);%преобразование в строку %fprintf(fw_bin,'%c ',MATR_Char);%запись файла в символах "0" и "1"(вариант, пробел после с имеет значение) fprintf(fw_bin,'%c\n',MATR_Char);%запись файла в символах "0" и "1" (вариант) fclose(fw_bin);% закрытие файла Файлы FilGL1.txt (восстановленный из шестнадцетиричного) и raw_bits_SV18J.dat (экспериментальный) одинаковы и применяя любой из них можно проводить декодирова- ние данных спутников ГЛОНАСС. Прежде чем приступить к обработке данных файла FilGL1.txt рассмотрим процесс формирования и декодирования в навигационном приемнике информационного навигаци- онного сигнала. Этот процесс может быть представлен в виде следующей блок – схемы (рис. 5.7). Информация, передаваемая с каждого навигационного спутника представляет собой непрерывный поток информационных символов. Начала отсчетов передачи данных нуль часов Московского времени (начало передачи первого суперкадра каждые сутки). Структурно поток цифровой информации (ЦИ) формируется в виде непрерывно по- вторяющихся суперкадров (рис. 5.8). 206 207Рис. 5.7. Формирование и декодирование радионавигационного сигнала 1. Формирование навигационной информации 2. Форматирование навигационной информации 3. Преобразование навигационной информации в двоичный код 4. Формирование информационных сигналов, соот- ветствующих двоичному коду методом двухпози- ционной фазовой манипуляции с частотой 50 Гц 5. Преобразование информационных сигналов в сигналы с относительной фазовой манипуляцией 6. Сложение по модулю 2 информационных сигналов с меандровым колебанием частоты 100 Гц и выделе- ние 170 символов длительностью 10 миллисекунд 7. Дополнение 170 символов тридцатью символами кода метки времени (символ кода метки времени име- ет длитальность 10 миллисеконд). Окончание формирования одной строки навигацион- ногокадра Итог: 200 символов по 10 миллисекунд (2 секунды- время передачи 1 строки) А 208Рис. 5.7. Формирование и декодирование радионавигационного сигнала (продолжение) А 9. Модуляция несущих частот спутника полученным сигналом 10. Демодуляция радионавигационного сигнала 11.Корреляционная обработка демодулированного сигнала. Итог: поток символов длительностью 10 миллисекунд 12. Определение границ строк Выделяются строки длиной по 170 бит 13. Выделение информационных символов каждой строки Выделяются 77 символов информации и 8 провероч- ных символов (код Хемминга) 14. Проверка правильности считанной информации и исправление ошибок 15. Дешифровка двоичных символом 15. Декодирование двоичных символом 16. Перевод декодированных сигналов в десятичный эквивалент 17. Конец декодирования 8. Сложение по модулю 2 информационных сигналов с псевдослучайным кодом стандартной точности с частотой 0.511 МГц 209Рис. 5.8. Изображение потока данных навигационного спутника В приведенной далее программе декодирования навигационных данных спутников ГЛОНАСС выполняются процедуры п. п. 12… 17 (рис. 5.7) по данным п. 11 в соответст- вии с ИКД ГЛОНАСС. 5. Задание 2. Откройте файл Decod_GL3.m. Изучите программные процедуры файла, ис- пользуя его текст и расширенный комментарий. Выполните файл и в командном окне MatLab прочитайте результат декодирования и занесите его в отчет.5.5.3 Задания и вопросы для самоподготовки 1. Изучите форматы сообщений спутников ГЛОНАСС? 2. Что понимается под суперкадром, кадром, строкой и словом в сообщениях спутников ГЛОНАСС? 3. Что понимается под оперативной и неоперативной информацией, передаваемой со спутников ГЛОНАСС? 4. Как контролируется достоверность информации, передаваемой со спутников ГЛО-НАСС? 5. Как выполняется преобразование между системами счисления: переход из двоичной системы счисления в десятичную и наоборот; переход из шестнадцатеричной системы счисления в десятичную и наоборот; переход из двоичной системы счисления в шест- надцатеричную и наоборот? 6. Что такой бит, байт, младший разряд, старший разряд? …Метеки времени Информационные символы строк (строка передается в течение 1.7 секунды) Начало суток Окончание суток 5.5.4 Листинг файла:Decod_GL3.m %Имя m- файла:Decod_GL3.m %Программа декодирования данных спутников ГЛОНАСС %Символы метки времени Met=[1 1 1 1 1 0 0 0 1 1 0 1 1 1 0 1 0 1 0 0 0 0 1 0 0 1 0 1 1 0]; %Открытие файла для считывания данных %fid = fopen('raw_bits_SV18J.dat','rb') %fid = fopen( 'SV18fw_bin.txt','rb') fid = fopen('FilGL1_bin.txt','rb') if fid=-1 kol = 0;%Порядковый номер считываемых данных while (kol < 24) %Поиск и обнаружение меток времени в массиве данных for i = 1:30 MAS(i) = fscanf(fid,'%d',1); end C = xor(MAS, Met); D = any(C); kk = 0; while ( (feof(fid)==0) & (D == 1) ) MAS(1:29) = MAS(2:30); ch = fscanf(fid,'%d',1); MAS(30) = ch; C = xor(MAS, Met); D = any(C); %0 - all 0 kk = kk + 1; end % while %Считывание данных строк kol = kol + 1 for i = 1:170 danDDD(i)= (fscanf(fid,'%d',1)); end; str_danDDD=sprintf('%d',danDDD(1:170)); %Сложение данных по модулю 2 с меандром for i = 1:170 meandr = rem((i-1), 2); %остаток от деления на 2 ((i-1)/2) danDD(i) = xor(danDDD(i), meandr); end; % for i = 1:170 210 % DD = danDD %Выделение информационных символов dd1 = danDD(1:2:169); dd2 = danDD(2:2:170); ddxor = xor(dd1, dd2); Dany = any(ddxor); str1(1) = dd1(1); %Перевод данных из относительного кода for i = 2:85 str1(i) = xor(dd1(i-1), dd1(i)); end; % Контроля четности (правильности) принятых данных Skod(1,1:85)=[0 0 0 0 0 0 0 0 1 1 0 1 1 0 1 0 1 0 1 1 0 1 0 1 0 1 ... 0 1 0 1 0 1 0 1 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1 ... 0 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1 1 0 1 0 1 0 1 ... 0 1 0 1 0 1 0 1 0 1 0 1 0]; Skod(2,1:85)=[0 0 0 0 0 0 0 0 1 0 1 1 0 1 1 0 0 1 1 0 1 1 0 0 1 1 ... 0 0 1 1 0 0 1 1 0 1 1 0 0 1 1 0 0 1 1 0 0 1 1 ... 0 0 1 1 0 0 1 1 0 0 1 1 0 0 1 1 0 1 1 0 0 1 1 ... 0 0 1 1 0 0 1 1 0 0 1 1 0]; Skod(3,1:85)=[0 0 0 0 0 0 0 0 0 1 1 1 0 0 0 1 1 1 1 0 0 0 1 1 1 1 ... 0 0 0 0 1 1 1 1 0 0 0 1 1 1 1 0 0 0 0 1 1 1 1 ... 0 0 0 0 1 1 1 1 0 0 0 0 1 1 1 1 0 0 0 1 1 1 1 ... 0 0 0 0 1 1 1 1 0 0 0 0 1]; Skod(4,1:85)=[0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 0 0 0 0 0 0 0 ... 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 ... 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 ... 1 1 1 1 1 1 1 1 0 0 0 0 0]; Skod(5,1:85)=[0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 ... 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ... 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 ... 0 0 0 0 0 0 0 0 1 1 1 1 0]; Skod(6,1:85)=[0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ... 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 ... 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 ... 0 0 0 0 0 0 0 0 0 0 0 0 0]; 211 Skod(7,1:85)=[0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ... 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ... 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 ... 1 1 1 1 1 1 1 1 1 1 1 1 1]; Skod(8,1:85)=[1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 ... 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 ... 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 ... 1 1 1 1 1 1 1 1 1 1 1 1 1]; E = rot90(eye(85)); for k=1:8 Vi1 = Skod(k,1:85) * E; STR1Vi = and(str1,Vi1); a1 = sum(STR1Vi); n = floor(a1 / 2); c1 = a1 - n * 2; if (k <8) cc_all(kol,k)=xor(c1,str1(86 - k)); else cc_all(kol,k)=c1; end; end; %Декодирование номеров строк num4 = str1(2:5); NN1 = sprintf('%d',num4); num(kol)=bin2dec(NN1); %Вывод данных строк str_danDDD1=sprintf('%d',num(kol),danDDD(1:85)); str_danDDD2=sprintf('%d',num(kol),danDDD(86:170)); s_danDDD=[str_danDDD1;str_danDDD2]; %++++++++++++++++++++++++++++++++++++++++++ %Вывод номера считываемых данных (kol), номера строки %[num(kol)] и коэффициентов четности cc_all(kol,1:8)= %C8 C7 C6 C5 C4 C3 C2 C1 out1 = [kol num(kol) cc_all(kol,8:-1:1)] %Вывод считываемых данных строки SS1=sprintf('%d',str1); %Декодирование данных строк согласно их номеров 212 % в кадре [например, num(kol) == 1- cтрока номер 1] if num(kol) == 1 zx=1; zxdot=1; zx2dot=1; Nstr=sprintf('%d',str1(2:5)); R1= sprintf('%d',str1(6:7)); P1=sprintf('%d',str1(8:9)); tk5=sprintf('%d',str1(10:14)); tk6=sprintf('%d',str1(15:20)); tk1=sprintf('%d',str1(21)); if sprintf('%d',str1(22))=='1' zxdot = -1; end; Xdot=sprintf('%d',str1(23:45)); if sprintf('%d',str1(46))=='1' zx2dot = -1; end; X2dot=sprintf('%d',str1(47:50)); if sprintf('%d',str1(51))=='1' zx = -1; end; X=sprintf('%d',str1(52:77)); St1=bin2dec({Nstr tk5 tk6 tk1 X Xdot X2dot}); Mng=[1 1 1 30 zx*2^(-11) zxdot*2^(-20) zx2dot*2^(-30)]; Str1=(St1'.*Mng)'%вывод декодированных данных end; if num(kol) == 2 zy=1; zydot=1; zy2dot=1; Nstr=sprintf('%d',str1(2:5)); Bn= sprintf('%d',str1(6:8)); P2=sprintf('%d',str1(9)); tb=sprintf('%d',str1(10:16)); R2=sprintf('%d',str1(17:21)); if sprintf('%d',str1(22))=='1' zydot = -1; end; Ydot=sprintf('%d',str1(23:45)); if sprintf('%d',str1(46))=='1' zy2dot = -1; end; 213 Y2dot=sprintf('%d',str1(47:50)); if sprintf('%d',str1(51))=='1' zy = -1; end; Y=sprintf('%d',str1(52:77)); St2=bin2dec({Nstr Bn P2 tb Y Ydot Y2dot}); Mng=[1 1 1 15*60 zy*2^(-11) zydot*2^(-20) zy2dot*2^(-30)]; Str2=(St2'.*Mng)'%вывод декодированных данных end; if num(kol) == 3 zz=1; zzdot=1; zz2dot=1; zgamma=1; Nstr=sprintf('%d',str1(2:5)); P3= sprintf('%d',str1(6)); if sprintf('%d',str1(7))=='1' zgamma = -1; end; Gamma=sprintf('%d',str1(8:17)); R3=sprintf('%d',str1(18:19)); P=sprintf('%d',str1(20)); ln=sprintf('%d',str1(21)); if sprintf('%d',str1(22))=='1' zzdot = -1; end; Zdot=sprintf('%d',str1(23:45)); if sprintf('%d',str1(46))=='1' zz2dot = -1; end; Z2dot=sprintf('%d',str1(47:50)); if sprintf('%d',str1(51))=='1' zz = -1; end; Z=sprintf('%d',str1(52:77)); St3=bin2dec({Nstr P3 ln Gamma Z Zdot Z2dot}); Mng=[1 1 1 zgamma*2^(-40) zz*2^(-11) zzdot*2^(-20) zz2dot*2^(-30)]; Str3=(St3'.*Mng)'%вывод декодированных данных end; 214 if num(kol) == 4 ztaun=1; zdtaun=1; Nstr=sprintf('%d',str1(2:5)); if sprintf('%d',str1(6))=='1' ztaun = -1; end; taun= sprintf('%d',str1(7:27)); if sprintf('%d',str1(28))=='1' zdtaun = -1; end; deltatau=sprintf('%d',str1(29:32)); En=sprintf('%d',str1(33:37)); R4=sprintf('%d',str1(38:51)); P4=sprintf('%d',str1(52)); FT=sprintf('%d',str1(53:56)); R4i4=sprintf('%d',str1(57:59)); NT=sprintf('%d',str1(60:70)); nSV=sprintf('%d',str1(71:75)); modifSV=sprintf('%d',str1(76:77)); St4=bin2dec({Nstr taun deltatau En FT NT nSV modifSV}); Mng=[1 ztaun*2^(-30) zdtaun*2^(-30) 1 1 1 1 1]; Str4=(St4'.*Mng)'%вывод декодированных данных end; if num(kol) == 5 ztauc=1; ztauGPS=1; Nstr=sprintf('%d',str1(2:5)); NA= sprintf('%d',str1(6:16)); if sprintf('%d',str1(17))=='1' ztauc = -1; end; tauc=sprintf('%d',str1(18:48)); R5=sprintf('%d',str1(49)); N4=sprintf('%d',str1(50:54)); if sprintf('%d',str1(55))=='1' ztauGPS = -1; end; tauGPS=sprintf('%d',str1(56:76)); ln=sprintf('%d',str1(77)); St5=bin2dec({Nstr NA tauc N4 tauGPS ln Nstr Nstr}); Mng=[1 1 ztauc*2^(-27) 1 ztauGPS*2^(-30) 1 0 0]; 215 Str5=(St5'.*Mng)'%вывод декодированных данных end; if num(kol) == 6 ztaun=1; zlambda=1;zdeltai=1; Nstr=sprintf('%d',str1(2:5)); CNA= sprintf('%d',str1(6)); Mn=sprintf('%d',str1(7:8)); nSV=sprintf('%d',str1(9:13)); if sprintf('%d',str1(14))=='1' ztaun = -1; end; taun=sprintf('%d',str1(15:23)); if sprintf('%d',str1(24))=='1' zlambda = -1; end; lambda=sprintf('%d',str1(25:44)); if sprintf('%d',str1(45))=='1' zdeltai= -1; end; deltai=sprintf('%d',str1(46:62)); ecs=sprintf('%d',str1(63:77)); St6=bin2dec({Nstr Mn nSV taun lambda deltai ecs}); Mng=[1 1 1 ztaun*2^(-18) zlambda*2^(-20) zdeltai*2^(-20) 2^(-20)]; Str6=(St6'.*Mng)'%вывод декодированных данных end; if num(kol) == 7 zomegan=1;zdeltaT=1; zdeltaTdot=1; Nstr=sprintf('%d',str1(2:5)); if sprintf('%d',str1(6))=='1' zomegan = -1; end; omegan= sprintf('%d',str1(7:21)); tlambda=sprintf('%d',str1(22:44)); if sprintf('%d',str1(43))=='1' zdeltaT = -1; end; deltaT=sprintf('%d',str1(44:64)); if sprintf('%d',str1(65))=='1' zdeltaTdot= -1; end; 216 deltaTdot=sprintf('%d',str1(66:71)); Hn=sprintf('%d',str1(72:76)); ln=sprintf('%d',str1(77)); St7=bin2dec({Nstr omegan tlambda deltaT deltaTdot Hn ln}); Mng=[1 zomegan*2^(-15) 2^(-5) zdeltaT*2^(-9) zdeltaTdot*2^(-14) 1 1]; Str7=(St7'.*Mng)'%вывод декодированных данных end; if num(kol) == 8 ztaun=1; zlambda=1;zdeltai=1; Nstr=sprintf('%d',str1(2:5)); CNA= sprintf('%d',str1(6)); Mn=sprintf('%d',str1(7:8)); nSV=sprintf('%d',str1(9:13)); if sprintf('%d',str1(14))=='1' ztaun = -1; end; taun=sprintf('%d',str1(15:23)); if sprintf('%d',str1(24))=='1' zlambda = -1; end; lambda=sprintf('%d',str1(25:44)); if sprintf('%d',str1(45))=='1' zdeltai= -1; end; deltai=sprintf('%d',str1(46:62)); ecs=sprintf('%d',str1(63:77)); St8=bin2dec({Nstr Mn nSV taun lambda deltai ecs}); Mng=[1 1 1 ztaun*2^(-18) zlambda*2^(-20) zdeltai*2^(-20) 2^(-20)]; Str8=(St8'.*Mng)' %вывод декодированных данных end; if num(kol) == 9 zomegan=1;zdeltaT=1; zdeltaTdot=1; Nstr=sprintf('%d',str1(2:5)); if sprintf('%d',str1(6))=='1' zomegan = -1; end; omegan= sprintf('%d',str1(7:21)); tlambda=sprintf('%d',str1(22:44)); if sprintf('%d',str1(43))=='1' zdeltaT = -1; 217 end; deltaT=sprintf('%d',str1(44:64)); if sprintf('%d',str1(65))=='1' zdeltaTdot= -1; end; deltaTdot=sprintf('%d',str1(66:71)); Hn=sprintf('%d',str1(72:76)); ln=sprintf('%d',str1(77)); St9=bin2dec({Nstr omegan tlambda deltaT deltaTdot Hn ln}); Mng=[1 zomegan*2^(-15) 2^(-5) zdeltaT*2^(-9) zdeltaTdot*2^(-14) 1 1]; Str9=(St9'.*Mng)'%вывод декодированных данных end; if num(kol) == 10 ztaun=1; zlambda=1;zdeltai=1; Nstr=sprintf('%d',str1(2:5)); CNA= sprintf('%d',str1(6)); Mn=sprintf('%d',str1(7:8)); nSV=sprintf('%d',str1(9:13)); if sprintf('%d',str1(14))=='1' ztaun = -1; end; taun=sprintf('%d',str1(15:23)); if sprintf('%d',str1(24))=='1' zlambda = -1; end; lambda=sprintf('%d',str1(25:44)); if sprintf('%d',str1(45))=='1' zdeltai= -1; end; deltai=sprintf('%d',str1(46:62)); ecs=sprintf('%d',str1(63:77)); St10=bin2dec({Nstr Mn nSV taun lambda deltai ecs}); Mng=[1 1 1 ztaun*2^(-18) zlambda*2^(-20) zdeltai*2^(-20) 2^(-20)]; Str10=(St10'.*Mng)'%вывод декодированных данных end; if num(kol) == 11 zomegan=1;zdeltaT=1; zdeltaTdot=1; Nstr=sprintf('%d',str1(2:5)); if sprintf('%d',str1(6))=='1' zomegan = -1; 218 end; omegan= sprintf('%d',str1(7:21)); tlambda=sprintf('%d',str1(22:44)); if sprintf('%d',str1(43))=='1' zdeltaT = -1; end; deltaT=sprintf('%d',str1(44:64)); if sprintf('%d',str1(65))=='1' zdeltaTdot= -1; end; deltaTdot=sprintf('%d',str1(66:71)); Hn=sprintf('%d',str1(72:76)); ln=sprintf('%d',str1(77)); St11=bin2dec({Nstr omegan tlambda deltaT deltaTdot Hn ln}); Mng=[1 zomegan*2^(-15) 2^(-5) zdeltaT*2^(-9) zdeltaTdot*2^(-14) 1 1]; Str11=(St11'.*Mng)'%вывод декодированных данных end; if num(kol) == 12 ztaun=1; zlambda=1;zdeltai=1; Nstr=sprintf('%d',str1(2:5)); CNA= sprintf('%d',str1(6)); Mn=sprintf('%d',str1(7:8)); nSV=sprintf('%d',str1(9:13)); if sprintf('%d',str1(14))=='1' ztaun = -1; end; taun=sprintf('%d',str1(15:23)); if sprintf('%d',str1(24))=='1' zlambda = -1; end; lambda=sprintf('%d',str1(25:44)); if sprintf('%d',str1(45))=='1' zdeltai= -1; end; deltai=sprintf('%d',str1(46:62)); ecs=sprintf('%d',str1(63:77)); St12=bin2dec({Nstr Mn nSV taun lambda deltai ecs}); Mng=[1 1 1 ztaun*2^(-18) zlambda*2^(-20) zdeltai*2^(-20) 2^(-20)]; Str12=(St12'.*Mng)'%вывод декодированных данных end; 219 if num(kol) == 13 zomegan=1;zdeltaT=1; zdeltaTdot=1; Nstr=sprintf('%d',str1(2:5)); if sprintf('%d',str1(6))=='1' zomegan = -1; end; omegan= sprintf('%d',str1(7:21)); tlambda=sprintf('%d',str1(22:44)); if sprintf('%d',str1(43))=='1' zdeltaT = -1; end; deltaT=sprintf('%d',str1(44:64)); if sprintf('%d',str1(65))=='1' zdeltaTdot= -1; end; deltaTdot=sprintf('%d',str1(66:71)); Hn=sprintf('%d',str1(72:76)); ln=sprintf('%d',str1(77)); St13=bin2dec({Nstr omegan tlambda deltaT deltaTdot Hn ln}); Mng=[1 zomegan*2^(-15) 2^(-5) zdeltaT*2^(-9) zdeltaTdot*2^(-14) 1 1]; Str13=(St13'.*Mng)'%вывод декодированных данных end; if num(kol) == 14 ztaun=1; zlambda=1;zdeltai=1; Nstr=sprintf('%d',str1(2:5)); CNA= sprintf('%d',str1(6)); Mn=sprintf('%d',str1(7:8)); nSV=sprintf('%d',str1(9:13)); if sprintf('%d',str1(14))=='1' ztaun = -1; end; taun=sprintf('%d',str1(15:23)); if sprintf('%d',str1(24))=='1' zlambda = -1; end; lambda=sprintf('%d',str1(25:44)); if sprintf('%d',str1(45))=='1' zdeltai= -1; end; deltai=sprintf('%d',str1(46:62)); ecs=sprintf('%d',str1(63:77)); 220 St14=bin2dec({Nstr Mn nSV taun lambda deltai ecs}); Mng=[1 1 1 ztaun*2^(-18) zlambda*2^(-20) zdeltai*2^(-20) 2^(-20)]; Str14=(St14'.*Mng)'%вывод декодированных данных end; %Строка 14 для 5 кадра if num(kol) == 14 zB1=1; zB2=1; Nstr=sprintf('%d',str1(2:5)); if sprintf('%d',str1(6))=='1' zB1 = -1; end; B1=sprintf('%d',str1(7:16)); if sprintf('%d',str1(17))=='1' zB2 = -1; end; B2=sprintf('%d',str1(18:26)); Kp=sprintf('%d',str1(27:28)); St14_5=bin2dec({Nstr B1 B2 Kp}); Mng=[1 zB1*2^(-10) zB2*2^(-16) 1]; Str14_5=(St14_5'.*Mng)'%вывод декодированных данных end; if num(kol) == 15 zomegan=1;zdeltaT=1; zdeltaTdot=1; Nstr=sprintf('%d',str1(2:5)); if sprintf('%d',str1(6))=='1' zomegan = -1; end; omegan= sprintf('%d',str1(7:21)); tlambda=sprintf('%d',str1(22:44)); if sprintf('%d',str1(43))=='1' zdeltaT = -1; end; deltaT=sprintf('%d',str1(44:64)); if sprintf('%d',str1(65))=='1' zdeltaTdot= -1; end; deltaTdot=sprintf('%d',str1(66:71)); Hn=sprintf('%d',str1(72:76)); ln=sprintf('%d',str1(77)); 221 St15=bin2dec({Nstr omegan tlambda deltaT deltaTdot Hn ln}); Mng=[1 zomegan*2^(-15) 2^(-5) zdeltaT*2^(-9) zdeltaTdot*2^(-14) 1 1]; Str15=(St15'.*Mng)'%вывод декодированных данных end; str1; end; %while (feof(fid)==0) | % while (kol < 14) end % if fid=-1 fclose(fid); %Технологические данные Qout=[num' cc_all] %Номер строки и проверка на четность %Qoutnum1=Qout(1,2:9); %SS1=sprintf('%d',str1); %Qout123=[Str1 Str2 Str3] %Qout45=[Str4 Str5] %Qout6789=[Str6 Str7 Str8 Str9] %Qout10_11_12=[Str10 Str11 Str12] %Qout13_14_15=[Str13 Str14 Str15] clear; 5.5.5 Пример выполнения файла:Decod_GL3.m В результате выполнения программы в командном окне Command Window MatLab отображаются декодированные данные по каждой из строк в последовательности и соответствии ИКД ГЛОНАСС. Так, например, структура строки 1согласно ИКД ГЛОНАСС имеет вид 1 2 - 5 6 - 7 8 - 9 10 - 21 22 - 45 85 84 -81 80 - 79 78 - 77 76 - 65 64 - 410 m P1 t k)t(nXb.Номер строки в навигацион- ном кадре Резерв Признак смены оператив- ной инфо- рмации Время начала кадра внутри текущих суток, исчисляемое в шкале бортового времени. Начало суток по бортовому времени спут- ника совпадает с началом очередного су- перкадра. Составляющие вектора скорости n-го спутника в системе коорди- нат ПЗ-90 на момент времени tb46 - 50 51 - 77 1 - 8 1 - 30 40 - 36 35 - 9 8 - 1 1 - 30 )t(nXb..)t(nXbСоставляющие ускорения n-го спутника на момент времени tb, Координаты n-го спутника в системе координат ПЗ-90 на момент времени t b; Код Хемминга Метка времени Рис. 5.9. - Структура 1 строки для (1- 4) кадров суперкадра 222 Декодированная информация из командного окна MatLab приведена на рис. 5.10. Рис. 5.10 Декодированные данные спутника ГЛОНАСС Данные, изображенные на рис. 5.10 соответствуют рис. 5.9, но записаны в столбец. 223 32>1   ...   6   7   8   9   10   11   12   13   14

Mean Anom(rad): -0.1107067493E+001
Af0(s): 0.5626678467E-004
Af1(s/s): 0.0000000000E+000 week: 371
******** Week 371 almanac for PRN-22 ********
ID: 22
Health: 000
Eccentricity: 0.5007266998E-002
Time of Applicability(s): 589824.0000
Orbital Inclination(rad): 0.9536770543
Rate of Right Ascen(r/s): -0.8343204671E-008
SQRT(A) (m 1/2): 5153.625000
Right Ascen at Week(rad): -0.8549179031E+000
Argument of Perigee(rad): -1.675495384
Mean Anom(rad): 0.2358868096E+001
Af0(s): 0.1401901245E-003
Af1(s/s): 0.3637978807E-011 week: 371
******** Week 371 almanac for PRN-23 ********
ID: 23
Health: 000
Eccentricity: 0.4817008972E-002
Time of Applicability(s): 589824.0000
Orbital Inclination(rad): 0.9696759945
Rate of Right Ascen(r/s): -0.7600316584E-008
SQRT(A) (m 1/2): 5153.688965
Right Ascen at Week(rad): 0.1404933156E+000
Argument of Perigee(rad): 2.539610056
Mean Anom(rad): 0.1512449255E+001
Af0(s): 0.1382827759E-003
Af1(s/s): 0.0000000000E+000 week: 371
******** Week 371 almanac for PRN-24 ********
ID: 24
Health: 000
Eccentricity: 0.8838176727E-002
Time of Applicability(s): 589824.0000 238

Orbital Inclination(rad): 0.9582849887
Rate of Right Ascen(r/s): -0.8263201338E-008
SQRT(A) (m 1/2): 5155.323242
Right Ascen at Week(rad): -0.1901905855E+001
Argument of Perigee(rad): -0.829185525
Mean Anom(rad): -0.7264556257E+000
Af0(s): 0.4196166992E-004
Af1(s/s): 0.3637978807E-011 week: 371
******** Week 371 almanac for PRN-25 ********
ID: 25
Health: 000
Eccentricity: 0.1263189316E-001
Time of Applicability(s): 589824.0000
Orbital Inclination(rad): 0.9556604435
Rate of Right Ascen(r/s): -0.8057478483E-008
SQRT(A) (m 1/2): 5153.598633
Right Ascen at Week(rad): 0.1094536862E+001
Argument of Perigee(rad): -1.322550598
Mean Anom(rad): -0.5200636789E+000
Af0(s): 0.3700256348E-003
Af1(s/s): 0.1818989404E-010 week: 371
******** Week 371 almanac for PRN-26 ********
ID: 26
Health: 000
Eccentricity: 0.1748418808E-001
Time of Applicability(s): 589824.0000
Orbital Inclination(rad): 0.9914573233
Rate of Right Ascen(r/s): -0.7406022776E-008
SQRT(A) (m 1/2): 5153.694336
Right Ascen at Week(rad): 0.1706179119E+000
Argument of Perigee(rad): 0.818631543
Mean Anom(rad): 0.6507314281E+000
Af0(s): -0.2841949463E-003
Af1(s/s): -0.1455191523E-010 week: 371
******** Week 371 almanac for PRN-27 ********
ID: 27 239

Health: 000
Eccentricity: 0.1995992661E-001
Time of Applicability(s): 589824.0000
Orbital Inclination(rad): 0.9601964726
Rate of Right Ascen(r/s): -0.8080336578E-008
SQRT(A) (m 1/2): 5153.531250
Right Ascen at Week(rad): 0.1126961306E+001
Argument of Perigee(rad): -1.892119987
Mean Anom(rad): -0.1971266429E+001
Af0(s): 0.7152557373E-004
Af1(s/s): 0.3637978807E-011 week: 371
******** Week 371 almanac for PRN-28 ********
ID: 28
Health: 000
Eccentricity: 0.1171159744E-001
Time of Applicability(s): 589824.0000
Orbital Inclination(rad): 0.9606279047
Rate of Right Ascen(r/s): -0.8091765626E-008
SQRT(A) (m 1/2): 5153.618652
Right Ascen at Week(rad): 0.2256309346E+001
Argument of Perigee(rad): -2.241381122
Mean Anom(rad): 0.3110969588E+001
Af0(s): 0.1716613770E-004
Af1(s/s): 0.0000000000E+000 week: 371
******** Week 371 almanac for PRN-29 ********
ID: 29
Health: 000
Eccentricity: 0.9555339813E-002
Time of Applicability(s): 589824.0000
Orbital Inclination(rad): 0.9880538034
Rate of Right Ascen(r/s): -0.7440309919E-008
SQRT(A) (m 1/2): 5153.555664
Right Ascen at Week(rad): 0.1356359595E+000
Argument of Perigee(rad): -0.841310565
Mean Anom(rad): 0.2484434685E+001
Af0(s): 0.4196166992E-004
Af1(s/s): 0.5093170330E-010 week: 371 240

******** Week 371 almanac for PRN-30 ********
ID: 30
Health: 000
Eccentricity: 0.9448528290E-002
Time of Applicability(s): 589824.0000
Orbital Inclination(rad): 0.9449225780
Rate of Right Ascen(r/s): -0.8308917528E-008
SQRT(A) (m 1/2): 5153.618164
Right Ascen at Week(rad): 0.2199777136E+001
Argument of Perigee(rad): 1.316367861
Mean Anom(rad): -0.2894497036E+001
Af0(s): 0.1239776611E-004
Af1(s/s): 0.0000000000E+000 week: 371
******** Week 371 almanac for PRN-31 ********
ID: 31
Health: 063
Eccentricity: 0.4019737244E-002
Time of Applicability(s): 589824.0000
Orbital Inclination(rad): 0.9588722157
Rate of Right Ascen(r/s): -0.7931758961E-008
SQRT(A) (m 1/2): 5171.043457
Right Ascen at Week(rad): 0.1184157143E+001
Argument of Perigee(rad): -2.911533361
Mean Anom(rad): 0.9753739691E+000
Af0(s): -0.8773803711E-004
Af1(s/s): 0.2546585165E-010 week: 371
Приложение 3 Альманах 002.txt
GPS: id=1 health=0 eccentricity=0.00668621063232421875 corr. to inclination=0.9891563653945922852 rate of right ascention=-8.011762260562438342e-12 square root of semi-major axis=5153.6787109375 right ascention parameter=-1.190934181213378906 argument of perigee=-1.770861387252807617 mean anomaly=1.673161864280700684
Af0m=0.09632110595703125 241

Af1=3.637978807091712952e-12
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=2 health=0 eccentricity=0.00882053375244140625 corr. to inclination=0.948260188102722168 rate of right ascention=-7.714606731212825252e-12 square root of semi-major axis=5153.7060546875 right ascention parameter=2.949104547500610352 argument of perigee=2.182616949081420898 mean anomaly=2.895400762557983398
Af0m=0.06580352783203125
Af1=3.637978807091712952e-12
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=3 health=0 eccentricity=0.0090236663818359375 corr. to inclination=0.9254841804504394531 rate of right ascention=-8.263201754787896647e-12 square root of semi-major axis=5153.63916015625 right ascention parameter=1.81049799919128418 argument of perigee=0.726136922836303711 mean anomaly=2.211142301559448242
Af0m=0.47397613525390625
Af1=1.818989403545856476e-11
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=4 health=0 eccentricity=0.00775909423828125 corr. to inclination=0.9477388858795166016 rate of right ascention=-7.748894408077244834e-12 square root of semi-major axis=5153.6650390625 right ascention parameter=2.96854710578918457 argument of perigee=0.2202932685613632202 mean anomaly=-0.8276365399360656738
Af0m=0.34809112548828125 242

Af1=3.637978807091712952e-12
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=5 health=0 eccentricity=0.007638454437255859375 corr. to inclination=0.9389963746070861816 rate of right ascention=-7.806039668822872812e-12 square root of semi-major axis=5153.82275390625 right ascention parameter=0.7738637328147888184 argument of perigee=1.15076756477355957 mean anomaly=-1.178214788436889648
Af0m=0.438690185546875
Af1=1.818989403545856476e-11
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=6 health=0 eccentricity=0.00571727752685546875 corr. to inclination=0.9330402016639709473 rate of right ascention=-8.18319804279932228e-12 square root of semi-major axis=5153.53857421875 right ascention parameter=1.870720624923706055 argument of perigee=-1.767207980155944824 mean anomaly=0.1111147403717041016
Af0m=0.49877166748046875
Af1=-1.818989403545856476e-11
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=7 health=0 eccentricity=0.010822296142578125 corr. to inclination=0.9345921874046325684 rate of right ascention=-8.18319804279932228e-12 square root of semi-major axis=5152.01513671875 right ascention parameter=1.84484410285949707 argument of perigee=-1.705631256103515625 mean anomaly=-0.2914162874221801758
Af0m=0.46825408935546875 243

Af1=-1.091393642127513886e-11
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=8 health=255 eccentricity=0.009907245635986328125 corr. to inclination=0.9763332605361938477 rate of right ascention=-7.611745435343042488e-12 square root of semi-major axis=5153.5458984375 right ascention parameter=-0.1443144083023071289 argument of perigee=2.730259895324707031 mean anomaly=0.2543621659278869629
Af0m=-0.09059906005859375
Af1=0
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=9 health=0 eccentricity=0.01850032806396484375 corr. to inclination=0.9640074372291564941 rate of right ascention=-7.783181217579926426e-12 square root of semi-major axis=5153.68310546875 right ascention parameter=-0.2266040891408920288 argument of perigee=1.299516916275024414 mean anomaly=-0.1355632990598678589
Af0m=0.049591064453125
Af1=0
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=10 health=0 eccentricity=0.00734424591064453125 corr. to inclination=0.9713417887687683105 rate of right ascention=-7.771751991958453232e-12 square root of semi-major axis=5153.69580078125 right ascention parameter=-2.259819507598876953 argument of perigee=0.4371939003467559814 mean anomaly=-2.844457387924194336
Af0m=0.09441375732421875 244

Af1=0
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=11 health=0 eccentricity=0.006408214569091796875 corr. to inclination=0.8953858017921447754 rate of right ascention=-8.274630113047631852e-12 square root of semi-major axis=5153.63232421875 right ascention parameter=2.79404592514038086 argument of perigee=0.440089195966720581 mean anomaly=0.6644263863563537598
Af0m=0.00095367431640625
Af1=0
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=12 health=0 eccentricity=0.003298282623291015625 corr. to inclination=0.9608616232872009277 rate of right ascention=-7.600316209721569294e-12 square root of semi-major axis=5153.68505859375 right ascention parameter=0.8550205230712890625 argument of perigee=-1.343068361282348633 mean anomaly=1.30158078670501709
Af0m=-0.08106231689453125
Af1=1.091393642127513886e-11
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=13 health=0 eccentricity=0.002939701080322265625 corr. to inclination=0.9933208823204040527 rate of right ascention=-7.97747545105975675e-12 square root of semi-major axis=5153.63232421875 right ascention parameter=-1.204549789428710938 argument of perigee=1.263086795806884766 mean anomaly=-3.038789987564086914
Af0m=0.12493133544921875 245

Af1=3.637978807091712952e-12
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=14 health=0 eccentricity=0.002999782562255859375 corr. to inclination=0.9873707294464111328 rate of right ascention=-8.03462071180538473e-12 square root of semi-major axis=5153.67041015625 right ascention parameter=-1.218866467475891113 argument of perigee=-2.083027839660644531 mean anomaly=2.695198297500610352
Af0m=0.0019073486328125
Af1=0
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=16 health=0 eccentricity=0.003540515899658203125 corr. to inclination=0.9628809094429016113 rate of right ascention=-7.634603886585988874e-12 square root of semi-major axis=5153.51513671875 right ascention parameter=0.8703708052635192871 argument of perigee=-0.7588751912117004394 mean anomaly=-1.416219711303710938
Af0m=0.10395050048828125
Af1=3.637978807091712952e-12
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=17 health=0 eccentricity=0.002378940582275390625 corr. to inclination=0.9593995213508605957 rate of right ascention=-7.88604251344970919e-12 square root of semi-major axis=5153.6953125 right ascention parameter=1.909842729568481445 argument of perigee=3.054111957550048828 mean anomaly=-2.438158750534057617
Af0m=0.0858306884765625 246

Af1=0
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=18 health=0 eccentricity=0.007780551910400390625 corr. to inclination=0.9554147720336914062 rate of right ascention=-7.908900964692655578e-12 square root of semi-major axis=5153.69287109375 right ascention parameter=-2.241117000579833984 argument of perigee=-2.594098567962646484 mean anomaly=-1.445593714714050293
Af0m=-0.2536773681640625
Af1=0
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=19 health=0 eccentricity=0.00347232818603515625 corr. to inclination=0.9573382139205932617 rate of right ascention=-7.908900964692655578e-12 square root of semi-major axis=5153.6171875 right ascention parameter=1.968618631362915039 argument of perigee=-1.0020391941070556641 mean anomaly=-2.867439746856689453
Af0m=0.0476837158203125
Af1=3.637978807091712952e-12
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=20 health=0 eccentricity=0.002883434295654296875 corr. to inclination=0.9549533724784851074 rate of right ascention=-7.897471739071182384e-12 square root of semi-major axis=5153.5830078125 right ascention parameter=-2.294228076934814453 argument of perigee=1.280813694000244141 mean anomaly=-1.452993631362915039
Af0m=-0.02193450927734375 247

Af1=0
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=21 health=0 eccentricity=0.01182270050048828125 corr. to inclination=0.9413033127784729004 rate of right ascention=-7.817468894444346006e-12 square root of semi-major axis=5153.6162109375 right ascention parameter=2.987459897994995117 argument of perigee=-2.899925470352172852 mean anomaly=-0.1171304434537887573
Af0m=0.06580352783203125
Af1=0
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=22 health=0 eccentricity=0.00497913360595703125 corr. to inclination=0.9529759883880615234 rate of right ascention=-7.94318777419533717e-12 square root of semi-major axis=5153.6376953125 right ascention parameter=-2.233473777770996094 argument of perigee=-1.647951126098632812 mean anomaly=-2.951039314270019531
Af0m=0.15735626220703125
Af1=3.637978807091712952e-12
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=23 health=0 eccentricity=0.00469684600830078125 corr. to inclination=0.9696160554885864258 rate of right ascention=-8.206056494042268668e-12 square root of semi-major axis=5153.6533203125 right ascention parameter=-1.235992670059204102 argument of perigee=2.586992979049682617 mean anomaly=2.443128108978271484
Af0m=0.1430511474609375 248

Af1=0
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=24 health=0 eccentricity=0.008758544921875 corr. to inclination=0.9579614400863647461 rate of right ascention=-7.566029400218887702e-12 square root of semi-major axis=5155.3740234375 right ascention parameter=3.00449061393737793 argument of perigee=-0.772795259952545166 mean anomaly=-0.7267537117004394531
Af0m=0.06198883056640625
Af1=3.637978807091712952e-12
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=25 health=0 eccentricity=0.01288080215454101562 corr. to inclination=0.9573562145233154297 rate of right ascention=-7.874614155189973986e-12 square root of semi-major axis=5153.671875 right ascention parameter=-0.2823831737041473389 argument of perigee=-1.308804750442504883 mean anomaly=0.4633183777332305908
Af0m=0.400543212890625
Af1=-7.275957614183425904e-12
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=26 health=0 eccentricity=0.017627716064453125 corr. to inclination=0.9914153814315795898 rate of right ascention=-8.000333034940965149e-12 square root of semi-major axis=5153.564453125 right ascention parameter=-1.204336762428283691 argument of perigee=0.8254280686378479004 mean anomaly=1.649268627166748047
Af0m=-0.057220458984375 249

Af1=-1.818989403545856476e-11
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=27 health=0 eccentricity=0.02038526535034179688 corr. to inclination=0.9619401693344116211 rate of right ascention=-7.794610443201399619e-12 square root of semi-major axis=5153.75048828125 right ascention parameter=-0.249639272689819336 argument of perigee=-1.867654919624328613 mean anomaly=-0.9617273211479187012
Af0m=0.08392333984375
Af1=3.637978807091712952e-12
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=28 health=0 eccentricity=0.0119800567626953125 corr. to inclination=0.9609574675559997559 rate of right ascention=-7.588887851461834089e-12 square root of semi-major axis=5153.6337890625 right ascention parameter=0.8796859383583068848 argument of perigee=-2.232036352157592774 mean anomaly=-2.177708625793457031
Af0m=0.011444091796875
Af1=0
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=29 health=0 eccentricity=0.00990962982177734375 corr. to inclination=0.9879399538040161133 rate of right ascention=-8.011762260562438342e-12 square root of semi-major axis=5153.56396484375 right ascention parameter=-1.239597320556640625 argument of perigee=-0.8064963817596435547 mean anomaly=-2.79503321647644043
Af0m=0.27370452880859375 250

Af1=2.182787284255027771e-11
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=30 health=0 eccentricity=0.009554386138916015625 corr. to inclination=0.9454199075698852539 rate of right ascention=-7.748894408077244834e-12 square root of semi-major axis=5153.60595703125 right ascention parameter=0.821842193603515625 argument of perigee=1.328736305236816406 mean anomaly=-1.895180821418762207
Af0m=0.02002716064453125
Af1=0
Af0l=0 week almanac=382 time of week almanac=405504000
GPS: id=31 health=0 eccentricity=0.00612545013427734375 corr. to inclination=0.9606039524078369141 rate of right ascention=-7.794610443201399619e-12 square root of semi-major axis=5153.56884765625 right ascention parameter=-0.1925415545701980591 argument of perigee=-1.421732783317565918 mean anomaly=0.668113410472869873
Af0m=-0.0171661376953125
Af1=7.275957614183425904e-12
Af0l=0 week almanac=382 time of week almanac=405504000
GLONASS: id=1 здоровье=0 номер литера=7 уход шкалы времени спутника=-0.057220458984375 долгота первого внутри суток восходящего узла=3.13790440559387207 наклонение орбиты=1.125113844871520996 эксцентриситет=0.00095844268798828125 аргумент перигея орбиты спутника=2.598275899887084961 время прохождения восходящего узла орбиты=21426938 драконический период=40543929.6875 251
скорость изменения драконического периода=-0.18310546875 номер дня=1086
GLONASS: id=2 здоровье=0 номер литера=1 уход шкалы времени спутника=0.030517578125 долгота первого внутри суток восходящего узла=2.796126365661621094 наклонение орбиты=1.127121090888977051 эксцентриситет=0.00434970855712890625 аргумент перигея орбиты спутника=1.660821795463562012 время прохождения восходящего узла орбиты=26365062 драконический период=40543832.03125 скорость изменения драконического периода=-0.18310546875 номер дня=1086
GLONASS: id=3 здоровье=255 номер литера=12 уход шкалы времени спутника=-0.1220703125 долгота первого внутри суток восходящего узла=2.377643346786499024 наклонение орбиты=1.12163841724395752 эксцентриситет=0.002330780029296875 аргумент перигея орбиты спутника=-2.618601083755493164 время прохождения восходящего узла орбиты=31583032 драконический период=40543917.96875 скорость изменения драконического периода=-0.18310546875 номер дня=1086
GLONASS: id=4 здоровье=0 номер литера=6 уход шкалы времени спутника=-0.209808349609375 долгота первого внутри суток восходящего узла=2.064510345458984375 наклонение орбиты=1.12721097469329834 эксцентриситет=0.00396633148193359375 аргумент перигея орбиты спутника=1.437052369117736816 время прохождения восходящего узла орбиты=36400220 драконический период=40543822.265625 скорость изменения драконического периода=-0.18310546875 номер дня=1086
GLONASS: id=5 здоровье=255 номер литера=7 уход шкалы времени спутника=0 252
долгота первого внутри суток восходящего узла=0 наклонение орбиты=1.099557399749755859 эксцентриситет=0 аргумент перигея орбиты спутника=0 время прохождения восходящего узла орбиты=0 драконический период=43200000 скорость изменения драконического периода=0 номер дня=1086
GLONASS: id=6 здоровье=0 номер литера=1 уход шкалы времени спутника=0.034332275390625 долгота первого внутри суток восходящего узла=3.117007017135620117 наклонение орбиты=1.126872420310974121 эксцентриситет=0.004566192626953125 аргумент перигея орбиты спутника=1.749025702476501465 время прохождения восходящего узла орбиты=22698750 драконический период=40543734.375 скорость изменения драконического периода=-0.244140625 номер дня=1083
GLONASS: id=7 здоровье=0 номер литера=4 уход шкалы времени спутника=-0.20599365234375 долгота первого внутри суток восходящего узла=2.752881288528442383 наклонение орбиты=1.124948978424072266 эксцентриситет=0.00047016143798828125 аргумент перигея орбиты спутника=2.485048770904541016 время прохождения восходящего узла орбиты=27434032 драконический период=40544041.015625 скорость изменения драконического периода=-0.244140625 номер дня=1083
GLONASS: id=8 здоровье=0 номер литера=6 уход шкалы времени спутника=0.0457763671875 долгота первого внутри суток восходящего узла=2.38481283187866211 наклонение орбиты=1.12516176700592041 эксцентриситет=0.0001850128173828125 аргумент перигея орбиты спутника=2.661936044692993164 время прохождения восходящего узла орбиты=32490438 драконический период=40544046.875 253
скорость изменения драконического периода=-0.30517578125 номер дня=1083
GLONASS: id=9 здоровье=255 номер литера=0 уход шкалы времени спутника=0 долгота первого внутри суток восходящего узла=0 наклонение орбиты=1.099557399749755859 эксцентриситет=0 аргумент перигея орбиты спутника=0 время прохождения восходящего узла орбиты=0 драконический период=43200000 скорость изменения драконического периода=0 номер дня=1083
GLONASS: id=10 здоровье=255 номер литера=0 уход шкалы времени спутника=0 долгота первого внутри суток восходящего узла=0 наклонение орбиты=1.099557399749755859 эксцентриситет=0 аргумент перигея орбиты спутника=0 время прохождения восходящего узла орбиты=0 драконический период=43200000 скорость изменения драконического периода=0 номер дня=1083
GLONASS: id=11 здоровье=255 номер литера=0 уход шкалы времени спутника=0 долгота первого внутри суток восходящего узла=0 наклонение орбиты=1.099557399749755859 эксцентриситет=0 аргумент перигея орбиты спутника=0 время прохождения восходящего узла орбиты=0 драконический период=43200000 скорость изменения драконического периода=0 номер дня=1083
GLONASS: id=12 здоровье=255 номер литера=0 уход шкалы времени спутника=0 254

Смотрите также файлы