title: quickest way to calculate neighbourhood of pixels
body: <p>I have a matrix which consists of positive and negative values. I need to do these things.</p>
  
  <p>Let <code>u(i,j)</code> denote the pixels of the matrix <code>u</code>.</p>
  
  <ol>
  <li>Calculate the zero crossing pixels. These are the pixels in the grid if <code>u(i-1,j)</code> and <code>u(i+1,j)</code> are of opposite signs or <code>u(i,j-1)</code> and <code>u(i,j+1)</code> are of opposite signs.</li>
  <li>Then I need to calculate the narrow band around these zero crossing pixels. The width of the narrow band is <code>(2r+1)X(2r+1)</code> for each pixel. I am taking <code>r=1</code> so for it I have to actually get the 8 neighborhood pixels of each zero crossing pixels. </li>
  </ol>
  
  <p>I have done this in a program. Please see it below. </p>
  
  <pre><code>%// calculate the zero crossing pixels  
  front = isfront(u);
  %// calculate the narrow band of around the zero crossing pixels
  band  = isband(u,front,1);
  </code></pre>
  
  <p>I am also attaching the <code>isfront</code> and <code>isband</code> functions.</p>
  
  <pre><code>function front = isfront( phi )
  %// grab the size of phi
  [n, m] = size(phi);
  
  %// create an boolean matrix whose value at each pixel is 0 or 1
  %// depending on whether that pixel is a front point or not
  front = zeros( size(phi) );
  
  %// A piecewise(Segmentation) linear approximation to the front is contructed by
  %// checking each pixels neighbour. Do not check pixels on border.
  for i = 2 : n - 1;
    for j = 2 : m - 1;
      if (phi(i-1,j)*phi(i+1,j)&lt;0) || (phi(i,j-1)*phi(i,j+1)&lt;0)
          front(i,j) = 100;
      else
          front(i,j) = 0;
      end
    end
  end
  
  function band = isband(phi, front, width)
  %// grab size of phi
  [m, n] = size(phi);
  
  %// width=r=1;
  width = 1;
  
  [x,y] = find(front==100);
  
  %// create an boolean matrix whose value at each pixel is 0 or 1
  %// depending on whether that pixel is a band point or not
  band = zeros(m, n);
  
  %// for each pixel in phi
  for ii = 1:m
    for jj = 1:n
      for k = 1:size(x,1)
          if (ii==x(k)) &amp;&amp; (jj==y(k))
              band(ii-1,jj-1) = 100;  band(ii-1,jj) = 100; band(ii-1,jj+1) = 100;
              band(ii  ,jj-1) = 100;  band(ii  ,jj) = 100; band(ii,jj+1) = 100;
              band(ii+1,jj-1) = 100;  band(ii+1,jj) = 100; band(ii+1,jj+1) = 100;
          end
      end
    end
  end
  </code></pre>
  
  <p>The outputs are given below:, as well as the computation time: </p>
  
  <p><img src="https://i.imgur.com/iX3WANw.jpg" alt="Figure"></p>
  
  <pre><code>%// Computation time
  
  %// for isfront function
  Elapsed time is 0.003413 seconds.
  
  %// for isband function
  Elapsed time is 0.026188 seconds.
  </code></pre>
  
  <p>When I run the code I do get the correct answers but the computation for the tasks is too much to my liking. Is there a better way to do it ? Especially the <code>isband</code> function? How can I optimize my code further ? </p>
  
  <p>Thanks in advance. </p>