From b394d3401750c3d1b1a43aa262abb0a6b95f5824 Mon Sep 17 00:00:00 2001 From: Igor Pisarev Date: Thu, 1 Nov 2018 11:54:25 -0400 Subject: [PATCH 1/2] add fit() overload to affine models with floating point precision arg --- .../java/mpicbg/models/AffineModel1D.java | 54 ++++++++++++++++--- .../java/mpicbg/models/AffineModel2D.java | 46 ++++++++++++++-- .../java/mpicbg/models/AffineModel3D.java | 51 +++++++++++++++--- mpicbg/src/main/java/mpicbg/util/Util.java | 36 +++++++++++++ 4 files changed, 170 insertions(+), 17 deletions(-) diff --git a/mpicbg/src/main/java/mpicbg/models/AffineModel1D.java b/mpicbg/src/main/java/mpicbg/models/AffineModel1D.java index 88f7138a..c9c984cd 100644 --- a/mpicbg/src/main/java/mpicbg/models/AffineModel1D.java +++ b/mpicbg/src/main/java/mpicbg/models/AffineModel1D.java @@ -18,6 +18,8 @@ import java.util.Collection; +import mpicbg.util.Util; + /** * * @author Stephan Saalfeld <saalfelds@janelia.hhmi.org> @@ -99,6 +101,20 @@ final public void fit( final double[][] q, final double[] w ) throws NotEnoughDataPointsException, IllDefinedDataPointsException + { + fit( p, q, w, 0 ); + } + + /** + * Closed form weighted least squares solution as described by + * \citet{SchaeferAl06}. + */ + final public void fit( + final double[][] p, + final double[][] q, + final double[] w, + final double epsilon ) + throws NotEnoughDataPointsException, IllDefinedDataPointsException { assert p.length >= 1 && @@ -113,7 +129,7 @@ final public void fit( final int l = p[ 0 ].length; if ( l < MIN_NUM_MATCHES ) - throw new NotEnoughDataPointsException( l + " data points are not enough to estimate a 2d affine model, at least " + MIN_NUM_MATCHES + " data points required." ); + throw new NotEnoughDataPointsException( l + " data points are not enough to estimate a 1d affine model, at least " + MIN_NUM_MATCHES + " data points required." ); double pcx = 0; double qcx = 0; @@ -143,7 +159,7 @@ final public void fit( b += wwpx * qx; } - if ( a == 0 ) + if ( Util.isApproxEqual( a, 0, epsilon ) ) throw new IllDefinedDataPointsException(); m00 = b / a; @@ -163,6 +179,20 @@ final public void fit( final float[][] q, final float[] w ) throws NotEnoughDataPointsException, IllDefinedDataPointsException + { + fit( p, q, w, 0 ); + } + + /** + * Closed form weighted least squares solution as described by + * \citet{SchaeferAl06}. + */ + final public void fit( + final float[][] p, + final float[][] q, + final float[] w, + final float epsilon) + throws NotEnoughDataPointsException, IllDefinedDataPointsException { assert p.length >= 1 && @@ -177,7 +207,7 @@ final public void fit( final int l = p[ 0 ].length; if ( l < MIN_NUM_MATCHES ) - throw new NotEnoughDataPointsException( l + " data points are not enough to estimate a 2d affine model, at least " + MIN_NUM_MATCHES + " data points required." ); + throw new NotEnoughDataPointsException( l + " data points are not enough to estimate a 1d affine model, at least " + MIN_NUM_MATCHES + " data points required." ); double pcx = 0; double qcx = 0; @@ -207,7 +237,7 @@ final public void fit( b += wwpx * qx; } - if ( a == 0 ) + if ( Util.isApproxEqual( a, 0, epsilon ) ) throw new IllDefinedDataPointsException(); m00 = b / a; @@ -219,15 +249,23 @@ final public void fit( /** * Closed form weighted least squares solution as described by * \citet{SchaeferAl06}. - * - * TODO */ @Override final public < P extends PointMatch >void fit( final Collection< P > matches ) throws NotEnoughDataPointsException, IllDefinedDataPointsException + { + fit( matches, 0 ); + } + + /** + * Closed form weighted least squares solution as described by + * \citet{SchaeferAl06}. + */ + final public < P extends PointMatch >void fit( final Collection< P > matches, final double epsilon ) + throws NotEnoughDataPointsException, IllDefinedDataPointsException { if ( matches.size() < MIN_NUM_MATCHES ) - throw new NotEnoughDataPointsException( matches.size() + " data points are not enough to estimate a 2d affine model, at least " + MIN_NUM_MATCHES + " data points required." ); + throw new NotEnoughDataPointsException( matches.size() + " data points are not enough to estimate a 1d affine model, at least " + MIN_NUM_MATCHES + " data points required." ); double pcx = 0; double qcx = 0; @@ -262,7 +300,7 @@ final public < P extends PointMatch >void fit( final Collection< P > matches ) b += wwpx * qx; } - if ( a == 0 ) + if ( Util.isApproxEqual( a, 0, epsilon ) ) throw new IllDefinedDataPointsException(); m00 = b / a; diff --git a/mpicbg/src/main/java/mpicbg/models/AffineModel2D.java b/mpicbg/src/main/java/mpicbg/models/AffineModel2D.java index d6ded2f6..1822db71 100644 --- a/mpicbg/src/main/java/mpicbg/models/AffineModel2D.java +++ b/mpicbg/src/main/java/mpicbg/models/AffineModel2D.java @@ -19,6 +19,8 @@ import java.awt.geom.AffineTransform; import java.util.Collection; +import mpicbg.util.Util; + /** * 2d-affine transformation models to be applied to points in 2d-space. * This model includes the closed form weighted least squares solution as @@ -117,6 +119,20 @@ final public void fit( final double[][] q, final double[] w ) throws NotEnoughDataPointsException, IllDefinedDataPointsException + { + fit( p, q, w, 0 ); + } + + /** + * Closed form weighted least squares solution as described by + * \citet{SchaeferAl06}. + */ + final public void fit( + final double[][] p, + final double[][] q, + final double[] w, + final double epsilon ) + throws NotEnoughDataPointsException, IllDefinedDataPointsException { assert p.length >= 2 && @@ -183,7 +199,7 @@ final public void fit( final double det = a00 * a11 - a01 * a01; - if ( det == 0 ) + if ( Util.isApproxEqual( det, 0, epsilon ) ) throw new IllDefinedDataPointsException(); m00 = ( a11 * b00 - a01 * b10 ) / det; @@ -207,6 +223,20 @@ final public void fit( final float[][] q, final float[] w ) throws NotEnoughDataPointsException, IllDefinedDataPointsException + { + fit( p, q, w, 0 ); + } + + /** + * Closed form weighted least squares solution as described by + * \citet{SchaeferAl06}. + */ + final public void fit( + final float[][] p, + final float[][] q, + final float[] w, + final float epsilon ) + throws NotEnoughDataPointsException, IllDefinedDataPointsException { assert p.length >= 2 && @@ -273,7 +303,7 @@ final public void fit( final double det = a00 * a11 - a01 * a01; - if ( det == 0 ) + if ( Util.isApproxEqual( det, 0, epsilon ) ) throw new IllDefinedDataPointsException(); m00 = ( a11 * b00 - a01 * b10 ) / det; @@ -294,6 +324,16 @@ final public void fit( @Override final public < P extends PointMatch >void fit( final Collection< P > matches ) throws NotEnoughDataPointsException, IllDefinedDataPointsException + { + fit( matches, 0 ); + } + + /** + * Closed form weighted least squares solution as described by + * \citet{SchaeferAl06}. + */ + final public < P extends PointMatch >void fit( final Collection< P > matches, final double epsilon ) + throws NotEnoughDataPointsException, IllDefinedDataPointsException { if ( matches.size() < MIN_NUM_MATCHES ) throw new NotEnoughDataPointsException( matches.size() + " data points are not enough to estimate a 2d affine model, at least " + MIN_NUM_MATCHES + " data points required." ); @@ -343,7 +383,7 @@ final public < P extends PointMatch >void fit( final Collection< P > matches ) final double det = a00 * a11 - a01 * a01; - if ( det == 0 ) + if ( Util.isApproxEqual( det, 0, epsilon ) ) throw new IllDefinedDataPointsException(); m00 = ( a11 * b00 - a01 * b10 ) / det; diff --git a/mpicbg/src/main/java/mpicbg/models/AffineModel3D.java b/mpicbg/src/main/java/mpicbg/models/AffineModel3D.java index b5210a59..1dc62b6a 100644 --- a/mpicbg/src/main/java/mpicbg/models/AffineModel3D.java +++ b/mpicbg/src/main/java/mpicbg/models/AffineModel3D.java @@ -19,6 +19,7 @@ import java.util.Collection; import mpicbg.util.Matrix3x3; +import mpicbg.util.Util; /** * 3d-affine transformation models to be applied to points in 3d-space. @@ -148,6 +149,20 @@ final public void fit( final double[][] q, final double[] w ) throws NotEnoughDataPointsException, IllDefinedDataPointsException + { + fit( p, q, w, 0 ); + } + + /** + * Closed form weighted least squares solution as described by + * \citet{SchaeferAl06}. + */ + final public void fit( + final double[][] p, + final double[][] q, + final double[] w, + final double epsilon ) + throws NotEnoughDataPointsException, IllDefinedDataPointsException { assert p.length >= 3 && @@ -162,7 +177,7 @@ final public void fit( final int l = p[ 0 ].length; if ( l < MIN_NUM_MATCHES ) - throw new NotEnoughDataPointsException( l + " data points are not enough to estimate a 2d affine model, at least " + MIN_NUM_MATCHES + " data points required." ); + throw new NotEnoughDataPointsException( l + " data points are not enough to estimate a 3d affine model, at least " + MIN_NUM_MATCHES + " data points required." ); double pcx = 0, pcy = 0, pcz = 0; double qcx = 0, qcy = 0, qcz = 0; @@ -244,7 +259,7 @@ final public void fit( a12 * a12 * a00 - a22 * a01 * a01; - if ( det == 0 ) + if ( Util.isApproxEqual( det, 0, epsilon ) ) throw new IllDefinedDataPointsException(); final double idet = 1.0 / det; @@ -285,6 +300,20 @@ final public void fit( final float[][] q, final float[] w ) throws NotEnoughDataPointsException, IllDefinedDataPointsException + { + fit( p, q, w, 0 ); + } + + /** + * Closed form weighted least squares solution as described by + * \citet{SchaeferAl06}. + */ + final public void fit( + final float[][] p, + final float[][] q, + final float[] w, + final float epsilon ) + throws NotEnoughDataPointsException, IllDefinedDataPointsException { assert p.length >= 3 && @@ -299,7 +328,7 @@ final public void fit( final int l = p[ 0 ].length; if ( l < MIN_NUM_MATCHES ) - throw new NotEnoughDataPointsException( l + " data points are not enough to estimate a 2d affine model, at least " + MIN_NUM_MATCHES + " data points required." ); + throw new NotEnoughDataPointsException( l + " data points are not enough to estimate a 3d affine model, at least " + MIN_NUM_MATCHES + " data points required." ); double pcx = 0, pcy = 0, pcz = 0; double qcx = 0, qcy = 0, qcz = 0; @@ -381,7 +410,7 @@ final public void fit( a12 * a12 * a00 - a22 * a01 * a01; - if ( det == 0 ) + if ( Util.isApproxEqual( det, 0, epsilon ) ) throw new IllDefinedDataPointsException(); final double idet = 1.0 / det; @@ -419,9 +448,19 @@ final public void fit( @Override final public < P extends PointMatch >void fit( final Collection< P > matches ) throws NotEnoughDataPointsException, IllDefinedDataPointsException + { + fit( matches, 0 ); + } + + /** + * Closed form weighted least squares solution as described by + * \citet{SchaeferAl06}. + */ + final public < P extends PointMatch >void fit( final Collection< P > matches, final double epsilon ) + throws NotEnoughDataPointsException, IllDefinedDataPointsException { if ( matches.size() < MIN_NUM_MATCHES ) - throw new NotEnoughDataPointsException( matches.size() + " data points are not enough to estimate a 2d affine model, at least " + MIN_NUM_MATCHES + " data points required." ); + throw new NotEnoughDataPointsException( matches.size() + " data points are not enough to estimate a 3d affine model, at least " + MIN_NUM_MATCHES + " data points required." ); double pcx = 0, pcy = 0, pcz = 0; double qcx = 0, qcy = 0, qcz = 0; @@ -494,7 +533,7 @@ final public < P extends PointMatch >void fit( final Collection< P > matches ) a12 * a12 * a00 - a22 * a01 * a01; - if ( det == 0 ) + if ( Util.isApproxEqual( det, 0, epsilon ) ) throw new IllDefinedDataPointsException(); final double idet = 1.0 / det; diff --git a/mpicbg/src/main/java/mpicbg/util/Util.java b/mpicbg/src/main/java/mpicbg/util/Util.java index 01d11cd9..2361bd0d 100644 --- a/mpicbg/src/main/java/mpicbg/util/Util.java +++ b/mpicbg/src/main/java/mpicbg/util/Util.java @@ -366,4 +366,40 @@ final public static void memset( final char[] array, final char value ) for ( int i = 1; i < len; i += i ) System.arraycopy( array, 0, array, i, ( ( len - i ) < i ) ? ( len - i ) : i ); } + + /** + * Checks if two real values are approximately equal. + * + * @param a + * @param b + * @param threshold + * @return + */ + final public static boolean isApproxEqual( final float a, final float b, final float threshold ) + { + if ( a == b ) + return true; + else if ( a + threshold > b && a - threshold < b ) + return true; + else + return false; + } + + /** + * Checks if two real values are approximately equal. + * + * @param a + * @param b + * @param threshold + * @return + */ + final public static boolean isApproxEqual( final double a, final double b, final double threshold ) + { + if ( a == b ) + return true; + else if ( a + threshold > b && a - threshold < b ) + return true; + else + return false; + } } From 8df5b0bebc3dc64218e2e63168ee48dd9580e946 Mon Sep 17 00:00:00 2001 From: Igor Pisarev Date: Thu, 1 Nov 2018 12:16:49 -0400 Subject: [PATCH 2/2] test detection of ill-defined point match sets with user-specified precision --- mpicbg/pom.xml | 7 ++ .../models/IllDefinedConfigurationTest.java | 78 +++++++++++++++++++ 2 files changed, 85 insertions(+) create mode 100644 mpicbg/src/test/java/mpicbg/models/IllDefinedConfigurationTest.java diff --git a/mpicbg/pom.xml b/mpicbg/pom.xml index 44986984..d63cde68 100644 --- a/mpicbg/pom.xml +++ b/mpicbg/pom.xml @@ -157,5 +157,12 @@ gov.nist.math jama + + + + junit + junit + test + diff --git a/mpicbg/src/test/java/mpicbg/models/IllDefinedConfigurationTest.java b/mpicbg/src/test/java/mpicbg/models/IllDefinedConfigurationTest.java new file mode 100644 index 00000000..8a8f1d7b --- /dev/null +++ b/mpicbg/src/test/java/mpicbg/models/IllDefinedConfigurationTest.java @@ -0,0 +1,78 @@ +/** + * License: GPL + * + * This program is free software; you can redistribute it and/or + * modify it under the terms of the GNU General Public License 2 + * as published by the Free Software Foundation. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. + */ +package mpicbg.models; + +import java.util.ArrayList; +import java.util.Collections; +import java.util.List; + +import org.junit.Assert; +import org.junit.Test; + +/** +* @author Igor Pisarev +*/ +public class IllDefinedConfigurationTest +{ + @Test( expected = IllDefinedDataPointsException.class ) + public void testIllDefined() throws NotEnoughDataPointsException, IllDefinedDataPointsException + { + final List< PointMatch > pointMatches = new ArrayList< PointMatch >(); + pointMatches.add( new PointMatch( new Point( new double[] { 11391.245522, -9557.644515, 3214.912133 } ), new Point( new double[] { 2334.168005583636, 3518.967770862933, 107.03217906932194 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 11070.939867, -8226.053258, 3213.153702 } ), new Point( new double[] { 2015.8052608873058, 4847.946467837462, 117.1739487701094 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 11070.928153, -8891.92001, 3213.146659 } ), new Point( new double[] { 2013.7415143287162, 4180.5686283211535, 114.917495945355 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 11391.257236, -8891.777763, 3214.919176 } ), new Point( new double[] { 2335.582386241431, 4185.54358979955, 111.22513943121862 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 11711.574604, -9557.502269, 3216.68465 } ), new Point( new double[] { 2652.076926238893, 3520.681050219568, 102.93656776245956 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 11711.586318, -8891.635516, 3216.691693 } ), new Point( new double[] { 2653.9875667012084, 4181.60078651819, 107.51853905527568 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 11070.91644, -9557.786762, 3213.139616 } ), new Point( new double[] { 2011.6923045944154, 3513.598776666881, 110.35896771552083 } ) ) ); + + final double weight = 1. / pointMatches.size(); + for ( final PointMatch pointMatch : pointMatches ) + pointMatch.setWeights( new double[] { weight } ); + + Collections.shuffle( pointMatches ); + + final AffineModel3D model = new AffineModel3D(); + model.fit( pointMatches, 1e-3 ); + } + + @Test + public void testWellDefined() throws NotEnoughDataPointsException, IllDefinedDataPointsException + { + final List< PointMatch > pointMatches = new ArrayList< PointMatch >(); + pointMatches.add( new PointMatch( new Point( new double[] { 10109.940901, -8892.346749, 3207.829103 } ), new Point( new double[] { 1053.264137949094, 4178.499556821791, 123.9601076750317 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 10750.610779, -8226.195504, 3211.38118 } ), new Point( new double[] { 1697.794268507087, 4847.981272290586, 119.63945933786255 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 9789.623532, -8226.622244, 3206.063629 } ), new Point( new double[] { 735.531633561916, 4840.987725209846, 127.4219850438575 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 10430.269983, -8892.204503, 3209.60162 } ), new Point( new double[] { 1374.5419728598247, 4176.975761289115, 121.43534016162218 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 10430.281697, -8226.337751, 3209.608663 } ), new Point( new double[] { 1375.6180691601485, 4841.696951542491, 120.5498694360961 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 10750.599065, -8892.062256, 3211.374137 } ), new Point( new double[] { 1696.3724170447979, 4180.403160693752, 117.54197364214336 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 9469.294444, -8226.764491, 3204.291107 } ), new Point( new double[] { 419.10957349136515, 4837.719300430445, 136.8057825553509 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 9469.28273, -8892.631243, 3204.284064 } ), new Point( new double[] { 414.84339554709277, 4172.678354893637, 131.63173672791973 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 10109.952614, -8226.479997, 3207.836146 } ), new Point( new double[] { 1053.9922003467032, 4843.204481890159, 123.80925600454901 } ) ) ); + pointMatches.add( new PointMatch( new Point( new double[] { 9789.611818, -8892.488996, 3206.056586 } ), new Point( new double[] { 733.8905377640995, 4173.987433223819, 126.7735051753456 } ) ) ); + + final double weight = 1. / pointMatches.size(); + for ( final PointMatch pointMatch : pointMatches ) + pointMatch.setWeights( new double[] { weight } ); + + Collections.shuffle( pointMatches ); + + final AffineModel3D model = new AffineModel3D(); + model.fit( pointMatches, 1e-3 ); + Assert.assertTrue( model.isInvertible ); + } +}