There is in fact no real reason to partition the sphere into a regular non-overlapping mesh, try this:
-
partition your sphere into semi-overlapping circles
see here for generating uniformly distributed points (your circle centers)
Dispersing n points uniformly on a sphere
-
you can identify the points in each circle very fast by a simple dot product..it really doesn’t matter if some points are double counted, the circle with the most points still represents the highest density
mathematica implementation
this takes 12 seconds to analyze 5000 points. (and took about 10 minutes to write )
testcircles = { RandomReal[ {0, 1}, {3}] // Normalize};
Do[While[ (test = RandomReal[ {-1, 1}, {3}] // Normalize ;
Select[testcircles , #.test > .9 & , 1] ) == {} ];
AppendTo[testcircles, test];, {2000}];
vmax = testcircles[[First@
Ordering[-Table[
Count[ (testcircles[[i]].#) & /@ points , x_ /; x > .98 ] ,
{i, Length[testcircles]}], 1]]];