filling up efficiently a sparse matrix based on 4-pixel or 8-pixel neighborhoods

Question

Given an image of size [hh,ww], I would like to create efficiently a sparse matrix of size [hh*ww, hh*ww]. For each 4- or 8-neighbor of a given pixel, the sparse matrix should be filled with a constant value (say -1) at the proper row (which represents the pixel under evaluation) and columns (the corresponding 4-neighbors or 8-neighbors of the pixel).

For instance, considering a 4-pixel neighborhood and a [2,2] matrix, the resulting sparse matrix of size [4,4] looks like:

 0  -1  -1   0
-1   0   0  -1
-1   0   0  -1
 0  -1  -1   0

The first row was filled with -1 at columns #2 and #3 since those pixels are in the 4-neighborhood of pixel 1 (linear indexing given below ):

1   3
2   4

The code below does the work, but becomes very slow whenever the matrices become too large, for instance, for a 2000x2000 matrix.

hh=2;
ww=2;
%hh=2000;
%ww=2000;
sG     = sparse(hh*ww,hh*ww);
linIdx = reshape(1:hh*ww, [hh ww]);
sG( sub2ind([hh*ww hh*ww], linIdx(:,1:end-1),linIdx(:,2:end)) ) = -1;
sG( sub2ind([hh*ww hh*ww], linIdx(1:end-1,:),linIdx(2:end,:)) ) = -1;
sG = max(sG, sG'); 

Any ideas how to make the code efficient when having large matrices? Ideally, it should work for 4-neighborhoods or 8-neighborhoods.

Answer 1

I wrote a function that computes efficiently the sparse adjacency matrix, as you described. See sparse_adj_matrix.

Usage:

[ii jj] = sparse_adj_matrix( [hh ww], 1, 1 ); % p is L1 for 4-connect 
sG = sparse( ii, jj, -1, hh*ww, hh*ww ); % construct the sparse matrix

For 8-connect

[ii jj] = sparse_adj_matrix( [hh ww], 1, inf ); % p is inf for 8-connect 
sG = sparse( ii, jj, -1, hh*ww, hh*ww ); % construct the sparse matrix

This function can handle arbitrary dimension (more than 2) regular grids with neighborhoods different than 4 or 8 (radius larger than 1, metric L1, L2 or Loo).