Convert latitude / longitude to pixels (x, y) on the projection of the Mercator

Question

I am trying to convert a lat / long point to a 2d point so that I can display it on a world image - a projection of a mercator.

I saw various ways to do this and a few questions about the stack overflow - I tried different pieces of code, and although I get the correct longitude to the pixel, the latitude is always off - it seems more reasonable.

I need a formula to account for image size, width, etc.

I tried this piece of code:

double minLat = -85.05112878; double minLong = -180; double maxLat = 85.05112878; double maxLong = 180; // Map image size (in points) Double mapHeight = 768.0; Double mapWidth = 991.0; // Determine the map scale (points per degree) double xScale = mapWidth/ (maxLong - minLong); double yScale = mapHeight / (maxLat - minLat); // position of map image for point double x = (lon - minLong) * xScale; double y = - (lat + minLat) * yScale; System.out.println("final coords: " + x + " " + y); 

In the example I'm trying, latitude seems to be off by about 30 pixels. Any help or advice?

Update

Based on this question: Lat / lon to xy

I tried using the provided code, but I still have problems converting the latitude, longitude is ok.

  int mapWidth = 991; int mapHeight = 768; double mapLonLeft = -180; double mapLonRight = 180; double mapLonDelta = mapLonRight - mapLonLeft; double mapLatBottom = -85.05112878; double mapLatBottomDegree = mapLatBottom * Math.PI / 180; double worldMapWidth = ((mapWidth / mapLonDelta) * 360) / (2 * Math.PI); double mapOffsetY = (worldMapWidth / 2 * Math.log((1 + Math.sin(mapLatBottomDegree)) / (1 - Math.sin(mapLatBottomDegree)))); double x = (lon - mapLonLeft) * (mapWidth / mapLonDelta); double y = 0.1; if (lat < 0) { lat = lat * Math.PI / 180; y = mapHeight - ((worldMapWidth / 2 * Math.log((1 + Math.sin(lat)) / (1 - Math.sin(lat)))) - mapOffsetY); } else if (lat > 0) { lat = lat * Math.PI / 180; lat = lat * -1; y = mapHeight - ((worldMapWidth / 2 * Math.log((1 + Math.sin(lat)) / (1 - Math.sin(lat)))) - mapOffsetY); System.out.println("y before minus: " + y); y = mapHeight - y; } else { y = mapHeight / 2; } System.out.println(x); System.out.println(y); 

When using the source code, if the latitude value is positive, it returned a negative point, so I changed it a bit and tested it with extreme latitudes, which should be point 0 and point 766, it works fine. However, when I try a different ex latitude value: 58.07 (just north of the UK), it displays as northern Spain.

Answer 1

The projection of the Mercator map is a special limiting case of the projection of a Lambert-Conic conformal mapping projector with equator as a single standard parallel. All other latitude parallels - straight lines and meridians are also straight lines at right angles to the equator, evenly distributed from each other. This is the basis for transverse and oblique projection shapes. It is little used for land mapping purposes, but is practically universal for navigational charts. Like conformal, it has the property that the straight lines drawn on it are lines of a permanent bearing. In this way, navigators can infer their course from the angle of a straight line the course does with the meridians. ^[one].

The formulas for deriving the predicted East and North coordinates from spherical latitude φ and longitude λ are:

 E = FE + R (λ – λₒ) N = FN + R ln[tan(π/4 + φ/2)] 

where λ _O is the longitude of natural origin, and FE and FN are false east and false north. In the spherical Mercator, these values ​​are not actually used, so you can simplify the formula to

An example of pseudocode, so this can be adapted to each programming language.

 latitude = 41.145556; // (φ) longitude = -73.995; // (λ) mapWidth = 200; mapHeight = 100; // get x value x = (longitude+180)*(mapWidth/360) // convert from degrees to radians latRad = latitude*PI/180; // get y value mercN = log(tan((PI/4)+(latRad/2))); y = (mapHeight/2)-(mapWidth*mercN/(2*PI)); 

Sources:

• OGP Geomatics Committee, Handbook No. 7, Part 2: Coordinate Transformation and Transformation
• Mercator projection output
• National Atlas: Map Projections
• Mercator Map Forecast

EDIT Created a working example in PHP (because I suck in Java)

https://github.com/mfeldheim/mapStuff.git

Answer 2

You cannot just transfer from longitude / latitude to x / y like this because the world is not flat. Are you looking at this post? Convert Longitude / Latitude to X / Y Coordinate

UPDATE - 1/18/13

I decided to give him a punch, and here's how I do it: -

 public class MapService { // CHANGE THIS: the output path of the image to be created private static final String IMAGE_FILE_PATH = "/some/user/path/map.png"; // CHANGE THIS: image width in pixel private static final int IMAGE_WIDTH_IN_PX = 300; // CHANGE THIS: image height in pixel private static final int IMAGE_HEIGHT_IN_PX = 500; // CHANGE THIS: minimum padding in pixel private static final int MINIMUM_IMAGE_PADDING_IN_PX = 50; // formula for quarter PI private final static double QUARTERPI = Math.PI / 4.0; // some service that provides the county boundaries data in longitude and latitude private CountyService countyService; public void run() throws Exception { // configuring the buffered image and graphics to draw the map BufferedImage bufferedImage = new BufferedImage(IMAGE_WIDTH_IN_PX, IMAGE_HEIGHT_IN_PX, BufferedImage.TYPE_INT_RGB); Graphics2D g = bufferedImage.createGraphics(); Map<RenderingHints.Key, Object> map = new HashMap<RenderingHints.Key, Object>(); map.put(RenderingHints.KEY_INTERPOLATION, RenderingHints.VALUE_INTERPOLATION_BICUBIC); map.put(RenderingHints.KEY_RENDERING, RenderingHints.VALUE_RENDER_QUALITY); map.put(RenderingHints.KEY_ANTIALIASING, RenderingHints.VALUE_ANTIALIAS_ON); RenderingHints renderHints = new RenderingHints(map); g.setRenderingHints(renderHints); // min and max coordinates, used in the computation below Point2D.Double minXY = new Point2D.Double(-1, -1); Point2D.Double maxXY = new Point2D.Double(-1, -1); // a list of counties where each county contains a list of coordinates that form the county boundary Collection<Collection<Point2D.Double>> countyBoundaries = new ArrayList<Collection<Point2D.Double>>(); // for every county, convert the longitude/latitude to X/Y using Mercator projection formula for (County county : countyService.getAllCounties()) { Collection<Point2D.Double> lonLat = new ArrayList<Point2D.Double>(); for (CountyBoundary countyBoundary : county.getCountyBoundaries()) { // convert to radian double longitude = countyBoundary.getLongitude() * Math.PI / 180; double latitude = countyBoundary.getLatitude() * Math.PI / 180; Point2D.Double xy = new Point2D.Double(); xy.x = longitude; xy.y = Math.log(Math.tan(QUARTERPI + 0.5 * latitude)); // The reason we need to determine the min X and Y values is because in order to draw the map, // we need to offset the position so that there will be no negative X and Y values minXY.x = (minXY.x == -1) ? xy.x : Math.min(minXY.x, xy.x); minXY.y = (minXY.y == -1) ? xy.y : Math.min(minXY.y, xy.y); lonLat.add(xy); } countyBoundaries.add(lonLat); } // readjust coordinate to ensure there are no negative values for (Collection<Point2D.Double> points : countyBoundaries) { for (Point2D.Double point : points) { point.x = point.x - minXY.x; point.y = point.y - minXY.y; // now, we need to keep track the max X and Y values maxXY.x = (maxXY.x == -1) ? point.x : Math.max(maxXY.x, point.x); maxXY.y = (maxXY.y == -1) ? point.y : Math.max(maxXY.y, point.y); } } int paddingBothSides = MINIMUM_IMAGE_PADDING_IN_PX * 2; // the actual drawing space for the map on the image int mapWidth = IMAGE_WIDTH_IN_PX - paddingBothSides; int mapHeight = IMAGE_HEIGHT_IN_PX - paddingBothSides; // determine the width and height ratio because we need to magnify the map to fit into the given image dimension double mapWidthRatio = mapWidth / maxXY.x; double mapHeightRatio = mapHeight / maxXY.y; // using different ratios for width and height will cause the map to be stretched. So, we have to determine // the global ratio that will perfectly fit into the given image dimension double globalRatio = Math.min(mapWidthRatio, mapHeightRatio); // now we need to readjust the padding to ensure the map is always drawn on the center of the given image dimension double heightPadding = (IMAGE_HEIGHT_IN_PX - (globalRatio * maxXY.y)) / 2; double widthPadding = (IMAGE_WIDTH_IN_PX - (globalRatio * maxXY.x)) / 2; // for each country, draw the boundary using polygon for (Collection<Point2D.Double> points : countyBoundaries) { Polygon polygon = new Polygon(); for (Point2D.Double point : points) { int adjustedX = (int) (widthPadding + (point.getX() * globalRatio)); // need to invert the Y since 0,0 starts at top left int adjustedY = (int) (IMAGE_HEIGHT_IN_PX - heightPadding - (point.getY() * globalRatio)); polygon.addPoint(adjustedX, adjustedY); } g.drawPolygon(polygon); } // create the image file ImageIO.write(bufferedImage, "PNG", new File(IMAGE_FILE_PATH)); } } 

RESULT: image width = 600 pixels, image height = 600 pixels, image fill = 50 pixels

RESULT: image width = 300 pixels, image height = 500 pixels, image fill = 50 pixels

Answer 3

The Java version of the original Google Maps API JavaScript API version 3 java script is as follows: it works without problems

 public final class GoogleMapsProjection2 { private final int TILE_SIZE = 256; private PointF _pixelOrigin; private double _pixelsPerLonDegree; private double _pixelsPerLonRadian; public GoogleMapsProjection2() { this._pixelOrigin = new PointF(TILE_SIZE / 2.0,TILE_SIZE / 2.0); this._pixelsPerLonDegree = TILE_SIZE / 360.0; this._pixelsPerLonRadian = TILE_SIZE / (2 * Math.PI); } double bound(double val, double valMin, double valMax) { double res; res = Math.max(val, valMin); res = Math.min(res, valMax); return res; } double degreesToRadians(double deg) { return deg * (Math.PI / 180); } double radiansToDegrees(double rad) { return rad / (Math.PI / 180); } PointF fromLatLngToPoint(double lat, double lng, int zoom) { PointF point = new PointF(0, 0); point.x = _pixelOrigin.x + lng * _pixelsPerLonDegree; // Truncating to 0.9999 effectively limits latitude to 89.189. This is // about a third of a tile past the edge of the world tile. double siny = bound(Math.sin(degreesToRadians(lat)), -0.9999,0.9999); point.y = _pixelOrigin.y + 0.5 * Math.log((1 + siny) / (1 - siny)) *- _pixelsPerLonRadian; int numTiles = 1 << zoom; point.x = point.x * numTiles; point.y = point.y * numTiles; return point; } PointF fromPointToLatLng(PointF point, int zoom) { int numTiles = 1 << zoom; point.x = point.x / numTiles; point.y = point.y / numTiles; double lng = (point.x - _pixelOrigin.x) / _pixelsPerLonDegree; double latRadians = (point.y - _pixelOrigin.y) / - _pixelsPerLonRadian; double lat = radiansToDegrees(2 * Math.atan(Math.exp(latRadians)) - Math.PI / 2); return new PointF(lat, lng); } public static void main(String []args) { GoogleMapsProjection2 gmap2 = new GoogleMapsProjection2(); PointF point1 = gmap2.fromLatLngToPoint(41.850033, -87.6500523, 15); System.out.println(point1.x+" "+point1.y); PointF point2 = gmap2.fromPointToLatLng(point1,15); System.out.println(point2.x+" "+point2.y); } } public final class PointF { public double x; public double y; public PointF(double x, double y) { this.x = x; this.y = y; } } 

Answer 4

I would like to indicate that the code within the procedure should read

  double bound(double val, double valMin, double valMax) { double res; res = Math.max(val, valMin); res = Math.min(res, valMax); return res; } 

Answer 5

  public static String getTileNumber(final double lat, final double lon, final int zoom) { int xtile = (int)Math.floor( (lon + 180) / 360 * (1<<zoom) ) ; int ytile = (int)Math.floor( (1 - Math.log(Math.tan(Math.toRadians(lat)) + 1 / Math.cos(Math.toRadians(lat))) / Math.PI) / 2 * (1<<zoom) ) ; if (xtile < 0) xtile=0; if (xtile >= (1<<zoom)) xtile=((1<<zoom)-1); if (ytile < 0) ytile=0; if (ytile >= (1<<zoom)) ytile=((1<<zoom)-1); return("" + zoom + "/" + xtile + "/" + ytile); } }