Finding distance between endpoints of image segment: PixelList vs 'endpoints'

New Code:

      %%Calculate Tortuosity
      disp('Calculating tortuosity of large segments');
      V = thin8-imdilate(crosslinks,strel('disk',2)); 
      V(V<=0) = 0; V(V>0) = 1; % return to binary
      V = bwareaopen(V,3); % get rid of small segments
      [K n] = bwlabel(V,8);
      props = regionprops(K,'Perimeter');
      tort = zeros(1,n);
      endpts = bwmorph(V,'endpoints');
      endpts_labeled = immultiply(K,endpts);
      for m = 1:n
          [rows,cols] = find(endpts_labeled==m);
          if numel(rows) ~= 2
              tort(m) = NaN;
          else
              p1 = [rows(1),cols(1)];
              p2 = [rows(2),cols(2)];
              d = sqrt((p1(1)-p2(1))^2+(p1(2)-p2(2))^2);
              c = props(m,1).Perimeter/2;
              tort(m) = c/d;
          end
      end
      tort(isnan(tort)) = [];
      tort_avg = mean(tort);
      tort_std = std(tort);

Original Post:

I am working to find the tortuosity of about 500-1000 segments like this one:

They are randomly spaced throughout a larger image and have varying lengths and orientation. I am using a simple definition of tortuosity (tort = distance between endpoints/length of curve). I am getting an accurate reading of the curve length verified via Matlab Blog Post, but I am not sure if I am getting an accurate distance between the endpoints.

I am trying two different methods (using PixelList and the slower method of Endpoints), but neither method is giving me the exact number expected for some of the segments I spot checked and I was wondering which is the more accurate method for determining the euclidean distance of the segments endpoints. My code is as follows:

Answer 1

[Image Analyst]

It's slow because you're using it in a loop. You don't need to call it in a loop. You can just call it on your labeled image once. But it gives you a binary image so you'd have to multiply that by your labeled image to get the endpoints associated with particular objects. Then you'd use that masked labeled image in the loop instead of calling bwmorph each time.

If you know for a fact that PixelList has the first element as one endpoint, and the last element as the last element, then you can use your method. However I doubt this is the case. What is you had a blob instead of a line? See what I mean? If not, run this simple example that illustrates why:

m = [0 0 2 0;0 2 0 0; 0 0 2 0]
measurements = regionprops(m, 'PixelList')
measurements.PixelList

ans =
     2     2
     3     1
     3     3

The first pixel is the middle pixel (location 2,2) of the arc! Now you see why you need to do it like I said at the beginning. The first pixel at 3,1 is actually the second element of PixelList and the last pixel is that last element (so only that one is correct). You need to understand that connected components works by finding the upper left most pixel and then proceeding down rows first, then across columns after there are no more pixels in the current column. I know that's not particularly easy to find or well explained anywhere in the documentation that I know of, but if there is, Sean, Jeff, or Steve will direct you to it.

[Ryan]

I did check the PixelList documentation to try to figure out how they were assigned, and my assumption that it was not along the path of the line appears to be correct. I am not sure why I didn't think of using the 'endpoints' as a mask on the labels. Thanks.

[Ryan]

So in applying this method today, I ended up with an odd exception about 1 in 3 images where a segment labeled would have either just 1 or 0 points corresponding to the 'endpoints.' I threw in the if statement to catch it and ignore them, but I am not sure why I am getting them. I thought maybe extremely small segments were doing it, but I get rid of those < 3 pixels in length. I've added the code to the original post for formatting purposes.

Answer 2

[Thomas]

You could lockup the Medical Image processing webinar which deals with toruosity.. You could also try bwdistgeodesic command..

<http://www.mathworks.com/company/events/webinars/wbnr66293.html?id=66293&p1=1075061617&p2=1075061622>

 doc bwdist
 doc bwdistgeodesic

[Ryan]

are the webinars free?

[Thomas]

yes. they are..

[Jay Cheng]

Can I download Tortuosity.m in the webinar?