0

I have a bounding box of:

Left -122.27671
Bottom 37.80445
Right -122.26673
Top 37.81449

It could also be converted into NE Lat/Long and SW Lat/Long

Within that bounding box, I'd like to find the X,Y position of a specific Lat/Long. This would be using the Mercator projection.

I've seen answers that find the X,Y of a position on a world map using Mercator, but not within a specific lat/lon.

Any help appreciated!

UPDATE Put this together from another question I saw. Can anyone validate if this seems legit?

map_width = 1240
map_height = 1279

map_lon_left = -122.296916
map_lon_right = -122.243380
map_lon_delta = map_lon_right - map_lon_left

map_lat_bottom = 37.782368
map_lat_bottom_degree = map_lat_bottom * Math::PI / 180

def convert_geo_to_pixel(lat, long)
  x = (long - map_lon_left) * (map_width / map_lon_delta)

  lat = lat * Math::PI / 180
  world_map_width = ((map_width / map_lon_delta) * 360) / (2 * Math::PI)
  map_offset_y = (world_map_width / 2 * Math.log((1 + Math.sin(map_lat_bottom_degree)) / (1 - Math.sin(map_lat_bottom_degree))))
  y = map_height - ((world_map_width / 2 * Math.log((1 + Math.sin(lat)) / (1 - Math.sin(lat)))) - map_offset_y)

  return [x, y]
end
hummmingbear
  • 2,294
  • 5
  • 25
  • 42

2 Answers2

1

Found a better solution that I've test and validated. Posting this for anyone else who might find it useful. It's written in Ruby but easy to convert to any other language

@north = to_radians(37.81449)
@south = to_radians(37.80445)
@east = to_radians(-122.26673)
@west = to_radians(-122.27671)
# Coordinates above are a subsection of Oakland, CA

@map_width = map_width
@map_height = map_height

def location_to_pixel(lat:, lon:)
  lat = to_radians(lat)
  lon = to_radians(lon)
  ymin = mercator_y(@south)
  ymax = mercator_y(@north)
  x_factor = @map_width/(@east - @west)
  y_factor = @map_height/(ymax - ymin)

  y = mercator_y(lat);
  x = (lon - @west) * x_factor
  y = (ymax - y) * y_factor
  [x, y]
end

def to_radians(deg)
  deg * Math::PI/180
end

def mercator_y(lat)
    Math.log(
      Math.tan(lat/2 + Math::PI/4)
    )
end
hummmingbear
  • 2,294
  • 5
  • 25
  • 42
0

Let's s is shift of map in world space, bottom latitude in radians B, top latitude T. (I assume y=0 is bottom)

enter image description here

C * Sin(B) = 0 + s
C * Sin(T) = map_height + s
=>
C = map_height / (Sin(T) - Sin(B))
s = C * Sin(B)
y = C * Sin(Lat) - s = 
    C * Sin(Lat) - C * Sin(B) = 
    C * (Sin(Lat) - Sin(B)) = 
    map_height * (Sin(Lat) - Sin(B) / (Sin(T) - Sin(B))

        // note - resembles linear interpolation is sine space
MBo
  • 77,366
  • 5
  • 53
  • 86
  • Is the code I used above incorrect? Having trouble understanding where your example outputs `x, y` when given a lat/lon within a bounding box? – hummmingbear Jan 10 '17 at 20:57
  • Your x calculation is correct. y seems wrong. My Lat corresponds to your lat, B=map_lat_bottom, T=map_lat_top (don't see it in your code) – MBo Jan 11 '17 at 05:23
  • I copied my solution from http://stackoverflow.com/questions/2103924/mercator-longitude-and-latitude-calculations-to-x-and-y-on-a-cropped-map-of-the/10401734#10401734 I believe map_lat left and right are the top left and top right corners. Since we know where the corner is, all we need is the lat of the bottom to finish the square. That's how I interpreted this solution? – hummmingbear Jan 11 '17 at 18:18
  • Posted a new answer with a refactored version above. Both pieces of code I posted work, but the answer I selected is more clear as it specifies NSEW – hummmingbear Jan 12 '17 at 18:50