Product Support

1212 - Using MESHGRID and GRIDDATA to Fit Vector Data and Plot Unevenly Spaced Data

1. Introduction
2. Plotting Evenly Spaced Data
3. Plotting Unevenly Spaced Data
4. Notes on GRIDDATA
5. Reconstructing Data by Interpolating Polled Data

In creating three-dimensional representations of data, it is often necessary to have data defined over evenly spaced, two-dimensional grids. The SURF, MESH, CONTOUR, and QUIVER functions are examples of graphics functions that use data in this form.

You can map data in the form of vectors or matrices containing X, Y, and Z coordinate points onto grids for use with the above plotting routines. MESHGRID is the most common way to create grids of X and Y point data.

Each (X,Y) point in the grid has a corresponding Z value, typically the height of a surface, Z = f(X,Y). After creating grids of X and Y points, you can use the GRIDDATA function to fit that surface over the grid domain using the original X, Y, and Z point data.

Mapping is often required when users perform field studies, such as measuring Z values at specific X-Y coordinates. For example, an oceanographer may be measuring the ocean's depth at specific latitude-longitude locations.

Evenly spaced data is best illustrated in the case of three spatial dimensions as depicted in Figure 1. In this case, each (X,Y) coordinate pair lies on an evenly spaced grid. For every coordinate pair in the grid, there is a corresponding Z coordinate.

Figure 1:   Illustration of data points on an evenly spaced X-Y grid

When measuring field study data, you generally end up with a table of (X,Y,Z) triplets. In the example above, that table could read as follows:

As an example, consider altitudes of a mountain. For each (X,Y) location on the surface of the earth, there is a measured altitude, Z. Each data point can then be represented by its (X,Y,Z) coordinate. When connected, all of these points form a three dimensional surface. For this to happen, it must be clear which points are adjacent with one another. Looking at the table of data above, it is difficult to visualize which points border each other. This section discusses how you can reorganize the data, and how this can be a powerful format for MATLAB to use.

Consider Figure 1 again. Imagine each point's (X,Y,Z) coordinate written at the base of the X-Y plane. There would be a two dimensional matrix of (X,Y,Z) coordinates. Now imagine making three matrices of that same size, representing each of the dimensions, X, Y, and Z. The matrices would be as follows:

1
1
1
1
1

2
2
2
2
2

3
3
3
3
3

4
4
4
4
4

5
5
5
5
5

1
2
3
4
5

1
2
3
4
5

1
2
3
4
5

1
2
3
4
5

1
2
3
4
5

152
103
 89
115
100

 89
  0
 13
100
 85

100
100
100
187
111

100
100
100
200
 97

100
100
100
111
 48

For three dimensional data, MATLAB requires each of the three dimensions to fit the form of a two dimensional matrix. (More generally, for N dimensional data, each dimension's data must fit the form of an (N-1) dimensional matrix. However, for the sake of example, this technical note refers only to two dimensions of grid space; the third dimension is the data.) To assist in converting your data from the form of Table 1 into the form of Table 2, MATLAB includes the functions MESHGRID and GRIDDATA . Together, these functions can define a grid space (the X and Y matrices), and reorganize your data (the Z matrix) into the corresponding matrix that fits the grid space.

MESHGRID can take specified vector ranges for your grid and return the grid matrices X and Y, such as those in Table 2. Using these grid matrices and your original data, GRIDDATA returns the data matrix, Z.

Below is an M-file that demonstrates this by using a basic example. In this example, the data is already defined at points in an evenly spaced X-Y grid. The grid represents a discretized set of points over which a function will be calculated. While this example introduces the use of MESHGRID, the use of GRIDDATA is unnecessary. An example using the two functions together to plot unevenly spaced data is described in Section 3.

(The following example also uses the LINSPACE function.)

% This example will create a surface describing the function:
% z = f(x,y) = x^2-y^2, over the range -5<x<5, -5<y<5 
% Define the data: 
% First define the ranges and resolution of the grid space:
xmin = -5;
xmax = 5;
ymin = -5;
ymax = 5;
gridresolution = 25; 


% Use linspace to define a vector of points representing
% the ranges of the x- and y-dimensions:
x = linspace(xmin, xmax, gridresolution);
y = linspace(ymin, ymax, gridresolution); 


% Use meshgrid to calculate the grid matrices:
% for the x- and y-coordinates
[X,Y] = meshgrid(x,y); 


% Now calculate the function using array operations to manipulate
% the matrices element by element. The result is a matrix of Z-data,
% the same size as the X and Y grid matrices:
Z = X.^2-Y.^2; 


% Plot the results by drawing the surface:
figure
surf(X,Y,Z)
xlabel X; ylabel Y; zlabel Z;
colormap(copper)

Figure 2:  Plot of evenly spaced data, calculated using array
operations on grid matrices defined using MESHGRID

While the data in the previous section was conveniently spaced evenly on a grid, often times field data is recorded at points that do not make a grid at all. The following is a more complicated example. In this example, X and Y are unevenly spaced vectors, and are not vertices of a rectangular array. The example uses MESHGRID to create an evenly spaced grid around the range of the data. This can be considered the X-Y "interpolation space". Then, using the original data and the X-Y "interpolation space", GRIDDATA calculates the interpolated Z data. To run the example, download the sample.mat file and place it on your MATLAB path.

0.1
0.1
0.1
0.1
0.2
0.2
0.2
0.2
0.2
0.2
0.3
0.3
0.3
0.3
0.3
0.4
0.4
0.4
0.4
0.5
0.5

0.5
1.0
1.5
5.0
0.5
1.0
1.5
2.5
4.0
5.0
0.5
2.0
3.5
4.0
5.0
0.5
1.0
4.0
4.5
1.5
2.0

0.0
2.0
8.0
8.0
0.5
1.0
4.0
6.0
5.0
7.0
3.0
6.0
4.0
2.0
0.0
5.0
5.0
0.0
1.0
4.0
4.0

0.5
0.5
0.5
0.6
0.6
0.6
0.6
0.6
0.7
0.7
0.7
0.7
0.7
0.8
0.8
0.8
0.8
0.9
0.9
0.9
0.9

3.0
4.5
5.0
1.0
2.5
3.5
4.0
5.0
0.5
1.5
3.0
3.5
4.5
1.0
2.0
3.0
5.0
0.5
2.5
4.0
5.0

5.0
0.0
0.0
3.0
6.0
5.0
3.0
1.0
4.0
4.0
7.0
6.0
3.0
6.0
5.0
7.0
4.0
6.0
7.0
7.0
4.0

The following M-file demonstrates how to fit this unevenly spaced data.

% Load the data and extract the (x,y,z) information:
load sample.mat 


% Determine the minimum and the maximum x and y values:
xmin = min(x); ymin = min(y);
xmax = max(x); ymax = max(y); 


% Define the resolution of the grid:
xres=10;
yres=15;


% Define the range and spacing of the x- and y-coordinates,
% and then fit them into X and Y
xv = linspace(xmin, xmax, xres);
yv = linspace(ymin, ymax, yres);
[Xinterp,Yinterp] = meshgrid(xv,yv); 


% Calculate Z in the X-Y interpolation space, which is an 
% evenly spaced grid:
Zinterp = griddata(x,y,z,Xinterp,Yinterp); 


% Generate the mesh plot (CONTOUR can also be used):
figure
mesh(Xinterp,Yinterp,Zinterp)
colormap(cool(8))
xlabel X; ylabel Y; zlabel Z;

The plot for this example is shown in Figure 3.

Figure 3:  Interpolated Z using GRIDDATA on unevenly spaced data

Mapping the original data onto a 10-by-15 grid produced a fairly smooth plot in Figure 3. To create a smoother surface plot, simply increase the resolution of the grid. For example, the following 20-by-20 interpolation space is used to produce the plot in Figure 4.

% Define the resolution of the grid:
xres=20;
yres=20;


% Define the range and spacing of the x- and y-coordinates,
% and then fit them into X and Y:
xv = linspace(xmin, xmax, xres);
yv = linspace(ymin, ymax, yres);
[Xinterp,Yinterp] = meshgrid(xv,yv); 


% Calculate Z in the X-Y interpolation space, which is an 
% evenly spaced grid:
Zinterp = griddata(x,y,z,Xinterp,Yinterp); 


% Generate the mesh plot (CONTOUR can also be used):
hf4 = figure;
mesh(Xinterp,Yinterp,Zinterp)
colormap(cool(8))
xlabel X; 
ylabel Y; 
zlabel Z;
        

Figure 4:  Interpolated Z using GRIDDATA on unevenly spaced data
This interpolation space has higher resolution than in Figure 3.

To check your work and to better visualize how the original data fits with the generated surface, overlay the original points on the generated surface. We do this by using the PLOT3 command. We can also show the points that are behind the mesh by using the HIDDEN command.

% Generate data point on the same axes with specified properties:
figure(hf4)
hold on
plot3(x,y,z,'marker','o','markerfacecolor','r','linestyle','none')
hidden off 


% Try rotating with the camera toolbar for more details.
% Choose "View-->Camera Toolbar" from the menu bar, and select
% the first toolbar button on the ensuing toolbar.
        

The plot for this example is continued in Figure 5.

Figure 5:  Interpolated Z data overlaid by original Z data

Several interpolation methods are available for use with GRIDDATA:

• Linear: triangle-based linear interpolation (default)
• Cubic: triangle-based cubic interpolation
• Nearest: Nearest Neighbor Interpolation
• v4: MATLAB 4 GRIDDATA method

More details can be found in the documentation on GRIDDATA.

Possible Warnings:

• Matrix is singular to working precision. This error means that there are multiple occurrences of the same [X,Y] pair in the data such that the output from GRIDDATA may contain NaNs and Infs.
• All the (X,Y) points must be distinct. This occurs in some versions of GRIDDATA. It means the same as the previously listed message.
• Data cannot be triangulated. This error generally means that the data forms an extremely small region when plotted.

To avoid either of the first two warnings, you can average the Z values at the repeated (X,Y) points before calling GRIDDATA.

In this example we reconstruct an image using random data points polled from a base image. The image is a matrix loaded from the durer.mat file, a data file shipped with MATLAB.

This is a two-part example. The first part shows the preparation of the random data from the image itself. The second part reconstructs the image. The sampfrac variable controls the number of points selected from the original image. It is the ratio of the number of sample points to the number of pixels in the original image. Try changing this ratio to see the effect on the output. If memory problems are encountered, try reducing the value of sampfrac.

% Data preparation
% Loading data from the durer.mat file
load durer
imdat=X;
figure
imagesc(imdat)
colormap(gray) 


% Ratio of sample points to number of pixels
sampfrac = 0.15;
samps = round(sampfrac*numel(imdat)); 


% Reset random seed to produce consistent output:
rand('state',0); 


% Randomly select sample points in the image:
x = ((size(imdat,2))*rand(samps,1));
y = ((size(imdat,1))*rand(samps,1)); 


% Find interpolated intensity value to each random point:
z = interp2([1:size(imdat,2)],[1:size(imdat,1)],imdat,x,y); 


% Image reconstruction 
% We now have unevenly spaced x, y, and z data. 
% Create interpolation plan - this will be the original image 
% size in pixels:
xd = 1:1:size(imdat,2);
yd = 1:1:size(imdat,1);
[Xd,Yd] = meshgrid(xd,yd);


% Form data to the interpolated plane - recreate image:
Z = griddata(x,y,z,Xd,Yd); 


% Plot results:
figure
imagesc(Z)
colormap(gray)
        

Figure 6:  Original Durer image

Figure 7:  Sampled Durer image (15%)