Converting points to a geodesic grid (subdivided icosahedron)

Recommended Posts

Numsgil    501
I have what's known as a geodesic grid (a subdivided icosahedron) used to represent the sphere of a planet(link 1 and link 2). Basically, you take an icosahedron and subdivide it an arbitrary number of times to build more and more verticies. Each vertex represents a single cell (the voronoi region around each vertex is the cell). What I want to do is "project" an arbitrary 3 Dimensional point (or more specifically, the ray between the origin of the sphere toward the point) onto the sphere. This involves finding which triangle (3 verticies of the geodesic grid) the point is "above", and projecting the point onto that triangle. Basically I'm working on how to perform "terrain following" on the surface of a sphere. The problem's a little more complex than the usual heightmaped terrain. At present, my basic algorithm consists of: 1. Convert the point to spherical coordinates 2. Convert the base 12 verticies of the subdivided icosahedron into spherical coordinates 3. Extract just the angles of the coordinates, leaving behind the distance term. 4. Using the two angles, create a transformation into R^2 space for the point and the base icosahedron. 5. Taking each triangle in the base icosahedron, find the barycentric coordinates of the point. If the barycentric coordinates are in the triangle, all three will be 0 <= u,v,w <= 1. 6. Once the base triangle that the point is above is found, I recursively subdivide the triangle an arbitrary number of times (however much my geodesic grid was subdivided), repeating the process. 7. I return the smallest triangle the point is in, along with the barycentric coordinates for the point. 8. I convert the smallest triangle's verticies back into the original R^3, and use the barycentric coordinates to find the R^3 position of the point on that triangle. I do take into account the fact that the spherical coordinates' theta term (that is, the term in the range [0, 2 pi]) has a seam at 0 and 2pi. Or I try to anyway. My implemenation of the above doesn't work. More specifically, it finds a triangle in the base icosahedron, but when it tries to subdivide it, none of the subdivided triangles appear to contain the point. I think my algorithm is sound, though if anyone notices any issues or incorrect assumptions I'm making, I'd appreciate it. It's probably just my implementation that's buggy, though I'd like a second opinion. I'm curious if anyone sees a simpler way to do what I'm trying to do, one with less transformations to different coordinate spaces. I've already implemented a brute force method that just tests every possible triangle in the subdivision to see if the ray between the sphere's origin along the direction of the point passes through it. But it quickly becomes computationally unwieldly as the number of subdivisions exceeds about 3 or 4.

Share on other sites
oliii    2196
Assuming you found the base triangle, and the corresponding (u, v) pair of that point in respect to a triangle of the icosahedron.

Then surely there must be a solution to finding the subdivision triangle that contains that [u, v] pair. I don't see how that would fail. Say your original triangle's vertices mappings have the (u, v) coordinates (0, 0), (1, 0), (0, 1). Then the 3 midpoints would be (0.5, 0), (0, 0.5), (0.5, 0.5). And so on, and so forth. So all matters of (u, v) pairs in the original range [0, 1] should be covered.

Share on other sites
oliii    2196
but if you just want a nice fly-by of the planet, I would work backwards for that. I would work from the 2D map directly instead. Pickup the point on the 2D map, moving it along in 2D, and then convert it to 3D.

I wouldn't mess with spherical coordinates and geodesic spheres myself. I tried that, it's not pretty.

Share on other sites
Numsgil    501
My interest is to create a 3D planet that you can explore on foot. Something similar to the Little Prince's asteroid in the French story The Little Prince.

Quote:
 I wouldn't mess with spherical coordinates and geodesic spheres myself. I tried that, it's not pretty.

No, they're not ;) But I've got most of the difficulties figured out (building the mesh, texturing it, subdividing it, linking it to a data structure, etc.), after several months of work. I just have terrain following and LOD issues left to solve.

You're indication then is that my algorithm should work? No major flaws?

Share on other sites
oliii    2196
I don't see how it would fail to pickup a sub-triangle, if that point is in the base triangle. That's my take on it. There must be something wrong in your algo. It's possible you get some weird artifacts at the sub-triangle edges and vertices.

In fact, it might even be possible to get the triangle directly from the base (u, v) pair.

Say, you have a subdivision depth of 2 (base, LOD1, LOD2, LOD3). Say the length of the base triangle edge is 1.

the vertices of your triangles will have pairs (u, v) as a multiple of 1/8 (1/2^3) at LOD3.

say, you have u = 0.3425, and v = 0.76556.

u is in range of vertices of LOD3 with u = [0.25, 0.375].
v is in range of vertices of LOD3 with v = [0.75, 0.875].

The choice is limited to two triangles. you know your point is inside the quad defined by vertices of LOD3 [(0.25, 0.75), (0.25, 0.875), (0.375, 0.75), (0.375, 0.875)]. Which one it is should be easy to find out.

Share on other sites
oliii    2196
BTW, I will stick with my idea of using the 2D map as your point of view, and converting it to 3D, since you have the texturing sussed out. I think it would be far easier to move a point on a 2D map (even with a non-continuous map), convert to 3D, than try to move the point in 3D and look up the triangle. I'm not sure you will get a nice smooth movement out of it. But I could be wrong :)

Share on other sites
Bob Janova    769
How I would approach this (and like you I use a subdivided icosahedron as a planetary mesh, though I use the triangles themselves as cells not the vertices):

• Perform a ray-triangle test on each of the 20 major faces in standard 3D cartesian space. This will give you a point in the plane of one of the faces.
• Each face has a uniform grid of cells within it, so find the intersection point in terms of the face's coordinates and you're done. This depends how you've done your coordinates but there will be an X and a Y direction on each face and a (0,0) point; project (intersection point - face origin) onto X and Y vectors.

Share on other sites
Numsgil    501
I think I've gotten it to work now. As you suggested Bob, I performed a ray test in 3D cartesian space against the 20 major faces. Then I continue recursively raytesting against the 4 inner triangles from the subdivided triangle, until I arrive at the proper depth.

Much easier than my original method. The trick is to realize that since the ray points towards the center of gravity (center of the sphere), you're implicitly handling the distance of the verticies from the center of the sphere.

Share on other sites
daftasbrush    252
I think there's an optimization for the initial face testing.

Since all the triangles are equilateral
And your ray "origin" is always the centre of the sphere.

If you pre-computed and stored the normalized direction vector to the "centre" of each face from the centre of the sphere. (? the outward facing normal of the face)
Then a simple dot-product between that direction and your ray direction, gives you a "fast rule-out" for each face....the test would be simply

if a . b < 0.8 (or thereabouts)

The limit would a fixed value, and would be the dot product of any centre point direction (?face normal) with the normalized direction vector from the centre to any vertex of that face.
Basically, the cosine of the angle, at the cente of the sphere, between the line to the "centre" of a face, and the line to any vertex.

I don't have the time (right now) to work out the exact value... but that would save you actually "ray testing" most of the faces... mostly only 1 face would pass that test.
And you'd only need to do the "full" ray-triangle intersection if it did.

You could even extend it...
There's another (higher) limit value on the dot product, where you can "instantly" say, the line must pass through this face... so you don't need to test whether the ray intersects the triangle, it DOES/MUST, you just need to work out "where" (u,v)
That value is the cosine of the angle between the face normal, and the line from the centre to the mid-point of any side of that face.
Note : These limit values will be the same for all the major faces, so just a constant really.

Basically you are testing whether your ray direction lies within a cone with it's "point" at the sphere origin, and it's base on the surface.
The "instant yes" test, is testing that the ray lies within a cone, whose base is the inscribed circle of the face triangle.
The "definitely not" test, is testing that the ray lies outside a cone, whose base is the exscribed circle of the face triangle.
And "in between" those 2, is the only place you actually need to test "whether" the ray and triangle intersect.

Also, I'm with Bob and Oliii...
If you calculate (u,v) for the intersected "major" face, you can use those values to determine which sub face you hit in "no time"...

[Edited by - daftasbrush on March 7, 2007 6:19:27 PM]

Share on other sites
daftasbrush    252
From the r,p,R table about halfway down the page.

For an Icosahedron with unit side length
r = radius of the in-sphere = 0.75576
p = radius of the mid-sphere = 0.80902
R = radius of the circum-sphere = 0.95106

The dot product "limits" can be easily determined from these values.

If a is your normalized ray direction (from centre to the point outside sphere)
and b is the "outward facing" face normal

Rule IN (Must hit this face)
a.b > r / p = 0.934167

Rule OUT (Cannot hit this face)
a.b < r / R = 0.79465

so only if the dot product is within this "small" range (>~0.8 and <~0.93) do you actually need to test whether the ray intersects the triangle.

You might want to use the original "analytical" derivations (also on the Mathworld page), rather than the numeric values above.

Share on other sites
Numsgil    501
Thanks, I'll look it over ;)