function [x, t] = uen(T, Fs, Lp, cal)
    // Generation of uniform excitation noise
    // Uniform excitation noise is a noise which has the
    // same intensity in each critical band.
    //
    // Usage:
    //       [x, t] = uen(T, Fs, Lp, cal)
    // where
    //       x:   vector containing uniform exciting noise 
    //       t:   time vector in s
    //       T:   duration in s
    //       Fs:  sampling rate in Hz
    //       Lp:  sound pressure level
    //       cal: sound pressure level corresponding to |x| = 1
    //
    // Output t is optional. Inputs fs, Lp and cal are 
    // optional. Their default values are 44100 Hz, 74 dB
    // and 94 dB 
    //
    // Critical band intensity of white noise is proportional 
    // to critical bandwidth. Hence, uniform excitation noise
    // is obtained by filtering white noise with a frequency 
    // response that is proportional to the square root of the
    // inverse of critical bandwidth
    // 
    // Reference
    // Zwicker & Fastl, "Psychoacoustics, Facts and Models"
    //
    // --------------------------
    // Author: Federico Miyara
    // Date: 2008-06-30
    //       2020-03-07 

    exec fftfilter.sci;

    // Default values if not enough arguments
    if nargin < 4
        cal = 94;
    end
    if nargin < 3
        Lp = 74;
    end
    if nargin < 2
        Fs = 44100;
    end

    // Window size 
    N  = 4096;

    // Frequency vector
    f = Fs * [0:N/2-1]/N;

    // Crítical bandwidth according to Zwicker-Fastl, page 164 
    cbw = 25 + 75*(1 + 1.4*(f/1000).^2).^0.69;

    // Size of signal vector
    siz = floor(T*Fs) + 1;

    // White noise
    wn = grand(1, siz, "nor", 0, 1);

    // Uniform exciting noise
    x = real(fftfilter(wn, sqrt(100 ./ cbw), N));
    x = x(1:siz);

    // Calibration
    x = x / sqrt(mean(x.^2)) * 10^((Lp - cal)/20);

    if nargout == 2
        // Time vector
        t = [0:siz-1]/Fs;
    end

endfunction
