algorithm - Given a line on the earth's surface, how do I plot a line perpendicular to it?

Given a line on the earth's surface, how do I plot a line perpendicular to it?

Apologies if this is a very simple question. I thought this would be a straightforward task, but it is proving counter-intuitive.

I start with the blue line in the figure below (see link -- I'm unable to upload figures yet). I found a perpendicular line by calculating the gradient of the blue line (m), then plotting another line (green) with gradient -1/m. When I plot the lines in Matlab (using 'plot' and 'axis equal'), they look perpendicular, as expected.

http://imgur.com/7qMkx

However, when I export these lines to Google Earth (using the KML Toolbox), they no longer look perpendicular (see link below; the shorter line is the blue line from the previous figure).

http://imgur.com/ncJQ7

I understand that strange things happen on curved surfaces, but I thought the lines should at least look perpendicular locally. I suspect this has something to do with the projection in Google Earth -- in particular, the fact that the grid cells appear to have roughly similar side lengths, yet the longitudinal edge has length = 1 degree, whereas the latitudinal edge has length = 0.5 degree.

So, in summary:

• is my method of finding a perpendicular line valid on a curved surface? (i.e., plotting a line with gradient -1/m)


• in the Google Earth image, do the perpendicular lines look as expected, or is something strange going on?

UPDATE:

To provide more context: I am looking at radar data taken from a plane. The multi-colored area is the 'swath', where observations have been recorded. The blue line I started with in the explanation above is parallel to the swath: this is the aircraft's flight line (the plane was moving in the roughly south-west direction). The radar looks in the direction orthogonal to the flight line, on the left. I am trying to draw a line perpendicular to the flight line; this should be the direction the radar is looking, and should cut the swath neatly. As you can see, this is not the case.

Answer

An elegant principle provides a simple answer:

All points on a smooth curved surface are flat at a sufficiently large scale.

This means that after affine change of coordinates (usually involving just a rescaling of one of them), we can use formulas of Euclidean geometry, such as the Pythagorean Theorem for computing distances and the negative-reciprocal-slope formula for finding perpendiculars.

With latitude and longitude coordinates on the sphere (away from the poles, where longitude becomes singular), all we need to do is rescale the east-west direction to reflect the shrinking length of a degree of longitude as one approaches the poles. With a spherical model of the earth, that shrinkage is given by the cosine of the latitude. This is merely a change in the plot's aspect ratio, nothing more.

This works for regions that extend no more than a few degrees of latitude north-south and do not approach either pole.

Therefore, all you have to do is:

1. 
   

   Multiply all longitudes by the cosine of a typical latitude.

   
   
   


2. 
   

   Compute the perpendicular line.

   
   


3. 
   

   Undo the coordinate adjustment.

   
   

For example, suppose the plane's track took it from (lon, lat) = (-78, 40) to (-79, 41). We may take a typical latitude to lie between 40 and 41, such as 40.5.

Step 1 The adjusted coordinates are (-78 * cos(40.5), 40) = (-59.31167, 40) and (-79 * cos(40.5), 41) = (-60.07207, 41).

Step 2 The question proposes doing this using a negative-reciprocal-slope method. That would be correct, but it will fail in some cases (where the slope is infinite). It is more general and more powerful to use vector arithmetic. Here's how the calculation goes.

The direction vector for the flight's path is the displacement from its beginning to its end,

v =  (-60.07207, 41) - (-59.31167, 40)
  =  (-0.7604, 1.0).

Turning any vector (x,y) at right angles clockwise produces (y,-x), whence a perpendicular direction to the right is

w = (1.0, 0.7604).

According to the Pythagorean Theorem, the length of this vector is the square root of the sum of the squares of its coefficients,

|w| = sqrt(1^2 + 0.7604^2) = 1.256268

Let's move, say, 0.2 degrees along this vector from the starting point of the plane's flight. The start is at (-59.31167, 40) and the displacement is 0.2 / |w| times w, ending up at

(-59.31167, 40) + 0.2 / 1.256268 * (1.0, 0.7604) = (-59.15247  40.12106).

Step 3 To undo the adjustment, divide the first coordinates of any resulting points by the same cosine used in Step 1:

(-59.15247/cos(40.5), 40.12106) = (-77.79064, 40.12106)

If you plot these points using a 1:1 aspect ratio, the angle will appear to be obtuse rather than a right angle. But if you change the aspect ratio to 1:cos(40.5) (about 4:3), the angle will correctly appear to be 90 degrees. When you plot the points using any conformal projection--including Google's Mercator--the angle will also be correct.