wolffd@0
|
1 function G = mk_2D_lattice(nrows, ncols, wrap_around)
|
wolffd@0
|
2 % MK_2D_LATTICE Return adjacency matrix for 4-nearest neighbor connected 2D lattice
|
wolffd@0
|
3 % G = mk_2D_lattice(nrows, ncols, wrap_around)
|
wolffd@0
|
4 % G(k1, k2) = 1 iff k1=(i1,j1) is connected to k2=(i2,j2)
|
wolffd@0
|
5 %
|
wolffd@0
|
6 % If wrap_around = 1, we use toroidal boundary conditions (default = 0)
|
wolffd@0
|
7 %
|
wolffd@0
|
8 % Nodes are assumed numbered as in the following 3x3 lattice
|
wolffd@0
|
9 % 1 4 7
|
wolffd@0
|
10 % 2 5 8
|
wolffd@0
|
11 % 3 6 9
|
wolffd@0
|
12 %
|
wolffd@0
|
13 % e.g., G = mk_2D_lattice(3, 3, 0) returns
|
wolffd@0
|
14 % 0 1 0 1 0 0 0 0 0
|
wolffd@0
|
15 % 1 0 1 0 1 0 0 0 0
|
wolffd@0
|
16 % 0 1 0 0 0 1 0 0 0
|
wolffd@0
|
17 % 1 0 0 0 1 0 1 0 0
|
wolffd@0
|
18 % 0 1 0 1 0 1 0 1 0
|
wolffd@0
|
19 % 0 0 1 0 1 0 0 0 1
|
wolffd@0
|
20 % 0 0 0 1 0 0 0 1 0
|
wolffd@0
|
21 % 0 0 0 0 1 0 1 0 1
|
wolffd@0
|
22 % 0 0 0 0 0 1 0 1 0
|
wolffd@0
|
23 % so find(G(5,:)) = [2 4 6 8]
|
wolffd@0
|
24 % but find(G(1,:)) = [2 4]
|
wolffd@0
|
25 %
|
wolffd@0
|
26 % Using wrap around, G = mk_2D_lattice(3, 3, 1), we get
|
wolffd@0
|
27 % 0 1 1 1 0 0 1 0 0
|
wolffd@0
|
28 % 1 0 1 0 1 0 0 1 0
|
wolffd@0
|
29 % 1 1 0 0 0 1 0 0 1
|
wolffd@0
|
30 % 1 0 0 0 1 1 1 0 0
|
wolffd@0
|
31 % 0 1 0 1 0 1 0 1 0
|
wolffd@0
|
32 % 0 0 1 1 1 0 0 0 1
|
wolffd@0
|
33 % 1 0 0 1 0 0 0 1 1
|
wolffd@0
|
34 % 0 1 0 0 1 0 1 0 1
|
wolffd@0
|
35 % 0 0 1 0 0 1 1 1 0
|
wolffd@0
|
36 % so find(G(5,:)) = [2 4 6 8]
|
wolffd@0
|
37 % and find(G(1,:)) = [2 3 4 7]
|
wolffd@0
|
38
|
wolffd@0
|
39 if nargin < 3, wrap_around = 0; end
|
wolffd@0
|
40
|
wolffd@0
|
41 % M contains the number of each cell e.g.
|
wolffd@0
|
42 % 1 4 7
|
wolffd@0
|
43 % 2 5 8
|
wolffd@0
|
44 % 3 6 9
|
wolffd@0
|
45 % North neighbors (assuming wrap around) are
|
wolffd@0
|
46 % 3 6 9
|
wolffd@0
|
47 % 1 4 7
|
wolffd@0
|
48 % 2 5 8
|
wolffd@0
|
49 % Without wrap around, they are
|
wolffd@0
|
50 % 1 4 7
|
wolffd@0
|
51 % 1 4 7
|
wolffd@0
|
52 % 2 5 8
|
wolffd@0
|
53 % The first row is arbitrary, since pixels at the top have no north neighbor.
|
wolffd@0
|
54
|
wolffd@0
|
55 npixels = nrows*ncols;
|
wolffd@0
|
56
|
wolffd@0
|
57 N = 1; E = 2; S = 3; W = 4;
|
wolffd@0
|
58 if wrap_around
|
wolffd@0
|
59 rows{N} = [nrows 1:nrows-1]; cols{N} = 1:ncols;
|
wolffd@0
|
60 rows{E} = 1:nrows; cols{E} = [2:ncols 1];
|
wolffd@0
|
61 rows{S} = [2:nrows 1]; cols{S} = 1:ncols;
|
wolffd@0
|
62 rows{W} = 1:nrows; cols{W} = [ncols 1:ncols-1];
|
wolffd@0
|
63 else
|
wolffd@0
|
64 rows{N} = [1 1:nrows-1]; cols{N} = 1:ncols;
|
wolffd@0
|
65 rows{E} = 1:nrows; cols{E} = [2:ncols 2];
|
wolffd@0
|
66 rows{S} = [2:nrows 2]; cols{S} = 1:ncols;
|
wolffd@0
|
67 rows{W} = 1:nrows; cols{W} = [1 1:ncols-1];
|
wolffd@0
|
68 end
|
wolffd@0
|
69
|
wolffd@0
|
70 M = reshape(1:npixels, [nrows ncols]);
|
wolffd@0
|
71 nbrs = cell(1, 4);
|
wolffd@0
|
72 for i=1:4
|
wolffd@0
|
73 nbrs{i} = M(rows{i}, cols{i});
|
wolffd@0
|
74 end
|
wolffd@0
|
75
|
wolffd@0
|
76
|
wolffd@0
|
77 G = zeros(npixels, npixels);
|
wolffd@0
|
78 if wrap_around
|
wolffd@0
|
79 for i=1:4
|
wolffd@0
|
80 if 0
|
wolffd@0
|
81 % naive
|
wolffd@0
|
82 for p=1:npixels
|
wolffd@0
|
83 G(p, nbrs{i}(p)) = 1;
|
wolffd@0
|
84 end
|
wolffd@0
|
85 else
|
wolffd@0
|
86 % vectorized
|
wolffd@0
|
87 ndx2 = sub2ind([npixels npixels], 1:npixels, nbrs{i}(:)');
|
wolffd@0
|
88 G(ndx2) = 1;
|
wolffd@0
|
89 end
|
wolffd@0
|
90 end
|
wolffd@0
|
91 else
|
wolffd@0
|
92 i = N;
|
wolffd@0
|
93 mask = ones(nrows, ncols);
|
wolffd@0
|
94 mask(1,:) = 0; % pixels in row 1 have no nbr to the north
|
wolffd@0
|
95 ndx = find(mask);
|
wolffd@0
|
96 ndx2 = sub2ind([npixels npixels], ndx, nbrs{i}(ndx));
|
wolffd@0
|
97 G(ndx2) = 1;
|
wolffd@0
|
98
|
wolffd@0
|
99 i = E;
|
wolffd@0
|
100 mask = ones(nrows, ncols);
|
wolffd@0
|
101 mask(:,ncols) = 0;
|
wolffd@0
|
102 ndx = find(mask);
|
wolffd@0
|
103 ndx2 = sub2ind([npixels npixels], ndx, nbrs{i}(ndx));
|
wolffd@0
|
104 G(ndx2) = 1;
|
wolffd@0
|
105
|
wolffd@0
|
106 i = S;
|
wolffd@0
|
107 mask = ones(nrows, ncols);
|
wolffd@0
|
108 mask(nrows,:)=0;
|
wolffd@0
|
109 ndx = find(mask);
|
wolffd@0
|
110 ndx2 = sub2ind([npixels npixels], ndx, nbrs{i}(ndx));
|
wolffd@0
|
111 G(ndx2) = 1;
|
wolffd@0
|
112
|
wolffd@0
|
113 i = W;
|
wolffd@0
|
114 mask = ones(nrows, ncols);
|
wolffd@0
|
115 mask(:,1)=0;
|
wolffd@0
|
116 ndx = find(mask);
|
wolffd@0
|
117 ndx2 = sub2ind([npixels npixels], ndx, nbrs{i}(ndx));
|
wolffd@0
|
118 G(ndx2) = 1;
|
wolffd@0
|
119 end
|
wolffd@0
|
120
|
wolffd@0
|
121 G = setdiag(G, 0);
|