Importing Faces Into CGAL Arrangements

This problems occurs repeatedly in the X-Plane scenery tools code: we need to import a series of polygons into a single arrangement_2 structure, and we want to tag the faces contained by these polygons as they appear in the arrangement. This isn't entirely trivial for two reasons:

• The polygons may collide with each other and thus there may be a 1:many relationship between the original data and the final map.
• The polygons may be self-intersected or in other ways "hosed", so we need a particular strategy for handling this situation.

CGAL provides a number of built-in tools to deal with these situations. Here are 3 basic and useful building blocks:

1. If you are using the general polygon set code that is built on top of arrangements, you can simply perform unions and intersections of a large number of faces. For merging a large number of areas, this code is faster than anything else you might code, because it can do an N-way divide and conquer, where N is larger than 2.
2. For custom merging and handling of multiple polygons, you can use the overlay free function, which lets you specify how the combinations of each set of face from two maps are handled.
3. You can simply insert a set of curves into an empty arrangement and they will be "swept" together. This is a useful way to turn a messy polygon into something useful - it finds intersections, builds topology, runs quickly, and handles input no matter how degenerate.

For example, if the goal is to have any area contained by any piece of polygon as "inside" you can simply insert all polygon sides into an empty arrangement and then tag every bounded face.

Finding Polygon Internal Areas

When building arrangements out of polygons of dubious origin (or simply building an arrangement out of a large number of unrelated polygons) I use a bulk insert to "sweep" the curves into the arrangement. How do I then find the faces? Here are three techniques:

1. The contained area of a polygon can be found by simply checking whether the face is bounded or not. This is not useful though when importing multiple polygons at the same time. (When I need this technique, each polygon is individually imported into its own arrangement, then all arrangements are merged later, typically with general-polygon-set code.)

2. We can implement a "toggle" policy (e.g. each line toggles interior vs. exterior) by doing a search from the outside to the inside of the arrangement, toggling whether we are "inside" or "outside" each time we cross a halfedge. The halfedges can retain curve-based properties; typically I use a consolidated data curve so that halfedges retain every property attached to them.

   One danger: an antenna will produce incorrect results in this technique because it won't toggle the data property twice. This can be hard to work around because data from the insert is maintained per edge, not half-edge

   Unfortunately, topology of the final arrangement doesn't help us resolve this either. Imagine an antenna (in the original polygon that crashes into another polygon, thus becoming part of a real partition. Technically the original face is not split by the antenna, but the faces in the final arrangement are, so saying face()==twin()->face() doesn't tell us we have an antenna in the original.

   Two ways to work around this: don't insert antennas, or don't tag known antennas with any data. Both cases require knowing that we have an antenna ahead of time.

3. Sometimes we want to use a "winding rule" - that is, contain areas inside closed left turning contours. This is, for example, useful when calculating offset buffers and minkowski sums; the artifacts from strange shapes being offset too much turn out not to be left turning contours and get thrown out.

   To find the winding rule areas, we have to look at the direction of the curves inserted.

   The way my code does this is to look at the direction of the underlying curve vs. the direction of the half-edge and then mark the half-edge as being on the "inside" or "outside" of the winding with a dat a field on the half-edge itself. This is reconstructed after bulk insert, and then we can traverse the whole arrangement, counting windings.

   The limitation here is similar to above: if we have an antenna, the underlying curve can have only one direction, not two, and one half-edge will be incorrectly tagged. Fortunately antennas are not typically necessary to produce offset buffers.
   

It should be noted that if you insert curves incrementally (insert one curve into the arrangement) an observer of the arrangement returns all generated and overlapped half-edges, which gives you the contained bounds of the contour inserted. I use this technique when inserting a low-side-count face into a very complex arrangement, to avoid re-sweeping a huge amount of data. Overlays and bulk inserts do not produce "per curve" announcements via the observer mechanism.

CGAL: It's All About the Mantissa

In a past post I described CGAL as having no rounding errors. It does this by using number types of variable size (using dynamically allocated memory per number!) so that it never runs out of digits. (It also maintains the numerator and denominator of fractions separately to avoid problems with repeating decimals.)

The advantage of this is that geometric algorithms that rely on precise calculations never go haywire due to rounding errors. For example, when using fixed-precision math (e.g. IEEE floats) the intersection of two near-parallel lines will be calculated inaccurately - sometimes with the intersection showing up miles from the original lines. CGAL always has more precision, so it avoids this problem.

But there is one down-side: when you perform a series of intersections, the result is exact numbers whose mantissas (the number of actual digits) have grown very long. And CGAL won't blink about making them even longer as you do more calculations.

Instead CGAL will become insanely slow.

I hit this case the other day. The first piece of processing I do is to combine a whole pile of vector data from OSM into one integrated map. While OSM is not particularly high precision (from a bits standpoint) the resulting intersecting points are calculated "perfectly", sometimes with very large mantissas.

I then wrote a piece of code to take a city block from that OSM map and perform some calculations to find the sidewalk calculation. The problem: the four corners of the city block were already very long numbers since they were the result of a CGAL calculation. Thus a long calculation on a long calculation becomes very slow.

The original algorithm took about 36 minutes for a fully optimized build to find all sidewalks in San Diego. That is way too slow, and unusable for our project.

I the put a rounding stage in: fore each corner of the block, I would convert it to a regular 64-bit IEEE float and then back to CGAL, throwing out any "extra" precision that CGAL was saving. Note that the 64-bit float already gives me better than 1 millimeter precision, which is more than overkill for a road. The algorithm run on the "simplified" data ran in 67 seconds.

Now there is one danger: if, due to mismatched road locations in OSM or conflicting edits, some of the "blocks" were really tiny (less than 1 mm) CGAL would have correctly built that block using infinite precision, and my "rounding" would have incorrectly reshaped those blocks, perhaps turning them inside out or in some other way damaging them.

So a necessary step to productizing this 'resolution reduction' is to do a sanity check on each resulting block. Fortunately most of the time if the block contains too-small-to-use data, we don't need the data in the first place.

Constitutional Opposition

One part of this post by Daring Fireball on the iPhone SDK licensing agreement made me chuckle:

If you are constitutionally opposed to developing for a platform where you’re expected to follow the advice of the platform vendor, the iPhone OS is not the platform for you. It never was. It never will be.

It inspired me to come up with a new quotable:

If you are constitutionally opposed to developing for a platform where you’re expected to follow the advice of the platform vendor, you should not be a computer programmer.

See also basically every post by Raymond Chen: "just because, in Win98SE2, you could call SomeRandomWin32API with a combination of NULL, -1, and Bill Gate's IQ and get an undocumented behavior that violates all of Microsoft's guidelines for applications development doesn't mean it will continue to work in Windows 7."

Thank You Jeeves, That Will Be All

The other day I went in to discover why a new piece of scenery code had mysteriously stopped working. Eventually I came to this:

(p,path.size()/2,def,degree,inExtrudeFunc,
inObjectFunc,inChecker,ag_mode_draped_obj);

Ah! Now it all makes sense. The code should have read:

AG_extrude_string(p,path.size()/2,def,degree,inExtrudeFunc,
inObjectFunc,inChecker,ag_mode_draped_obj);

After having done a global search, clearly I had hit the space bar by accident, nuking my function call. The charming thing is that C++ doesn't question why I have a giant list of paranthetical "stuff", it just blissfully compiles it into an expression that does...well, pretty much nothing.

Some of my other favorite C++ isms:

case a: do_it(); break;
b: do_x(); break; // no case, not illegal - now "b'' is a label!
defaultl: do_more(); break; // typo in default? That's a label too!

Of course we are all familiar with the fun that emerges from swapping = and ==. And having a stray semi-colon never hurt anything.

Propsman had an apt characterization: C++ is like an overly polite butler. "A...hamburger on the rocks, Sir? Certainly, Sir, I'll bring you one directly..."

Santa to Ben: You're An Idiot

A few months ago I posted a request (to Santa) for parallel command dispatch. The idea is simple: if I am going to render several CSM shadow map levels and the scene graph contained in each does not overlap, then each one is both (1) independent in render target and (2) independent in the actual work being done. Because the stuff being rendered to each shadow map is different, using geometry shaders and multi-layer FBOs doesn't help.* My idea was: well I have 8 cores on the CPU - if the GPU could slurp down and run 8 streams, it'd be like having 8 independent GLs running at once and I'd get through my shadow prep 8x as fast.

I was talking with a CUDA developer and finally I got a clue. The question at hand was whether CUDA actually runs parallel kernels. Her comment was that while you can queue multiple kernels asynchronously, the goal of such a technique is to keep the GPU busy - that is, to keep the GPU from going idle between batches of kernel processing. The technique of multiple kernels isn't necessary to keep the GPU fully busy, because even iwth hundreds of shader units, the kernel is going to run over thousands or tens of thousands of data points. That is, CUDA is intended for wildly parallel processing, so the entire swarm of "cores" (or "shaders"?) is still smaller than the number of units of work in a batch.

If you submit a tiny batch (only 50 items to work over) there's a much bigger problem than keeping the GPU hardware busy - the overhead of talking to the GPU at all is going to be worse than the benefit of using the GPU. For small numbers of items, the CPU is a better bet - it has better locality to the rest of your program!

So I thought about that, then turned around to OpenGL and promptly went "man am I an idiot". Consider a really trivial case: we're preparing an environment map, it's small (256 x 256) and the shaders have been radically reduced in complexity because the environment map is going to be only indirectly shown to the user.

That's still at least 65,536 pixels to get worked over (assuming we don't have over-draw, which we do). Even on our insane 500-shader modern day cards, the number of shaders is still much smaller than the amount of fill we have to do. The entire card will be busy - just for a very short time. (In other words, graphics are still embarrassingly parallel.)

So, at least on a one GPU card, there's really no need for parallel dispatch - serial dispatch will still keep the hardware busy.

So...parallel command dispatch? Um...never mind.

This does beg the question (which I have not been able to answer with experimentation): if I use multiple contexts to queue up multiple command queues to the GPU using multiple cores (thus "threading the driver myself") will I get faster command-buffer fill and thus help keep the card busy? This assumes that the card is going idle when it performs trivially simple batches that require a fair amount of setup.

To be determined: is the cost of a batch in driver overhead (time spent deciding whether we need to change the card configuration or real overhead (e.g. we have to switch programs and GPU isn't that fast at it). It can be very hard to tell from an app standpoint where the real cost of a batch lives.

Thanks to those who haunt the OpenGL forums for smacking me around^H^H^H^H^H^H^H^H^Hsetting me straight re: parallel dispatch.

* geometry shaders and multilayer FBO help, in theory, when the batches and geometry for each rendering layer are the same. But for a cube map if most of the scene is not visible from each cube face, then the work for each cube face is disjoint and we are simply running our scene graph, except now we're going through the slower geometry shader vertex path.

The Value Of Granularity

OpenGL is a very leaky abstraction. It promises to draw in 3-d. And it does! But it doesn't say a lot about how long that drawing will take, yet performance is central to GL-based games and apps. Filling in this gap is transient information about OpenGL and its current dominant implementations that isn't easy to come by - it comes from a mix of insight from other developers, connecting the dots, reading many disparate documents, and direct experimentation. This isn't easy for someone who isn't working full time as an OpenGL developer, so I figure there may be some value to blogging things I have learned the hard way about OpenGL while working on X-Plane.

OpenGL presents new functionality via extensions. (It also presents new functionality via version numbers, but the extensions tend to range ahead of the version numbers because the version number can only be bumped when allrequired extensions are available.) When building an OpenGL game you need a strategy for coping with different hardware with different capabilities. X-Plane dates back well over a decade, and has been using OpenGL for a while, so the app has had to cope with pretty much every extension being not available at one point or another.

Our overall strategy is to categorize hardware into "buckets". For X-Plane 9 we have 2.5 buckets:

• Pre-shader hardware, running on a fixed function pipeline.
• Modern shader enabled hardware, using shaders whenever possible.
• We have a few shaders that get cased off into a special bucket for the first-gen shader hardware (R300, NV25), since that hardware has some performance and capability limitations.

These buckets then get sliced up by features the user select, but these don't complicate the buckets - we simply make sure we can shade without per pixel lighting, for example, if the user wants higher framerate.

So here is what has turned out to be surprising: we were basically forced to allow X-Plane to run with a very granular set of extensions for debugging purposes. An example will ilWhat lustrate.

Using the buckets strategy you might say: "The shader bucket uses GLSL, FBOs, and VBOs. Any hardware in that category has all three, so don't write any code that uses GLSL but not FBOs, or GLSL but not VBOs." The idea is to save coding by reducing the combination of all possible OpenGL hardware (we have eight combos of these three extensions) to only two combinations (have them all, don't have them all).

What we found in practice was that being able to run in a semi-useful state without FBOs but with GLSL was immensely useful for in-field debugging. This is not a configuration we'd ever want to really support or use, but at least during the time period that we started using FBOs heavily, the driver support for them was spotty on the configurations we hit in-field. Being able to tell a user to run with --no_fbos was an invaluable differential to demonstrate that a crash or corrupt screen was related specifically to FBOs and not some other part of OpenGL.

As a result, X-Plane 9 can run with any of these "core" extensions in an optional mode: FBOs, GLSL, VBOs (!), PBOs, point sprites, occlusion queries, and threaded OpenGL. That list matches a series of driver problems we ran across pretty directly.

Maintaining a code base that supports virtually every combination is not sustainable indefinitely, and in fact we've started to "roll up" some of these extensions. For example, X-Plane 9.45 requires a threaded OpenGL driver, whereas X-Plane 9.0 would run without it. We remove support for individual extensions going missing when tech support calls indicate that "in field" the extension is now reliable.

At this point it looks like FBOs, threaded OpenGL, and VBOs are pretty much stable. But I believe that as we move forward into newer, weirder OpenGL extensions, we will need to keep another set of extensions optional on a per-feature basis as we find out the hard way what isn't stable in-field.

One More On VBOs - glBufferSubData

So if you survived the timing of VBO updates (or rather, my speculations on what is possible with VBO updates), now you're in a position to ask the question: how fast might glBufferSubData be? In particular, developers like myself are often astonished when glBufferSubData does things like block.

In a world before manual synchronizing of VBOs (via the 3.0 buffer management APIs or Apple's buffer range extensions) we can now see why a sub-data buffer on a streamed VBO might perform quite badly.

The naive code goes something like this:

1. Fill half the buffer with buffer sub-data.
2. Issue a draw call to that half of the buffer.
3. Flip which half of the buffer we are using and go back to step 1.

In other words, double buffering by dividing the buffer in half, or treating it like a ring buffer.

This implementation is going to perform terribly. T sub-data call is going to block until the previous draw call has completed, even though they use opposite halves of the buffer, and we'll lose all of our concurrency. Let's see if we can understand why.

If we go to respecify a VBO in AGP memory using glBufferSubData while that VBO is in progress, glBufferSubData must block; it can't rewrite the buffer until the last draw finishes because we would see the new vertices, not the old, or maybe half and half. In order for the "fill" to complete, the driver would have to be able to determine that the pending draws and the new fill are completely disjoint.

There are two reasons why the driver might not be able to figure this out:

1. You've drawn using glDrawElements, and thus the actual part of the vertex VBO you draw from is determined by the index table. The cost of figuring out the "extent" of this draw is to process all of the indices. The cure is worse than the disease. Any sane driver is going to simply assume that any part of the VBO could be used.
2. Let's assume you use glDrawRangeElements to tell the driver that you're really only going to use half the VBO. Even then, the structure to mark "locked" regions would be a complex one - a series of draws over overlapping regions would require a complex data structure. For this one special case, you're asking the drivers to replace a simple time-stamp based lock (e.g. this VBO is locked until this many commands have executed) with a dynamic range marking structure. If I were a driver writer I'd say "let's keep it simple and not eat this cost on all VBOs."

I think it's safe to assume that some implementations (and all if you use glDrawElements) are simply going to mark the entire VBO as in use until the draw happens, and thus the partial rewrite is going to block as if there was a conflict, even if there was not.

Can we do anything about this? Besides falling back to an "orphaned" approach where we get a fresh buffer each time, our alternative is to use the more exact APIs from ARB_map_buffer_range or APPLE_flush_buffer_range. With these APIs we can map only the part of the VBO we know is not in use, with the unsynchronized bit set to avoid blocking because the other half is pending. We can use flush explicit to then flush only the areas we modified. (With the 3.0 APIs we can also use the discard range option to simply say "we are rewriting what we map".)

Of course, this technique isn't without peril - all synchronization is up to the client. The main danger is an over-run: your app is so fast that it needs to modify a range that the GL isn't done with - we made it all the way around our ring buffer. Probably the safest way to cope with this is to put explicit fences in place to wait until the last dependent draw call that we issued is finished.