I set myself a task. I have a Garmin GPSMap 62st: a great little device (if a little limited), but the utility I can derive from it is tightly coupled to the availability of maps. What's the point in knowing you're at

`13°09′47″S 72°32′44″W`

if the map around is just empty black space? Armed with a 1024 x 1024 pixel map (the maximum tile size for the 62st) in Mercator projection, I isolated two points, in diagonally opposite corners, to get the most distance between them. Google Earth allowed me to drop markers and get the corresponding latitude and longitude co-ordinates for each pair of (x,y) pixels.Degrees are easily converted into radians (and back again), and radians of latitude can easily be converted into the unitless measures of height used in the Mercator projection, while degrees of longitude don't need any conversion at all. If you look at a complete map of the earth, you'll see the areas surrounding the poles are most distorted. Anyway, long story short:

Given a pair of (N,E) coordinates and their respective (X,Y) positions on the raster map, it's pretty simple to work out what the left, top, right and bottom boundaries of the entire map are, and hence easy to create the KML and KMZ file required to correctly position a raster map (of the Mercator projection) on the 62st.

Longitude first because it's easiest:

`X2-X1`

gives the pixel width between the two selected points on the raster map; `X1-X0`

gives the pixel padding to the left edge, and X3 - X2 gives the pixel padding to the right. We know the angle of longitude at E1 and E2, so we can work out a ratio of degrees to pixels using `(X2-X1)/(E2-E1)`

. Substituting the ratio into `X3-X2`

and `X1-X0`

, and adding the correct starting point, we get the longitudinal boundaries.Latitude has some tricky hyperbolic and trigonmetric functions (as seen on Wolfram): an angle of latitude N (in radians) is represented at height H (no units) above the equator. N1's height would be

`H1 = atan(sinh(N1))`

and N2's height `H2 = atan(sinh(N2))`

. The formula `(Y2-Y1)/(N2-N1)`

gives us the ratio of pixels to "heights". It's simple to get to Y0 and Y3 from Y2 (and/or Y1) using just addition and subtraction. Also, we know Y1 and Y2's latitudes so it's simple to get the "height" of Y0 and Y3: H0 and H3. Getting back to radians from "height" is as simple as `Y0 = log(tan(H0) + 1/cos(H0))`

. All that's left is to convert the answer from radians back into degrees. Simples.
## No comments:

Post a Comment