distancelatitude-longitudegeographysqlgeographynettopologysuite

Using NetTopologySuite, how can I convert meters to degrees?


I have no language preference for an answer; either VB.NET or C# will be fine.

I've seen this problem discussed elsewhere here on SO, for example here and here, but I've not found solutions that employ the NetToplogySuite and ProjNet libraries. Further, these other Q&As deal with two known points—I only have one. Therefore, I believe this question to not be a duplicate.

We can easily convert degrees to meters using the sample code found in Microsoft's documentation here. This works well.

However, in my case I'm querying on a database column. One location is known, but the others aren't. I need to be able to make a call like this:

oZipCodes = oDb.ZipCodes.Where(Function(Unknown) Unknown.Location.Distance(oKnown.Location) < 22000)

EF Core successfully translates this into workable T-SQL. The resulting query, run against the lat/lon coordinates of 33,000 U.S. Zip Codes, is quite fast. The trouble, though, is that the response returned by Distance() is in degrees. The ProjectTo() calculation converts degrees to meters as expected, but EF Core fails when it tries to translate this call:

oZipCodes = oDb.ZipCodes.Where(Function(Unknown) Unknown.Location.ProjectTo(2855).Distance(oKnown.Location.ProjectTo(2855)) < 22000)

The LINQ expression 'DbSet<ZipCode>.Where(z => (Geometry)z.Location.ProjectTo(_2855).Distance(__ProjectTo_0) < 22000)' could not be translated.

I took a brief look at user-defined function mapping, but that's not a good fit. Any such internalized function would be called 33,000 times for each query, which of course would amount to a significant performance hit.

So it seems that the only solution is to convert the maximum number of meters to degrees before making our .Location.Distance() call. In other words, the reverse of what's described in Microsoft's documentation linked above. Best to convert it once for each call, rather than 33,000 times.

How can I use NetTopologySuite and ProjNet to 'reverse' convert degrees into meters, based on geographical lat/lon coordinates?

(All that said, if there's an easier—and still fast—way to get a list of lat/lon coordinates within a given distance from a given point I'd love to hear about it!)


Solution

  • I ended up using a stored procedure (with miles instead of meters):

    Imports Microsoft.EntityFrameworkCore.Migrations
    
    Namespace Global.Persistence.Migrations.MsSql
      ''' <inheritdoc />
      Partial Public Class _002
        Inherits Migration
    
        ''' <inheritdoc />
        Protected Overrides Sub Up(migrationBuilder As MigrationBuilder)
          migrationBuilder.Sql("
            CREATE PROCEDURE GetZipCodesWithinRadius
                @Longitude FLOAT,
                @Latitude FLOAT,
                @Radius FLOAT
            AS
            BEGIN
                DECLARE @KnownPoint geography = geography::Point(@Latitude, @Longitude, 4326);
                SELECT *
                FROM ZipCodes
                WHERE @KnownPoint.STDistance(Location) / 1609.34 <= @Radius;
            END")
        End Sub
    
        ''' <inheritdoc />
        Protected Overrides Sub Down(migrationBuilder As MigrationBuilder)
          migrationBuilder.Sql("DROP PROCEDURE GetZipCodesWithinRadius")
        End Sub
      End Class
    End Namespace
    

    It's called like this:

    Dim iMaxMiles As Integer
    Dim oLongitude As SqlParameter
    Dim oLatitude As SqlParameter
    Dim oRadius As SqlParameter
    
    iMaxMiles = 100
    oLongitude = New SqlParameter("@Longitude", oKnown.Location.Longitude)
    oLatitude = New SqlParameter("@Latitude", oKnown.Location.Latitude)
    oRadius = New SqlParameter("@Radius", iMaxMiles)
    
    oZipCodes = oDb.ZipCodes.FromSql($"
      GetZipCodesWithinRadius
        @Longitude = {oLongitude},
        @Latitude = {oLatitude},
        @Radius = {oRadius}
      ")
    
    oZipCodes.ForEach(Sub(ZipCode)
                        Console.WriteLine($"{ZipCode.Code}: {ZipCode.City.Name}")
                      End Sub)
    

    This neatly bypasses all of the ProjNet and SRID complexities (which ultimately are unnecessary for this simple task).

    It works perfectly.