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

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

[Walter Roberson]

http://www.mathworks.com/matlabcentral/answers/6200-tutorial-how-to-ask-a-question-on-answers-and-get-a-fast-answer

Answer 1

[ANSHU MALINI]

      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

[MABROUK ABAZA]

I'm not understand what is m and n ? this code was not worked in my case

you can see the attached files

Any help ! thank you

[Muhammad Usman Saleem]

@any expert of MATLAB elaborate what are the m and n used in this question?

Please

[Walter Roberson]

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.

Answer 2

[Pai-Feng Teng]

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?