CODEWITHSUNDEEP
arraysmatlabvectorization
noir
noir

Reputation: 185

Efficiently vectorize an element-wise operation in matlab

I have an nx4 matrix A representing n spheres, and an mx3 matrix B representing m points. I need to test whether these m points are inside any of the spheres. I can do this using a for loop, but with large n and m this method is very inefficient. How can I vectorize this operation? My current method is

A = [0.8622    1.1594    0.7457    0.6925;
     1.4325    0.2559    0.0520    0.4687;
     1.8465    0.3979    0.2850    0.4259;
     1.4387    0.8713    1.6585    0.4616;
     0.2383    1.5208    0.5415    0.9417;
     1.6812    0.2045    0.1290    0.1972];

B = [0.5689    0.9696    0.8196;
     0.5211    0.4462    0.6254;
     0.9000    0.4894    0.2202;
     0.4192    0.9229    0.4639];

for i=1:size(B,1)

    mask = vecnorm(A(:, 1:3) - B(i,:), 2, 2) < A(:, 4);

    if sum(mask) > 0
        C(i) = true;
    else
        C(i) = false;
    end %if

end %for

I tested the method suggested by @LuisMendo, and it seems it only speeds up the calculation for quite small m and n, but for large m and n, say, around 10000 for my problem, the improvement is very limited. But @NickyMattsson gave me some hint. Because logical operation in matlab is faster than vecnorm, I first use a rough check to find the spheres near the point, and then do a fine check: 

A = [0.8622    1.1594    0.7457    0.6925;
     1.4325    0.2559    0.0520    0.4687;
     1.8465    0.3979    0.2850    0.4259;
     1.4387    0.8713    1.6585    0.4616;
     0.2383    1.5208    0.5415    0.9417;
     1.6812    0.2045    0.1290    0.1972];

B = [0.5689    0.9696    0.8196;
     0.5211    0.4462    0.6254;
     0.9000    0.4894    0.2202;
     0.4192    0.9229    0.4639];

ids = 1:size(A, 1);

for i=1:size(B,1)

    % first a rough check
    xbound = abs(A(:, 1) - B(i, 1)) < A(:, 4);
    ybound = abs(A(:, 2) - B(i, 2)) < A(:, 4);
    zbound = abs(A(:, 3) - B(i, 3)) < A(:, 4);
    nears = ids(xbound & ybound & zbound);
    if isempty(nears)
        C(i) = false;
    else 
        % then a fine check
        mask = vecnorm(A(nears, 1:3) - B(i,:), 2, 2) < A(nears, 4);

        if sum(mask) > 0
            C(i) = true;
        else
            C(i) = false;
        end 
    end

end

This may reduce the time to 1/2 or 1/3, which is acceptable, and if I divide m and n into batches it may be even faster without too heavy memory burden. @CrisLuengo mentioned the R*-tree method, but it seems that the implementation is quite complicated XD

Upvotes: 3

Views: 111

Answers (1)

Luis Mendo
Luis Mendo

Reputation: 112679

This uses implicit expansion to compute all distances between points and sphere centers, and then to compare those with the sphere radii:

C = any(vecnorm(permute(B, [1 3 2]) - permute(A(:,1:3), [3 1 2]), 2, 3) < A(:,4).', 2);

This is probably faster than the loop approach, but also more memory-intensive, because an intermediate m×n×3 array is computed.

Upvotes: 3

Related Questions