2

I'm trying to create a bounding box at a specific scale from a center coordinate. I'm trying to keep it within the aspect ratio of a 8.5x11inch piece of paper (612x792 pixels @ 72dpi).

The code I'm using below mostly works, but the heigh seems a bit too tall for the aspect ratio of a letter. Am I not accounting for mercator projection? What am I missing here?

def bounding_box_from_point(center:, size:, scale_denominator:)
  dpi = 72
  inches_per_unit = 4374754
  resolution = 1 / (scale_denominator * inches_per_unit * dpi)
  half_width_deg = (size.width * resolution) / 2
  half_height_deg = (size.height * resolution) / 2

  BoundingBox.new(
    north: center.lat + half_height_deg,
    south: center.lat - half_height_deg,
    east: center.lon + half_width_deg,
    west: center.lon - half_width_deg
  )
end

Calling bounding_box_from_point(center: center, size: size, scale_denominator: scale) with:

scale = 0.0008861342166177423 (i.e. 1/18055.955520)
center = Geometry::Location.new(lat: 37.806336, lon: -122.270625)
size.width = 612,
size.height = 792

It returns:

west: -122.27172131608657,
east: -122.26952868391342,
south: 37.804917238005615
north: 37.80775476199439

If you go to http://www.openstreetmap.org/export and enter those bounding box coordinates, you can see that the ratio does not match that of a 8.5x11in piece of paper...it's slightly too tall. What am I doing wrong here or not understanding?

hummmingbear
  • 2,294
  • 5
  • 25
  • 42
  • Its not so much the Mercator projection as it is the fact that longitude lines get closer together as you move away from the equator. Your map isn't too tall, its too skinny, – Ron Jensen Feb 13 '17 at 17:50
  • See http://stackoverflow.com/a/2913249/4504607 – Ron Jensen Feb 13 '17 at 18:17
  • @RonJensen so how do you adjust for this? – hummmingbear Feb 13 '17 at 18:17
  • @ronjensen I read through that post, but doesn't quite answer it. The article that's linked doesn't directly/clearly answer it – hummmingbear Feb 14 '17 at 17:47
  • What language are you using in your code example? I believe you need to do something like: double deltaLon = Math.asin(Math.sin(radDist) / Math.cos(radLat)); but use your desired width for radDist – Ron Jensen Feb 14 '17 at 18:46
  • @ronjensen what units does radDist take? – hummmingbear Feb 14 '17 at 19:35
  • The fragment I posted was from the link. He says there // angular distance in radians on a great circle double radDist = distance / radius; For your case you could have a vertical and horizontal radDist – Ron Jensen Feb 14 '17 at 19:56

1 Answers1

0

This was my implementation that solved it!

def self.bounding_box_from_location(location:, scale:)
  scale_meters_per_pixel = 0.0003

  # scale.zoom is an integer between 0-19
  # See http://wiki.openstreetmap.org/wiki/Zoom_levels for a description on `zoom`
  # location.lat_rad is the latitude in radians, same for lon_rad

  meters_per_pixel = EARTH_CIR * Math.cos(location.lat_rad) / (2 ** (scale.zoom + 8))
  earth_meters_per_pictureMeters = meters_per_pixel / scale_meters_per_pixel


  # height/width in meters is the height/width of the box you're applying this to
  # In my case I'm using the heigh/width of a Letter of paper
  # You can convert Pixels to meters (@ 72dpi) with
  # pixels * 0.00035277777777777776

  meters_north = earth_meters_per_pictureMeters * scale.height_in_meters / 2
  meters_east = earth_meters_per_pictureMeters * scale.width_in_meters / 2

  meters_per_degree_lat = 111111.0
  meters_per_degree_long = 111111.0 * Math.cos(location.lat_rad)

  degrees_north = meters_north / meters_per_degree_lat
  degrees_east = meters_east / meters_per_degree_long

  BoundingBox.new(
    north: (location.lat + degrees_north),
    south: location.lat - degrees_north,
    east: location.lon + degrees_east,
    west: location.lon - degrees_east,
    width: scale.width_in_pixels,
    height: scale.height_in_pixels
  )
end
hummmingbear
  • 2,294
  • 5
  • 25
  • 42