function dist = soft_identity(n, p)
% SOFT_IDENTITY Make a noisy identity matrix
% 
% DIST = SOFT_IDENTITY(N, P) creates an N*N matrix with P on the diagonal,
% and 1-P spread over the remaining columns. When N=2, this is a binary symmetric noisy channel.

prest = (1-p)/(n-1);
dist = prest*ones(n,n);
dist = setdiag(dist, p);
