]> Creatis software - CreaPhase.git/blobdiff - octave_packages/communications-1.1.1/qaskenco.m
Add a useful package (from Source forge) for octave
[CreaPhase.git] / octave_packages / communications-1.1.1 / qaskenco.m
diff --git a/octave_packages/communications-1.1.1/qaskenco.m b/octave_packages/communications-1.1.1/qaskenco.m
new file mode 100644 (file)
index 0000000..760d453
--- /dev/null
@@ -0,0 +1,180 @@
+## Copyright (C) 2003 David Bateman
+##
+## This program is free software; you can redistribute it and/or modify it under
+## the terms of the GNU General Public License as published by the Free Software
+## Foundation; either version 3 of the License, or (at your option) any later
+## version.
+##
+## This program is distributed in the hope that it will be useful, but WITHOUT
+## ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+## FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+## details.
+##
+## You should have received a copy of the GNU General Public License along with
+## this program; if not, see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Function File}  {} qaskenco (@var{m})
+## @deftypefnx {Function File} {} qaskenco (@var{msg},@var{m})
+## @deftypefnx {Function File} {@var{y} = } qaskenco (@var{...})
+## @deftypefnx {Function File} {[@var{inphase}, @var{quadr}] =} qaskenco (@var{...})
+##
+## Map a digital signal using a square QASK constellation. The argument
+## @var{m} must be a positive integer power of 2. With two input arguments
+## the variable @var{msg} represents the message to be encoded. The values
+## of @var{msg} must be between 0 and @code{@var{m}-1}. In all cases
+## @code{qaskenco(@var{M})} is equivalent to @code{qaskenco(1:@var{m},@var{m})}
+##
+## Three types of outputs can be created depending on the number of output 
+## arguments. That is
+##
+## @table @asis
+## @item No output arguments
+## In this case @dfn{qaskenco} plots the constellation. Only the
+## points in @var{msg} are plotted, which in the case of a single input
+## argument is all constellation points.
+## @item A single output argument
+## The returned variable is a complex variable representing the in-phase 
+## and quadrature components of the mapped  message @var{msg}. With, a 
+## single input argument this effectively gives the mapping from symbols
+## to constellation points
+## @item Two output arguments
+## This is the same as two ouput arguments, expect that the in-phase
+## and quadrature components are returned explicitly. That is
+##
+## @example
+## octave> c = qaskenco(msg, m);
+## octave> [a, b] = qaskenco(msg, m);
+## octave> all(c == a + 1i*b)
+## ans = 1
+## @end example
+## @end table
+##
+## If @code{sqrt(@var{m})} is an integer, then @dfn{qaskenco} uses a Gray
+## mapping. Otherwise, an attempt is made to create a nearly square mapping 
+## with a minimum Hamming distance between adjacent constellation points.
+## @end deftypefn
+## @seealso{qaskdeco}
+
+## 2005-04-23 Dmitri A. Sergatskov <dasergatskov@gmail.com>
+##     * modified for new gnuplot interface (octave > 2.9.0)
+
+function [a, b] = qaskenco(msg, M)
+
+  if (nargin == 1)
+    M = msg;
+  elseif (nargin == 2)
+    if ((min(msg(:)) < 0) || (max(msg(:)) > M-1))
+      error ("qaskenco: message invalid");
+    endif
+  else
+    error ("qaskenco: incorrect number of arguments");
+  endif
+
+  if (!isscalar(M) || (M != ceil(M)) || (M < 2))
+    error ("qaskenco: order of modulation must be a positive integer greater than 2");
+  endif
+
+  if (log2(M) != ceil(log2(M)))
+    error ("qaskenco: the order must be a power of two");
+  endif
+  
+  if (M == 2)
+    inphase = [-1,  1];
+    quadr =   [ 0,  0];
+  elseif (M == 4)
+    inphase = [-1, -1,  1,  1];
+    quadr =   [-1,  1, -1,  1];
+  elseif (M == 8)
+    inphase = [-1, -1,  1,  1, -3, -3,  3,  3];
+    quadr =   [-1,  1, -1,  1, -1,  1, -1,  1];
+  else
+    NC =2^floor(log2(sqrt(M)));
+    MM = NC * NC;
+    Gray = [0, 1];
+    for i=2:ceil(log2(NC))
+      Gray = [Gray 2^(i-1) + fliplr(Gray)];
+    end
+    Gray = fliplr(de2bi(shift(Gray,length(Gray)/2 - 1)));
+    Gray2 = zeros(MM,log2(MM));
+    Gray2(:,1:2:log2(MM)) = repmat(Gray,NC,1);
+    for i=1:NC
+      Gray2(i:NC:MM,2:2:log2(MM)) = Gray;
+    end
+    layout = reshape(bi2de(fliplr(Gray2)),NC,NC);
+
+    if (MM != M)
+      ## Not sure this is the best that can be done for these mappings. If 
+      ## anyone wants to improve this, go ahead, but do it in qaskdeco too.
+      OFF = sqrt(M/32);
+      NR = NC + 2*OFF;
+      layout2 = NaN * ones(NR);
+      layout2(1+OFF:OFF+NC,1+OFF:OFF+NC) = layout;
+      
+      layout2(1:OFF,1+OFF:OFF+NC) = MM + layout(OFF:-1:1,:);
+      layout2(NR-OFF+1:NR,1+OFF:OFF+NC) = MM + layout(NC:-1:NC-OFF+1,:);
+
+      layout2(1+2*OFF:NC,1:OFF) = MM + layout(OFF+1:NC-OFF,OFF:-1:1);
+      layout2(1+2*OFF:NC,NR-OFF+1:NR) = MM + ...
+         layout(OFF+1:NC-OFF,NC:-1:NC-OFF+1);
+
+      layout2(1+OFF:2*OFF,1:OFF) = MM + ...
+         layout(NC/2:-1:NC/2-OFF+1,NC/2:-1:OFF+1);
+      layout2(NC+1:OFF+NC,1:OFF) = MM + ...
+         layout(NC-OFF:-1:NC/2+1,NC/2:-1:OFF+1);
+
+      layout2(1+OFF:2*OFF,NR-OFF+1:NR) = MM + ...
+         layout(NC/2:-1:NC/2-OFF+1,NC-OFF:-1:NC/2+1);
+      layout2(NC+1:OFF+NC,NR-OFF+1:NR) = MM + ...
+         layout(NC-OFF:-1:NC/2+1,NC-OFF:-1:NC/2+1);
+      NC = NR;
+      layout = layout2;
+    endif
+
+    inphase = repmat([0:NC-1]*2 - NC + 1,1,NC);
+    for i=1:NC
+      quadr(i:NC:NC*NC) = [0:NC-1]*2 - NC + 1;
+    end
+    [dummy, indx] = sort(layout(:));
+    indx = indx(1:M);                  ## Get rid of remaining NaN's
+    inphase = inphase(indx);
+    quadr = quadr(indx);
+  endif
+
+  if (nargin == 2)
+    inphase = inphase(msg+1);
+    quadr = quadr(msg+1);
+    ## Fix up indexing if using column vector
+    if (size(msg,2) == 1)
+      inphase = inphase';
+      quadr = quadr';
+    endif
+  endif
+
+  if (nargout == 0)
+    inphase = inphase(:);
+    quadr = quadr(:);
+    plot (inphase, quadr, "r+");
+    title("QASK Constellation");
+    xlabel("In-phase");
+    ylabel("Quadrature");
+    axis([min(inphase)-1, max(inphase)+1, min(quadr)-1, max(quadr)+1]);
+    xd = 0.02 * max(inphase);
+    if (nargin == 2)
+      msg = msg(:);
+      for i=1:length(inphase)
+       text(inphase(i)+xd,quadr(i),num2str(msg(i)));
+      end
+    else
+      for i=1:length(inphase)
+       text(inphase(i)+xd,quadr(i),num2str(i-1));
+      end
+    endif
+  elseif (nargout == 1)
+    a = inphase + 1i * quadr;
+  else
+    a = inphase;
+    b = quadr;
+  endif
+
+endfunction