50

I want to convert GPS location (latitude, longitude) into x,y coordinates. I found many links about this topic and applied it, but it doesn't give me the correct answer!

I am following these steps to test the answer: (1) firstly, i take two positions and calculate the distance between them using maps. (2) then convert the two positions into x,y coordinates. (3) then again calculate distance between the two points in the x,y coordinates and see if it give me the same result in point(1) or not.

one of the solution i found the following, but it doesn't give me correct answer!

latitude = Math.PI * latitude / 180;
longitude = Math.PI * longitude / 180;

// adjust position by radians
latitude -= 1.570795765134; // subtract 90 degrees (in radians)

// and switch z and y 
xPos = (app.radius) * Math.sin(latitude) * Math.cos(longitude);
zPos = (app.radius) * Math.sin(latitude) * Math.sin(longitude);
yPos = (app.radius) * Math.cos(latitude);

also i tried this link but still not work with me well!

any help how to convert from(latitude, longitude) to (x,y) ?

Thanks,

Pedro77
  • 5,176
  • 7
  • 61
  • 91
userInThisWorld
  • 1,361
  • 4
  • 18
  • 35

6 Answers6

74

No exact solution exists

There is no isometric map from the sphere to the plane. When you convert lat/lon coordinates from the sphere to x/y coordinates in the plane, you cannot hope that all lengths will be preserved by this operation. You have to accept some kind of deformation. Many different map projections do exist, which can achieve different compromises between preservations of lengths, angles and areas. For smallish parts of earth's surface, transverse Mercator is quite common. You might have heard about UTM. But there are many more.

The formulas you quote compute x/y/z, i.e. a point in 3D space. But even there you'd not get correct distances automatically. The shortest distance between two points on the surface of the sphere would go through that sphere, whereas distances on the earth are mostly geodesic lengths following the surface. So they will be longer.

Approximation for small areas

If the part of the surface of the earth which you want to draw is relatively small, then you can use a very simple approximation. You can simply use the horizontal axis x to denote longitude λ, the vertical axis y to denote latitude φ. The ratio between these should not be 1:1, though. Instead you should use cos(φ0) as the aspect ratio, where φ0 denotes a latitude close to the center of your map. Furthermore, to convert from angles (measured in radians) to lengths, you multiply by the radius of the earth (which in this model is assumed to be a sphere).

  • x = r λ cos(φ0)
  • y = r φ

This is simple equirectangular projection. In most cases, you'll be able to compute cos(φ0) only once, which makes subsequent computations of large numbers of points really cheap.

MvG
  • 57,380
  • 22
  • 148
  • 276
  • Thanks @MvG, I am not involved with these topics, i need just to draw points on a map, these points the distances between them don't exceed 50 meter, is there an easy way for that? – userInThisWorld Apr 29 '13 at 13:03
  • @alsadi90: I updated my answer to describe a simple approximation for your use case. – MvG Apr 29 '13 at 14:27
  • Thanks, but i have these questions: here λ and φ in radian? also about φ0, what do you suggest to me put it, if my area nearly 50mx50m – userInThisWorld Apr 29 '13 at 20:07
  • 2
    @alsadi90: Whether you measure your angles in degrees or radians only affects the size but not the shape of the map. You have to make sure to compute the cosine correctly, which in most programming languages means using radians. For φ₀ I'd suggest you simply take the average between the lower bound and the upper bound of your latitudes. But for such small maps, taking any latitude from within that area will likely work well enough. – MvG Apr 30 '13 at 04:48
  • Thanks for your help, but i applied the steps on this example:: point1(32.49194144785378,35.99031005054712) ,, point2(32.492285325127966,35.99031541496515) ,, the actual distance between them = 38.284 meter _ using this sight http://www.daftlogic.com/projects-google-maps-distance-calculator.htm _ then i applied the steps assuming φ₀ = 32.492113386490873 which is the avg of lat of the two points ,, then in x,y coordination i got point1(0.3559,0.5668) and ppoint2(0.3559,0.5668), then distance between them = 5.9990e-06 !! which is wrong value! – userInThisWorld Apr 30 '13 at 07:16
  • @alsadi90: If you want the result in meters, then you have to multiply radians by the radius of the earth. Updated my answer accordingly. – MvG Apr 30 '13 at 08:20
18

I want to share with you how I managed the problem. I've used the equirectangular projection just like @MvG said, but this method gives you X and Y positions related to the globe (or the entire map), this means that you get global positions. In my case, I wanted to convert coordinates in a small area (about 500m square), so I related the projection point to another 2 points, getting the global positions and relating to local (on screen) positions, just like this:

First, I choose 2 points (top-left and bottom-right) around the area where I want to project, just like this picture:

enter image description here

Once I have the global reference area in lat and lng, I do the same for screen positions. The objects containing this data are shown below.

//top-left reference point
var p0 = {
    scrX: 23.69,        // Minimum X position on screen
    scrY: -0.5,         // Minimum Y position on screen
    lat: -22.814895,    // Latitude
    lng: -47.072892     // Longitude
}

//bottom-right reference point
var p1 = {
    scrX: 276,          // Maximum X position on screen
    scrY: 178.9,        // Maximum Y position on screen
    lat: -22.816419,    // Latitude
    lng: -47.070563     // Longitude
}

var radius = 6371;      //Earth Radius in Km

//## Now I can calculate the global X and Y for each reference point ##\\

// This function converts lat and lng coordinates to GLOBAL X and Y positions
function latlngToGlobalXY(lat, lng){
    //Calculates x based on cos of average of the latitudes
    let x = radius*lng*Math.cos((p0.lat + p1.lat)/2);
    //Calculates y based on latitude
    let y = radius*lat;
    return {x: x, y: y}
}

// Calculate global X and Y for top-left reference point
p0.pos = latlngToGlobalXY(p0.lat, p0.lng);
// Calculate global X and Y for bottom-right reference point
p1.pos = latlngToGlobalXY(p1.lat, p1.lng);

/*
* This gives me the X and Y in relation to map for the 2 reference points.
* Now we have the global AND screen areas and then we can relate both for the projection point.
*/

// This function converts lat and lng coordinates to SCREEN X and Y positions
function latlngToScreenXY(lat, lng){
    //Calculate global X and Y for projection point
    let pos = latlngToGlobalXY(lat, lng);
    //Calculate the percentage of Global X position in relation to total global width
    pos.perX = ((pos.x-p0.pos.x)/(p1.pos.x - p0.pos.x));
    //Calculate the percentage of Global Y position in relation to total global height
    pos.perY = ((pos.y-p0.pos.y)/(p1.pos.y - p0.pos.y));

    //Returns the screen position based on reference points
    return {
        x: p0.scrX + (p1.scrX - p0.scrX)*pos.perX,
        y: p0.scrY + (p1.scrY - p0.scrY)*pos.perY
    }
}

//# The usage is like this #\\

var pos = latlngToScreenXY(-22.815319, -47.071718);
$point = $("#point-to-project");
$point.css("left", pos.x+"em");
$point.css("top", pos.y+"em");

As you can see, I made this in javascript, but the calculations can be translated to any language.

P.S. I'm applying the converted positions to an HTML element whose id is "point-to-project". To use this piece of code on your project, you shall create this element (styled as position absolute) or change the "usage" block.

allexrm
  • 181
  • 1
  • 2
  • 2
    There's an error: the radius of the Earth in km is `6371`, not `6.371` (which in JS is a little more than 6). Of course over such a small area it doesn't change much, but future readers should take note of it. – Giulio Muscarello Sep 08 '19 at 10:32
4

Since this page shows up on top of google while i searched for this same problem, I would like to provide a more practical answers. The answer by MVG is correct but rather theoratical.

I have made a track plotting app for the fitbit ionic in javascript. The code below is how I tackled the problem.

//LOCATION PROVIDER
index.js
var gpsFix = false;
var circumferenceAtLat = 0;
function locationSuccess(pos){
  if(!gpsFix){
    gpsFix = true;
    circumferenceAtLat = Math.cos(pos.coords.latitude*0.01745329251)*111305;
  }
  pos.x:Math.round(pos.coords.longitude*circumferenceAtLat),
  pos.y:Math.round(pos.coords.latitude*110919), 
  plotTrack(pos);
}

plotting.js

plotTrack(position){

let x = Math.round((this.segments[i].start.x - this.bounds.minX)*this.scale);
let y = Math.round(this.bounds.maxY - this.segments[i].start.y)*this.scale; //heights needs to be inverted

//redraw?
let redraw = false;

//x or y bounds?
 if(position.x>this.bounds.maxX){
   this.bounds.maxX = (position.x-this.bounds.minX)*1.1+this.bounds.minX; //increase by 10%
   redraw = true;
 }
 if(position.x<this.bounds.minX){
   this.bounds.minX = this.bounds.maxX-(this.bounds.maxX-position.x)*1.1;
    redraw = true;
 };
 if(position.y>this.bounds.maxY){
   this.bounds.maxY = (position.y-this.bounds.minY)*1.1+this.bounds.minY; //increase by 10%
    redraw = true;
 }
 if(position.y<this.bounds.minY){
   this.bounds.minY = this.bounds.maxY-(this.bounds.maxY-position.y)*1.1;
    redraw = true;
 }
 if(redraw){
   reDraw();
 }
}


function reDraw(){

let xScale = device.screen.width / (this.bounds.maxX-this.bounds.minX);
let yScale = device.screen.height / (this.bounds.maxY-this.bounds.minY); 
if(xScale<yScale) this.scale = xScale; 
else this.scale = yScale;

//Loop trough your object to redraw all of them
}
  • Are you converting to SI (km) or Imperial (mi) unit? – Cloud Cho Jan 11 '20 at 09:22
  • 1
    Hi @CloudCho. I use metric, as any programming language does. I do a conversion only for visual reasons in the very end of my code.111305 is the difference in meters between 1 degree of latitude. – Pieter Oskam Jan 15 '22 at 22:15
3

For completeness I like to add my python adaption of @allexrm code which worked really well. Thanks again!

radius = 6371    #Earth Radius in KM

class referencePoint:
    def __init__(self, scrX, scrY, lat, lng):
        self.scrX = scrX
        self.scrY = scrY
        self.lat = lat
        self.lng = lng


# Calculate global X and Y for top-left reference point        
p0 = referencePoint(0, 0, 52.526470, 13.403215)
# Calculate global X and Y for bottom-right reference point
p1 = referencePoint(2244, 2060, 52.525035, 13.405809) 


# This function converts lat and lng coordinates to GLOBAL X and Y positions
def latlngToGlobalXY(lat, lng):
    # Calculates x based on cos of average of the latitudes
    x = radius*lng*math.cos((p0.lat + p1.lat)/2)
    # Calculates y based on latitude
    y = radius*lat
    return {'x': x, 'y': y}


# This function converts lat and lng coordinates to SCREEN X and Y positions
def latlngToScreenXY(lat, lng):
    # Calculate global X and Y for projection point
    pos = latlngToGlobalXY(lat, lng)
    # Calculate the percentage of Global X position in relation to total global width
    perX = ((pos['x']-p0.pos['x'])/(p1.pos['x'] - p0.pos['x']))
    # Calculate the percentage of Global Y position in relation to total global height
    perY = ((pos['y']-p0.pos['y'])/(p1.pos['y'] - p0.pos['y']))

    # Returns the screen position based on reference points
    return {
        'x': p0.scrX + (p1.scrX - p0.scrX)*perX,
        'y': p0.scrY + (p1.scrY - p0.scrY)*perY
    }


pos = latlngToScreenXY(52.525607, 13.404572);

pos['x] and pos['y] contain the translated x & y coordinates of the lat & lng (52.525607, 13.404572)

I hope this is helpful for anyone looking like me for the proper solution to the problem of translating lat lng into a local reference coordinate system.

Best

Alexander
  • 1,422
  • 11
  • 15
2

Its better to convert to utm coordinates, and treat that as x and y.

import utm
u = utm.from_latlon(12.917091, 77.573586)

The result will be (779260.623156606, 1429369.8665238516, 43, 'P') The first two can be treated as x,y coordinates, the 43P is the UTM Zone, which can be ignored for small areas (width upto 668 km).

  • 1
    I would assume that just converting spherical coordinates to UTM and then ignoring the zone might cause problems for regions that cross zone boundaries. You'd get results from different zones that don't fit together. There should be a way to force all points to use the same zone, then they should fit better. – MvG May 01 '22 at 08:24
0

The easiest way to do it is to use Mercator Projection. Be careful that there might be some stretching around the poles.

And the reference code (Kotlin):

fun latLngToXY(point: LatLng): Pair<Double, Double> {
    // Converts to radians.
    val latRad = Math.toRadians(point.latitude)
    val longRad = Math.toRadians(point.longitude)
    
    // Keep the x axis.
    val x = longRad

    // Mercator
    val y = ln(tan(latRad) - (1 / cos(latRad)))

    return Pair(x, y)
}

Note that it returns values in [0;pi/2] so you might want to scale it.

DatGeoudon
  • 342
  • 3
  • 15