diff --git a/Changelog.md b/Changelog.md index 84fd62a30..2c5a1bd88 100644 --- a/Changelog.md +++ b/Changelog.md @@ -25,6 +25,7 @@ All notable changes to this project will be documented in this file. - Removed obsolete methods\properties (#797) ### Fixed +- inconsistent use of affine coefficients (#822) - Fixed the shift in x-coordinate when reprojecting from WGS84 to LAEA (#815) - Fixed LAEA reprojected y coordinate that resulted in n.def (#813) - ShapeReader skipping one entry when switching the page (#774) diff --git a/Source/DotSpatial.Data.Rasters.GdalExtension/GdalImage.cs b/Source/DotSpatial.Data.Rasters.GdalExtension/GdalImage.cs index 5c99e9213..044c0f672 100644 --- a/Source/DotSpatial.Data.Rasters.GdalExtension/GdalImage.cs +++ b/Source/DotSpatial.Data.Rasters.GdalExtension/GdalImage.cs @@ -124,6 +124,7 @@ private void ReadHeader() NumBands = _dataset.RasterCount; var test = new double[6]; _dataset.GetGeoTransform(test); + test = (new AffineTransform(test)).TransfromToCorner(0.5, 0.5);//shift origin by half a cell ProjectionString = _dataset.GetProjection(); Bounds = new RasterBounds(Height, Width, test); WorldFile.Affine = test; diff --git a/Source/DotSpatial.Data.Rasters.GdalExtension/GdalImageSource.cs b/Source/DotSpatial.Data.Rasters.GdalExtension/GdalImageSource.cs index 6e688e1c6..47d82a45a 100644 --- a/Source/DotSpatial.Data.Rasters.GdalExtension/GdalImageSource.cs +++ b/Source/DotSpatial.Data.Rasters.GdalExtension/GdalImageSource.cs @@ -184,6 +184,7 @@ private void ReadHeader() double[] test = new double[6]; _dataset.GetGeoTransform(test); + test = (new AffineTransform(test)).TransfromToCorner(0.5, 0.5);//shift origin by half a cell Bounds = new RasterBounds(_dataset.RasterYSize, _dataset.RasterXSize, test); _red = _dataset.GetRasterBand(1); _numOverviews = _red.GetOverviewCount(); diff --git a/Source/DotSpatial.Data.Rasters.GdalExtension/GdalRaster.cs b/Source/DotSpatial.Data.Rasters.GdalExtension/GdalRaster.cs index 72e9303b3..bec76b016 100644 --- a/Source/DotSpatial.Data.Rasters.GdalExtension/GdalRaster.cs +++ b/Source/DotSpatial.Data.Rasters.GdalExtension/GdalRaster.cs @@ -582,6 +582,9 @@ private void ReadHeader() } double[] affine = new double[6]; _dataset.GetGeoTransform(affine); + // in gdal (row,col) coordinates are defined relative to the top-left corner of the top-left cell + // shift them by half a cell to give coordinates relative to the center of the top-left cell + affine = (new AffineTransform(affine)).TransfromToCorner(0.5, 0.5); ProjectionString = projString; Bounds = new RasterBounds(base.NumRows, base.NumColumns, affine); PixelSpace = Marshal.SizeOf(typeof(T)); diff --git a/Source/DotSpatial.Data.Rasters.GdalExtension/GdalTiledImage.cs b/Source/DotSpatial.Data.Rasters.GdalExtension/GdalTiledImage.cs index 727f2ae25..a3906d09b 100644 --- a/Source/DotSpatial.Data.Rasters.GdalExtension/GdalTiledImage.cs +++ b/Source/DotSpatial.Data.Rasters.GdalExtension/GdalTiledImage.cs @@ -68,6 +68,7 @@ private void ReadHeader() WorldFile = new WorldFile { Affine = new double[6] }; double[] test = new double[6]; _dataset.GetGeoTransform(test); + test = (new AffineTransform(test)).TransfromToCorner(0.5, 0.5);//shift origin by half a cell Bounds = new RasterBounds(Height, Width, test); WorldFile.Affine = test; Close(); diff --git a/Source/DotSpatial.Data.Tests/AffineCoefficientsTests.cs b/Source/DotSpatial.Data.Tests/AffineCoefficientsTests.cs new file mode 100644 index 000000000..d30a1d68b --- /dev/null +++ b/Source/DotSpatial.Data.Tests/AffineCoefficientsTests.cs @@ -0,0 +1,71 @@ +using System; +using System.IO; +using GeoAPI.Geometries; +using DotSpatial.Data.Rasters.GdalExtension; +using DotSpatial.Tests.Common; +using NUnit.Framework; +using TestClass = NUnit.Framework.TestFixtureAttribute; +using TestMethod = NUnit.Framework.TestAttribute; +using TestCleanup = NUnit.Framework.TearDownAttribute; +using TestInitialize = NUnit.Framework.SetUpAttribute; +using ClassCleanup = NUnit.Framework.TestFixtureTearDownAttribute; +using ClassInitialize = NUnit.Framework.TestFixtureSetUpAttribute; + +namespace DotSpatial.Data.Tests +{ + /// + ///This is a class testing the implementation of affine coefficients. + /// + [TestClass()] + public class AffineCoefficientsTests + { + /// + ///A test for affine coefficients in AffineTransform: generate random points with + ///known (row, column) and check whether AffineTransform.ProjToCell returns the correct values. + ///This test fails when using Math.Floor in AffineTransform.ProjToCell. + /// + [TestMethod] + public void AffineTransformTest() + { + var c = new double[6]; + c[0] = -179.75; // x-coordinate of the *center* of the top-left cell + c[1] = 0.5; // column width + c[2] = 0.0; // rotation/skew term + c[3] = 89.75; // y-coordinate of the *center* of the top-left cell + c[4] = 0.0; // rotation/skew term + c[5] = -0.5; // negative row height + var at = new AffineTransform(c); + + var rand = new Random(); + + for (int row = 0; row < 9; row++) + { + for (int col = 0; col < 9; col++) + { + Coordinate corner = at.CellTopLeft_ToProj(row, col); + double frac_col = rand.NextDouble(); + double frac_row = rand.NextDouble(); + double dX = c[1] * frac_col + c[2] * frac_row; + double dY = c[4] * frac_col + c[5] * frac_row; + Coordinate point = new Coordinate(corner.X + dX, corner.Y + dY); + Assert.AreEqual(at.ProjToCell(point), new RcIndex(row, col)); + } + } + } + + /// + ///A test for affine coefficients in GdalRaster: use gdal to load a raster file with + ///known geolocation and test the location of the center of the first grid cell. + ///Without a half-cell shift applied to the origin of the grid in GdalRaster.ReadHeader, + ///this test fails because affine coefficients are defined differently in gdal and dotspatial. + /// + [TestMethod] + public void GdalRasterTest() + { + var rp = new GdalRasterProvider(); + var raster = rp.Open(@"Data\Grids\sample_geotiff.tif"); + var at = new AffineTransform(raster.Bounds.AffineCoefficients); + Assert.AreEqual(at.CellCenter_ToProj(0, 0), new Coordinate(-179.9499969, 89.9499969));// correct location from sample_geotiff.tfw + } + } +} \ No newline at end of file diff --git a/Source/DotSpatial.Data.Tests/Data/Grids/sample_geotiff.tfw b/Source/DotSpatial.Data.Tests/Data/Grids/sample_geotiff.tfw new file mode 100644 index 000000000..35e75cf52 --- /dev/null +++ b/Source/DotSpatial.Data.Tests/Data/Grids/sample_geotiff.tfw @@ -0,0 +1,6 @@ + 0.1000000 + 0.0000000 + 0.0000000 + -0.1000000 + -179.9499969 + 89.9499969 diff --git a/Source/DotSpatial.Data.Tests/Data/Grids/sample_geotiff.tif b/Source/DotSpatial.Data.Tests/Data/Grids/sample_geotiff.tif new file mode 100644 index 000000000..617a73b0a Binary files /dev/null and b/Source/DotSpatial.Data.Tests/Data/Grids/sample_geotiff.tif differ diff --git a/Source/DotSpatial.Data.Tests/DotSpatial.Data.Tests.csproj b/Source/DotSpatial.Data.Tests/DotSpatial.Data.Tests.csproj index 7db32c6f9..c22e07041 100644 --- a/Source/DotSpatial.Data.Tests/DotSpatial.Data.Tests.csproj +++ b/Source/DotSpatial.Data.Tests/DotSpatial.Data.Tests.csproj @@ -86,6 +86,7 @@ + diff --git a/Source/DotSpatial.Data/AffineTransform.cs b/Source/DotSpatial.Data/AffineTransform.cs index 37e5d8a3f..484a310c0 100644 --- a/Source/DotSpatial.Data/AffineTransform.cs +++ b/Source/DotSpatial.Data/AffineTransform.cs @@ -119,10 +119,10 @@ public Coordinate CellBottomRight_ToProj(int row, int column) /// /// Given the row and column, this returns new affine coefficients transformed to that cell /// - /// The integer row index from 0 to numRows - 1 - /// The integer column index from 0 to numColumns - 1 - /// Transoformed affine coefficients - public double[] TransfromToCorner(int startColumn, int startRow) + /// The row index from 0 to numRows - 1 + /// The column index from 0 to numColumns - 1 + /// Transformed affine coefficients + public double[] TransfromToCorner(double startColumn, double startRow) { // X = [0] + [1] * column + [2] * row; // Y = [3] + [4] * column + [5] * row; @@ -165,8 +165,8 @@ public RcIndex ProjToCell(Coordinate location) rw = (c[3] + c[4] * location.X / c[1] - c[4] * c[0] / c[1] - location.Y) / div; cl = (location.X - c[2] * rw - c[0]) / c[1]; } - var iRow = (int)Math.Floor(rw); - var iCol = (int)Math.Floor(cl); + var iRow = (int)Math.Round(rw); + var iCol = (int)Math.Round(cl); return new RcIndex(iRow, iCol); } diff --git a/Source/DotSpatial.Plugins.MapWindowProjectFileCompatibility.Tests/DotSpatial.Plugins.MapWindowProjectFileCompatibility.Tests.csproj b/Source/DotSpatial.Plugins.MapWindowProjectFileCompatibility.Tests/DotSpatial.Plugins.MapWindowProjectFileCompatibility.Tests.csproj index ad04a4bb4..f8c3e87c6 100644 --- a/Source/DotSpatial.Plugins.MapWindowProjectFileCompatibility.Tests/DotSpatial.Plugins.MapWindowProjectFileCompatibility.Tests.csproj +++ b/Source/DotSpatial.Plugins.MapWindowProjectFileCompatibility.Tests/DotSpatial.Plugins.MapWindowProjectFileCompatibility.Tests.csproj @@ -73,6 +73,9 @@ + + + diff --git a/Source/DotSpatial.Positioning.Tests/DotSpatial.Positioning.Tests.csproj b/Source/DotSpatial.Positioning.Tests/DotSpatial.Positioning.Tests.csproj index b8a82a43c..4c4cfa6ad 100644 --- a/Source/DotSpatial.Positioning.Tests/DotSpatial.Positioning.Tests.csproj +++ b/Source/DotSpatial.Positioning.Tests/DotSpatial.Positioning.Tests.csproj @@ -62,6 +62,9 @@ DotSpatial.Positioning + + +