Code covered by the BSD License  

Highlights from
GridSphere

4.66667

4.7 | 3 ratings Rate this file 24 Downloads (last 30 days) File Size: 16.9 KB File ID: #28842

GridSphere

by

 

26 Sep 2010 (Updated )

Produces a nearly even grid over the surface of a sphere.

| Watch this File

File Information
Description

The grid may have a total of exactly 12, 42, 162, 642, ... points. Mathematically speaking, the grid may have either 12 or 2 + ( 10 * (4 ^ k) ) points, where k is a positive integer. The user may request any number of points, and the closest attainable value will be produced. All code is compatible with GNU Octave. The algorithm was developed by Nick A. Teanby of Oxford University. Refer to his website for the publication describing the approach and more elaborate geodesic grid software written in IDL: http://www.atm.ox.ac.uk/user/teanby/software.html#icos. Use the FindNearestNeighbors function, available on the MATLAB file exchange at http://www.mathworks.com/matlabcentral/fileexchange/28844-find-nearest-neighbors-on-sphere, to find the grid points closest to arbitrary query points.

GridSphere and FindNearestNeighbors share some functions in common. Each package contains a copy of these functions so that both can stand alone. To eliminate duplicates, simply move all the files into a single folder and replace the shared files when prompted.

Sample usage:

[latGridInDegrees, longGridInDegrees] = GridSphere(12)

latGridInDegrees =

  -58.28253
  -58.28253
  -31.71747
  -31.71747
    0.00000
    0.00000
    0.00000
    0.00000
   31.71747
   31.71747
   58.28253
   58.28253

longGridInDegrees =

     0.00000
   180.00000
   -90.00000
    90.00000
  -148.28253
   -31.71747
    31.71747
   148.28253
   -90.00000
    90.00000
     0.00000
   180.00000

Required Products MATLAB
MATLAB release MATLAB 7.9 (R2009b)
Tags for This File   Please login to tag files.
Please login to add a comment or rating.
Comments and Ratings (14)
02 Jul 2013 Andrew

Fantastic file.

10 Jun 2013 Kurt

Hi, Ben. Thanks for rating. I haven't tested it, but see the first comment on this file (from Paul on February 7th, 2012). I believe this will work if you're using MATLAB but possibly not in GNU Octave.

06 Jun 2013 Ben

Hi Kurt, thanks a lot for the sharing! I have a naive question. How to visualize the result? Could you also kindly provide some examples on using this code?

14 Feb 2013 Kurt

The fix has been approved and is now available for download. Thanks again for catching the bug, Francisco.

14 Feb 2013 Kurt

Finally had time to come up with a more elegant fix and submitted it for review. Will post again once the submission is approved.

25 Jan 2013 Kurt

Sweet; glad to hear it. Yeah, my fix is pretty brutish too.

25 Jan 2013 Francisco

Dont worry !, y already fixed it (in a quite brute way, but it's ok) just wanted to report the bug.

24 Jan 2013 Kurt

Please disregard my previous post.

Thanks! Glad you liked the comments.

Thank you for the clear bug report, Francisco. I don't have MATLAB on hand at the moment, but I have an untested fix ready. If you like I can email it straight to you. Either way, I'll make sure I didn't introduce any additional bugs and then go through the normal MATLAB File Exchange upload process. (It usually takes them a few days to review submissions.)

24 Jan 2013 Kurt

comments*

23 Jan 2013 Francisco

Hi Kurt ! awesome comments in the code, i could understand it entirely.

just one thing, i found something i believe is an error, if i calculate de size of the results (i.e. length(GridSphere(n))

im getting this:
>> length(GridSphere(42)) = 41

>> length(GridSphere(162)) = 161

>> length(GridSphere(600)) = 641

the ammount of points are not correct, so i entered the code and found this at the 'RemoveEqualPairs' function:

uniquePairIndices = DistinctFromPrevious(sortedTwos);
% Assuming the greatest two in each cluster does not equal the least two in
% the next cluster, the pairs for which the two does not equal the previous
% two are the unique pairs.

the asumption there is wrong i believe, for instance for the 42 point grid, when i get the sorted points (with repeated points), the first ones are:


>> [sortedOnes, sortedTwos]

ans =

-90.0000 0
-90.0000 0
-58.2825 0
-58.2825 0
-58.2825 180.0000
-54.0000 -121.7175
... ....

so the 'cleaning' function 'DistinctFromPrevious' should be changed, in particular it must consider the 'sortedOnes' vector.

That's all,
Regards,
Francisco.

19 May 2012 Kurt

Thank you, LU! That's a great question. How concerned are you concerned with performance? If you can stomach an O(n ^ 2) algorithm, where n is the number of points output by GridSphere, then I think the following approach will be easiest.

1. Find the distance between all the pairs of points output by GridSphere.
2. Sort by these distances to find all the pairs of points that are closest to one another. There are at least three gotchas here: (1) you'll need to make sure that either you're not computing the distance between a point and itself or that you're filtering out the pairs with distance 0, (2) you'll probably want to avoid computing both the distance from point a to point b and the distance from point b to point a, and (3) you'll need to account for the imperfect precision of floating point arithmetic, which will result in clumps of distances that are very close but not quite equal. Some of the functions that come with GridSphere may be helpful to you, such as DistinctFromPrevious, which will find the boundaries between the aforementioned clumps.
3. Draw great circles (or just straight lines) connecting the pairs of points you found in step 2. Each of these pairs corresponds to one edge of a triangle. As Paul mentioned, you may need to call conversion functions like deg2rad and sph2cart on your great circle (or line) vectors to project them appropriately.

16 May 2012 LU

Hi Kurt, Wonderful program. I was just wondering if it is possible to connect the vertices with lines? Using Paul's example plots the vertiecs but i would like to display it as a mesh grid showing the edges of the triangles.

Thanks!

08 Feb 2012 Kurt

Thank you for the feedback, Paul. I went ahead and put in a request for the description to be modified to include a usage example. It should appear within five business days.

07 Feb 2012 Paul

Blazingly fast and well documented. Only an example would not have been unbecoming, e.g. to illustrate the fact that returned values are in degrees (whereas Matlab's sph2cart requires radians).

[lat,long] = GridSphere(1000);
latrad = deg2rad(lat);
longrad = deg2rad(long);
[x,y,z] = sph2cart(longrad, latrad, 1);
scatter3(x, y, z, 22);
axis equal vis3d;

Updates
27 Sep 2010

Added a missing function file.

27 Sep 2010

Under the description, I added a link to the FindNearestNeighbors function.

10 Oct 2010

Replaced the ElementWiseMax function with an equivalent built-in function.

08 Feb 2012

Added usage example per user request.

21 May 2012

Replaced tabs with spaces so that the source code displays consistently in all text editors.

14 Feb 2013
14 Feb 2013

Fixes a serious bug in RemoveEqualPairs that resulted in a small number of latitude, longitude pairs being omitted from the final output of GridSphere. Big thanks to Francisco for the bug report.

01 Apr 2013

A MATLAB user reported a crash on line 23 of RemoveEqualPairs.m. This change fixes the type mismatch error by converting a double to a uint32. GNU Octave users will be unaffected by this change as GNU Octave handles the type conversion automatically.

Contact us