How to get bounds of a google static map?
Asked Answered
C

7

19

How to get bounds in degrees of google static map which has been returned, for example, for following request

http://maps.googleapis.com/maps/api/staticmap?center=0.0,0.0&zoom=10&size=640x640&sensor=false

As I know, full Earth map is 256x256 image. This means that n vertical pixels contain x degrees, but n horizontal pixels contain 2x degrees. Right?

As google says center defines the center of the map, equidistant from all edges of the map. As I understood equidistant in pixels (or in degrees?). And each succeeding zoom level doubles the precision in both horizontal and vertical dimensions. So, I can find delta value of Longitude of map for each zoom value as:

dLongitude = (HorizontalMapSizeInPixels / 256 ) * ( 360 / pow(2, zoom) );

Same calculations for Latitude:

dLatitude = (VerticalMapSizeInPixels / 256 ) * ( 180 / pow(2, zoom) );

VerticalMapSizeInPixels and HorizontalMapSizeInPixels are parameters of map size in URL.

It's good to calculate delta value of Longitude, but for Latitude it is wrong. I cannot find delta value of Latitude, there is some delta error.

Chromomere answered 20/9, 2012 at 6:43 Comment(1)
Possible duplicate of #4731385Center
I
34

As I know, full Earth map is 256x256 image.

Yes.

This means that n vertical pixels contain x degrees, but n horizontal pixels contain 2x degrees. Right?

No. One pixel will represent varying amounts of latitude depending on the latitude. One pixel at the Equator represents less latitude than one pixel near the poles.

The corners of the map will depend on center, zoom level and map size, and you'd need to use the Mercator projection to calculate them. If you don't want to load the full API, here's a MercatorProjection object:

var MERCATOR_RANGE = 256;

function bound(value, opt_min, opt_max) {
  if (opt_min != null) value = Math.max(value, opt_min);
  if (opt_max != null) value = Math.min(value, opt_max);
  return value;
}

function degreesToRadians(deg) {
  return deg * (Math.PI / 180);
}

function radiansToDegrees(rad) {
  return rad / (Math.PI / 180);
}

function MercatorProjection() {
  this.pixelOrigin_ = new google.maps.Point( MERCATOR_RANGE / 2, MERCATOR_RANGE / 2);
  this.pixelsPerLonDegree_ = MERCATOR_RANGE / 360;
  this.pixelsPerLonRadian_ = MERCATOR_RANGE / (2 * Math.PI);
};

MercatorProjection.prototype.fromLatLngToPoint = function(latLng, opt_point) {
  var me = this;

  var point = opt_point || new google.maps.Point(0, 0);

  var origin = me.pixelOrigin_;
  point.x = origin.x + latLng.lng() * me.pixelsPerLonDegree_;
  // NOTE(appleton): 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.
  var siny = bound(Math.sin(degreesToRadians(latLng.lat())), -0.9999, 0.9999);
  point.y = origin.y + 0.5 * Math.log((1 + siny) / (1 - siny)) * -me.pixelsPerLonRadian_;
  return point;
};

MercatorProjection.prototype.fromPointToLatLng = function(point) {
  var me = this;

  var origin = me.pixelOrigin_;
  var lng = (point.x - origin.x) / me.pixelsPerLonDegree_;
  var latRadians = (point.y - origin.y) / -me.pixelsPerLonRadian_;
  var lat = radiansToDegrees(2 * Math.atan(Math.exp(latRadians)) - Math.PI / 2);
  return new google.maps.LatLng(lat, lng);
};

//pixelCoordinate = worldCoordinate * Math.pow(2,zoomLevel)

You can save that to a separate file, for example "MercatorProjection.js", and then include it in your application.

<script src="MercatorProjection.js"></script>

With the above file loaded, the following function calculates the SW and NE corners of a map of a given size and at a given zoom.

function getCorners(center,zoom,mapWidth,mapHeight){
    var scale = Math.pow(2,zoom);
    var centerPx = proj.fromLatLngToPoint(center);
    var SWPoint = {x: (centerPx.x -(mapWidth/2)/ scale) , y: (centerPx.y + (mapHeight/2)/ scale)};
    var SWLatLon = proj.fromPointToLatLng(SWPoint);
    alert('SW: ' + SWLatLon);
    var NEPoint = {x: (centerPx.x +(mapWidth/2)/ scale) , y: (centerPx.y - (mapHeight/2)/ scale)};
    var NELatLon = proj.fromPointToLatLng(NEPoint);
    alert(' NE: '+ NELatLon);
}

and you'd call it like this:

var proj = new MercatorProjection();
var G = google.maps;
var centerPoint = new G.LatLng(49.141404, -121.960988);
var zoom = 10;
getCorners(centerPoint,zoom,640,640);
Indoiranian answered 20/9, 2012 at 11:38 Comment(4)
How precise is this? I ported it to java, I must say I'm really happy with it! I'm only a few meters of. I wonder if that is the algorithm or me making a mistake. For some reason I had to multiple the scale in getCorners with 2.0 else the bounds where way to big.Kab
GOT IT I had a request with size=1024x1024. It gave me an image of 1080x1080 which I was not aware of! So I did bounds checking on a to small region. Resulting in a to small region...Kab
Hey @marcelo, I have a problem with your code translated to C#. More info on the problem hereDartmoor
How would I use this to go from a lat long within the view to a pixel x, y?Adermin
F
16

Thanks Marcelo for your answer. It has been quite helpful. In case anyone would be interested, here the Python version of the code (a rough translation of the PHP code, probably not as pythonic as it could be):

from __future__ import division
import math
MERCATOR_RANGE = 256

def  bound(value, opt_min, opt_max):
  if (opt_min != None): 
    value = max(value, opt_min)
  if (opt_max != None): 
    value = min(value, opt_max)
  return value


def  degreesToRadians(deg) :
  return deg * (math.pi / 180)


def  radiansToDegrees(rad) :
  return rad / (math.pi / 180)


class G_Point :
    def __init__(self,x=0, y=0):
        self.x = x
        self.y = y



class G_LatLng :
    def __init__(self,lt, ln):
        self.lat = lt
        self.lng = ln


class MercatorProjection :


    def __init__(self) :
      self.pixelOrigin_ =  G_Point( MERCATOR_RANGE / 2, MERCATOR_RANGE / 2)
      self.pixelsPerLonDegree_ = MERCATOR_RANGE / 360
      self.pixelsPerLonRadian_ = MERCATOR_RANGE / (2 * math.pi)


    def fromLatLngToPoint(self, latLng, opt_point=None) :
      point = opt_point if opt_point is not None else G_Point(0,0)
      origin = self.pixelOrigin_
      point.x = origin.x + latLng.lng * self.pixelsPerLonDegree_
      # NOTE(appleton): 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.
      siny = bound(math.sin(degreesToRadians(latLng.lat)), -0.9999, 0.9999)
      point.y = origin.y + 0.5 * math.log((1 + siny) / (1 - siny)) * -     self.pixelsPerLonRadian_
      return point


def fromPointToLatLng(self,point) :
      origin = self.pixelOrigin_
      lng = (point.x - origin.x) / self.pixelsPerLonDegree_
      latRadians = (point.y - origin.y) / -self.pixelsPerLonRadian_
      lat = radiansToDegrees(2 * math.atan(math.exp(latRadians)) - math.pi / 2)
      return G_LatLng(lat, lng)

#pixelCoordinate = worldCoordinate * pow(2,zoomLevel)

def getCorners(center, zoom, mapWidth, mapHeight):
    scale = 2**zoom
    proj = MercatorProjection()
    centerPx = proj.fromLatLngToPoint(center)
    SWPoint = G_Point(centerPx.x-(mapWidth/2)/scale, centerPx.y+(mapHeight/2)/scale)
    SWLatLon = proj.fromPointToLatLng(SWPoint)
    NEPoint = G_Point(centerPx.x+(mapWidth/2)/scale, centerPx.y-(mapHeight/2)/scale)
    NELatLon = proj.fromPointToLatLng(NEPoint)
    return {
        'N' : NELatLon.lat,
        'E' : NELatLon.lng,
        'S' : SWLatLon.lat,
        'W' : SWLatLon.lng,
    }

Usage :

>>> import MercatorProjection
>>> centerLat = 49.141404
>>> centerLon = -121.960988
>>> zoom = 10
>>> mapWidth = 640
>>> mapHeight = 640
>>> centerPoint = MercatorProjection.G_LatLng(centerLat, centerLon)
>>> corners = MercatorProjection.getCorners(centerPoint, zoom, mapWidth, mapHeight)
>>> corners
{'E': -65.710988,
'N': 74.11120692972199,
'S': 0.333879313530149,
'W': -178.210988}
>>> mapURL = "http://maps.googleapis.com/maps/api/staticmap?center=%f,%f&zoom=%d&size=%dx%d&scale=2&maptype=roadmap&sensor=false"%(centerLat,centerLon,zoom,mapWidth,mapHeight)
>>> mapURL
http://maps.googleapis.com/maps/api/staticmap?center=49.141404,-121.960988&zoom=10&size=640x640&scale=2&maptype=roadmap&sensor=false'
Fitful answered 14/1, 2014 at 14:49 Comment(1)
Buggy : zoom should be defined as 2**scale, not 2^scale which is a bitwise XOR between 2 and scaleFitful
R
5

Simple Python Version

Having spent a long time on a similar problem and using this thread for help, here is a simplified Python version of Marcelo (and Jmague)'s code:

import math
import requests

def latLngToPoint(mapWidth, mapHeight, lat, lng):

    x = (lng + 180) * (mapWidth/360)
    y = ((1 - math.log(math.tan(lat * math.pi / 180) + 1 / math.cos(lat * math.pi / 180)) / math.pi) / 2) * mapHeight

    return(x, y)

def pointToLatLng(mapWidth, mapHeight, x, y):

    lng = x / mapWidth * 360 - 180

    n = math.pi - 2 * math.pi * y / mapHeight
    lat = (180 / math.pi * math. atan(0.5 * (math.exp(n) - math.exp(-n))))

    return(lat, lng)

def getImageBounds(mapWidth, mapHeight, xScale, yScale, lat, lng):

    centreX, centreY = latLngToPoint(mapWidth, mapHeight, lat, lng)

    southWestX = centreX - (mapWidth/2)/ xScale
    southWestY = centreY + (mapHeight/2)/ yScale
    SWlat, SWlng = pointToLatLng(mapWidth, mapHeight, southWestX, southWestY)

    northEastX = centreX + (mapWidth/2)/ xScale
    northEastY = centreY - (mapHeight/2)/ yScale
    NElat, NElng = pointToLatLng(mapWidth, mapHeight, northEastX, northEastY)

    return[SWlat, SWlng, NElat, NElng]

lat = 37.806716
lng = -122.477702
zoom = 16
picHeight = 640 #Resulting image height in pixels (x2 if scale parameter is set to 2)
picWidth = 640

mapHeight = 256 #Original map size - specific to Google Maps
mapWidth = 256

xScale = math.pow(2, zoom) / (picWidth/mapWidth)
yScale = math.pow(2, zoom) / (picHeight/mapWidth)

corners = getImageBounds(mapWidth, mapHeight, xScale, yScale, lat, lng)

Here I have used x and y to represent pixel values and lat lng as, Latitude and Longitude. lat, lng, zoom, picHeight and picWidth can all be changed to your specific use case. Changing the scale/ maptype etc will not affect this calculation.

I used this code to tile Static Maps images with no gaps/ overlap. If you want to see more of it in use/ how it can work in that sense there is more information on my GitHub.

Repugnant answered 6/2, 2021 at 14:21 Comment(0)
C
2

Here is a line by line translation of Marcelo's code in PHP, which can probably be cleaned up a bit. Works great! Thanks to Marcelo for doing the hard part.

define("MERCATOR_RANGE", 256);

function degreesToRadians($deg) {
  return $deg * (M_PI / 180);
}

function radiansToDegrees($rad) {
  return $rad / (M_PI / 180);
}

function bound($value, $opt_min, $opt_max) {
  if ($opt_min != null) $value = max($value, $opt_min);
  if ($opt_max != null) $value = min($value, $opt_max);
  return $value;
}

class G_Point {
    public $x,$y;
    function G_Point($x=0, $y=0){
        $this->x = $x;
        $this->y = $y;
    }
}

class G_LatLng {
    public $lat,$lng;
    function G_LatLng($lt, $ln){
        $this->lat = $lt;
        $this->lng = $ln;
    }
}

class MercatorProjection {

    private $pixelOrigin_, $pixelsPerLonDegree_, $pixelsPerLonRadian_;

    function MercatorProjection() {
      $this->pixelOrigin_ = new G_Point( MERCATOR_RANGE / 2, MERCATOR_RANGE / 2);
      $this->pixelsPerLonDegree_ = MERCATOR_RANGE / 360;
      $this->pixelsPerLonRadian_ = MERCATOR_RANGE / (2 * M_PI);
    }

    public function fromLatLngToPoint($latLng, $opt_point=null) {
      $me = $this;

      $point = $opt_point ? $opt_point : new G_Point(0,0);

      $origin = $me->pixelOrigin_;
      $point->x = $origin->x + $latLng->lng * $me->pixelsPerLonDegree_;
      // NOTE(appleton): 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.
      $siny = bound(sin(degreesToRadians($latLng->lat)), -0.9999, 0.9999);
      $point->y = $origin->y + 0.5 * log((1 + $siny) / (1 - $siny)) * -$me->pixelsPerLonRadian_;
      return $point;
    }

    public function fromPointToLatLng($point) {
      $me = $this;

      $origin = $me->pixelOrigin_;
      $lng = ($point->x - $origin->x) / $me->pixelsPerLonDegree_;
      $latRadians = ($point->y - $origin->y) / -$me->pixelsPerLonRadian_;
      $lat = radiansToDegrees(2 * atan(exp($latRadians)) - M_PI / 2);
      return new G_LatLng($lat, $lng);
    }

    //pixelCoordinate = worldCoordinate * pow(2,zoomLevel)
}

function getCorners($center, $zoom, $mapWidth, $mapHeight){
    $scale = pow(2, $zoom);
    $proj = new MercatorProjection();
    $centerPx = $proj->fromLatLngToPoint($center);
    $SWPoint = new G_Point($centerPx->x-($mapWidth/2)/$scale, $centerPx->y+($mapHeight/2)/$scale);
    $SWLatLon = $proj->fromPointToLatLng($SWPoint);
    $NEPoint = new G_Point($centerPx->x+($mapWidth/2)/$scale, $centerPx->y-($mapHeight/2)/$scale);
    $NELatLon = $proj->fromPointToLatLng($NEPoint);
    return array(
        'N' => $NELatLon->lat,
        'E' => $NELatLon->lng,
        'S' => $SWLatLon->lat,
        'W' => $SWLatLon->lng,
    );
}

Usage:

$centerLat = 49.141404;
$centerLon = -121.960988;
$zoom = 10;
$mapWidth = 640;
$mapHeight = 640;
$centerPoint = new G_LatLng($centerLat, $centerLon);
$corners = getCorners($centerPoint, $zoom, $mapWidth, $mapHeight);
$mapURL = "http://maps.googleapis.com/maps/api/staticmap?center={$centerLat},{$centerLon}&zoom={$zoom}&size={$mapWidth}x{$mapHeight}&scale=2&maptype=roadmap&sensor=false";
Cowpox answered 31/10, 2013 at 17:1 Comment(0)
T
2

For big zoom factors (>=8), where non-uniformity of map scale on y axis can be neglected, there is much easier method, where we just takes into accout 1/cos(latitude) correction for pixels/(degrees of latitude) resolution. Initial resolution for zoom=0 is 256 pixels per 360 degrees both for x and y at 0 lattitude.

def get_static_map_bounds(lat, lng, zoom, sx, sy):
    # lat, lng - center
    # sx, sy - map size in pixels

    # 256 pixels - initial map size for zoom factor 0
    sz = 256 * 2 ** zoom

    #resolution in degrees per pixel
    res_lat = cos(lat * pi / 180.) * 360. / sz
    res_lng = 360./sz

    d_lat = res_lat * sy / 2
    d_lng = res_lng * sx / 2

    return ((lat-d_lat, lng-d_lng), (lat+d_lat, lng+d_lng))
Thereabouts answered 30/11, 2017 at 17:32 Comment(0)
B
1

Here is translation to Delphi/Pascal with some optimizations to correspond more strict Pascal language and memory management.

unit Mercator.Google.Maps;

interface

uses System.Math;

type TG_Point = class(TObject)
    private
        Fx: integer;
        Fy: integer;
    public
        property x: integer read Fx write Fx;
        property y: integer read Fy write Fy;

        constructor Create(Ax: integer = 0; Ay: integer = 0);
end;

type TG_LatLng = class(TObject)
    private
        FLat: double;
        FLng: double;
    public
        property Lat: double read FLat write FLat;
        property Lng: double read FLng write FLng;

        constructor Create(ALat: double; ALng: double);
end;

type TMercatorProjection = class(TObject)
    private
        pixelOrigin_: TG_Point;
        pixelsPerLonDegree_, pixelsPerLonRadian_: double;

        function degreesToRadians(deg: double): double;
        function radiansToDegrees(rad: double): double;
        function bound(value: double; opt_min: double; opt_max: double): double;
    public
        constructor Create;
        procedure fromLatLngToPoint(latLng: TG_LatLng; var point: TG_Point);
        procedure fromPointToLatLng(point: TG_point; var latLng: TG_LatLng);
        procedure getCorners(center: TG_LatLng; zoom: integer; mapWidth: integer; mapHeight: integer;
                             var NELatLon: TG_LatLng; var SWLatLon: TG_LatLng);
end;

implementation

const MERCATOR_RANGE = 256;

constructor TG_Point.Create(Ax: Integer = 0; Ay: Integer = 0);
begin
     inherited Create;
     Fx := Ax;
     Fy := Ay
end;

// **************

constructor TG_LatLng.Create(ALat: double; ALng: double);
begin
     inherited Create;
     FLat := ALat;
     FLng := ALng
end;

// **************

constructor TMercatorProjection.Create;
begin
    inherited Create;

    pixelOrigin_ := TG_Point.Create( Round(MERCATOR_RANGE / 2), Round(MERCATOR_RANGE / 2));
    pixelsPerLonDegree_ := MERCATOR_RANGE / 360;
    pixelsPerLonRadian_ := MERCATOR_RANGE / (2 * PI);
end;

// Translate degrees to radians
function TMercatorProjection.degreesToRadians(deg: double): double;
begin
  Result := deg * (PI / 180);
end;

// Translate radians to degrees
function TMercatorProjection.radiansToDegrees(rad: double): double;
begin
  Result := rad / (PI / 180);
end;

// keep value insid defined bounds
function TMercatorProjection.bound(value: double; opt_min: double; opt_max: double): double;
begin
  if Value < opt_min then Result := opt_min
  else if Value > opt_max then Result := opt_max
  else Result := Value;
end;

procedure TMercatorProjection.fromLatLngToPoint(latLng: TG_LatLng; var point: TG_Point);
var
    siny: double;
begin
   if Assigned(point) then
   begin
      point.x := Round(pixelOrigin_.x + latLng.lng * pixelsPerLonDegree_);
      // NOTE(appleton): 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.
      siny := bound(sin(degreesToRadians(latLng.lat)), -0.9999, 0.9999);
      point.y := Round(pixelOrigin_.y + 0.5 * ln((1 + siny) / (1 - siny)) * -pixelsPerLonRadian_);
   end;
end;

procedure TMercatorProjection.fromPointToLatLng(point: TG_point; var latLng: TG_LatLng);
var
    latRadians: double;
begin
    if Assigned(latLng) then
    begin
      latLng.lng := (point.x - pixelOrigin_.x) / pixelsPerLonDegree_;
      latRadians := (point.y - pixelOrigin_.y) / -pixelsPerLonRadian_;
      latLng.lat := radiansToDegrees(2 * arctan(exp(latRadians)) - PI / 2);
    end;
end;

//pixelCoordinate = worldCoordinate * pow(2,zoomLevel)

procedure TMercatorProjection.getCorners(center: TG_LatLng; zoom: integer; mapWidth: integer; mapHeight: integer;
                     var NELatLon: TG_LatLng; var SWLatLon: TG_LatLng);
var
    scale: double;
    centerPx, SWPoint, NEPoint: TG_Point;
begin
    scale := power(2, zoom);

    centerPx := TG_Point.Create(0, 0);
    try
        fromLatLngToPoint(center, centerPx);
        SWPoint := TG_Point.Create(Round(centerPx.x-(mapWidth/2)/scale), Round(centerPx.y+(mapHeight/2)/scale));
        NEPoint := TG_Point.Create(Round(centerPx.x+(mapWidth/2)/scale), Round(centerPx.y-(mapHeight/2)/scale));
        try
            fromPointToLatLng(SWPoint, SWLatLon);
            fromPointToLatLng(NEPoint, NELatLon);
        finally
            SWPoint.Free;
            NEPoint.Free;
        end;
    finally
        centerPx.Free;
    end;
end;

end.

Usage example:

            with TMercatorProjection.Create do
            try
                CLatLon := TG_LatLng.Create(Latitude, Longitude);
                SWLatLon := TG_LatLng.Create(0,0);
                NELatLon := TG_LatLng.Create(0,0);
                try
                    getCorners(CLatLon, Zoom,
                                MapWidth, MapHeight,
                                SWLatLon, NELatLon);
                finally
                    ShowMessage('SWLat='+FloatToStr(SWLatLon.Lat)+' | SWLon='+FloatToStr(SWLatLon.Lng));
                    ShowMessage('NELat='+FloatToStr(NELatLon.Lat)+' | NELon='+FloatToStr(NELatLon.Lng));

                    SWLatLon.Free;
                    NELatLon.Free;
                    CLatLon.Free;
                end;
            finally
                Free;
            end;
Barnie answered 11/5, 2016 at 14:17 Comment(0)
S
0

Here is the NodeJS version, inspired from the "Simple Python Version" above.

import { create, all } from "mathjs";

// create a new instance
const math = create(all, {});


function latLngToPoint(mapWidth, mapHeight, lat, lng) {
  const x = (lng + 180) * (mapWidth / 360);
  const y = ((1 - math.log(math.tan(lat * math.pi / 180) + 1 / math.cos(lat * math.pi / 180)) / math.pi) / 2) * mapHeight;
  return [x, y];
}

function pointToLatLng(mapWidth, mapHeight, x, y) {
  const lng = x / mapWidth * 360 - 180;
  const n = math.pi - 2 * math.pi * y / mapHeight;
  const lat = (180 / math.pi * math.atan(0.5 * (math.exp(n) - math.exp(-n))));
  return [lat, lng];
}

function getImageBounds(mapWidth, mapHeight, xScale, yScale, lat, lng) {
  const [centreX, centreY] = latLngToPoint(mapWidth, mapHeight, lat, lng);

  const southWestX = centreX - (mapWidth / 2) / xScale;
  const southWestY = centreY + (mapHeight / 2) / yScale;
  const [SWlat, SWlng] = pointToLatLng(mapWidth, mapHeight, southWestX, southWestY);

  const northEastX = centreX + (mapWidth / 2) / xScale;
  const northEastY = centreY - (mapHeight / 2) / yScale;
  const [NElat, NElng] = pointToLatLng(mapWidth, mapHeight, northEastX, northEastY);

  return [SWlat, SWlng, NElat, NElng];
}
const lat = 37.806716
const lng = -122.477702

const zoom = 16;
const picHeight = 640; // Resulting image height in pixels (x2 if scale parameter is set to 2)
const picWidth = 640;

const mapHeight = 256; // Original map size - specific to Google Maps
const mapWidth = 256;

const xScale = math.pow(2, zoom) / (picWidth / mapWidth);
const yScale = math.pow(2, zoom) / (picHeight / mapWidth);

const corners = getImageBounds(mapWidth, mapHeight, xScale, yScale, lat, lng);

console.log(corners);


Shirting answered 15/4, 2023 at 8:44 Comment(0)

© 2022 - 2024 — McMap. All rights reserved.