Reputation: 102
I have a very large 3D image stored in a matrix (approx. 500x500x40
voxels) in Matlab. In this matrix about 30000 points are selected (suppose randomly). Suppose the selected voxels have a value of one and the non-selected points are zero. Now I need to calculate for every voxel in the entire 3D-image the Euclidian distance to the nearest selected point.
So for example in 2D, given a 4x4
matrix:
selection = 0 0 1 0
1 0 0 0
0 0 0 1
0 0 0 0
The corresponding distance matrix would be:
distance = 1 1 0 1
0 1 1 1
1 √2 1 0
2 √5 √2 1
Is there an efficient way to do this, in terms of both time and memory?
Upvotes: 1
Views: 337
Reputation: 46
If you don't need to know which combinations had which distances, then you could calculate the 3D cross-correlation.
To illustrate this in 2D, take the following matrix and calculate the 2D correlation as described in the reference
M =
0 0 0
1 0 0
1 0 0
convn(v,v(end:-1:1,end:-1:1)) =
0 0 0 0 0
0 0 1 0 0
0 0 2 0 0
0 0 1 0 0
0 0 0 0 0
One can read off the distances between the ones because the correlation matrix here can be understood as the differences in indices. The central column in convn means the horizontal distance is zero. Likewise, the middle row means zero vertical distance. As a result, the central value gives you the amount of zero distance values, which is the sum of ones in your matrix. The two ones correspond to a vertical distance of 1 between the ones in M. One combination has a positive distance and the other combination has a negative distance.
As a result, you now have all the distances, but they contain the horizontal and vertical directions too. But you can still process that as you wish.
A way to calculate all possible squared distances without a for loop would be
n = size(M,1) =
3
tmp = repmat([-(n-1):(n-1)].^2,2*n-1,1)
d2 = tmp+tmp' =
8 5 4 5 8
5 2 1 2 5
4 1 0 1 4
5 2 1 2 5
8 5 4 5 8
Both matrices together contain basically the histogram of distances.
Upvotes: 1
Reputation: 46
If your points are given as coordinates X = n x 3, you can use
D = pdist(X,'euclidean')
to efficiently calculate all combinations of distances.
Upvotes: 0