PX4-Autopilot icon indicating copy to clipboard operation
PX4-Autopilot copied to clipboard

Geocoordinate projection errors by assuming the earth is a sphere

Open Jaeyoung-Lim opened this issue 3 years ago • 13 comments
trafficstars

Describe the bug

The current Map projection assumes that the earth is a sphere,

https://github.com/PX4/PX4-Autopilot/blob/9246d38667311b8f53c3b38a13e35bf4a3969f70/src/lib/geo/geo.cpp#L71-L92

However, this can result in large errors when projecting WGS84 coordinates to a tangential local coordinate frame. The projection should be done through a geodetic to cartesian conversion as in https://proj.org/operations/conversions/cart.html#cart

To Reproduce

Steps to reproduce the behavior:

  1. Run SITL (make px4_sitl gazebo_plane) with the following branch: https://github.com/PX4/PX4-SITL_gazebo/tree/pr-geographic-lib
  2. Fly out 20km or more and compare the geoprojection
  3. Compare the projection errors using Geographiclib using topocentric projection vs projection using spherical earth assumptions

Expected behavior

A clear and concise description of what you expected to happen.

Log Files and Screenshots

The figure below shows the difference between the altitude differences when flying 20km out between the spherical assumption, and a conversion using geographic lib

The log was generated in SITL Gazebo with the plane model. sitl_gazebo uses the same projection: https://github.com/PX4/PX4-SITL_gazebo/blob/3162519888b825a167fc3fff49c32903f81f7bce/include/common.h#L260-L280

Logged topics show

  • /vehicle_global_position/alt : GeographicLib conversion
  • /vehicle_global_position_groundtruth/alt: spherical assumption (PX4)

It can be clearly seen that the altitude errors amount to ~30m at 20km distance image

Always provide a link to the flight log file:

  • Download the flight log file from the vehicle (tutorial).
  • Upload the log to the PX4 Flight Review
  • Share the link to the log (Copy and paste the URL of the log)

Jaeyoung-Lim avatar Nov 09 '22 09:11 Jaeyoung-Lim

Is there any follow up on this topic? This seems like an important issue to be resolved.

murundb avatar Feb 03 '24 20:02 murundb

Has the earth spheroid been updated since 1984?

mcsauder avatar Feb 03 '24 23:02 mcsauder

You know.... wgs 84?

mcsauder avatar Feb 03 '24 23:02 mcsauder

Once in a while you need an old guy to remind you that the earth is not burning down. Someone born before 1980.

mcsauder avatar Feb 03 '24 23:02 mcsauder

The magnetic fields, are, however, about to flip.

mcsauder avatar Feb 04 '24 00:02 mcsauder

@mcsauder I think the earth has been burning down in a lot of places because of climate change. The current projection to NED is not using WGS84 reference ellipsoid. The issue was further discussed here and it seems like the altitude error is reasonable.

murundb avatar Feb 04 '24 00:02 murundb

I don't agree, but that's fine. :)

mcsauder avatar Feb 04 '24 00:02 mcsauder

Which part you don't agree with?

murundb avatar Feb 04 '24 00:02 murundb

It's not new news that the earth is a spheroid, not a sphere. What do mean?

mcsauder avatar Feb 04 '24 00:02 mcsauder

It's a thing that happens because the earth has mass and spins on an axis. It's wider around the middle than it is from top to bottom.

mcsauder avatar Feb 04 '24 00:02 mcsauder

If it's not wgs84 and you need to fly more than 300miles, submit a pr. And yes,the earth is not a sphere. @Jaeyoung-Lim is quite capable of adding the math, he is as sharp as they get. More than 300miles is military and you can get yourself into trouble quick.

mcsauder avatar Feb 04 '24 01:02 mcsauder

I still don't understand what you don't agree with or what you are trying to tell me. MapProjection::project() is not doing standard curvilinear (with respect to WGS84 reference ellipsoid) to NED local tangent plane. Assuming you have established the reference origin of your local tangent plane, you would do transformations from curvilinear to ECEF to NED local tangent plane.

I am not saying that there is something wrong with MapProjection::project(). For all I know, it might have been the original intent to not use WGS84 datum. Another observation is that for the standard conversion, you'd use the altitude above the reference ellipsoid, but in MapProjection::project(), altitude above mean sea level is being used.

murundb avatar Feb 04 '24 02:02 murundb

Can you do it with wgs84?

mcsauder avatar Feb 04 '24 02:02 mcsauder