function missmatch=mm(P,A,B)

% function missmatch=mm(P,A,B) : if missmatch is not <<1 then P,A,B are not consistent
%
% Computes the singular values of
% [P-A*P*A' B; B' zeros(m,m)] beyond the first 2*m, normalized by the
%                                                   smallest of the first
%                                                   2*m.
% These ought to be small... otherwise P is not consistent with the
% assumption of being a state covariance of the system (A,B)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

[n,m]=size(B);
svdM=svd([P-A*P*A' B; B' zeros(m,m)]);
missmatch=svdM(2*m+1:length(svdM));
missmatch=missmatch/min(svdM(1:2*m));