How to create contour of water table at 0.2 m intervals?

Question

Asanka Subasinghe 2020-4-13

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/517450-how-to-create-contour-of-water-table-at-0-2-m-intervals

评论： Asanka Subasinghe 2020-4-14

I have some wells that are located by a lake. I have the ground watertable elevation at each well, as well as the elevation at the lake. I created x and y coordinates for these elevations and am now trying to combine it to create contour lines of the water table. How can I do this at 0.2 m intervals?

Not sure what the problem is with using the reshape function here as well.

clc;clear all;close all
% x and y coordinates of well locations and lake points in cm  (3.4cm/km)
w1=[5.9,5.6];  
w8=[8,7.6];
w17=[11.9,2];
w21=[9.7,5.3];
w27=[4.6,4.6];
w56=[5.9,1.7];
w58=[5.9,1];
w59=[6.7,4.1];
lakepoint1=[4,6.55];
lakepoint2=[5.75,6.9];
lakepoint3=[7.7,7.65];
% Combined coordinates into a single matrix 
wells=[w1;w8;w17;w21;w27;w56;w58;w59;lakepoint1;lakepoint2;lakepoint3];
% Converted well locations to metres
scaledwells=wells*0.441176*1000;
%Groundwater table elevations at each well and at 3 points of the lake. (in metres)
gwtablew1=236.07;
gwtablew8=236.65;
gwtablew17=240.37;
gwtablew21=238.13;
gwtablew27=235.96;
gwtablew56=236.66;
gwtablew58=237.70;
gwtablew59=237.90;;
gwtablelake1=234.19;
gwtablelake2=234.19;
gwtablelake3=234.19;
%Combined into single matrix below
gwtable=[gwtablew1;gwtablew8;gwtablew17;gwtablew21;gwtablew27;gwtablew56;gwtablew58;gwtablew59;gwtablelake1;gwtablelake2;gwtablelake3];
%created scatter plot of x and y coordinates of wells and lake points)
scatter(scaledwells(:,1),scaledwells(:,2))
xlabel('m')
ylabel('m')
grid on
%Trying to produce contours
scaledwellsx=scaledwells(:,1)
scaledwellsy=scaledwells(:,2)
nx=length(scaledwellsx);
ny=length(scaledwellsy);
z=reshape(gwtable,ny,nx)

I am receiving the error:

Error using reshape
To RESHAPE the number of elements must not change.
Error in gwass1 (line 45)
z=reshape(gwtable,ny,nx)

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Ameer Hamza 2020-4-13

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/517450-how-to-create-contour-of-water-table-at-0-2-m-intervals#answer_425697

在 MATLAB Online 中打开

You first need to interpolate your data to a grid before plotting contour. Try the following code

clc;clear all;close all
% x and y coordinates of well locations and lake points in cm  (3.4cm/km)
w1=[5.9,5.6];  
w8=[8,7.6];
w17=[11.9,2];
w21=[9.7,5.3];
w27=[4.6,4.6];
w56=[5.9,1.7];
w58=[5.9,1];
w59=[6.7,4.1];
lakepoint1=[4,6.55];
lakepoint2=[5.75,6.9];
lakepoint3=[7.7,7.65];
% Combined coordinates into a single matrix 
wells=[w1;w8;w17;w21;w27;w56;w58;w59;lakepoint1;lakepoint2;lakepoint3];
% Converted well locations to metres
scaledwells=wells*0.441176*1000;
%Groundwater table elevations at each well and at 3 points of the lake. (in metres)
gwtablew1=236.07;
gwtablew8=236.65;
gwtablew17=240.37;
gwtablew21=238.13;
gwtablew27=235.96;
gwtablew56=236.66;
gwtablew58=237.70;
gwtablew59=237.90;;
gwtablelake1=234.19;
gwtablelake2=234.19;
gwtablelake3=234.19;
%Combined into single matrix below
gwtable=[gwtablew1;gwtablew8;gwtablew17;gwtablew21;gwtablew27;gwtablew56;gwtablew58;gwtablew59;gwtablelake1;gwtablelake2;gwtablelake3];
%created scatter plot of x and y coordinates of wells and lake points)
scatter(scaledwells(:,1),scaledwells(:,2))
xlabel('m')
ylabel('m')
grid on
hold on
%Trying to produce contours
scaledwellsx=scaledwells(:,1);
scaledwellsy=scaledwells(:,2);
[X,Y] = meshgrid(sort(scaledwellsx), sort(scaledwellsy));
F = scatteredInterpolant(scaledwellsx,scaledwellsy,gwtable);
Z = F(X,Y);
contour_levels = min(Z,[],'all'):0.2:max(Z,[],'all');
contour(X, Y, Z, contour_levels);

9 个评论
显示 7更早的评论隐藏 7更早的评论

Tommy 2020-4-14

在 MATLAB Online 中打开

Alternatively, if you want to keep the contour lines at 0.2 m intervals as initially requested, you can use the clabel function to add your labels. clabel allows you to specify which contour lines should be labeled. For example, to label only every 4th contour line:

contour_levels = min(min(Z)):0.2:max(max(Z));
[M, c] = contour(X, Y, Z, contour_levels);
clabel(M, c, contour_levels(1:4:end));

or to label only the upper half of contour lines:

clabel(M,c, contour_levels(end/2:end));

You might have to play around with it to get the labels exactly where you want.

Below was my attempt to find and plot the flow directions using gradient:

[DX,DY] = gradient(Z);
zi = arrayfun(@(i) find(X(1,:)==scaledwellsx(i)&Y(:,1)==scaledwellsy(i),1), 1:numel(scaledwellsx));
% ^ linear indices within Z (and DX, DY) of each well/lake
quiver(scaledwellsx',scaledwellsy',-DX(zi),-DY(zi));