Quick Reference

3.1.4 plc

The plc command plots contours of a function z(x,y). If z is a two dimensional array of the same shape as the x and y mesh coordinate arrays, then

plc, z, y, x

plots contours of z. The levs= keyword specifies particular contour levels, that is, the values of z at which contours are drawn; by default you get eight equally spaced contours spanning the full range of z. Each contour is actually a set of polylines. Here is the algorithm Yorick uses to find contour polylines:

Each edge of the mesh has a z value specified at either end. For those edges with one end above and the other below the desired contour level, linearly interpolate z along the edge to find the (x,y) point on the edge where z has the level value.

Start at any such point, choose one of the two zones containing that edge, and move to the point on another edge of that zone, stepping across the new edge to the zone on its opposite side. Continue until you reach the boundary of the mesh or the point at which you started. If any points remain, repeat the process to get another disjoint contour. If you start at boundary points as long as any remain, you will first walk all open contours -- those which run from one point on the mesh boundary to another -- then all closed contours. Note that each contour level may consist of several polylines.

One additional complication can arise: Some zones might have two diagonal corners above the contour value and the two other corners below, so all four edges have points on the contour. In a sense, your mesh does not resolve this contour -- it could represent a true X-point where contour lines cross, in which case you want to move directly across the zone to the point on the opposite edge. Unfortunately, if you connect the opposite edges, you will make very ugly contour plots. (Even if you disagree on my taste in other places, trust me on this one.)

So the question is, which adjacent edge do you pick? If you don't make the same choice for every contour level you plot, contours at different levels cross (and so will your eyes when you try to interpret the picture). By default, plc will use a sort of minimum curvature algorithm: it turns the opposite direction that it turned crossing the previous zone. The triangulate= keyword to plc can be used both to force a particular triangulation of any or all zones in the mesh, and to return automatic triangulation decisions made during the course of the plc command. Again, if triangulation decisions become important to you, your mesh is probably not fine enough to resolve the contour you are trying to draw.

Making a good contour plot is an art, not an algorithm. Contours work best if used sparingly; more than a half dozen levels is likely to result in an unintelligible picture. I suggest you try drawing a very few levels on top of a filled mesh plot (draw with plf). You will need to select a contrasting color, and you may want to call plc once per level to get a different line type for each level.

If you need lots of levels, you should try the plfc function, which fills the regions between successive contour levels with different colors. Beyond about six levels, contour lines drawn with plc will be unintelligible, but with plfc you can get perhaps as many as twenty distinguishable levels.