how to write function for thessian polygon method in hydrology by using matlab ?

14 次查看（过去 30 天）
ANSHU MALINI
ANSHU MALINI 2012-1-28
回答： Pai-Feng Teng 2018-10-2
Consider a rectangular area with vertices ([0,0],[10,0],[10,10] & [0,10]). For the gages given in Table 1, estimate
(a) Thiessen weights,
(b) average rainfall using Thiessen weights, and
(c)average rainfall using arithmetic average.
Plot the area, gages and associated Thiessen polygons on a figure.
Table 1: Rain Gages
Gage Co-ordinates Rainfall (mm)
A (2,5) 10.0
B (8,4) 7.0
C (5,5) 5.0
D (5,9) 12.0
E (6,1) 2.0
F (9,9) 8.0

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

采纳的回答

ANSHU MALINI
ANSHU MALINI 2012-2-6
编辑：Walter Roberson 2017-1-12
function [theiss_avg_rainfall ] = hw2(X,Y,x,y,m,n,rainfall)
%HW2: theissan polygon analysis-calulation of thessian weights for calculating
%average rainfall
% catchment area coordinates
X = input('X =');
Y = input('Y =');
plot(X,Y,'r','LineWidth',3)
% gages
x = [ 2 8 5 5 6 9]'
y = [ 5 4 5 9 1 9]'
gages = [x ...
y];
hold on
scatter(x,y,'b','Marker','*')
text([2; 8; 5; 5; 6; 9],[5; 4; 5; 9; 1; 9],['A';'B'; 'C'; 'D'; 'E'; 'F'],...
'HorizontalAlignment','right',...
'FontSize',12,...
'FontName','Impact');
voronoi(x(:),y(:))
% calculations of grid points
m = input('m =');
n = input('n =');
for a=1:m
for b=1:n
gx(a,b)=b;
gy(b,a)=b;
scatter(a,b,...
'r','Marker','.')
end
end
p = cat(2,gx,gy);
% grid points
grdpts = (m*n);
r = reshape(p,grdpts,2);
index=1:grdpts;
for d1 = 1:6
for d2 = 1:grdpts
distn = pdist2(gages,r);
end
[~, pos] = min(distn);
indx= index(pos);
end
% Nearest gage to each grid point
numb =zeros(grdpts,1);
for g=1:6
numb = find(indx==g);
wt(g,1) = length(numb);
end
numb;
wt;
% theissen weights
theisswts = wt./grdpts
rainfall = input('rainfall =');
%Avg rainfall by theissen ploygon analysis
theiss_avg_rainfall=sum(rainfall.*theisswts)
% Average rainfall by arithmetic mean
aritm_avg_rainfall = mean(rainfall)
end
  3 个评论
显示 1更早的评论
Walter Roberson
Walter Roberson 2017-1-12
The code creates an m x n grid and finds the distances between each point in the grid and each of the gages and uses that to figure out the closest gage to each grid point. The number of grid points closest to each gage is calculated and used to create a weight matrix that is used in the remaining calculations.
I speculate that the idea is that water is assumed to run from each grid point into the closest gage, and so some gages would have a larger catch-basin than others, and so to calculate average rainfall you need to divide by the size of the catch basin. If I am correct then your m and n should only be as large as is required to account for the known catch basins.

请先登录，再进行评论。

更多回答（1 个）

Pai-Feng Teng
Pai-Feng Teng 2018-10-2
I tried something very similar, with 9 gages, but I couldn't find the answer. Here is the error message:
Index exceeds matrix dimensions.
Error in HW8 (line 43) indx= index(pos);
Why does that happen, and how can I fix it?

类别

在 Help Center 和 File Exchange 中查找有关 Environment and Clutter 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by