Broken feImage and lib2geom matrix inversion
Hi,
I've traced the problem with the broken feImage filter primitive (see https://bugs.launchpad.net/inkscape/+bug/382313 ) to a problem with lib2geom. More specifically, Geom::Matrix.inverse() returns an incorrect matrix. I run this simple test program:
#include <iostream> #include <2geom/matrix.h>
int main() {
Geom::Matrix test_matrix1( 0.003333, 0.000000, 0.000000, 0.003333, 0.083333, 0.083333 );
Geom::Matrix test_inverse1 = test_matrix1.inverse();
std::cout << "Test input" << std::endl; std::cout << test_matrix1 << std::endl;
std::cout << "Test output" << std::endl; std::cout << test_inverse1 << std::endl;
Geom::Matrix test_matrix2( 0.003077, 0.000000, 0.000000, 0.003077, 0.076923, 0.076923 );
Geom::Matrix test_inverse2 = test_matrix2.inverse();
std::cout << "Test input" << std::endl; std::cout << test_matrix2 << std::endl;
std::cout << "Test output" << std::endl; std::cout << test_inverse2 << std::endl; }
The output I see is:
Test input A: 0.003333 C: 0 E: 0.083333 B: 0 D: 0.003333 F: 0.083333
Test output A: 300.03 C: 0 E: 25.0024 B: 0 D: 300.03 F: 25.0024
Test input A: 0.003077 C: 0 E: 0.076923 B: 0 D: 0.003077 F: 0.076923
Test output A: 1 C: 0 E: 0 B: 0 D: 1 F: 0
You can see that the first matrix inversion is correct but the second is wrong. This was tested with SVN as of yesterday on Linux.
Tav
A bit more info.
On Sun, 20090830 at 13:10 +0200, Tavmjong Bah wrote:
Hi,
I've traced the problem with the broken feImage filter primitive (see https://bugs.launchpad.net/inkscape/+bug/382313 ) to a problem with lib2geom. More specifically, Geom::Matrix.inverse() returns an incorrect matrix.
In Geom::Matrix::inverse() the matrix determinate is calculated and then compared to the constant EPSILON which is defined as 10^5 (in coord.h). If it is less than EPSILON the matrix is set to unity. This EPSILON is clearly way too large for use in inverse(). For some reason a value of 10^18 in commented out in coord.h. My guess is that the value 10^5 was selected for comparing pixel positions to each other but it is being used for other purposes in the lib2geom library.
Tav
On Sun, Aug 30, 2009 at 03:20:44PM +0200, Tavmjong Bah wrote:
... https://bugs.launchpad.net/inkscape/+bug/382313 ... More specifically, Geom::Matrix.inverse() returns an incorrect matrix.
In Geom::Matrix::inverse() the matrix determinate is calculated and then compared to the constant EPSILON which is defined as 10^5 (in coord.h). If it is less than EPSILON the matrix is set to unity. This EPSILON is clearly way too large for use in inverse().
njh thinks that EPSILON shouldn't be smaller than 1e8 with the current code & callers.
Sorry for not checking myself (I'm just about to leave), but is changing EPSILON to 1.5e8 enough to make this test case work?
pjrm.
On Tue, 20090901 at 10:39 +1000, Peter Moulder wrote:
On Sun, Aug 30, 2009 at 03:20:44PM +0200, Tavmjong Bah wrote:
... https://bugs.launchpad.net/inkscape/+bug/382313 ... More specifically, Geom::Matrix.inverse() returns an incorrect matrix.
In Geom::Matrix::inverse() the matrix determinate is calculated and then compared to the constant EPSILON which is defined as 10^5 (in coord.h). If it is less than EPSILON the matrix is set to unity. This EPSILON is clearly way too large for use in inverse().
njh thinks that EPSILON shouldn't be smaller than 1e8 with the current code & callers.
Sorry for not checking myself (I'm just about to leave), but is changing EPSILON to 1.5e8 enough to make this test case work?
The inverse() function is used to find the bounding box of the image in drawing coordinates. Setting EPSILON to 1.5e8 would make this test case work but would cause trouble if anybody wanted to use the feImage filter with an image greater than 8000x8000 px (not common... but you never know). (There is the question of why a matrix inversion is necessary in the first place to find the bounding box but that is a separate issue for a filters expert.)
Fundamentally, I think having one value of EPSILON used everywhere in lib2geom is wrong as the way it is used is different in different functions. Not knowing much about lib2geom, I would guess that using EPISLON^2 in the matrix inverse function would be more reasonable.
Tav
On Tue, Sep 01, 2009 at 08:39:35AM +0200, Tavmjong Bah wrote:
The inverse() function is used to find the bounding box of the image in drawing coordinates. Setting EPSILON to 1.5e8 would make this test case work but would cause trouble if anybody wanted to use the feImage filter with an image greater than 8000x8000 px (not common... but you never know). (There is the question of why a matrix inversion is necessary in the first place to find the bounding box but that is a separate issue for a filters expert.)
If we're literally talking px in the CSS sense (1px = 1/2688 of expected viewing distance) then it is indeed OK not to work on images wider than about 4/3 of viewing distance (3584px) wide. The number 4/3 is one I just made up based on trying to look at a picture and deciding how far away I need to move from a picture for it to be comfortable to look at the whole picture; certainly 8000px (112°) is too big for me to take in the whole picture, even though my peripheral vision extends beyond that.
Of course it's still useful to have an error margin for e.g. zoomrelated uses, preferably an error margin bigger than a factor of 8000/3584, it's just a question of what costs we're willing to take to do so.
Fundamentally, I think having one value of EPSILON used everywhere in lib2geom is wrong as the way it is used is different in different functions.
Agreed. Presumably it arises because it's difficult to come up with the right number to use :) .
Anyway, note that njh has now changed the calculation a little (in upstream 2geom at least), so the best number to use may now be different again.
pjrm.
On Tue, 20090901 at 18:08 +1000, Peter Moulder wrote:
On Tue, Sep 01, 2009 at 08:39:35AM +0200, Tavmjong Bah wrote:
The inverse() function is used to find the bounding box of the image in drawing coordinates. Setting EPSILON to 1.5e8 would make this test case work but would cause trouble if anybody wanted to use the feImage filter with an image greater than 8000x8000 px (not common... but you never know). (There is the question of why a matrix inversion is necessary in the first place to find the bounding box but that is a separate issue for a filters expert.)
If we're literally talking px in the CSS sense (1px = 1/2688 of expected viewing distance) then it is indeed OK not to work on images wider than about 4/3 of viewing distance (3584px) wide. The number 4/3 is one I just made up based on trying to look at a picture and deciding how far away I need to move from a picture for it to be comfortable to look at the whole picture; certainly 8000px (112°) is too big for me to take in the whole picture, even though my peripheral vision extends beyond that.
While it's probably not common, one might use Inkscape to create a picture that is meant to be looked at piece by piece (Chinese scrolls, for example). Also, remember that one might want to use create a drawing with higher resolution to avoid sampling errors upon export at different resolutions.
Of course it's still useful to have an error margin for e.g. zoomrelated uses, preferably an error margin bigger than a factor of 8000/3584, it's just a question of what costs we're willing to take to do so.
Fundamentally, I think having one value of EPSILON used everywhere in lib2geom is wrong as the way it is used is different in different functions.
Agreed. Presumably it arises because it's difficult to come up with the right number to use :) .
Anyway, note that njh has now changed the calculation a little (in upstream 2geom at least), so the best number to use may now be different again.
Can someone please see that this gets patched into Inkscape 0.47 before release. I will be gone for a week with little if any Internet connection so I won't be able to help. Thanks.
pjrm.
Tav
Peter Moulder wrote:
On Tue, Sep 01, 2009 at 08:39:35AM +0200, Tavmjong Bah wrote:
Fundamentally, I think having one value of EPSILON used everywhere in lib2geom is wrong as the way it is used is different in different functions.
Agreed. Presumably it arises because it's difficult to come up with the right number to use :) .
Might I suggest to do it in a fundamentally different way in this particular case (matrix inverse). There are a number of different options that come to mind that would sound a lot more reasonable than any arbitrary constant:
 Use range arithmetic. If you using different rounding modes for the determinant can get you numbers in a range that straddles zero I'd say it's fair to say the matrix is essentially singular.  If changing rounding modes is too problematic, we could always do a poor man's version of range arithmetic using nextafter/nexttoward.  We compute 1.0/det. If det is so small that this (or the value in the resulting matrix with the largest magnitude) is inf(inite) there wouldn't be much sense in computing the inverse, as we can't represent it.
Jasper van de Gronde wrote:
 Use range arithmetic. If you using different rounding modes for the determinant can get you numbers in a range that straddles zero I'd say it's fair to say the matrix is essentially singular.
This is not going to work, because a small determinant will always round to 0 regardless of the rounding mode. Rounding modes only affect what is done to halfway cases like 0.5 when nearbyint() is called; it doesn't affect round() or trunc().
 If changing rounding modes is too problematic, we could always do a poor man's version of range arithmetic using nextafter/nexttoward.
Won't work for the same reason as above; rounding has little to do with this issue.
 We compute 1.0/det. If det is so small that this (or the value in the resulting matrix with the largest magnitude) is inf(inite) there wouldn't be much sense in computing the inverse, as we can't represent it.
I think this is the correct way to fix this (and of course review the way EPSILON is used). We could use some more unit tests (right now I don't even know how to run the existing ones...)
Regards, Krzysztof
participants (4)

Jasper van de Gronde

Krzysztof Kosiński

Peter Moulder

Tavmjong Bah