How do you indicate a line between a user and a GPS satellite?

The distance between the satellite and the user was obtained, but it was obtained with the ECEF coordinate system.

If you use ecef2enu, a function that changes the distance to the ENU coordinate system, errors continue to occur. I am curious about the solution, and furthermore, I want to use the unit vector to reduce the distance and then visualize the satellite-user distance.

receiverLat = 37.498759; 
receiverLon = 127.027487; 
receiverAlt = 0;         

t = datetime('now', 'TimeZone', 'Local');
[satPos, ~] = gnssconstellation(t); 

wgs84 = wgs84Ellipsoid('meters'); 
[recX, recY, recZ] = geodetic2ecef(wgs84, receiverLat, receiverLon, receiverAlt);

distances = zeros(size(satPos, 1), 1);

for i = 1:size(satPos, 1)
    satX = satPos(i, 1);
    satY = satPos(i, 2);
    satZ = satPos(i, 3);
    distances(i) = sqrt((satX - recX)^2 + (satY - recY)^2 + (satZ - recZ)^2);
end

disp('Satellite location:');
disp(satPos);

disp('Receiver location:');
disp([recX, recY, recZ]);

disp('Satellite-Receiver distances:');
disp(distances);

Answer 1

[William Rose]

@인혁

Why do you think there are errors?

I modified your script slightly to display the receiver radius from Earth center, and the distance from receiver to the nearest and farthest satellites, and the difference between the nearest and farthest distances. One would expect the difference betwen nearest and farthest to be on the order of 2 Earth radii (~2 x 6.37e6 m). It is.

receiverLat = 37.498759; 
receiverLon = 127.027487; 
receiverAlt = 0;         
t = datetime('now', 'TimeZone', 'Local');
[satPos, ~] = gnssconstellation(t); 
wgs84 = wgs84Ellipsoid('meters'); 
[recX, recY, recZ] = geodetic2ecef(wgs84, receiverLat, receiverLon, receiverAlt);
distances = zeros(size(satPos, 1), 1);
for i = 1:size(satPos, 1)
    satX = satPos(i, 1);
    satY = satPos(i, 2);
    satZ = satPos(i, 3);
    distances(i) = sqrt((satX - recX)^2 + (satY - recY)^2 + (satZ - recZ)^2);
end
%disp('Satellite location:');
%disp(satPos);
fprintf('Receiver location: %.3e, %.3e, %.3e m.\n',recX,recY,recZ);

Receiver location: -3.051e+06, 4.045e+06, 3.861e+06 m.

fprintf('Receiver radius from center = %.3e m.\n',sqrt(recX^2+recY^2+recZ^2));

Receiver radius from center = 6.370e+06 m.

%disp('Satellite-Receiver distances:');
fprintf('Closest=%.2e, farthest=%.2e, difference=%.2e.\n',...
    min(distances),max(distances),max(distances)-min(distances));

Closest=2.02e+07, farthest=3.27e+07, difference=1.25e+07.

[William Rose]

@인혁,

In the script below, I find the distances to the GPS satellites*, using ECEF coordinates, and I display the minimum and maximum distances. Then I find the closest satellite, and compute its ENU coordinates, using ecef2enu(). Then I compute and display the distance to that satellite, using the ENU coordinates. This distance matches the minimum distance displayed earlier, which was computed using the ECEF coordinates. There are not any obvious errors.

*I changed the distances calculation, so that it does not use a for loop.

receiverLat = 37.498759; 
receiverLon = 127.027487; 
receiverAlt = 0;         
t = datetime('now', 'TimeZone', 'Local');
[satPos, ~] = gnssconstellation(t); 
wgs84 = wgs84Ellipsoid('meters'); 
[recX, recY, recZ] = geodetic2ecef(wgs84, receiverLat, receiverLon, receiverAlt);
distances=((satPos(:,1)-recX).^2+(satPos(:,2)-recY).^2+(satPos(:,3)-recZ).^2).^0.5;
fprintf('Closest=%.2e, farthest=%.2e, difference=%.2e.\n',...
    min(distances),max(distances),max(distances)-min(distances));

Closest=2.05e+07, farthest=3.29e+07, difference=1.24e+07.

[~,idx]=min(distances);  %index of closest satellite 
% Next: find the ENU coords of closest satellite.
[satE,satN,satU]=ecef2enu(satPos(idx,1),satPos(idx,2),satPos(idx,3),...
    receiverLat,receiverLon,receiverAlt,wgs84);
satDist=sqrt(satE^2+satN^2+satU^2);  % distance to closest, using ENU
fprintf('Distance to satellite %d=%.3e m.',idx,satDist)

Distance to satellite 20=2.045e+07 m.

OK

[William Rose]

@인혁,

Here is an easier way to obtain the ranges.

llaSeoul = [37.498759,127.027487,0];  % Seoul lat lon alt         
[satPos,~] = gnssconstellation(datetime('now'));
distances=pseudoranges(llaSeoul,satPos);
[minDist,idx]=min(distances);  %index of closest satellite
fprintf('Satellite %d is closest. Distance= %.3e m.\n',...
    idx,minDist)

Satellite 8 is closest. Distance= 2.037e+07 m.

More later

[William Rose]

@인혁,

[Edit: Add comments in the code. Add markers at the ends of the line segments. Update screenshots to show the lines with markers.]

The code below displays the earth, with lines from a ground station to the four closest GPS satellites.

llaSeoul=[37.498759,127.027487,0];  % Seoul lat lon height
rEarth=6.37e6;                      % Earth radius (m)
[satPos,~] = gnssconstellation(datetime('now'));
distances=pseudoranges(llaSeoul,satPos);
[minDist,idx]=mink(distances,4);    % indices of 4 closest satellites
satPoslla=ecef2lla(satPos);         % lat,lon,alt for the satellites
line1=[];
% Loop below creates 4 line segments, from ground station to
% four satellites. Each segment is followed by a row of NaNs, 
% to disconnect the segments from each other.
for i=1:4
    line1=[line1;llaSeoul;satPoslla(idx(i),:);[NaN,NaN,NaN]];
end
uif = uifigure;
g = geoglobe(uif);
geoplot3(g,line1(:,1),line1(:,2),line1(:,3),"-ro",LineWidth=2,MarkerSize=3)
camheight(g,8*rEarth)               % set camera height

The code creates the image below. You can tell that it must be night time in Seoul, because the sun is visible behind the earth.

I can click and drag, or use my fingers, to rotate the image, or zoom in and out. For example, I can rotate to get the image below.

[inhyeok]

thank you so much!

[William Rose]

You’re welcome.