Converting Pixels to LatLong

[BaldPenguin]

I need to find either some CPAN modules or at least a good math page to help me convert pixels to latlong on a map image. Given a map image with the dimensions of 200 x 200 and a known latlong center point, I am trying to compute the next center point if I pan the map 100 pixels in any direction.

Does anyone know of some Perl modules or even some reference pages to help me on my quest?

Update: For purposes of documenting the quest. The Map projection is a Mercator Projection and I have scale lookups.

[Albannach]

This problem is not necessarily a simple one. Depending on the scale of the map, and the map projection used (the map is flat, the earth isn't), your task may be a relatively simple one or virtually impossible.

For starters, check out the Geo::Coordinates modules. You can convert the lat/lon to an appropriate UTM mapping given an ellipsoid. Then, as kyle points out, you will need to know how many metres to the pixel, or perhaps the lat/lon of another point on the map. With just one point and no scale, you cannot do what you want.

[akho]

Things are even worse: it's not enough to know the pixels/meters ratio. You have to know a) which projection was used and b) its exact parameters. It seems obvious that a pixel near the North Pole on a rectangular map represents a larger e-w distance than a pixel on the equator.

I'd read http://en.wikipedia.org/wiki/Map_projection if I were you.

BTW, you can compute the east-west distance represented by 1 pixel on the equator easily, since you know the length of the equator and its length in pixels on your map. That will give you scale if you find which projection was used.

[akho]

If it's a Mercator projection, then

long = x
lat = atan(sinhy)

Where lat, long, x, y of the center point == 0.

Most probably.

This should be adjusted to scale, of course. This can be done for long, as it's just linear. Lat is a bit harder, however. You should find a point with known lat_real on your map, then find it's lat_formula using the above formula. If the map was linearly stretched (which it probably was), you can use

lat = (lat_real / lat_formula) * atan(sinhy)

I strongly advise you to check these formulas for sanity before using them, though (just check several more points, like cities, peninsulas or whatever is visible on your map).

[blue_cowdawg]

Does anyone know of some Perl modules or even some reference pages to help me on my quest?

Google is your friend. I did a Google search and found this page which shows how to take two latitudes and longitudes and calculate distance.

Now the problem becomes one of calculating the distance represented by one pixel and applying some algebra to the formula pointed out above and you have your answer.

Hope this is a help...

[shotgunefx]

Geo::Forward and Geo::Distance?
That's what I used for Google Earth Cache Builder to compute the view points Lat/Longs when stepping through the area to be cached.

[moklevat]

If akho's solution does not get you what you want, you might check out PDL::Transform::Cartography to see if any of its methods would help.

[ikegami]

That depends on the projection used to obtain the map image.

This post is wrong from this point on.

If you draw a three points on the map such that
A_y == B_y && A_x != B_x
A_x == C_x && A_y != C_y

+=====================+
#                     #
#   A-------------B   #
#   |                 #
#   |                 #
#   C                 #
#                     #
+=====================+

Is A_lat == B_lat && A_long == C_long?

• If so, every pixel represents the same long and lat difference.

  Pixel_lat = (B_lat - A_lat) / (B_x - A_x)
  Pixel_long = (C_long - A_long) / (C_y - A_y)

• If not, the math gets a bit more complicated, but it's still pretty straightforward. You'll probably need a fourth point D such that
  D_x == B_x && D_y == C_y

  +=====================+
  #                     #
  #   A-------------B   #
  #   |             |   #
  #   |             |   #
  #   C-------------D   #
  #                     #
  +=====================+
  

  There surely exists modules to help you here.

[akho]

No.

If the Mercator projection was used, then A_lat == B_lat && A_long == C_long, but pixels near the poles and pixels on the equator represent very different long values.

It's hard to find out which projection was used if you only look at a finite number of points.

Finding a map with a known projection for use in the program seems a better way to do it.

[ikegami]

Shoot! I mixed up distances with long/lat somewhere along the way.

[kyle]

In addition to the data you've described, you'd also have to know the scale of the map. What real distance is represented by those 100 pixels?

[Moron]

IMO The most important module for this purpose is Math::Trig.

Mercator projects from the centre of the earth onto an infinite cylinder of equatorial radius running North-South, scaling being only applicable on the projection (degrees = angles can't have scale).

Longitudinal recentering means simply moving the arc difference between two arctangents whereas lateral recentering is linear.

Apart from having the required tan and atan for Mercator projection work, the module also has more specific functions for the lazier mind for converting between spherical and cylindrical co-ordinates.

It isn't possible for the recentering in any direction to traverse either pole but there are potential issues if suitable measures are not taken when crossing zero and the extreme lateral boundaries. This can be made transparent for latitude by converting to and from modulo 180, for example: (updated to include earth radius in mercator projection algorithm)

# untested
use Math::Trig;
my $Radius = 6378569.135 # in metres, so need $scale also in metres
sub ReCentre {
    my ( $lat, $long, $direction, $scale ) = @_;

    # scale should be passed as the change in degrees
    # of latitude afforded by a 100 pixel move
    # - this routine then handles the effects of
    # non-linear longitude automatically

    $lat %= 180; # transform from signed to unsigned ring element

    # handle lateral move
    if ( $direction eq 'W' ) {
        $direction = 'E';
        $scale = -$scale;
    }
    if ( $direction eq 'E' ) {
        return ( Ensign( $lat + 100 * $scale), $long );
            # i.e. move while in unsigned transformation
            # before converting back to signed
    }

    # handle logitudinal move
    $long = tan( $long*pi()/180.0 ) * $Radius; # Project Mercatorially

    $scale = -$scale if ( $direction eq 'S' );
 
    return ( $lat, 180.0 * atan( $long + 100*$scale )/pi() );
        # i.e. move scaled distance while on Mercator projection,
        # and then reverse projection to get degrees.
}

sub Ensign { # convert back from unsigned ring element to  degrees of 
+latitude
    ( $_[0] > 90 ) ? $_[0] - 180 : $_[0];
}

[punkish]

The short answer is, you can't, not without inaccuracies. That is because lat/lon is not a projection and is not cartesian. If your map is already Mercator (and you say it is), then you can use the following code --

my $h_img_dims   => {img_width => 200, img_height => 200};
my $h_map_bbox   => {map_minx => 2090491, map_miny => 482208, map_maxx
+ => 2224055, map_maxy => 587708};
my $h_new_center => {x => 200, y => 100}; # for pan right/East
#my $h_new_center => {x => 0, y => 100}; # for pan left/West
#my $h_new_center => {x => 100, y => 200}; # for pan bottom/South
#my $h_new_center => {x => 100, y => 0}; # for pan up/North


# Convert the mouse click to earth coordinates
my ($x, $y) = _pixel2earth(
  h_img_dims   => $h_img_dims,
  h_map_bbox   => $h_map_bbox,
  h_new_center => $h_new_center,
);
  
  
sub _pixel2earth {
  my (%args) = @_;
  
  my $h_img_dims   = $args{h_img_dims};
  my $h_map_bbox   = $args{h_map_bbox};
  my $h_new_center = $args{h_new_center};

  # Mouse clicks equivalent in earth coordinates
  return (
    $h_map_bbox->{map_minx} + ($h_new_center->{x} * ( ($h_map_bbox->{m
+ap_maxx} - $h_map_bbox->{map_minx}) / $h_img_dims->{img_width}  )),
    $h_map_bbox->{map_miny} + ($h_new_center->{y} * ( ($h_map_bbox->{m
+ap_maxy} - $h_map_bbox->{map_miny}) / $h_img_dims->{img_height} ))
  );
  
}

If you want to display a rectangular image, but measure non-rect coordinates such as lat/lon, you will have to project on the fly. Look at Proj.4, or better yet, look at MapServer which can do the image generation for you as well.