In image processing, images can have multiple-channels. RGB images have three channels and remote sensing data can have many bands. These images can be stored as a three dimensional array in MATLAB, with the first two dimensions representing raster coordinates and the third dimension representing the channel. These 3D arrays might look like this:

Channel 1 | Channel 2 | Channel 3 | |||
---|---|---|---|---|---|

A1 | A3 | B1 | B3 | C1 | C3 |

A2 | A4 | B2 | B4 | C2 | C4 |

Each channel is a two dimensional array but these are stacked together in a 3D array in MATLAB. The channels in these 3D arrays can be interleaved so that the entire 3D array is stored as a 2D array with the use of the `permute`

function (kind of but not exactly like band interleaving for remote sensing images). I came up with this when I was thinking of ways to make coding an image feature extraction algorithm easier by putting all the bands into the same 2D matrix.

- Interleave matrices in MATLAB with reshape: Similar idea but interleaves two regular 2-dimensional matrices.

Interleaving by row will create a single 2D matrix from the 3D matrix above and the entries will look like this:

A1 | A3 |

B1 | B3 |

C1 | C3 |

A2 | A4 |

B2 | B4 |

C2 | C4 |

The code to do this is shown below:

% create some data. a = [11 13; 12 14] b = [21 23; 22 24] c = [31 33; 32 34] % make a 3D array / 3 channel image: mat = cat(3,a,b,c); % shuffle the dimensions around. This switches the 2nd and 3rd dimensions t = permute(mat, [1 3 2]); % Reshape to a single 2D array. col_interleave = reshape(t, [], 3*size(mat,1))

If you are familiar with remote sensing file formats, the interleaved 2D array looks suspiciously like band interleaved by line (BIL). However, with the way MATLAB stores data, you cannot simply `fwrite`

the data and expect it to be stored in BIL format. MATLAB stores the entries in the matrix by iterating over rows, then columns (i.e. `fwrite`

would store them as A1, B1, … C2, A3, … C4).

Interleaving by column will create a single 2D matrix from the 3D matrix above with the following entries:

A1 | B1 | C1 | A3 | B3 | C3 |

A2 | B2 | C2 | A4 | B4 | C4 |

The code to do this is shown below:

% create some data. a = [11 13; 12 14] b = [21 23; 22 24] c = [31 33; 32 34] % make a 3D array / 3 channel image: mat = cat(3,a,b,c); % shuffle the dimensions around so reshape will produce desired result. % [3 1 2] indicates that the original dimension D3 is now the first dimension, % D2 is the last dimension and D1 is the second dimension. t = permute(mat, [3 1 2]); % Reshape to a single 2D array. row_interleave = reshape(t, 3*size(mat,1), [])

As with row interleaving, the interleaved 2D array looks suspiciously like band interleaved by pixel (BIP). However, with the way MATLAB stores data, you cannot simply `fwrite`

the data and expect it to be stored in BIP format. MATLAB stores the entries in the matrix by iterating over rows, then columns (i.e, when written with `fwrite`

, the order is: A1, A2, B1, B2, … C3, C4).

MATLAB array elements can be accessed by a single index number which is the array index. In a 3D array, MATLAB's indexing order increments over rows, then columns and finally channels. The element at (row, column, channel) = (1,1,1) has index 1, the element at (2,1,1) has index 2 and so on, until the last row is reached. The next element index is assigned to the element at (1,2,1). The rows are incremented until the last row and the indexing continues with the next column. When all columns of channel 1 are exhausted, then the indexing continues in channel 2 in the same manner. Multi-dimensional arrays simply extend this process. If you access the array elements with subscripts (D1, D2, D3, … DN), then the index of each element always increments over D1, then D2, then D3, … then DN.

The `permute`

command switches the order of the dimensions. The vector that you supply as the second argument into `permute`

lists the original dimensions in the new order. This has the effect of also changing the indexing of the individual elements. It is easier to see this in a 2D array:

>> a = [1 3 5; 2 4 6] a = 1 3 5 2 4 6

We can flatten the array:

>> a(:) ans = 1 2 3 4 5 6

The order that the elements appear in is the same as the order of the indices of each element.

An array can be permuted, switching dimension 2 with dimension 1 (completely equivalent to a transpose operation; `permute`

is a generalization of transpose to more than 3 dimensions):

>> b = permute(a, [2 1]) b = 1 2 3 4 5 6 >> b(:) ans = 1 3 5 2 4 6

As shown above, the array is transposed by the `permute`

operation but the order of the indices is completely changed. This is important for interleaving the stack of image channels.

The other side of the operation is the reshape operation. Reshape simply changes the number of rows, columns and channels of the array while keeping total number of elements the same. Reshape does not affect the index of each of the individual elements and the reshaped matrix will have elements appear in the original index order, just “re-flowed” to fit the new shape.

In order to get the desired interleaved arrays with `reshape`

, the `permute`

command is used to get the right index order by swapping the dimensions. For example, the following creates a 2-channel image as a 3D array:

a = [11 13; 12 14]; b = [21 23; 22 24]; mat = cat(3,a,b);

The matrix looks like this:

>> mat mat(:,:,1) = 11 13 12 14 mat(:,:,2) = 21 23 22 24

Flatten the array to a 1D vector to see the index ordering:

>> mat(:) ans = 11 12 13 14 21 22 23 24

If you `permute`

the matrix:

mat = permute(mat, [3 1 2]);

The matrix dimensions are shuffled and the index order is different:

>> mat mat(:,:,1) = 11 12 21 22 mat(:,:,2) = 13 14 23 24 >> mat(:) ans = 11 21 12 22 13 23 14 24

It is the permuted matrix that actually has the correct index order to be reshaped into an interleaved matrix (in this case, column interleave).

In order to work out the correct `permute`

command, I just worked backwards from the desired interleaved matrix and found the proper command to use.

## Discussion