bugGNU Octave - Bugs: bug #49613, image package - regionprops.m...

 
 

bug #49613: image package - regionprops.m majoraxislength broken

Submitter:  Hartmut <hardy>
Submitted:  Tue 15 Nov 2016 09:33:59 PM UTC
   
 
Category:  Octave Package Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Regression
Status:  Fixed Assigned to:  None
Originator Name:  Open/Closed:  * Closed
Release:  * dev Operating System:  * Any
Fixed Release:  None Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup
   

Jump to the original submission

Thu 13 Apr 2017 06:38:26 PM UTC, comment #23: 

I have pushed the last patch to the stable branch after splitting into two commits

http://hg.code.sf.net/p/octave/image/rev/ed4dbfec84cf

http://hg.code.sf.net/p/octave/image/rev/794d73956664

I have also simplified the tests.  There is no need to do:


bw = logical ([0 0 0 0; 0 1 1 0; 0 1 1 0]);
l = bwlabel (bw);
props = regionprops (l, "Orientation");


this is enough


bw = logical ([0 0 0 0; 0 1 1 0; 0 1 1 0]);
props = regionprops (bw,"Orientation");


Carnë Draug <carandraug>
Group Member
Mon 12 Dec 2016 01:02:47 PM UTC, comment #22: 

I have tested Avinoam's latest patch (file #39195) and it looks all fine to me:

  • it resolves all mentioned issues of this bug report
    • it solves the original issue (major/minoraxislength with multiple regions), see comment #0
    • it removes the division by zero in the calculation of the ellipse orientation property, see comment #1
  • it gives Matlab compatible results on all mentioned examples in this bug report
  • it adds new tests for those improvements
    • those new tests areMatlab compatible
    • the tests all pass


If noone else finds a problem wit this patch, then I would like to see this patch commited to the image repository.

Thanks for supplying this patch, Avinoam.

Hartmut <hardy>
Sat 10 Dec 2016 09:55:19 PM UTC, comment #21: 

Thanks for checking the patch. I fixed the typo.

I hope this time it will pass the tests.

(file #39195)

Avinoam Kalma <avinoam>
Group Member
Sat 10 Dec 2016 07:20:52 PM UTC, comment #20: 

There is a small typo on your last patch, Avinoam. (It also makes one of your newly added tests fail): The number "90" from the last verbatim-box in comment #17 is mistakenly written as a "0".

I will try to run all your new tests in Matlab during the next days. If they all pass, I would be happy to see this patch applied (with the number 90 corrected).

Hartmut <hardy>
Thu 08 Dec 2016 06:43:48 PM UTC, comment #19: 


Attached a new cs, according to @Hartmut suggestios.

Please review.


(file #39177)

Avinoam Kalma <avinoam>
Group Member
Tue 06 Dec 2016 10:22:44 PM UTC, comment #18: 

Yes, I think the last "patch" in comment #17 would be fine.

Could you add some relevant tests from the discussed example images?

Hartmut <hardy>
Tue 06 Dec 2016 09:22:50 PM UTC, comment #17: 

Indeed, my patch was wrong. I tried to be compatible, but in the wrong way. It should be:


   if (major_vec(2) == 0)
      orientation(idx) = 90;
    else
      orientation(idx) = -(180/pi) .* atan (major_vec(1) ./ major_vec(2));
    end


which works fine on your last example, and on your example from comment #0, except entry 12, which matlab  result is 0, and mine is 90.

The same problem is with


BW_C = logical([0 0 0 0; 0 1 1 0; 0 1 1 0]);
L_C = bwlabel(BW_C);
props = regionprops(L_C,'Orientation');
orientation_C = props.Orientation


In which Matlab result is 0, and mine is 90.

I think that the problem is that both eigenvalues are equal, so
min = max.

What about:


    if (major(idx) == minor(idx))
      orientation(idx) = 0;
    elseif (major_vec(2) == 0)
      orientation(idx) = 90;
    else
      orientation(idx) = -(180/pi) .* atan (major_vec(1) ./ major_vec(2));
    endif


Avinoam Kalma <avinoam>
Group Member
Tue 06 Dec 2016 03:56:04 PM UTC, comment #16: 

With Matlab I get for the test script from comment #15:


orientation_A = 0
orientation_B = 90


This is different compared to the results with the patch in comment #14.

@Avinoam: Could you change your if-statement to also make those results ML compatible?

Hartmut <hardy>
Mon 05 Dec 2016 09:39:50 PM UTC, comment #15: 

I like the fact that your "patch" from comment #14 gives fully Matlab compatible results on my test image (comment #0).

But I suspect that the if-statement of this "patch" will not be Matlab compatible. We should check the Matlab result for something like this

clear;

BW_A = logical([0 0 0 0; 0 1 1 0; 0 0 0 0]);
L_A = bwlabel(BW_A);
props = regionprops(L_A,'Orientation');
orientation_A = props.Orientation

BW_B = BW_A';
L_B = bwlabel(BW_B);
props = regionprops(L_B,'Orientation');
orientation_B = props.Orientation

and might need to adjust the if-statement in your "patch" accordingly.

Hartmut <hardy>
Mon 05 Dec 2016 09:12:45 PM UTC, comment #14: 

Current implementation of orientation


    orientation(idx) = -(180/pi) .* atan (major_vec(2) ./ major_vec(1));


gives:


orientation =

  -35.96725
   82.06034
  -56.90598
  -86.07882
  -62.03271
   42.99091
   69.42293
   37.99984
   84.71827
  -90.00000
   31.92248
   -0.00000
  -31.02163
   14.03068
   34.96791
   64.73045
  -13.02547
   23.02190
   12.62266
  -77.01055
  -20.73767


Changing to:


    if (major_vec(1) == 0 || major_vec(2) == 0)
      orientation(idx) = 0;
    else
      orientation(idx) = -(180/pi) .* atan (major_vec(1) ./ major_vec(2));
    end

 
gives a Matlab compatible results. What do you think?


Avinoam Kalma <avinoam>
Group Member
Thu 01 Dec 2016 03:57:48 PM UTC, comment #13: 

Yes, the answer of Matlab R2014b is the very same as in commment #12.

Hartmut <hardy>
Thu 01 Dec 2016 07:01:28 AM UTC, comment #12: 

Thanks @Hartmut, and sorry, I did not write the correct lines:

for the commands

props1 = regionprops(L,'Orientation');
orientation = [props1.Orientation]'


I get with Matlab 2010a


orientation =

  -54.0328
    7.9397
  -33.0940
   -3.9212
  -27.9673
   47.0091
   20.5771
   52.0002
    5.2817
         0
   58.0775
         0
  -58.9784
   75.9693
   55.0321
   25.2696
  -76.9745
   66.9781
   77.3773
  -12.9895
  -69.2623


Do you get the same results?

Thanks.

Avinoam Kalma <avinoam>
Group Member
Wed 30 Nov 2016 06:26:06 PM UTC, comment #11: 

Here is the Matlab result (R2014b) from the script in comment #0 and the final line as in comment #10 :


Orientation =

   15.0780
   19.6747
   26.5785
   26.6431
   12.7492
    8.1949
   15.0308
   12.8203
   10.4529
   12.7017
   19.5944
   17.3205
   15.0855
   17.3745
   21.9570
   10.4287
    9.1946
   35.8453
   34.7005
   28.9482
   33.4851


Hartmut <hardy>
Tue 29 Nov 2016 05:00:03 PM UTC, comment #10: 

To check Matlab compatability, can some check the result of the
script from comment #0 with the line


Orientation = [props.MinorAxisLength]'


Thanks.

Avinoam Kalma <avinoam>
Group Member
Sun 27 Nov 2016 07:22:12 PM UTC, comment #9: 

About the different results, I suspected it was due to
roundoff errors, so I have changed the tolerance.

Avinoam Kalma <avinoam>
Group Member
Fri 25 Nov 2016 09:34:35 PM UTC, comment #8: 

I can confirm that Avinoam's first patch (file #39060) fixes the original issue. And all tests of regionprops.m now pass on linux and on Windows (w32 tested).

@Avinoam: Do you understand why Windows gives different results (some eps only) compared to Linux for those tests you have changed?

Hartmut <hardy>
Wed 23 Nov 2016 09:52:23 PM UTC, comment #7: 

Attached a changeset for the original bug.
I also added a test for the change, and change tolerances for some
tests that failed in Windows.

Please review.

I will fix later the zero divided problem in Orientation calculation.


(file #39060)

Avinoam Kalma <avinoam>
Group Member
Thu 17 Nov 2016 06:08:00 PM UTC, comment #6: 

Sorry, the argumentation in my last post was wrong. A single pixel gets 1/12 in size, this leads to the guarantee that the majoraxislength and the minoraxislength will never become zero, I think. But this doesn't prove yet that the major_vec(1) will never be 0 in our code. Sorry.

The following image really gives me a value of major_vec(1) = 0, nevertheless the current regionprops.m DOES return the right value for the orientation:


BW = logical([0 0 0 0; 0 1 1 0; 0 0 0 0]);
L = bwlabel(BW);
props = regionprops(L,'Orientation');
props.Orientation


This is due to the result of "atan(1/0)" which does give a warning (no error!), but nevertheless returns the "correct" value of pi/2 (that results in the correct orientation value of -90 degrees).

I have also tried the simple change of atan to atan2, and I do get a failing test of the Orientation feature. So a quick substitute of atan into atan2 will not do.

To conclude: I have no objection about switching from atan to atan2, but the resulting orientation values SHOULD stay the same, please. (So this would mean some additional adjustments to this code line, and double checking that the orientation results stay the same.)

Hartmut <hardy>
Thu 17 Nov 2016 05:41:51 PM UTC, comment #5: 

Regarding Avinoam's second topic (the tan function in the calculation of the ellipse orientation):

I have implemented this calculation myself (see bug #44100) two years ago. When I remember right, I succeeded to make this orientation value also Matlab compatible. So chances are high that changes in the resulting orientation values would be a regression. Do the included test cases still pass when you switch to atan2 ? (But even if yes, I would hesitate to change this).

Here is the Matlab blog article, that I used at that time to work towards Matlab compatibility of the ellipse orientation: http://blogs.mathworks.com/steve/2010/07/30/visualizing-regionprops-ellipse-measurements/

The problem with tan would only occur if major_vec(1)=0. But can this happen at all? Even a region with only one single pixel shoud receive the 1/12 value for its width (it's second moment). So I would guess it is impossible that one of the eigenvalues can become zero here. That would mean that the ellipse degrades to a 0-thin line, but there are still some pixels in it, because without any pixels it wouldn't count as a region at all... (But you could test some corner cases and try to prove me wrong, that's fine!)

Regarding Avinoam's first topic (the missing +1 for the loop counter): This sounds reasonable to fix. Thanks's for spotting this. Could you add some test for calculation the Orientation two regions at once (that have different orientations)? Maybe BW = [1 0 0 0 1; 0 1 0 1 0] or something similar?

Hartmut <hardy>
Thu 17 Nov 2016 04:36:29 AM UTC, comment #4: 

I will make two csets with tests.

Note that the tan/atan2 should be discussed, because it will
not be back compatible, since atan range is [-90:90] and atan2
range is [-180:180]

Avinoam Kalma <avinoam>
Group Member
Wed 16 Nov 2016 08:46:29 PM UTC, comment #3: 

Avinoam, could you prepare two csets (one for each issue) with your proposed fixes and a test for each of them? I will only be able to look at this next week anyway.

Carnë Draug <carandraug>
Group Member
Wed 16 Nov 2016 07:53:11 AM UTC, comment #2: 

Adding package maintainer.

Avinoam Kalma <avinoam>
Group Member
Wed 16 Nov 2016 07:51:48 AM UTC, comment #1: 

The problem is in the loop


  for idx = 1:no
    sel = c_idx:(c_idx + area(idx) -1);
    X = pixellist(sel, 2);
    Y = pixellist(sel, 1);
  ....


c_idx is constant in the loop. Adding the line


    c_idx += area(idx);


fixes the problem.

Another problem is in the line:


orientation(idx) = -(180/pi) .* atan (major_vec(2) ./ major_vec(1));


if major_vec(1)=0 we get "division by zero" warning, and not the correct result. Maybe atan2 can fix this.

Avinoam Kalma <avinoam>
Group Member
Tue 15 Nov 2016 09:33:59 PM UTC, original submission:  

This problem happens with the current version 2.6.1 of the image package. It still worked fine with the older version of regionprops.m that is in the release image-2.4.1.

The property "majoraxislenght" and "minoraxislength" seem to be corrupted somewhere on their way through regionprops.m. If I observe this right, the results are fine when I only request those properties from a single region. But once I request those properties simultanously from several regions at once, I get corrupted results.

Here is a small script to reproduce the behavior. It uses the attached png-file "rechtecke.png":


clear, close all

I = imread('rechtecke.png');
schwelle = graythresh(I);
BW = im2bw(I, schwelle);
L = bwlabel(BW);
L_rgb = label2rgb(L, 'hsv');
imshow(L_rgb);
title('labelled regions');

props = regionprops(L,'MajorAxisLength','MinorAxisLength');
major=[props.MajorAxisLength]';
minor=[props.MinorAxisLength]';
figure, plot(major, minor, 'x');
xlabel('major axis length');
ylabel('minor axis length');


The second plot obviously does not show the correct major and minor axis length with image-2.6.1. If you use image-2.4.1 instead you can see the correct and expected result. I will also attach those two screenshots.

Hartmut <hardy>

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:
   
   
Comment:
   

Attached Files
file #39195:  regioprops2.diff added by avinoam (6KiB - application/octet-stream)
file #39177:  regioprops1.diff added by avinoam (6KiB - application/octet-stream)
file #39060:  regioprops.diff added by avinoam (4KiB - application/octet-stream)
file #38977:  rechtecke.png added by hardy (6KiB - image/png)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by avinoam
  • -email is unavailable- added by avinoam (Posted a comment)
  • -email is unavailable- added by hardy (Submitted the item)
  •  

    There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

    Only group members can vote.

     

    Follow 10 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2017-04-13 carandraug StatusPatch Submitted Fixed
        Open/ClosedOpen Closed
    2016-12-10 avinoam Attached File- Added regioprops2.diff, #39195
    2016-12-08 avinoam Attached File- Added regioprops1.diff, #39177
    2016-11-23 avinoam StatusNone Patch Submitted
    2016-11-23 avinoam Attached File- Added regioprops.diff, #39060
    2016-11-16 avinoam Carbon-Copy- Added -email is unavailable-
    2016-11-15 hardy Attached File- Added rechtecke.png, #38977
        Attached File- Added correct-result_image-2.4.1.png, #38978
        Attached File- Added wrong-result_image-2.6.1.png, #38979

    Back to the top

    Powered by Savane 3.13-758e.
    Corresponding source code