Untitled

function S = firesim( N, alpha, K, v, theta, M, T )
% FIRESIM simulates a fire in an area
%
% 
% Call:     firesim can be called in 4 different ways as shown below:
%           S = firesim(N, alpha, K) or
%           S = firesim(N, alpha, K, v, theta) or
%           S = firesim(N, alpha, K, v, theta, M) or
%           S = firesim(N, alpha, K, v, theta, M, T)
% 
% Inputs:   N, the size of the area of the fire
%           alpha, the spread of vegetation in the area.
%           K, the time of the simulation.
%           v, the speed of the wind. 0<= v < 1
%           theta, the angle of the wind
%           M, a NxN array with values that describe the flammability in
%              the area.
%           T, a vector with timestamps for changing the wind. The length
%              of the vector must be the same length as v and theta.
% 
% Outputs:  S, a NxN array with the simulated area after the forestfire

% Author:   Nicolas Borup - s083124
% Date:     May 2, 2011

% Initial values
midt = round(N/2);
S = zeros(N);
S(midt,midt) = 0.1;

if (alpha>1) || (alpha<=0)
    error('Alpha value does not meet the range restrictions');
end;

A = [3/8 3/4 3/8;
    3/4 0 3/4;
    3/8 3/4 3/8];
B = [0 0 0;
    0 3*sqrt(3) 0;
    0 0 0];
W = (1-alpha)*A+alpha*B;

% Checking for arguments and setting default values
if nargin == 3
    M = ones(N,N);
    T = 0;
    theta = 0;
    v = 0;
elseif nargin == 5
    if (0<=v) && (v<1) && (0<= theta) && (theta< 2*pi)
        M = ones(N,N);
        T = 0;
    else
        error('The specified values does not meet the range restrictions');
    end;
elseif nargin == 6
    if (0<=v) && (v<1) && (0<= theta) && (theta< 2*pi)
        T = 0;
    else
        error('The specified values does not meet the range restrictions');
    end;
else
    if (length(T) ~= length(v)) || (length(T) ~= length(theta))
        error('The specified vector is to large');
    end;
    if T(end) >= K
        error('Specified time(s) for change of wind exeeds simulation time');
    end;
end;


for i = 1:length(T)
    if (0<=v(i)) && (v(i)<1) && (0<= theta(i)) && (theta(i)< 2*pi)
        vx = v(i)*cos(theta(i));
       vy = v(i)*sin(theta(i));

       Wnew = W.*[(1+(-vx-vy)/sqrt(2)) (1-vy)  (1+(vx-vy)/sqrt(2));
              (1-vx) (1) (1+vx);
              (1+(-vx+vy)/sqrt(2)) (1+vy) (1+(vx+vy)/sqrt(2))];
          
        if T(i) == T(end)
          for j=T(i):K-1
            S = update(S,Wnew,N,M);
          end; 
          
        else
           for j=T(i):T(i+1)-1
               S = update(S,Wnew,N,M);
           end;
        end;
    else
        error('The specified values does not meet the range restrictions');
    end;
   
end;



end



function Snew = update(S,W,N,M)
% UPDATE updates an array S according to the spread of fire in the
% surrounding area.
% 
% Call:     Snew = update(S, W, N, M)
% 
% Inputs:   S, a NxN array needing to be updated.
%           W, a 3x3 array with weighted values.
%           N, the row and column length of S and M.
%           M, a NxN array with values that describe the flammability in
%              the area.
% Outputs:  Snew, a NxN array with the updated values of S

% Author:   Nicolas Borup - s083124
% Date: May 2, 2011

Snew = S;

for i=2:N-1
    for j=2:N-1
        I = [i-1 i i+1];
        J = [j-1 j j+1];
        if (S(i,j) < 1) && (S(i,j) > 0.01)
            Snew(I,J) = (Snew(I,J) + S(i,j).^2*(1-S(i,j))*(1-S(I,J)).*W.*M(I,J));
        end;
    end;

end;
end