Grid Sphere
- 6 files
- 31 downloads
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
Cite As
Kurt von Laven (2019). Grid Sphere (https://www.mathworks.com/matlabcentral/fileexchange/28842-grid-sphere), MATLAB Central File Exchange. Retrieved .
Comments and Ratings (20)
Updates
| 1.15 | Added a link to another FX submission that this submission uses code from. |
|
| 1.14.0.0 | Moved some words in comments onto previous lines where they fit without exceeding 80 columns. |
|
| 1.13.0.0 | Removed write and execute permissions from group/others on all *.m files. |
|
| 1.12.0.0 | Removed trailing whitespace and switched from 2 spaces after each period to 1 in the interest of readability. |
|
| 1.11.0.0 | Added space in title. |
|
| 1.10.0.0 | Made code available as a toolbox and removed some unnecessary hidden files. |
|
| 1.9.0.0 | 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. |
|
| 1.8.0.0 | 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. |
|
| 1.7.0.0 | n/a |
|
| 1.5.0.0 | Replaced tabs with spaces so that the source code displays consistently in all text editors. |
|
| 1.4.0.0 | Added usage example per user request. |
|
| 1.3.0.0 | Replaced the ElementWiseMax function with an equivalent built-in function. |
|
| 1.2.0.0 | Under the description, I added a link to the FindNearestNeighbors function. |
|
| 1.1.0.0 | Added a missing function file. |
MATLAB Release Compatibility
Platform Compatibility
Windows macOS LinuxCategories
Tags
Acknowledgements
Inspired by: Find Nearest Neighbors on Sphere, Golden Ratio, Spherical To Azimuthal Equidistant, Geodesic Midpoints
Discover Live Editor
Create scripts with code, output, and formatted text in a single executable document.
GridSphere/
Select a Web Site
Choose a web site to get translated content where available and see local events and offers. Based on your location, we recommend that you select: .
Select web siteYou can also select a web site from the following list:
Americas
- América Latina (Español)
- Canada (English)
- United States (English)
Europe
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)
Bobbie WU (view profile)
W Ryan Williamson (view profile)
Elapsed time: 0.0035 seconds for 642 points. Seems plenty fast. Suited my needs perfectly, along with Find Nearest Neighbors (http://www.mathworks.com/matlabcentral/fileexchange/28844-find-nearest-neighbors-on-sphere)
Peter Gagarinov (view profile)
spheretri is much faster than GridSphere, see http://www.mathworks.com/matlabcentral/fileexchange/58453-spheretri
sumana (view profile)
Kurt von Laven (view profile)
Hello, Peter! Sorry for the delayed response. Yes, it is intentional. Also, the formula you linked to assume the following coordinate system: http://en.wikipedia.org/wiki/File:3D_Spherical.svg, where r corresponds to rho. The theta computed by arccos(z / ((x ^ 2 + y ^ 2 + z ^ 2) ^ (1 / 2))) is not the latitude, but rather the angle between the positive z-axis and the point in question. Also, acosd (http://www.mathworks.com/help/matlab/ref/acosd.html) returns results in the range [0, 180] rather than the range [-90, 90] returned by http://www.mathworks.com/help/matlab/ref/asind.html that we want for longitude.
Peter Markowsky (view profile)
Hi Kurt,
Is the use of the arcus sinus instead of the tradition arcus cosinus (see
http://en.wikipedia.org/wiki/List_of_common_coordinate_transformations#From_Cartesian_coordinates)
for sphere coordinates inentional? It created some confusion for me :-)
Andrew (view profile)
Fantastic file.
Kurt von Laven (view profile)
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.
Ben (view profile)
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?
Kurt von Laven (view profile)
The fix has been approved and is now available for download. Thanks again for catching the bug, Francisco.
Kurt von Laven (view profile)
Finally had time to come up with a more elegant fix and submitted it for review. Will post again once the submission is approved.
Kurt von Laven (view profile)
Sweet; glad to hear it. Yeah, my fix is pretty brutish too.
Francisco (view profile)
Dont worry !, y already fixed it (in a quite brute way, but it's ok) just wanted to report the bug.
Kurt von Laven (view profile)
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.)
Kurt von Laven (view profile)
comments*
Francisco (view profile)
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.
Kurt von Laven (view profile)
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.
LU (view profile)
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!
Kurt von Laven (view profile)
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.
Paul (view profile)
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;