1

我正在尝试对以下代码进行矢量化:

% code before 
% code before 
% a lot of code before we got to the current comment 
% 
% houghMatrix holds some values 
for i=1:n
    for j=1:m
        for k = 1:maximalRadius            
            % get the maximal threshold 
            if houghMatrix(i,j,k) > getMaximalThreshold(k)                           
                lhs = [j i k];

                % verify that the new circle is not listed 
                isCircleExist = verifyCircleExists(circles,lhs,circleCounter);

                % not listed - then we put it in the circles vector 
                if isCircleExist == 0
                    circles(circleCounter,:) = [j i k];                    
                    fprintf('Circle % d: % d, % d, % d \n', circleCounter, j, i, k);
                    circleCounter = circleCounter + 1;                    
                end
            end
        end
    end
end

使用 tic tac 我得到以下输出:

>> x = findCircles(ii);
Circle  1:  38,  38,  35 
Circle  2:  89,  51,  34 
Circle  3:  72,  66,  11 
Circle  4:  33,  75,  30 
Circle  5:  90,  81,  31 
Circle  6:  54,  96,  26 

    Elapsed time is 3.111176 seconds.
>> x = findCircles(ii);
Circle  1:  38,  38,  35 
Circle  2:  89,  51,  34 
Circle  3:  72,  66,  11 
Circle  4:  33,  75,  30 
Circle  5:  90,  81,  31 
Circle  6:  54,  96,  26 

    Elapsed time is 3.105642 seconds.
>> x = findCircles(ii);
Circle  1:  38,  38,  35 
Circle  2:  89,  51,  34 
Circle  3:  72,  66,  11 
Circle  4:  33,  75,  30 
Circle  5:  90,  81,  31 
Circle  6:  54,  96,  26 

    Elapsed time is 3.135818 seconds.

意义——平均3.1秒。

我试图对代码进行矢量化,但问题是我需要使用i,j,k内部for(第 3 个for)主体中的索引。

任何有关如何矢量化代码的建议将不胜感激

谢谢

编辑 :

% -- function [circleExists] = verifyCircleExists(circles,lhs,total) --
%
%
function [circleExists] = verifyCircleExists(circles,lhs,total)

    MINIMUM_ALLOWED_THRESHOLD = 2;

    circleExists = 0;
    for index = 1:total-1                
        rhs = circles(index,:);
        absExpr = abs(lhs - rhs);
        maxValue = max( absExpr );
        if  maxValue <= MINIMUM_ALLOWED_THRESHOLD + 1
            circleExists = 1;
            break
        end
    end

end
4

1 回答 1

0

这是我认为您想要做的:对于每个有效的坐标三元组,您要检查附近是否已经存在三元组,否则,将其添加到列表中。如果没有“链接”的可能性,即如果每个可能的候选体素集群只能容纳一个中心,则该操作可以完全矢量化。在这种情况下,您只需使用:

%# create a vector of thresholds
maximalThreshold = getMaximalThreshold(1:maximalRadius);

%# make it 1-by-1-by-3
maximalThreshold = reshape(maximalThreshold,1,1,[]);

%# create a binary array the size of houghMatrix with 1's 
%# wherever we have a candidate circle center
validClusters = bsxfun(@gt, houghMatrix, maximalThreshold);

%# get the centroids of all valid clusters
stats = regionprops(validClusters,'Centroid');

%# collect centroids, round to get integer pixel values
circles = round(cat(1,stats.Centroid));

或者,如果您想遵循选择有效圆圈的方案,您可以从validClusters以下获取 ijk 索引:

[potentialCircles(:,1),potentialCircles(:,2), potentialCircles(:,3)]= ...  
       sub2ind(size(houghMatrix),find(validClusters));

nPotentialCircles = size(potentialCircles,1);


for iTest = 2:nPotentialCircles
    absDiff = abs(bsxfun(@minus,potentialCircles(1:iTest-1,:),potentialCircles(iTest,:)));

    if any(absDiff(:) <= MINIMUM_ALLOWED_THRESHOLD + 1)
       %# mask the potential circle
       potentialCircles(iTest,:) = NaN;
    end
end

circles = potentialCircles(isfinite(potentialCircles(:,1)),:);
于 2013-01-05T21:36:47.633 回答