Can operations on submatrices (and subvectors) be vectorized?

Question

I'm currently working on an edge detector in octave. Coming from other programming languages like Java and Python, I'm used to iterating in for loops, rather than performing operations on entire matrices. Now in octave, this causes a serious performance hit, and I'm having a bit of difficulty figuring out how to vectorize my code. I have the following two pieces of code:

1)

function zc = ZeroCrossings(img, T=0.9257)
  zc = zeros(size(img));

  # Iterate over central positions of all 3x3 submatrices
  for y = 2:rows(img) - 1
    for x = 2:columns(img) - 1
      ndiff = 0;

      # Check all necessary pairs of elements of the submatrix (W/E, N/S, NW/SE, NE/SW)
      for d = [1, 0; 0, 1; 1, 1; 1, -1]'
        p1 = img(y-d(2), x-d(1));
        p2 = img(y+d(2), x+d(1));
        if sign(p1) != sign(p2) && abs(p1 - p2) >= T
          ndiff++;
        end
      end

      # If at least two pairs fit the requirements, these coordinates are a zero crossing
      if ndiff >= 2
        zc(y, x) = 1;
      end
    end
  end
end

2)

function g = LinkGaps(img, k=5)
  g = zeros(size(img));

  for i = 1:rows(img)
    g(i, :) = link(img(i, :), k);
  end
end

function row = link(row, k)
  # Find first 1
  i = 1;
  while i <= length(row) && row(i) == 0
    i++;
  end

  # Iterate over gaps
  while true
    # Determine gap start
    while i <= length(row) && row(i) == 1
      i++;
    end
    start = i;

    # Determine gap stop
    while i <= length(row) && row(i) == 0
      i++;
    end

    # If stop wasn't reached, exit loop
    if i > length(row)
      break
    end

    # If gap is short enough, fill it with 1s
    if i - start <= k
      row(start:i-1) = 1;
    end
  end
end

Both of these functions iterate over submatrices (or rows and subrows in the second case), and particularly the first one seems to be slowing down my program quite a bit.

This function takes a matrix of pixels (img) and returns a binary (0/1) matrix, with 1s where zero crossings (pixels whose corresponding 3x3 neighbourhoods fit certain requirements) were found.

The outer 2 for loops seem like they should be possible to vectorize somehow. I can put the body into its own function (taking as an argument the necessary submatrix) but I can't figure out how to then call this function on all submatrices, setting their corresponding (central) positions to the returned value.

Bonus points if the inner for loop can also be vectorized.
This function takes in the binary matrix from the previous one's output, and fills in gaps in its rows (i.e. sets them to 1). A gap is defined as a series of 0s of length <= k, bounded on both sides by 1s.

Now I'm sure at least the outer loop (the one in LinkGaps) is vectorizable. However, the while loop in link again operates on subvectors, rather than single elements so I'm not sure how I'd go about vectorizing it.

Task 1 can likely be solved with convolution with a kernel designed to detect zero crossings. — Suever

Florian Florian · Accepted Answer · 2017-01-17T16:00:54

Not a full solution, but here is an idea how you could do the first without any loops:

% W/E
I1 = I(2:end-1,1:end-2);
I2 = I(2:end-1,3:end  );
C = (I1 .* I2 < 0) .* (abs(I1 - I2)>=T);

% N/S
I1 = I(1:end-2,2:end-1);
I2 = I(3:end,  2:end-1);
C = C + (I1 .* I2 < 0) .* (abs(I1 - I2)>=T);

% proceed similarly with NW/SE and NE/SW
% ...

% zero-crossings where count is at least 2
ZC = C>=2;

Idea: form two subimages that are appropriately shifted, check for the difference in sign (product negative) and threshold the difference. Both tests return a logical (0/1) matrix, the element-wise product does the logical and, result is a 0/1 matrix with 1 where both tests have succeeded. These matrices can be added to keep track of the counts (ndiff).

Can operations on submatrices (and subvectors) be vectorized?

1 Answers