Matlab Function for Joint Diagonalize for BSS

JointDiagonalize.m

function G = JointDiagonalize( S )

  [ ~, N ] = size( S );
  halfN = N / 2;

  % Make two Correlation matrices
  S1 = S(:,1:halfN); 
  S2 = S(:,halfN + 1:N);

  Ex1 = (S1 * S1')/halfN;
  [ P1, D1 ] = eig( Ex1 );
  rD1_inv =  D1^(-1/2);
  Gd = rD1_inv * P1';

  Ex2 = (S2 * S2')/halfN;
  Ex2d = Gd * Ex2 * Gd';
  [ P2, ~ ] = eig( Ex2d );

  G = P2' * Gd;

end