I've been trying to generate an evenly spaced grid centred at a given point on a sphere, such that the angular separation between any neighbouring pair of points is the same (e.g., 1 degree). The grid should be oriented with the central column of points along the same meridian as the given point. I have little experience with spherical trigonometry (aside from knowing about the Haversine formula and the spherical law of cosines) but after some reading I'm fairly certain that each cell will require four unique great circle arcs. The previous approach I followed using vectors only found great circles perpendicular to the central column (like the coloured grids shown in the first plot) and so the points were not evenly spaced, especially at the edges. Each point in the second plot (showing the result of my latest attempt) is plotted with a 1 degree-radius circle surrounding it, which ideally should coincide with the surrounding four points. The updated method described below finally generates a regular grid. However, for a 20-by-20 degree grid the angular separation of the corner point to the central column is ~4.5 arcseconds smaller than the ideal 10 degree separation. Does anyone know of a method that can determine the required grid points to greater precision?

Current method: The procedure I've followed so far to generate the n-by-n degree grid centred at a given point (e.g. longitude=0, latitude=0; let's call this point A) involves three stages:

Note: a 10-by-10 degree grid would consist of 11-by-11 points each separated from the nearest four points in the adjacent rows and columns by 1 degree. All calculations assume a sphere of radius 1.

  1. Find the points along the `central' column:
    Simply increment the latitude (i.e. following the meridian running through the centre point) and wrap the latitude and longitude if a pole is crossed (i.e. if the current point's latitude is less than -90 or greater than +90).

  2. Find the points along the `central' row (through the grid centre): Since the row of points through A follow a great circle that intersects the meridian at a right angle, we need to find the plane that this second great circle lies on. Consider two points at the given central longitude; let point A be the grid centre at 0 degrees longitude and 0 degrees latitude, and point B at 0 degrees longitude and 1 degree north. The two points can be expressed in Cartesian coordinates, and so the unit vector (N) normal to the great circle that contains A and B is found using the cross product, i.e. N = BxA. The cross product of the unit vector N with a unit vector in the direction of A gives another normal vector (M). The vector M can then be used as an axis to rotate N by (90-x) degrees^, where x is the offset (in this case 1 deg) from the central column. This is repeated outwards in steps of 1 degree for one half of the central row then the other half using a unit vector in the opposite direction to N. ^http://en.wikipedia.org/wiki/Rotation_matrix#Rotation_matrix_from_axis_and_angle

  3. Fill in the remaining points in each quadrant: Each point in the grid forms a rectangle lying on the same plane as three other points in the grid - let A be the point at the grid centre, B is a point 'above' or 'below' A in the same column and C is a point on the central row, and the unknown point D is in the same row as B and same column as C. For the upper left quadrant, the rectangle (sometimes a square) would have corners clockwise from the upper-left D, B, A, and C.

First compute vectors AB and AC then a unit vector normal to the triangle BAC can be found from their cross-product, i.e. t=ABxAC. To find D, the approach I've used is to rotate the plane on which A/B/C/D lie to align it with the xy-plane. The rotation matrix requires a rotation axis and rotation angle that I found here: https://gamedev.stackexchange.com/a/48102

The dot product of the rotation matrix and each of vectors A, B and C allows vector D' (the equivalent of vector D in the xy-plane) to be found from their components, i.e., $D' = (b_x-a_x+c_x, b_y-a_y+c_y, b_z-a_z+c_z)$

The dot product of the rotation matrix using the same rotation axis but rotation angle of opposite sign to that given above with vector D' gives the coordinates of vector D, which can be converted to the latitude and longitude for point D. This is repeated for every remaining point in the four quadrants.

I've spent a few days looking for a better method but without any luck.

Fig 1. Sphere with grids: Sphere with grids

Fig 2. Grid points - note circles shown are 1 degree in radius (i.e. the perimeter of a 1 degree spherical cap centred at the given grid point). I had used a gnomonic projection to present the grid in figure 2 (using the python plotting package APLpy) but replaced it with an orthographic projection since the former shows great circles as straight lines but in this case the rows of points are joined by different great circle arcs, except the (n/2+1)th row. Grid points with 1 degree radius circles


You can't do what you are hoping because you can't cover a sphere with squares. In fact, you can't even put one square on the sphere. If you look closely, the angles of your squares add up to more than $360^\circ$, by an amount that is determined by the solid angle taken up by the square.

A few choices that maintain some of what you want. 1) The black grid in your first figure comes from the classic lines of latitude and longitude. You get quadrilaterals except near the poles. You can choose what latitude to make the cells almost square, but they will be much larger near the equator than near the poles. 2) Use the latitude/longitude lines as in 1, but have fewer cells as you get near the poles. You can keep the cells close to square, but won't have a nice grid. 3) Use the square grid as you are, but recognize that you can't have one grid that covers the sphere. Maybe a grid that covers about half the sphere centered on a point of current interest meets your needs.

  • $\begingroup$ Thank you for your reply. When you say 'squares' regarding the cells of the grid, do you mean spherical rectangles with equal (angular, e.g. 1-degree) length sides? (That's the type of cell I was hoping to find - like four rays from the centre from the centre of the sphere each separated by 1 degree from the nearest two.) Also, the grid would only ever extend over a 10 by 10 or 20 by 20 degree area. I've update my current method with additional details. $\endgroup$ – Marcin Feb 24 '15 at 17:00
  • 1
    $\begingroup$ There are no spherical rectangles, ’cause the four angles of a quadrilateral have to add up to more than $360^\circ$. And don’t forget: the “parallels” of latitude are not “straight”, i.e. not geodesics, the way the lines of longitude are. They’re what are known as “small circles” on the surface of the sphere. $\endgroup$ – Lubin Feb 24 '15 at 17:05
  • $\begingroup$ I came across this diagram <math.stackexchange.com/a/1016129/217456> that describes a shape that is like the cells of the grid I had hoped construct. If alpha were to equal beta (degrees) then the left and right sides would be arcs along great circles (for this particular cell, also meridians) as would the top and bottom lines (which would be great circles pivoted upwards and downwards alpha degrees from the equator, respectively). As you pointed out these would not follow lines of latitude. If it's possible to tesselate this shape then one of the corners could be the centre of the grid $\endgroup$ – Marcin Feb 25 '15 at 7:22
  • $\begingroup$ You can't tesselate that shape exactly, because the angles are greater than $90^\circ$ and four won't fit together. If you make them small compared to the radius of the sphere, the errors are very small because the sphere looks like a plane at that scale. $\endgroup$ – Ross Millikan Feb 25 '15 at 14:12

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.