PolyMatAlign

[Shift,error] = PolyMatAlign(A,B);

   PolyMatAlign(A,B) takes input arguments A and B to represent polynomial 
   matrix A(z) and B(z) of the same dimension but potentially different order.
   It determines the best delay such that A(z) - B(z)z^{-Shift} is minised.
   
   [Shift,error] = PolyMatAlign(A,B) returns the determined shift and norm of
   of the mismatch between the aligned matrices A(z) and B(z)z^{-Shift}.

   [Shift,error,A2,B2] = PolyMatAlign(A,B) additionally returns the aligned 
   matrices of identical orders in A2 and B2.

   Input parameters:
      A           MxNxL1 matrix
      B           MxNxL2 matrix
   
   Output parameters:
      Shift       determined delay between A and B
      error       mismatch between the aligned matrices
      A2          MxBxL3 matrix A after alignment
      B2          MxBxL3 matrix B after alignment

0001 function [Shift,error,A2,B2] = PolyMatAlign(A,B);
0002 %[Shift,error] = PolyMatAlign(A,B);
0003 %
0004 %   PolyMatAlign(A,B) takes input arguments A and B to represent polynomial
0005 %   matrix A(z) and B(z) of the same dimension but potentially different order.
0006 %   It determines the best delay such that A(z) - B(z)z^{-Shift} is minised.
0007 %
0008 %   [Shift,error] = PolyMatAlign(A,B) returns the determined shift and norm of
0009 %   of the mismatch between the aligned matrices A(z) and B(z)z^{-Shift}.
0010 %
0011 %   [Shift,error,A2,B2] = PolyMatAlign(A,B) additionally returns the aligned
0012 %   matrices of identical orders in A2 and B2.
0013 %
0014 %   Input parameters:
0015 %      A           MxNxL1 matrix
0016 %      B           MxNxL2 matrix
0017 %
0018 %   Output parameters:
0019 %      Shift       determined delay between A and B
0020 %      error       mismatch between the aligned matrices
0021 %      A2          MxBxL3 matrix A after alignment
0022 %      B2          MxBxL3 matrix B after alignment
0023 
0024 % S. Weiss, 19/2/2023
0025 
0026 %----- check arguments
0027 [M1,N1,L1] = size(A); [M,N,L2] = size(B);
0028 if (M~=M1)||(N~=N1), error('dimension mismatch between input arguments'); end;
0029 
0030 %------------------------------------------------------------------------------
0031 %   determine shift
0032 %------------------------------------------------------------------------------
0033 %----- elementwise auto-correlation, add moduli
0034 r_ab = zeros(L1+L2-1,1);
0035 for m = 1:M, 
0036    for n = 1:N,
0037        r_ab = r_ab + abs(conv(squeeze(A(m,n,:)),flipud(conj(squeeze(B(m,n,:))))));
0038    end;    
0039 end;
0040 
0041 %----- determine delay
0042 TimeScale=(-L2+1:L1-1);
0043 [~,TIndex] = max(r_ab);
0044 Shift = -TimeScale(TIndex);
0045 
0046 %------------------------------------------------------------------------------
0047 %   determine mismatch
0048 %------------------------------------------------------------------------------
0049 % compensate for shift
0050 if Shift>0,    % a positive shift means that A needs zero-padding at the front
0051    dummy = zeros(M,N,L1+Shift);
0052    dummy(:,:,Shift+1:end) = A;
0053    A = dummy;
0054 else           % a negative shift means that B needs zeropadding at the front
0055    dummy = zeros(M,N,L2-Shift);
0056    dummy(:,:,1-Shift:end) = B;
0057    B = dummy;
0058 end;    
0059 % compensate for difference in orders
0060 L1 = size(A,3); L2 = size(B,3);
0061 if L1>L2,
0062    B2 = zeros(M,N,L1);
0063    B2(:,:,1:L2) = B;
0064    A2 = A;
0065 else
0066    A2 = zeros(M,N,L2);
0067    A2(:,:,1:L1) = A;
0068    B2 = B;
0069 end;
0070 % mismatch as sum of squared Frobenius norms
0071 error = PolyMatNorm(A2-B2);