views:

257

answers:

2

I need to find the corners in lat/lng of a ground overlay given in a kml-file either in php or javascript.

I.e. for a specific example I need to get from:

  <LatLonBox>
    <north>60.406505416667</north>
    <south>60.400570555556</south>
    <east>5.3351572222222</east>
    <west>5.3190577777778</west>
    <rotation>3.7088732260919</rotation>
  </LatLonBox>

to corner coordinates

SW: 60.400316388889;5.3194425
SE: 60.400824722222;5.3355405555556
NE: 60.406759444444;5.3347738888889
NW: 60.406251388889;5.3186730555556

I can get the other way (approximately at least, php code given) by

$w=($nw_lng+$sw_lng)/2;
$e=($ne_lng+$se_lng)/2;
$n=($ne_lat+$nw_lat)/2;
$s=($se_lat+$sw_lat)/2;
$rot= rad2deg (atan ( ( $nw_lng - $sw_lng ) / ($sw_lat - $nw_lat ) / 2  ) );

Should be easy to get back, but I've used hours for this without getting there. Any tips?

A: 

You need to use spherical trigonometry, part of spherical geometry for full accuracy. However, since you are dealing with only a small piece of the sphere, euclidian geometry will do if you remember one thing.

As latitude increases, the lines of longitude get closer together. For example, near the North Pole, the latitude lines are almost touching. So condition your latitude differences, diminishing them by mulitlying by a factor of cos(latitude). That will give you good enough accuracy for your app.

 $n = 60.406505416667;
 $s = 60.400570555556;
 $e = 5.3351572222222;
 $w = 5.3190577777778;
 $rotn = 3.7088732260919;

 $a = ($e + $w) / 2.0;
 $b = ($n + $s) / 2.0;
 $squish = cos(deg2rad($b));
 $x = $squish * ($e - $w) / 2.0;
 $y = ($n - $s) / 2.0;

 $ne = array(
   $a + ($x * cos(deg2rad($rotn)) - $y * sin(deg2rad($rotn))) /$squish,
   $b + $x * sin(deg2rad($rotn)) + $y *cos(deg2rad($rotn))
   );
 $nw = array(
   $a - ($x * cos(deg2rad($rotn)) + $y * sin(deg2rad($rotn))) /$squish,
   $b - $x * sin(deg2rad($rotn)) + $y *cos(deg2rad($rotn))
   );
 $sw = array(
   $a - ($x * cos(deg2rad($rotn)) - $y * sin(deg2rad($rotn))) /$squish,
   $b - $x * sin(deg2rad($rotn)) - $y *cos(deg2rad($rotn))
   );
 $se = array(
   $a + ($x * cos(deg2rad($rotn)) + $y * sin(deg2rad($rotn))) /$squish,
   $b + $x * sin(deg2rad($rotn)) - $y *cos(deg2rad($rotn))
   );
 print_r(array(
 'sw'=>$sw,
 'se'=>$se,
 'ne'=>$ne,
 'nw'=>$nw,
 ));

My $squish variable is the cos(lat) I mentioned. There is de-squishing for the relative part of horizontal lengths. The sine table looks like this:

NE: (a + x cos A - y sin A, b + x sin A + y cos A)
NW: (a - x cos A - y sin A, b - x sin A + y cos A)
SW: (a - x cos A + y sin A, b - x sin A - y cos A)
SE: (a + x cos A + y sin A, b + x sin A - y cos A)

Perhaps tttppp could account for the differences from tttppp's table.

Ewan Todd
the distances seem small to me, so maybe an euclidean approximation would be enough
Victor
I already know that I need to use some kind of trigonometry. I think an euclidean approximation would do (that's what the inverse formula uses which works fine for the maps I've considered does). However, I have not managed to get the trigonometry right for the forward case, even after hours of work. The inverse formula was easy to figure...
jankoc
This looks very good now! Probably my inverse formula would also need some fixing to get back the exact answer. This gave Array ( [sw] => Array ( [0] => 5.3194632953201 [1] => 60.400319597408 ) [se] => Array ( [0] => 5.3355290212874 [1] => 60.400833943604 ) [ne] => Array ( [0] => 5.3347517046799 [1] => 60.406756374815 ) [nw] => Array ( [0] => 5.3186859787126 [1] => 60.406242028619 ) )
jankoc
A: 

Find the width and height of the original rectangle and then define x = half the width, y = half the height. Let A be the angle rotated (anticlockwise) in radians. The origin (a,b) is the position of the centre of the rectangle.

The four corners end up at:

NE: (a + x cos A - y sin A, b + x sin A + y cos A)
NW: (a - x sin A - y cos A, b + x cos A - y sin A)
SW: (a - x cos A + y sin A, b - x sin A - y cos A)
SE: (a + x sin A + y cos A, b - x cos A + y sin A)
tttppp
This is along the lines I've been trying, but does not seem to work. I did x=(east-west)/2; y=(north-south)/2; a=(east+west)/2; b=(north+south)/2; A=-rot*pi/180; does give me quite far off results, e.g. ne_lng = 60.4059784, ne_lat=5.3353323. The definition of north, south, east, west and how the rotation is applied for ground overlays in the kml-format is what makes it tricky.
jankoc