function [Out,BlkSize] = mDIL(In,BlkSize,Ncpc)
% Deinterleaver of IEEE 802.16-2004 (WiMAX, OFDMA only)
%
%     [Out,BlkSize] = mDIL(In,BlkSize,Ncpc)
%
% In     : input column vector or matrix
%          one column per code block
% Out    : output column vector or matrix
%          one column per code block
% BlkSize: row vector, length = number of code blocks
%          number of data elements in each code block (or column)
%
% Ncpc   : number of coded bits per subcarrier (modulation mode)
%          1 = BPSK
%          2 = QPSK
%          4 = 16QAM
%          6 = 64QAM
%
% The output block size is the same as the input block size.
%
% Matlab 7 Release 14 SP2

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%  Property of Freescale
%  Freescale Confidential Proprietary
%  Freescale Copyright (C) 2005 All rights reserved
%  ----------------------------------------------------------------------------
%  $RCSfile: mDIL.m.rca $
%  $Revision: 1.6 $
%  $Date: Fri Oct 20 15:04:34 2006 $
%  Target: Matlab
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

[m,n]=size(In);

%check dimensions
[x,y]=size(BlkSize);
if ((x~=1)|(y~=n))
    error('Error: BlkSize must be row vector with the same number of columns as In.');
end

%-- Check dimension Ncpc value
if ~any([1 2 4 6]==Ncpc)
    error(['Error: Ncpc type "' int2str(Ncpc) '" not supported'])
end


%allocate output matrix
if islogical(In)
    Out = logical(zeros(m,n));
else
    Out = zeros(m,n,class(In));
end

%loop over all blocks
for b=1:n,

    %Stage1: calculation of Ncbps and S
    Ncbps  = BlkSize(b);
    s      = ceil(Ncpc/2);

    %check size
    if (Ncbps/16/s)~=floor(Ncbps/16/s)
        error('Error: length(input)/ceil(Ncbc/2)/16 must be an integer!');
    end

    %Stage2: First Permutation ... Alternately maps adjacent bits to more or less
    %        significant constellation bits.
    j      = 0:Ncbps-1;
    mj     = s*floor(j/s) + mod(j+floor(16*j/Ncbps),s);

    %Stage3: Second Permutation ... Maps adjacent bits to different modulated words
    kj     = 16*mj - (Ncbps-1)*floor(16*mj/Ncbps);

    %Stage4: Deinterleave data
    Out(kj+1,b) = In(1:Ncbps,b);

end