MATLAB®

Contour Plots

Functions for Creating Contour Displays

The contouring display functions compute, plot, and label isolines (contour lines) for one or more matrices. These displays vary according to whether they plot plain contour lines, filled contour lines, raised contours, or contours in concert with mesh or surface plots. Two of the functions listed below support contouring. The low-level contourc function computes isolines but does not plot them. The clabel function places elevation labels on previously generated contours.

Back to Top

Creating Simple Contour Plots

contour and contour3 display 2- and 3-D contours, respectively. They can be called with separate x, y, and z matrices, but need only one input argument — a z matrix interpreted as heights with respect to a plane. In this case, the contour functions determine the number of contours to display based on the minimum and maximum data values.

To explicitly set the number of contour levels displayed by the functions, you specify a second optional argument.

Contour Plot of the Peaks Function

The statements

[X,Y,Z] = peaks;
contour(X,Y,Z,20)

display 20 contours of the peaks function in a 2-D view.

The statements

[X,Y,Z] = peaks;
contour3(X,Y,Z,20)
h = findobj('Type','patch');
set(h,'LineWidth',2)
title('Twenty Contours of the peaks Function')

display 20 contours of the peaks function in a 3-D view and increase the line width to 2 points.

Back to Top

Labeling Contours

Each contour level has a value associated with it. clabel uses these values to display labels for 2-D contour lines. The contour matrix contains the values clabel uses for the labels. This matrix is returned by contour, contour3, and contourf and is described in The Contouring Algorithm.

clabel optionally returns the handles of the text objects used as labels. You can then use these handles to set the properties of the label string.

For example, display 10 contour levels of the peaks function,

Z = peaks;
[C,h] = contour(Z,10);

then label the contours and display a title.

clabel(C,h)
title({'Contour Labeled Using','clabel(C,h)'})

Note that clabel labels only those contour lines that are large enough to have an inline label inserted.

The 'manual' option enables you to add labels by selecting the contour you want to label with the mouse.

You can also use this option to label only those contours you select interactively.

For example,

clabel(C,h,'manual')

displays a crosshair cursor when your cursor is inside the figure. Pressing any mouse button labels the contour line closest to the center of the crosshair.

Back to Top

Filled Contours

The contourf displays a two-dimensional contour plot and fills the areas between contour lines. Use caxis to control the mapping of contour to color. For example, this filled contour plot of the peaks data uses caxis to map the fill colors into the center of the colormap.

Z = peaks;
[C,h] = contourf(Z,10);
caxis([-20 20])
title({'Filled Contour Plot Using','contourf(Z,10)'}) 

Back to Top

Specifying Contour Levels

The contouring functions permit you to specify the number of contour levels or the particular contour levels to draw. In the case of contour, the two forms of the function are contour(Z,n) and contour(Z,v). Z is the data matrix, n is the number of contour lines, and v is a vector of specific contour levels.

When you specify n (the number of contour levels to plot), you are setting the LevelStep property of the contourgroup. If you want to always draw n contour levels even when the data range of Z changes, you should obtain the contourgroup's handle and set its LevelStepMode to 'manual'.

When you specify v (a vector itemizing contour levels), you are setting the contourgroup's LevelList property. By default, the LevelList is recomputed whenever contours are redrawn. If you always want to depict the same contour levels even if the data changes, set the contourgroup's LevelListMode property to 'manual'.

Drawing a Single Contour

MATLAB^® functions do not differentiate between a scalar and a one-element vector. So, if v is a one-element vector specifying a single contour at that level, contour interprets it as the number of contour lines, not the contour level. Consequently, contour(Z,v) behaves in the same manner as contour(Z,n).

To display a single contour line, define v as a two-element vector with both elements equal to the desired contour level. For example, create a 3-D contour of the peaks function.

xrange = -3:.125:3;
yrange = xrange;
[X,Y] = meshgrid(xrange,yrange);
Z = peaks(X,Y);
contour3(X,Y,Z)

To display only one contour level at Z = 1, define v as [1 1].

v = [1 1]
contour3(X,Y,Z,v)

Example — Visualizing Contour Construction

You can think of a contour as the intersection of a 3–D surface with a horizontal plane. The intersection defines 0 or more level lines that trace contours. The level lines either form loops or terminate at the outer edges of the surface. Contour loops can intersect at saddle points, and therefore require special handling in their vicinity.

Run the following interactive code to visualize how contour lines are constructed. Use the slider to move the plane up or down through the range of z-values, and click the Plot Contour button to draw a contour line that delineates where the plane slices through the surface. Click the Plot Labels button to add a label to the contour you just plotted. Click Clear Contours to remove all the contours and labels.

% Create x, y, and z arrays for a parametric surface
[x y]=meshgrid(linspace(0,1,10),linspace(0,1,10));
z = .5*x + y - 1.5*x.*y;
% Display with the surface function in 3-D
fh = figure; colormap cool;
hpl = uipanel(fh,'Units','normalized','position',[.025 .025 .95 .95]);
s=surface('xdata',x,'ydata',y,'zdata',z,'cdata',z);
view(3);hold on;
% Display a second surface, a horixontal plane at z = 0
p=surface('xdata',[0 1;0 1],'ydata',[0 0; 1 1],...
          'zdata',[0 0; 0 0],'cdata',[0 0;0 0]);
set(p,'facealpha',.25,'facecolor','red'); % Make cut plane transparent
% Create a slider control for contour elevations
hs = uicontrol(hpl,'style','slider','min',0,'max',100,...
            'units','normalized','position',[.05 .05 .2 .05],...
            'sliderstep', [.01 .05]);
set(hs,'callback',...   % Tell the slider what it should do
    ['lvl=get(hs,''value'')/100;,' ...
     'set(p,''zdata'',[lvl lvl; lvl lvl]),' ...
     'set(hto,''string'',num2str(lvl)),' ...
     'set(hbc,''enable'',''on'')']);
 lvl = 0;       % Initialize the z-level of the cutting plane
% Create a label for the slider and a text box to show its value
hst = uicontrol(hpl,'Style','text', 'String','Z-level',...
                'units','normalized','Position',[.05 .10 .1 .05]);
hto = uicontrol(hpl,'Style','text', 'String','0',...
                'units','normalized','Position',[.13 .10 .1 .05]);
% Create a pushbutton control for drawing contours with CONTOUR3
hbc = uicontrol(hpl,'style','pushbutton','enable','off',...
    'string','Plot Contour',...
    'units','normalized','position',[.80 .05 .15 .05]);
set(hbc,'callback',['[C hc] = contour3(x,y,z, [lvl lvl],''r'');' ...
    'set(hbl,''enable'',''on''), set(hbe,''enable'',''on''),' ...
    'set(hbc,''enable'',''off'')']);
% Create a pushbutton control for labelling with CLABEL,
% which uses the "contour matrix" returned from CONTOUR3
hbl = uicontrol(hpl,'style','pushbutton','enable','off',...
    'string','Plot Labels',...
    'units','normalized','position',[.80 .90 .15 ,.05]);
set(hbl,'callback',['clabel(C, hc,''color'',''r'',' ...
   '''fontweight'',''bold'');' 'set(hbl,''enable'',''off''), '...
   'set(hbe,''enable'',''on'')']);
% Create a pushbutton to clear away the contours and labels
hbe = uicontrol(hpl,'style','pushbutton','enable','off',...
    'string','Clear Contours',...
    'units','normalized','position',[.05 .90 .15 .05]);
set(hbe,'callback',['delete(findall(gca,''color'',''r''));' ...
    'set(hbe,''enable'',''off'')']);

Here is what the figure and its controls look like with a contour plotted at the cut line:

See The Contouring Algorithm, below, for an explanation of how contour lines are computed.

Back to Top

Index Contours

You can index contours to visually emphasize certain contour levels. This technique, commonly used on topographic maps to highlight contours at set altitudes such as 25, 50, 75 ... meters above sea level, provides visual cues analogous to major ticks on a graph's axis. It is much easier to read a contour display that shows index contours because the heavier lines lessen the chance that one's eye will jump between adjacent contours in scanning across the plot.

Example — Specifying Index Contours

The following code example highlights contours at elevations of -6, -5, -4, ... 7 for the output of the peaks function.

1. Generate a data matrix to contour:

   z = peaks(100);

2. Compute 40 contour levels. Select contour levels so as to be round numbers; zlevs is the vector of contour levels to be plotted:

   zmin = floor(min(z(:))); zmax = ceil(max(z(:)));
   zinc = (zmax - zmin) / 40;
   zlevs = zmin:zinc:zmax;

3. Specify the vertical distance between index contours; here it is unity, but it can be any modulus of values in zlevs.

   zindex = 1;

4. Plot 2-D level lines with the CONTOUR function:

   [c2,hc2] = contour(z,zlevs);

5. Create Index Contours by thickening level lines every zindex units

   nc = get(hc2,'Children');
   for i = 1:length(nc)
      ud = get(nc(i),'UserData');
      if (mod(ud,zindex) == 0)
          set(nc(i),'LineWidth',2);
      end
   end

   A contour line will thicken with each call to set.

6. Annotate to identify the contouring parameters used

   s = sprintf('%s %g %s %g %s', 'Peaks Function Contoured at', ...
       zinc, 'Units, Indexed every', zindex, 'Units');
   title(s)

The loop of code in step 5 above works for contour but not contour3, because contour3 does not create contourgroup objects containing Children. To accomplish the same result with contour3 you must dereference the handle to the contours returned by contour3 (hc3, below) differently, as follows:

figure;
[c3,hc3] = contour3(z,zlevs);
for i = 1:length(hc3)
   ud = get(hc3(i),'UserData');
   if (mod(ud,zindex) == 0)
       set(hc3(i),'LineWidth',2);
   end
end
s = sprintf('%s %g %s %g %s',...
    'Peaks Function Contoured in 3-D at', ...
    zinc, 'Units, Indexed every', zindex, 'Units');
title(s)

Back to Top

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 either has associated X and Y matrices that locate each value in Z at the intersection of a row and a column, or these are inferred when they are unspecified. The width of rows and columns can vary, but typically is constant (i.e. Z is a regular grid).

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, then a contour at that level crosses the edge, and a linear interpolation is performed:

t=(c-Z0)/(Z1-Z0), where 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 next edge with z-values that bracket c is where the contour exits. Increment i, compute t for the edge, and then 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 below).

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 above 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 specifying all the contour lines. The format of the matrix is

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 value of the contour, 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

The circled values begin each definition of a contour line.

Back to Top

Changing the Offset of a Contour

The surfc and meshc functions display contours beneath a surface or a mesh plot. These functions draw the contour plot at the axes' minimum z-axis limit. To specify your own offset, you must change the ZData values of the contour lines. First, save the handles of the graphics objects created by meshc or surfc.

h = meshc(peaks(20));

The first handle belongs to the mesh or surface. The remaining handles belong to the contours you want to change. To raise the contour plane, increment the z coordinate of each contour line by some amount by resetting its Zdata value:

for i = 2:length(h);
    newz = get(h(i),'Zdata') + 5;
    set(h(i),'Zdata',newz)
end

Back to Top

Displaying Contours in Polar Coordinates

You can contour data defined in the polar coordinate system. As an example, set up a grid in polar coordinates and convert the coordinates to Cartesian coordinates.

[th,r] = meshgrid((0:5:360)*pi/180,0:.05:1);
[X,Y] = pol2cart(th,r);

Then generate the complex matrix Z on the interior of the unit circle.

Z = X+i*Y;

X, Y, and Z are points inside the circle.

Create and display a surface of the function .

f = (Z.^4-1).^(1/4);
surf(X,Y,abs(f))

Display the unit circle beneath the surface using the statements

hold on
surf(X,Y,zeros(size(X)))
hold off

Labeling the Graph

These statements add labels.

xlabel('Real','FontSize',14);
ylabel('Imaginary','FontSize',14);
zlabel('abs(f)','FontSize',14);

Contours in Cartesian Coordinates

These statements display a contour of the surface in Cartesian coordinates and label the x- and y-axis.

contour(X,Y,abs(f),30)
axis equal
xlabel('Real','FontSize',14);
ylabel('Imaginary','FontSize',14);

Contours on a Polar Axis

You can also display the contour within a polar axes. Create a polar axes using the polar function, and then delete the line specified with polar.

h = polar([0 2*pi], [0 1]);
delete(h)

With hold on, display the contour on the polar grid.

hold on
contour(X,Y,abs(f),30)

Back to Top

Preparing Data for Contouring

The various contour plotting functions, as well as the mesh and surface families of functions, accept 2-D matrices as inputs. For most applications, these input grids will represent continuous functions of two variables or relatively continuous fields of data. In many applications, source data might consist of z-values sampled over a two-dimensional domain in an irregular fashion, such as discrete spot elevations from GPS measurements (in the form of x, y, and z data vectors). To prepare such data for contour or mesh display, you need to interpolate it in some fashion.

There are several MATLAB methods for interpolating data into vectors, grids and triangulated (Delaunay) tessellations. Input observations can be one-, two-, three- or higher-dimensional. By choosing and using these functions carefully you can control parameters and constraints for interpolation to model your assumptions about the underlying nature of the raw data. Typically, you would use the interp2, meshgrid, and griddata functions to interpolate z-values for scattered x-y data points into a 2-D grid. See in the MATLAB Mathematics documentation for discussion and examples of data interpolation using these and other functions.

If the surface you are contouring is "noisy," contours depicting it will exhibit jaggedness. When you analyze and explore such data, you can filter it to attentuate high-frequency variations. One way to do this is with a convolution (with conv2 or filter2) filter, as the following example demonstrates:

Example — Smoothing a Matrix for Plotting Contours

The conv2 and filter functions can remove high-frequency components from a matrix representing a continuous surface or field to make the underlying data easier to visualize. This requires few assumptions

1. Create a function of two variables and plot contour lines at a specified, fixed interval.

   Z = peaks(100);
   figure; 
   set (gcf,'position',[400,100,600,600], 'color','w')
   subplot(2,2,1); 
   cl = [-7:1:10];         % Define contour levels for all plots
   contour(Z, cl)
   axis([0 100 0 100]); colormap autumn;
   set(gca,'Xtick',[0 100],'Ytick',[0 100]);
   title('Peaks Surface (underlying data)')

2. Add uniform random noise with mean of 0 to the surface and plot resulting contours. Irregularities in the contours tend to obscure the trend of the data:

   ZN = Z + rand(100) - .5;
   subplot (2,2,2)
   contour(ZN, cl)
   axis([0 100 0 100]);
   set(gca,'Xtick',[0 100],'Ytick',[0 100]);
   title('Peaks Surface (noise added)')

3. Specify a 3x3 convolution kernal, F, for smoothing the matrix and use the conv2 function to attenuate high spatial frequencies in the surface data:

   F = [.05 .1 .05; .1 .4 .1; .05 .1 .05];
   ZC = conv2(ZN,F,'same');

4. Visually compare the smoothed surface to the original and the noisy ones:

   subplot (2,2,3)
   contour(ZC, cl)
   axis([0 100 0 100]);
   set(gca,'Xtick',[0 100],'Ytick',[0 100]);
   title('Noisy Surface (smoothed once)')

5. Smooth the surface one more time using the same operator and compare (a larger or more uniform kernal could have achieved this in one pass):

   ZC2 = conv2(ZC,F,'same');
   subplot (2,2,4)
   contour(ZC2, cl)
   axis([0 100 0 100]);
   set(gca,'Xtick',[0 100],'Ytick',[0 100]);
   title('Noisy Surface (smoothed twice)')

   

Back to Top

Direction and Velocity Vector Graphs

Interactive Plotting