function X = normrnd_1D(rho,M,N);
%NORMRND_1D Random correlated sequences from the normal distribution
%   X = NORMRND_1D(RHO,N) returns a row or N random numbers from the normal
%   distrubution with a correlation of RHO between adjacent observations.
%
%   X = NORMRND_1D(RHO,M,N) returns M sequences of correlated random
%   variables of length N.  Each sequence is one row of X.
%
%   For accuracy, N should be large enough that RHO^N<.05.  If shorter N
%   are desired, please make N large enough and extract the first portion
%   after the fuction is executed.

if abs(rho)>1, error('RHO must be between -1 and 1'); end

if nargin == 2
    T = [0:M M-1:-1:1]';
    R = rho.^T;
    Z = randn(2*M,1);
 %   X = real(ifft((fft(Z).*sqrt(fft(R))));
    X = X(1:M,:)';
else
    T = [0:N N-1:-1:1]';
    R = rho.^T;
    Z = randn(2*N,M);
%    Z = randn(2*N,M)+sqrt(-1)*randn(2*N,M);
%    FR = fft(R);
%    FR = sqrt((2*N)^-1*FR);
%    Z = Z.*repmat(sqrt(FR/2/N),1,M);
%    X = fft(Z);
    X = real(ifft(fft(Z).*repmat(sqrt(fft(R(:))),1,M)));
    X = X(1:N,:)';
end