ddemod.m

数字通信第四版原书的例程

function x = ddemod(y, Fc, Fd, Fs, method, M, opt1, opt2, opt3, opt4)
%DDEMOD Demodulates the complex envelope of a digital modulated signal.
%       DDEMOD(Y, Fc, Fd, Fs) plots eye-pattern diagram for the signal Y.
%       The sampling frequency for Y is Fs (Hz) and the digital sampling
%       frequency is Fd (Hz). When Fd is a two element vector, the second
%       element is the offset value for the decision point. The offset
%       should be an integer. The offset timing in the plot is offset/Fs.
%       The parameter Fc is not used.
%
%       Z = DDEMOD(Y, Fc, Fd, Fs, METHOD, OPT1, OPT2, OPT3, OPT4, OPT5)
%       demodulates a digital modulated signal Y with carrier frequency 
%       Fc (Hz), digital signal frequency Fd (Hz), and computation sampling
%       frequency Fs (Hz). This function requires Fs/Fd to be a positive
%       integer. When Fd is a two element vector, the second element is the
%       offset value for the decision point. The offset timing in the plot is
%       offset/Fs. The default offset is 0. The received digital signal is 
%       sampled at time point i/Fd + offset/Fs. The sample value keeps until
%       time (i+1)/Fd+offset/Fs when the next sample point is taken place. 
%       The demapping process finds the distance from the sample value
%       to all of the possible digital symbols. The digital symbols with the
%       shortest distance to the current sampling point are the demodulated
%       output. When Fs is a two element vector, the second element specifies
%       the initial phase of the carrier signal in modulation.
%
%       METHOD is a string, which can be one of the following:
%       'ask'       M-ary Amplitude shift keying modulation.
%       'psk'       M-ary Phase shift keying modulation.
%       'qask'      M-ary Quadrature amplitude shift-keying modulation, with 
%                   choices of 'qask/sqr' (default), 'qask/cir', or 'qask/arb'
%                   for different mapping constellations.
%       'fsk'       M-ary Frequency shift keying modulation.
%       'msk'       Minimum shift keying modulation.
%       'sample'    Down sampling from a higher sampling frequency signal.
%
%       When there is a '/eye' or 'scat' string appended to the end of the
%       variable METHOD, the function plots eye-pattern diagram or the scatter 
%       plot.
%
%       The demodulation uses a lowpass filter. When it is not specified, a
%       default lowpass filter [NUM, DEN] = butter(5, Fc*2/Fs) is used.
%
%       Use DDEMOD(METHOD) to view the help for a specific method.
%
%       See also MODMAP, AMOD, ADEMOD, DMOD, DDEMODCE.
%

%       Wes Wang 1/9/95, 10/2/95
%       Copyright (c) 1995-96 by The MathWorks, Inc.
%       $Revision: 1.1 $  $Date: 1996/04/01 17:56:09 $

%position for optional parameters.
opt_pos = 7;
plot_const = 0;
if nargin < 1
    feval('help','ddemod');
    return;
elseif isstr(y)
    method = lower(deblank(y));
    if length(method) < 3
        error('Invalid method option for ddemod.')
    end
    if nargin == 1
        addition = 'See also DMOD, DMODCE, DDEMODCE, AMOD, ADEMOD.';
        if method(1:3) == 'qas'
          callhelp('ddemod.hlp',method(1:4),addition);
        else
          callhelp('ddemod.hlp',method(1:3),addition);
        end
    else
        disp('Warning: Worng number of input variable, use MODMAP for plotting constellations.');
        return
    end;
else

if length(Fs) > 1
    ini_phase = Fs(2);
    Fs = Fs(1);
else
    ini_phase = 0;
end;

if length(Fd) > 1
    offset = Fd(2);
    Fd = Fd(1);
else
    offset = 0;
end;

if (Fs == 0) | (Fd == 0) | isempty(Fd) | isempty(Fs)
    FsDFd = 1;
    offset = 0;
else
    FsDFd = Fs / Fd;
    if (ceil(FsDFd) ~= FsDFd) | (FsDFd <= 0)
        error('Fs / Fd must be a positive integer.')
    end;
end;

offset = rem(offset/FsDFd, 1)*FsDFd;
if offset < 0
    offset = rem(offset/FsDFd+1, 1)*FsDFd;
end;
%offset = rem(offset, FsDFd);
%if offset < 0
%    offset = rem(offset + FsDFd, FsDFd);
%end;

[r, c] = size(y);
if r * c == 0
    x = [];
    return;
end;
if r == 1
    y = y(:);
    len_y = c;
else
    len_y = r;
end;

if nargin < 3
    disp('Usage: Z=DDEMOD(Y, Fc, Fd, Fs, METHOD, OPT1, OPT2, OPT3) for modulation mapping');
    return;
elseif nargin < opt_pos  - 2
    if nargout < 1
        method = 'eye';
    else
        method = 'sample';
    end
end;

method = lower(method);

if strcmp(method(1:2), 'as')
    % M has to be specified in 'ask' option.
    if nargin < opt_pos - 1
        error('Not enough input varibale for DDEMOD.')
    end;
    if offset == 0
        offset = FsDFd;
    end;
    index = ([0 : M - 1] - (M - 1) / 2) * 2 / (M - 1);
    if nargin < opt_pos + 1
        [num, den] = butter(5, Fc * 2 / Fs);
    else
        num = opt1;
        den = opt2;
    end;
    if findstr(method, '/cos')
        y = ademod(y, Fc, [Fs, ini_phase], 'amdsb-sc/costas', num, den);
    else
        y = ademod(y, Fc, [Fs, ini_phase], 'amdsb-sc', num, den);
    end;
    if findstr(method, '/eye')
        ddemod(y, Fc, [Fd, offset], [Fs, ini_phase], 'eye')
    end;
    if findstr(method, '/sca')
        ddemod(y, Fc, [Fd, offset], [Fs, ini_phase], 'sca')
    end;
    yy = y([offset : FsDFd : len_y], :);
    x = [];
    for i = 1 : size(y, 2)
        [tmp, x_p] = min(abs(yy(:, i*ones(1, M)) - index(ones(1, size(yy, 1)), :))');
        x = [x x_p'-1];
    end;
elseif findstr(method, 'fsk')
    if nargin < opt_pos - 1
        Tone = 2 * Fd / M;
    else
        Tone = opt1;
    end;

    %calculate the correlation of fsk.
    z = [0: M-1] * Tone * 2 * pi / Fs;
    z = [0; ones(FsDFd-1, 1)] * z;
    t = [0 : 1/Fs : 1/Fd-1/Fs]';
    t = t(:, ones(1, M)) * 2 * pi * Fc;
    z = cos(t + cumsum(z) + ini_phase);
    zz = cos(t + cumsum(z) + ini_phase * 2);
    %leave space for x
    x = y([offset+1 : FsDFd : len_y], :);
    if findstr(method, '/eye')
        t = [0 : FsDFd-1]/Fs;
        t = t + offset/Fs;
        plot([min(t), max(t), max(t)], [-1/2, NaN, 1])
        axis([min(t) max(t), -1/2, 1]);
        hold on
    end;
    [len_x, wid_x] = size(x);
    for i = 1 : wid_x
        comp_low = 1;
        if offset <= 0
            comp_upp = FsDFd;
        else
            comp_upp = offset;
        end;
        for k = 1 : len_x
            corr = y(comp_low:comp_upp)*ones(1,M) .* z(1:comp_upp-comp_low+1,:);
            if findstr(method, '/eye')
                if findstr(method, '/nonc')
                    cor2 = y(comp_low:comp_upp)*ones(1,M) .* zz(1:comp_upp-comp_low + 1, :);
                    corr = cumsum(corr).^2 + cumsum(cor2).^2;
                else
                    corr = cumsum(corr);
                end;
                corr = corr / max(max(max(abs(corr))), eps);
                if k == 1
                    plot(t(FsDFd-size(corr, 1)+1 : FsDFd)', corr);
                else
                    plot(t(1:size(corr, 1))', corr)
                end
                corr = corr(size(corr, 1), :);
            else
                corr = sum(corr);
                if findstr(method, '/nonc')
                    cor2 = y(comp_low:comp_upp)*ones(1,M) .* zz(1:comp_upp-comp_low + 1, :);
                    corr = corr.^2 + sum(cor2).^2;
                end;
            end;
            [tmp_v, tmp_i] = max(corr);
            x(k, i) = tmp_i - 1;
            comp_low = min(comp_low + FsDFd, len_y);
            comp_upp = min(comp_upp + FsDFd, len_y);
        end
    end
    if findstr(method, '/eye')
        hold off
    end;
%psk has combined with qask.
elseif findstr(method, 'msk')
    %This is a special case of fsk call back to get fsk
    method(1) = 'f';
    M = 2;
    Tone = Fd;
    x = ddemod(y, Fc, [Fd, offset], [Fs, ini_phase], method, M, Tone);
elseif ~isempty(findstr(method, 'qask')) |...
         ~isempty(findstr(method, 'qam')) |...
         ~isempty(findstr(method, 'qsk')) |...
         ~isempty(findstr(method, 'psk'))
    if findstr(method, '/ar')
        % arbitraryly defined I, Q. 
        if nargin < opt_pos
            error('In correct format for METHOD=''qask/arbitrary''.');
        end;
        I = M;
        Q = opt1;
        if nargin < opt_pos + 2
            [num, den] = butter(5, Fc*2/Fs);
        else
            num = opt2;
            den = opt3;
        end;
        M = length(I);
    elseif ~isempty(findstr(method, '/ci')) | ~isempty(findstr(method, 'psk'))
        % circle defined NIC, AIC, PIC.
        if nargin < opt_pos - 1
            if findstr(method, 'psk')
                error('M-Ary number must be specified for psk demap.')
            else
                error('Incorrect format for METHOD=''qask/arbitrary''.');
            end;
        end;
        NIC = M;
        M = length(NIC);
        if nargin < opt_pos
            AIC = [1 : M];
        else
            AIC = opt1;
        end;
        if nargin < opt_pos + 1
            PIC = NIC * 0;
        else
            PIC = opt2;
        end;
        if nargin < opt_pos + 3
            [num, den] = butter(5, Fc*2/Fs);
        else
            num = opt3;
            den = opt4;
        end;
        inx = apkconst(NIC, AIC, PIC);
        I = real(inx);
        Q = imag(inx);
        M = sum(NIC);
    else
        %consider as square style.
        [I, Q] = qaskenco(M);
        if nargin < opt_pos + 1
            [num, den] = butter(5, Fc*2/Fs);
        else
            num = opt1;
            den = opt2;
        end;
    end;
    y = ademod(y, Fc, [Fs, ini_phase], 'qam', num, den);
    if findstr(method, '/eye')
        ddemod(y, Fc, [Fd, offset], [Fs, ini_phase], 'eye')
    end;
    if findstr(method, '/sca')
        ddemod(y, Fc, [Fd, offset], [Fs, ini_phase], 'sca')
    end;
    if offset <= 0
        offset = FsDFd;
    end;
    yy = y([offset : FsDFd : len_y], :);
    [len_y, wid_y] = size(yy);
    if (ceil(wid_y/2) ~= wid_y/2)
        error('qask demap requires input is a matrix with even number of columns.');
    end;
    x = []; I = I(:)'; Q = Q(:)';
    for i = 1 : 2 : wid_y
        [tmp, x_p] = min(...
            ((yy(:, i*ones(1, M)) - I(ones(1, len_y), :)).^2 + ...
            (yy(:, (i+1)*ones(1, M)) - Q(ones(1, len_y), :)).^2)');
        x = [x x_p'-1];
    end;
elseif findstr(method, 'samp')
    %This is made possible to convert an input signal from sampling frequency Fd
    %to sampling frequency Fs.
    if offset <= 0
        offset = FsDFd;
    end;
    x = y([offset : FsDFd : len_y], :);
elseif findstr(method, 'eye')
    % plot eye-pattern plot
    eyescat(y, Fd, Fs, offset, offset);
elseif findstr(method, 'scat')
    % plot scatter plot by using eyescat
    if nargin >= opt_pos
        if isstr(M)
            M = M(1);
        else
            M = '.';
        end;
    else
        M = '.';
    end;
    eyescat(y, Fd, Fs, offset, M);
else
    %The choice is not a valid one.
    disp('You have used an invalid method. The method should be one of the following string:')
    disp('  ''ask'' Amplitude shift keying modulation;')
    disp('  ''psk'' Phase shift keying modulation;')
    disp('  ''qask'' Quadrature amplitude shift-keying modulation, square constellation;')
    disp('  ''qask/cir'' Quadrature amplitude shift-keying modulation, circle constellation;')
    disp('  ''qask/arb'' Quadrature amplitude shift-keying modulation, user defined constellation;')
    disp('  ''fsk'' Frequency shift keying modulation;')
    disp('  ''msk'' Minimum shift keying modulation;')
    disp('  ''sample'' Convert sample frequency Fd input to sample frequency Fs output.')
end;

end;