%                                gmp1.m
%  Scope:   This MATLAB macro generates first order Gauss-Markov sequence.
%  Usage:   x = gmp1(nstep,beta,sigx,deltat,iseed)
%           x = gmp1(nstep,beta,sigx,deltat)      when the seed is not reset
%  Description of parameters:
%           nstep  - input, number of steps of the Gauss-Markov process
%           beta   - input, parameter of the Gauss-Markov process, 1/beta is
%                    the correlation time in seconds
%           sigx   - input, standard deviation of the Gauss-Markov process  x
%           deltat - input, time step in seconds
%           iseed  - input, initial value of the seed (optional); if it is
%                    specified the seed is reset to iseed value
%           x      - output, computed value of the Gauss-Markov process x,
%                    vector with nstep elements
%  External Matlab macros used:  genrn
%  Method: 
%           Direct implementation of the formula
%              x_new = exp(-beta*delta) * x_old + w
%           where
%              w  is a normally distributed random number with zero mean 
%                 and standard deviation equals to 
%                 sigma_x * sqrt(1 - exp(-2*beta*delta))
%  Last update:  01/07/01
%  Copyright (C) 1996-98 by LL Consulting. All Rights Reserved.

function  x = gmp1(nstep,beta,sigx,deltat,iseed)

if ( (nargin < 4) | (nargin > 5) )
   disp('Error - GMP1.m  - check the argument list');
   disp('  ');
   return
elseif  (nargin == 5)
   rand('seed',iseed);  
end

clear x
x = zeros(nstep,1);

for k = 1:nstep
   if k == 1
      temp =  1.0 - exp(-2.0 * beta * deltat);
      stdevw = sigx * sqrt(temp);
      expbdt = exp(-beta * deltat);
      xk = 0.;
   end
   x(k) = expbdt * xk + genrn(1,0.,stdevw);
   xk = x(k);
end