Conclusões Finais e Propostas de Continuidade 125
% Sérgio Luciano Avila – Sistema de Coordenadas – Posicionamento do Satélite
% Data base
load topo.mat % Load world bitmap load br % Load Brasil contour % longitude x latitude in radians
% Parameters
Ra = 6.37e6; % Average Earth radius (m) Rg = 42.2e6; % Geosynchronous orbit radius (m)
h = Rg-Ra; % Antenna height LONc = -51*pi/180; LATc = -13*pi/180; % Map center in LAT-LON
Pc = Ra*[cos(LATc)*cos(LONc) cos(LATc)*sin(LONc) sin(LATc)]; % Map center in cartesian
Laz = [-4 4]; % Azimute limits Lel = [-4 4]; % Elevation limits
% Transformation: earth coordinates Brasil contour LAT-LON to antenna coordinates (KEY OPERATION!!)
urc = [cos(LATc)*cos(LONc) cos(LATc)*sin(LONc) sin(LATc)]; % r axis at Pc
utetac = [sin(LATc)*cos(LONc) sin(LATc)*sin(LONc) -cos(LATc)]; % theta axis at Pc
uphic = [-sin(LONc) cos(LONc) 0]; % phi axis at Pc urhoc = [cos(LONc) sin(LONc) 0]; % rho axis at Pc
Pa = Rg*urhoc; % antenna position (translation) vector
uya = uphic; % antenna coordinates unit y axis (rotation matix)
uza = (Pc-Pa)/sqrt(sum((Pa-Pc).^2)); % antenna coordinates unit z axis
uxa = cross(uya,uza); % antenna coordinates unit x axis
% Apply transformation
brasile = Ra*[cos(brasil(:,2)).*cos(brasil(:,1)), ... % earth coordinates
cos(brasil(:,2)).*sin(brasil(:,1)), ...
sin(brasil(:,2))];
brasila = [brasile(:,1)-Pa(1), ... % antenna coordinates
brasile(:,2)-Pa(2), ...
brasile(:,3)-Pa(3)]; % translation
brasila = [dot(brasila,repmat(uxa,size(brasila,1),1),2), ...
dot(brasila,repmat(uya,size(brasila,1),1),2), ...
dot(brasila,repmat(uza,size(brasila,1),1),2)]; % rotation
AZbr = atan2(brasila(:,2), brasila(:,3)); ELbr = atan2(brasila(:,1), brasila(:,3)); % cartesian to Azimute-Elevation
% Map grid
[theta,phi] = ndgrid((-90:10:90)*pi/180,(-180:10:180)*pi/180); % Grid: latitude and longitude
x = Ra*cos(theta).*cos(phi); y = Ra*cos(theta).*sin(phi); z = Ra*sin(theta);
%--- 3D view % Plot globe
hfg = figure('color',[1 1 1],'NumberTitle','off','Name','3D View'); hold on
hg = surface(x,y,z,'FaceColor','texture','CData',[topo(:,181:360) topo(:,1:180)]); %earth surface
colormap(topomap1)
hl = line([0 Pa(1) Pc(1)],[0 Pa(2) Pc(2)],[0 Pa(3) Pc(3)],'LineStyle',':','Color',0*[1 1 1]); % quide lines
hc = plot3(Pc(1),Pc(2),Pc(3),'r.');
L = 2*Ra; % Plot earth axis
hl = line([0 0 0; L 0 0],[0 0 0; 0 L 0],[0 0 0; 0 0 L]);
text([L,0,0],[0,L,0],[0,0,L],['x'; 'y'; 'z'], 'Color', [0 0 0])
L = .5*Ra; % Plot antenna axis
hl = line([0 0 0; L*[uxa(1) uya(1) uza(1)]]+Pa(1),[0 0 0; L*[uxa(2) uya(2) uza(2)]]+Pa(2),[0 0 0; L*[uxa(3) uya(3)
uza(3)]]+Pa(3)); % antenna axis
text(L*[uxa(1) uya(1) uza(1)]+Pa(1), L*[uxa(2) uya(2) uza(2)]+Pa(2), L*[uxa(3) uya(3) uza(3)]+Pa(3),['xa'; 'ya'; 'za'], 'Color',
[0 0 0])
cam3d(gca),cam3d('reset'),view(180,0),set(gca,'units','norm'),alpha(.5),axis off
%--- Map view % Plot world map
grid on
[x,y] = ndgrid(0:10:360,-90:10:90);
hfm = figure('color',[1 1 1],'NumberTitle','off','Name','Flat map'); hold on
hm = surface(phi*180/pi,theta*180/pi,zeros(size(phi)),'FaceColor','texture','CData',[topo(:,181:360) topo(:,1:180)]);
colormap(topomap1)
title('World Map'),xlabel('longitude'),ylabel('latitude')
axis image
hc = plot(LONc*180/pi,LATc*180/pi,'r.'); % Plot center
hp = plot(brasil(:,1)*180/pi,brasil(:,2)*180/pi,'k'); % Plot Brasil contour
%--- CONBRA
hfc = figure('color',[1 1 1],'NumberTitle','off','Name','CONBRA'); hold on
hbrc = plot(AZbr*180/pi,ELbr*180/pi,'k','LineWidth',2);
hmjc = plot(AZmj*180/pi,ELmj*180/pi,'k','LineWidth',2);
axis equal, axis([Laz Lel]), grid on
title('CONBRA'),xlabel('AZ (º)'),ylabel('EL (º)','Rotation',0)
Figura A1.3. Pseudocódigo: Transformação de coordenadas.