GenerateMIMOSources

PURPOSE ^

[H,D,F] = GenerateMIMOSources(L,P,M,K,gamma)

SYNOPSIS ^

function [H,D,F] = GenerateMIMOSources(L,P,M,K,gamma)

DESCRIPTION ^

[H,D,F] = GenerateMIMOSources(L,P,M,K,gamma)

  GenerateMIMOSources(L,P,M,K) produces a source model for M-array data emanating
  from L independent Gaussian sources spectrally shaped by Pth order moving 
  average (MA, i.e. finite impulse response) innovation filters. The source 
  signals are then mixed by a random MxL paraunitary matrix of order K.

  GenerateMIMOSources(L,P,M,K,gamma) operates the same but the additional input
  parameter gamma limits the radii of zeros in the MA innovation filters. This 
  restricts the dynamic range of the generated spectra.

  A flow graph of this source model is given in [1]. In the context of multichannel
  coding, M=L is required to avoid a division by zero, as the geometric mean of 
  variances is zero for M>L. The purpose of limiting the dynamic range of source
  PSDs is driven by the same aim, and keeps the problem numerically tractable.

  The source model does currently not permit the overdetermined case where the 
  number of sources exceeds the number of sensors.

  [H,D,F]=GenerateMIMOSources(L,P,M,K,gamma) returns the paraunitary mixing matrix
  H(z) of dimension MxL and order K in H; D represents a diagonal parahermitian 
  matrix D(z) containing the power spectra of the L sources prior to convolutive
  mixing. Minimum phase filters of order P that generate such PSDs from unit 
  variance complex Gaussian white noise are contained in the columns of the matrix 
  F.

  Input parameters
     L       number of independent sources
     P       order of MA innovation filter
     M       dimension of array generated (M>=L)
     K       order of paraunitary mixing matrix
     gamma   max radius for zeros in the innovation filters
           (default value is 1)

  Output parameters
     H       paraunitary matrix
     D       diagonal, parahermitian matrix
     F       columns are the innovation filters, whose 
             ACFs generate the diagonal entries of D.

  Reference:

  [1] S. Redif, S. Weiss and J.G. McWhirter, "Sequential Matrix Diagonalisation 
      Algorithms for Polynomial EVD of Parahermitian Matrices," to appear in 
      IEEE Transactions on Signal Processing.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [H,D,F] = GenerateMIMOSources(L,P,M,K,gamma)
0002 %[H,D,F] = GenerateMIMOSources(L,P,M,K,gamma)
0003 %
0004 %  GenerateMIMOSources(L,P,M,K) produces a source model for M-array data emanating
0005 %  from L independent Gaussian sources spectrally shaped by Pth order moving
0006 %  average (MA, i.e. finite impulse response) innovation filters. The source
0007 %  signals are then mixed by a random MxL paraunitary matrix of order K.
0008 %
0009 %  GenerateMIMOSources(L,P,M,K,gamma) operates the same but the additional input
0010 %  parameter gamma limits the radii of zeros in the MA innovation filters. This
0011 %  restricts the dynamic range of the generated spectra.
0012 %
0013 %  A flow graph of this source model is given in [1]. In the context of multichannel
0014 %  coding, M=L is required to avoid a division by zero, as the geometric mean of
0015 %  variances is zero for M>L. The purpose of limiting the dynamic range of source
0016 %  PSDs is driven by the same aim, and keeps the problem numerically tractable.
0017 %
0018 %  The source model does currently not permit the overdetermined case where the
0019 %  number of sources exceeds the number of sensors.
0020 %
0021 %  [H,D,F]=GenerateMIMOSources(L,P,M,K,gamma) returns the paraunitary mixing matrix
0022 %  H(z) of dimension MxL and order K in H; D represents a diagonal parahermitian
0023 %  matrix D(z) containing the power spectra of the L sources prior to convolutive
0024 %  mixing. Minimum phase filters of order P that generate such PSDs from unit
0025 %  variance complex Gaussian white noise are contained in the columns of the matrix
0026 %  F.
0027 %
0028 %  Input parameters
0029 %     L       number of independent sources
0030 %     P       order of MA innovation filter
0031 %     M       dimension of array generated (M>=L)
0032 %     K       order of paraunitary mixing matrix
0033 %     gamma   max radius for zeros in the innovation filters
0034 %           (default value is 1)
0035 %
0036 %  Output parameters
0037 %     H       paraunitary matrix
0038 %     D       diagonal, parahermitian matrix
0039 %     F       columns are the innovation filters, whose
0040 %             ACFs generate the diagonal entries of D.
0041 %
0042 %  Reference:
0043 %
0044 %  [1] S. Redif, S. Weiss and J.G. McWhirter, "Sequential Matrix Diagonalisation
0045 %      Algorithms for Polynomial EVD of Parahermitian Matrices," to appear in
0046 %      IEEE Transactions on Signal Processing.
0047 
0048 % Stephan Weiss, August 2012
0049 
0050 %--------------------------------------------------------
0051 %   Input Parameter Check
0052 %--------------------------------------------------------
0053 if nargin==4,
0054    gamma=1;
0055 end;
0056 if M < L,
0057   error('more sources than array elements');
0058 end;
0059 
0060 %--------------------------------------------------------
0061 %   Paraunitary Matrix
0062 %--------------------------------------------------------
0063 [H,dummy] = qr(randn(M,M) + sqrt(-1)*randn(M,M));  
0064 H = eye(M);  
0065 for k = 1:K,
0066    v = randn(M,1) + sqrt(-1)*randn(M,1);
0067    v = v./norm(v,2);
0068    Hnew = zeros(M,M,2);
0069    [Q,dummy] = qr(randn(M,M) + sqrt(-1)*randn(M,M));  
0070    Hnew(:,:,1) = v*(v'*Q);
0071    Hnew(:,:,2) = Q - Hnew(:,:,1);
0072       
0073    H = PolyMatConv(H,Hnew);
0074 
0075 end;
0076 
0077 %--------------------------------------------------------
0078 %   Source Innovation Filters
0079 %--------------------------------------------------------
0080 Fz = (rand(P-1,L)*gamma).*exp(sqrt(-1)*2*pi*rand(P-1,L));  % zeros
0081 F = zeros(P,L);
0082 for l = 1:L,
0083    F(:,l) = poly( Fz(:,l)).';
0084 end;
0085 F(:,1) = F(:,1)./norm(F(:,1),2);
0086 Ffd = abs(fft(F,1024,1));
0087 for l = 1:L-1,
0088    alpha = min(Ffd(:,l)./Ffd(:,l+1));
0089    Ffd(:,l+1) = alpha*Ffd(:,l+1);
0090    F(:,l+1) = alpha*F(:,l+1);
0091 end;
0092 
0093 %--------------------------------------------------------
0094 %   Output Parameters
0095 %--------------------------------------------------------
0096 D1 = zeros(L,L,P);
0097 for l = 1:L,
0098   for p = 1:P,
0099     D1(l,l,p) = F(p,l);
0100   end;
0101 end;
0102 D = PolyMatConv(D1,ParaHerm(D1));
0103 H = H(:,1:L,:);
0104 
0105 
0106

Generated on Wed 03-Dec-2014 19:38:11 by m2html © 2005