Extract interpolated data from contourf?

Hello,

I have an elevation matrix (xq,yq,zq) that I have refined using meshgrid. I am hoping to extract a continuous 1 m contour, but my data is a little interrupted.

contourf provides a really nice output plot of continuous filled data between the extracted 1 m contour vertices. That's exactly what I am after. However, when I extract the contourline x,y data, I only obtain the vertices and not the data interpolated in between the vertices.

How do I extract the interpolated data shape that contourf plots? My goal is to fill in the gaps between the black and red lines shown in the image, tracking along the border of the blue data.

V = [1,1]
[M, h] = contourf(xq,yq,zq, V);  %% really nice plot of contour (black) + filled/interpolated data (blue). See Fig. 1/
contourTable = getContourLineCoordinates(M);
plot(contourTable_20170918_1.X, contourTable_20170918_1.Y, 'r.', 'MarkerSize', 10)  %% only contour vertices (red). See Fig. 2.

---> How do I extract the data between the black (Fig. 1) and red (Fig. 2) data, following the blue data border from contourf (Fig. 1)?

Answer 1

[Star Strider]

The contour functions return (x,y) matrices for each contour, and in all likelihood, the ‘x’ values are not the same for both contours, so it may be necessary to use interp1 to interpolate the ‘y’ vectors to a common reference. (This is actually straightforward, although you may want to do this with the ‘y’ vectors and interpolate the ‘x’ vectors instead in this instance.) After that, use either trapz or cumtrapz to calculate the areas, depending on what you want for the result.

The actual data and code would be necessary to provide a specific solution.

EDIT — (12 Sep 2023 at 15:35)

That might go something like this with your contours —

x = linspace(0,10);
y = linspace(0, 20);
[X,Y] = ndgrid(x,y);
Z = exp(-((X-2).^2+0.1*(Y-10).^2)/50);

figure
surfc(X,Y,Z)
xlabel('X')
ylabel('Y')
zlabel('Z')

figure
[M,h] = contourf(X,Y,Z, 'ShowText',1);
xlabel('X')
ylabel('Y')

Lvls = h.LevelList

Lvls = 1×8

0.2276 0.3000 0.4000 0.5000 0.6000 0.7000 0.8000 0.9000

Extract = Lvls([4 5])

Extract = 1×2

0.5000 0.6000

for k = 1:numel(Extract)
    idx = find(M(1,:)==Extract(k));
    Len = M(2,idx);
    rngv = idx+1:idx+Len;
    xv{k} = M(1,rngv);
    yv{k} = M(2,rngv);
end

vi = linspace(yv{k}(1), yv{k}(end), 250);                   % Common Interpolation Vector
for k = 1:numel(Extract)
    xi(k,:) = interp1(yv{k}, xv{k}, vi);                    % Interpolated 'x' Vector
    Area(k) = trapz(vi, xi(k,:));
end

AreaBetweenContours = diff(Area);

fprintf('\nThe Area between contour %.2f and %.2f is %.2f\n', Extract, AreaBetweenContours)

The Area between contour 0.50 and 0.60 is 17.78

.

[Laura Szczyrba]

Hi Star Strider,

Thank you for your response. I apologize, the wording in my question was misleading. I am not trying to calculate the area value of the blue zone, though I see how my wording indicated that. I updated my question to clarify.

I am hoping to extract the shape produced by contourf (the blue shaded zone). Specfficially, I want the border of the blue shaded zone that fills in the gaps of the extracted contour vertices.

In Figure 1, you can see the black vertices and the fill in between that contourf interpolated based on the elevation data. I would like to extract that continuous border. Any idea?

[Star Strider]

I still do not understand exactly what you want. If you only want to extract the border coordinats and not do any integrations or other processing, the (x,y) coordinates of those lines aer givne in the ‘xv’ and‘yv’ cell arrays in my code. They contain the ‘x’ and ‘y’ vectors respectively for the chosen contours. The edges of the region between them are then given by:

xedge1 = [xv{1}(1) xv{2}(1)] % Upper Edge

yedge1 = [yv{1}(1) yv{2}(1)]

xedge2 = [xv{1}(end) xv{2}(end)] % Lower Edge

yedge2 = [yv{1}(end) yv{2}(end)]

Those will connect them with straight lines if you need to do that.

x = linspace(0,10);
y = linspace(0, 20);
[X,Y] = ndgrid(x,y);
Z = exp(-((X-2).^2+0.1*(Y-10).^2)/50);

figure
surfc(X,Y,Z)
xlabel('X')
ylabel('Y')
zlabel('Z')

figure
[M,h] = contourf(X,Y,Z, 'ShowText',1);
xlabel('X')
ylabel('Y')

Lvls = h.LevelList

Lvls = 1×8

0.2276 0.3000 0.4000 0.5000 0.6000 0.7000 0.8000 0.9000

Extract = Lvls([4 5])

Extract = 1×2

0.5000 0.6000

for k = 1:numel(Extract)
    idx = find(M(1,:)==Extract(k));
    Len = M(2,idx);
    rngv = idx+1:idx+Len;
    xv{k} = M(1,rngv);
    yv{k} = M(2,rngv);
    Q = [xv{k}(1) xv{k}(end); yv{1}(1) yv{k}(end)]
end

Q = 2×2

6.9656 6.9656 20.0000 0

Q = 2×2

5.9422 5.9422 20.0000 0

vi = linspace(yv{k}(1), yv{k}(end), 250);                   % Common Interpolation Vector
for k = 1:numel(Extract)
    xi(k,:) = interp1(yv{k}, xv{k}, vi);                    % Interpolated 'x' Vector
    Area(k) = trapz(vi, xi(k,:));
end

AreaBetweenContours = diff(Area);

fprintf('\nThe Area between contour %.2f and %.2f is %.2f\n', Extract, AreaBetweenContours)

The Area between contour 0.50 and 0.60 is 17.78

xedge1 = [xv{1}(1) xv{2}(1)];                               % Upper Edge
yedge1 = [yv{1}(1) yv{2}(1)];
xedge2 = [xv{1}(end) xv{2}(end)];                           % Lower Edge
yedge2 = [yv{1}(end) yv{2}(end)]; 

figure
plot(xv{1}, yv{1})
hold on
plot(xv{2}, yv{2})
plot(xedge1, yedge1)
plot(xedge2, yedge2)
hold off
axis([-1 11  -1 21])
title('Extracted Contour')

.

[Laura Szczyrba]

Shoot, I apologize for the confusion.

Yes, I have extracted the vertices in the 'xv' and 'yv' cell arrays for one contour line (see the red points in my Figure 2 and the black points in Figure 1). Especially in Figure 2, you can see how the contour vertices are discontinuous.

If you look at Figure 1, notice how the blue shaded region fills in the gaps in between the extracted vertices (black points). In order to produce this fill, the function contourf must interpolate between the contour vertices based on the zq data. It looks great. I really want to extract the matrix that represents the shape of the fill with the ultimate goal of extracting the data along the border of the blue shape to fill in the gaps in the extracted vertices (black points).

For the life of me, I can't find the data representing that blue fill in neither the contour matrix nor contour object. Do you know if this data is associated with the contourf outputs or is there a better way for me to accomplish this manually?

[Star Strider]

I am not certain what you want.

My latest guess —

x = linspace(0,10);
y = linspace(0, 20);
[X,Y] = ndgrid(x,y);
Z = exp(-((X-2).^2+0.1*(Y-10).^2)/50);

figure
surfc(X,Y,Z)
xlabel('X')
ylabel('Y')
zlabel('Z')

figure
[M,h] = contourf(X,Y,Z, 'ShowText',1);
xlabel('X')
ylabel('Y')

Lvls = h.LevelList

Lvls = 1×8

0.2276 0.3000 0.4000 0.5000 0.6000 0.7000 0.8000 0.9000

Extract = Lvls([4 5])

Extract = 1×2

0.5000 0.6000

for k = 1:numel(Extract)
    idx = find(M(1,:)==Extract(k));
    Len = M(2,idx);
    rngv = idx+1:idx+Len;
    xv{k} = M(1,rngv);
    yv{k} = M(2,rngv);
    Q = [xv{k}(1) xv{k}(end); yv{1}(1) yv{k}(end)]
end

Q = 2×2

6.9656 6.9656 20.0000 0

Q = 2×2

5.9422 5.9422 20.0000 0

vi = linspace(yv{k}(1), yv{k}(end), 250);                   % Common Interpolation Vector
for k = 1:numel(Extract)
    xi(k,:) = interp1(yv{k}, xv{k}, vi);                    % Interpolated 'x' Vector
    Area(k) = trapz(vi, xi(k,:));
end

AreaBetweenContours = diff(Area);

fprintf('\nThe Area between contour %.2f and %.2f is %.2f\n', Extract, AreaBetweenContours)

The Area between contour 0.50 and 0.60 is 17.78

xedge1 = [xv{1}(1) xv{2}(1)];                               % Upper Edge
yedge1 = [yv{1}(1) yv{2}(1)];
xedge2 = [xv{1}(end) xv{2}(end)];                           % Lower Edge
yedge2 = [yv{1}(end) yv{2}(end)]; 

figure
patch([xv{1} flip(xv{2})], [yv{1} flip(yv{2})], 'b', 'FaceAlpha',0.5)
hold on
plot(xv{1}, yv{1})
plot(xv{2}, yv{2})
plot(xedge1, yedge1)
plot(xedge2, yedge2)
hold off
axis([-1 11  -1 21])
title('Filled Extracted Contour')

.