Matlab如何计算等高线？

Question

Matlab使用哪种算法生成等高线？换句话说，它如何将网格上的等值数据转换为一组线条？

我想要的是：

• 获取位于等高线上点的局部标准？
• 捕获所有等高线的全局过程？

我不需要关于底层代码的详细说明，但了解一般原理对我解释输出结果会有帮助。我在研究中使用contour（及其衍生品），希望了解此步骤引入的数值误差。

这看起来是一个非常简单的问题，但我在Matlab的文档中找不到任何解释，在SO或其他地方也没有找到任何信息。如果最终证明很容易找到，请原谅我的无知。

- tvo

---

contourc 是 contour、contour3 和 contourf 使用的底层函数。这是一个内置函数，因此它是编译代码，用户无法看到。在较新版本的MATLAB中，内联文档中并没有对算法进行真正的讨论。旧版本的MATLAB可能会有更多的内联细节。 - sco1

@excaza 感谢您的评论并列出这些相关函数。我知道它们的存在，但对其他人可能有帮助。我没有考虑过寻找代码，但由于它已经编译了，所以也不会起作用。 - tvo

@tvo 对我来说找到答案并不容易 ;-) 最新版本在 http://www.mathworks.com/help/pdf_doc/matlab/graphg.pdf 中省略了算法的细节。 - Alessandro Jacopson

谢谢。干得好！+1并接受。 - tvo

Answer 1

来自MATLAB®-图形-R2012a，从第5-73页到第5-76页：

The Contouring Algorithm

The contourc function calculates the contour matrix for the other contour functions. It is a low-level function that is not called from the command line. The contouring algorithm first determines which contour levels to draw. If you specified the input vector v, the elements of v are the contour level values, and length(v) determines the number of contour levels generated. If you do not specify v, the algorithm chooses no more than 20 contour levels that are divisible by 2 or 5.

The height matrix Z has associated X and Y matrices that locate each value in Z at the intersection of a row and a column, or these matrices are inferred when they are unspecified. The row and column widths can vary, but typically they are constant (i.e., Z is a regular grid). Before calling contourc to interpolate contours, contourf pads the height matrix with an extra row or column on each edge. It assigns z-values to the added grid cells that are well below the minimum value of the matrix. The padded values enable contours to close at the matrix boundary so that they can be filled with color. When contourc creates the contour matrix, it replaces the x,y coordinates containing the low z-values with NaNs to prevent contour lines that pass along matrix edges from being displayed. This is why contour matrices returned by contourf sometimes contain NaN values. Set the current level, c, equal to the lowest contour level to be plotted within the range [min(Z) max(Z)]. The contouring algorithm checks each edge of every square in the grid to see if c is between the two z values for the edge points. If so, a contour at that level crosses the edge, and a linear interpolation is performed:

t=(c-Z0)/(Z1-Z0)

Z0 is the z value at one edge point, and Z1 is the z value at the other edge point.

Start indexing a new contour line (i=1) for level c by interpolating x and y:

cx(i) = X0+t*(X1-X0)
cy(i) = Y0+t*(Y1-Y0)

Walk around the edges of the square just entered; the contour exits at the next edge with z values that bracket c. Increment i, compute t for the edge, and then compute cx(i) and cy(i), as above. Mark the square as having been visited. Keep checking the edges of each square entered to determine the exit edge until the line(cx,cy) closes on its initial point or exits the grid. If the square being entered is already marked, the contour line closes there. Copy cx, cy, c, and i to the contour line data structure (the matrix returned by contouring functions, described shortly).

Reinitialize cx, cy, and i. Move to an unmarked square and test its edges for intersections; when you find one at level c, repeat the preceding operations. Any number of contour lines can exist for a given level. Clear all the markers, increment the contour level, and repeat until c exceeds max(Z). Extra logic is needed for squares where a contour passes through all four edges (saddle points) to determine which pairs of edges to connect. contour, contour3, and contourf return a two-row matrix that specifies all the contour lines:

C = [ value1 xdata(1) xdata(2)...
        numv ydata(1) ydata(2)...]

The first row of the column that begins each definition of a contour line contains the contour value, as specified by v and used by clabel. Beneath that value is the number of (x,y) vertices in the contour line. Remaining columns contain the data for the (x,y) pairs. For example, the contour matrix calculated by C = contour(peaks(3)) is as follows.

The circled values begin each definition of a contour line.

- Alessandro Jacopson

Answer 2

您可以了解“Marching Squares”（https://en.wikipedia.org/wiki/Marching_squares）。

通常，沿着网格单元的边缘进行线性插值，得到轮廓点，您可以将它们连接起来形成折线。

如果单元格过于粗糙，则可以进行初步插值（例如双三次插值），以精细化网格。

- user1196549

---

但是不知道Matlab实际使用了什么。 - user1196549

+1 这是一个很好的参考。如果我需要开发一个算法，我基本上会采用相同的方法。这似乎是Matlab使用的方法。我可以验证从数据到等高线的线性插值是否精确（最多四舍五入）。嗯，Matlab的“linear”插值实际上是双线性的，但这在边缘上不应该有影响。你觉得呢？ - tvo

@tvo：关于什么？ - user1196549

事后看来，我才意识到在之前的工作中已经隐式地进行了这个验证（使用双线性插值到轮廓线上）。 - tvo

@tvo：这没有意义。沿着网格边缘没有双线性插值这样的东西。 - user1196549