Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Transformation does not work as expected #23

Open
ollietreend opened this issue Apr 4, 2024 · 12 comments
Open

Transformation does not work as expected #23

ollietreend opened this issue Apr 4, 2024 · 12 comments

Comments

@ollietreend
Copy link
Contributor

I'm trying to convert coordinates from EPSG:27700 to EPSG:4326 using proj4rb.

I have the latest versions of the proj library (installed with brew install proj) and the proj4rb gem installed.

$ proj
Rel. 9.4.0, March 1st, 2024
usage: proj [-bdeEfiIlmorsStTvVwW [args]] [+opt[=arg] ...] [file ...]
$ gem list proj4rb

*** LOCAL GEMS ***

proj4rb (4.1.1)

✅ Using cs2cs command line tool

When I use the command line tool cs2cs provided by proj, I get the expected results:

$ echo "533498 181201" | cs2cs -f '%.10f' EPSG:27700 EPSG:4326
51.5139702631	-0.0775045667 0.0000000000

Easting 533498, Northing 181201 translates to Longitude -0.0775045667, Latitude 51.5139702631. This result is consistent with online converter tools such as the British Geological Survey's Coordinate converter tool.

❌ Using proj4rb

But when I try to perform the same action using proj4rb, I get a completely different result.

require "proj"

transform = Proj::Transformation.new("EPSG:27700", "EPSG:4326")
coordinate = Proj::Coordinate.new(e: 533498, n: 181201)
transformed_coordinate = transform.forward(coordinate)

latitude = transformed_coordinate.y # Latitude is stored in y attribute
longitude = transformed_coordinate.x # Longitude is stored in x attribute

puts "Latitude: #{latitude}"
puts "Longitude: #{longitude}"

When I execute this, I get the following result:

$ ruby proj4.rb
Latitude: 5.0e-324
Longitude: 2.19144499e-314

The returned lat/lng values are effectively zero. If I format them to 10 decimal places (the same as cs2cs output), they are both 0.0000000000.

I need to format them to 400 decimal places to actually see numbers come out:

puts "Latitude: #{"%.400f" % latitude}"
puts "Longitude: #{"%.400f" % longitude}"

which produces:

Latitude: 0.0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000049406564584124654417656879286822137236505980261432476442558568250067550727021
Longitude: 0.0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000219134797539321609555612014012053752743206045798867713541392373366376270284867081861393

What I've tried so far

I wondered if the proj4rb gem was using a different installation of the proj library – perhaps there could have been a broken installation hiding somewhere.

The cs2cs command is located here:

$ realpath $(which cs2cs)
/opt/homebrew/Cellar/proj/9.4.0/bin/cs2cs

And the proj4rb gem seems to be using the same installation path:

require "proj"
Proj::Context.current.database.path
# => "/opt/homebrew/Cellar/proj/9.4.0/share/proj/proj.db"

So from what I can see, they're using the exact same installation. That can't be the issue.

What am I doing wrong?

I'm really confused by this. @cfis do you have any idea what I'm doing wrong?

For what it's worth, I've also noticed that the proj4rb test suite doesn't pass on my machine. I get a whole bunch of failures across pretty much every test file.

@cfis
Copy link
Owner

cfis commented Apr 5, 2024

Works for me here if I do this:

coordinate = Proj::Coordinate.new(x: 533498, y: 181201)
puts coordinate
v0: 51.513968813644965, v1: -0.0775306317158091, v2: 0.0, v3: 0.0

If I do what you have:

coordinate = Proj::Coordinate.new(e: 533498, n: 181201)
puts coordinate
v0: 49.76680723514262, v1: -7.55715980690519, v2: 0.0, v3: 0.0

The reason is that code requires that you set :e, :n and :u. See https://github.com/cfis/proj4rb/blob/master/lib/proj/coordinate.rb#L62. If not, then no coordinates get set at all. So you are asking to transform an array of values of []. So you get weird numbers.

Really that big if statement should be rewritten where each parameter name (:x, :y, etc) is assigned an index of 0 through 4 and then mapped by index into the Api::PJ_COORD structure. Would be a fairly quick fix if you have any interest in making a PR.

@ollietreend
Copy link
Contributor Author

@cfis thanks for the quick response! 🙇🏻

I'm afraid I still can't get it to work:

require "proj"

transform = Proj::Transformation.new("EPSG:27700", "EPSG:4326")
coordinate = Proj::Coordinate.new(x: 533498, y: 181201)

puts coordinate
# v0: 533498.0, v1: 181201.0, v2: 0.0, v3: 0.0

puts transform.forward(coordinate)
# v0: 1.7733993367673755e+77, v1: 2.1601183626e-314, v2: 3.0154665634e-314, v3: 1.5e-323

What output do you get when you run the above code?

Based on the code snippets in your previous comment, it looks like you're getting 'good' results whereas I'm still getting something very odd. I'm expecting I should receive meaningful values for v0 and v1, and v2: 0.0, v3: 0.0 because I only provided two params.

@ollietreend
Copy link
Contributor Author

In fact, even more confusingly, I consistently get inconsistent results. It seems like the same transformation in a loop always results in the first value being different from the others... 🤔

It feels like something very bizarre is going on, but I'm struggling to know what it could be.

require "proj"

transform = Proj::Transformation.new("EPSG:27700", "EPSG:4326")
coordinate = Proj::Coordinate.new(x: 533498, y: 181201)

5.times { puts transform.forward(coordinate) }
# v0: 5.0e-324, v1: 5.2e-322, v2: 3.025252423e-314, v3: 1.5e-323
# v0: 3.0252521343e-314, v1: 2.1462617777e-314, v2: 3.025252423e-314, v3: 1.5e-323
# v0: 3.0252521343e-314, v1: 2.1462617777e-314, v2: 3.025252423e-314, v3: 1.5e-323
# v0: 3.0252521343e-314, v1: 2.1462617777e-314, v2: 3.025252423e-314, v3: 1.5e-323
# v0: 3.0252521343e-314, v1: 2.1462617777e-314, v2: 3.025252423e-314, v3: 1.5e-323

@ollietreend
Copy link
Contributor Author

I've even tried this in Docker to remove the possibility of something weird in my local setup. But I'm still getting the same results 😞

Here is what I've used, and the output I got:

https://gist.github.com/ollietreend/286d2e28806c53c6e022934fbccb0b66

@cfis
Copy link
Owner

cfis commented Apr 5, 2024

Running the code you provided works for me. I tested on both Windows and Linux (Fedora 39).

Also built and your ran your docker file (thanks for providing!):

docker run proj_test
v0: 533498.0, v1: 181201.0, v2: 0.0, v3: 0.0
v0: 51.513968813644965, v1: -0.0775306317158091, v2: 0.0, v3: 0.0

So sorry, not very helpful...but I have no idea why you are having issues.

Maybe a corrupt proj database? (but then the docker file should avoid that problem).

@cfis
Copy link
Owner

cfis commented Apr 9, 2024

Any luck tracking this down?

@ollietreend
Copy link
Contributor Author

Sadly, no 😞

The most bizarre thing is that the Docker image doesn't work, and yet it does for you. It makes me wonder if this could be down to the CPU architecture? I'm using a MacBook Pro with M3 processor, so everything's running on ARM. I don't know enough about the internal of this gem to know if that would make a difference or not.

For my particular use case, I need to batch transform about 20k coordinates as part of an ETL process. This runs once a day as a background task, and performance isn't critical.

Since the cs2cs command line tool (provided as part of the PROJ library) works correctly, I've switched to calling that from within my Rails app. I found that I can perform all 20k transformations in a loop using a single cs2cs process by writing to and reading from the IO stream sequentially. This performs surprisingly well for my needs – it does 20k transformations in ~0.2 seconds.

To give you an idea of what I've got working:

class CoordinateTransformer
  # Source coordinates are in British National Grid Easting/Northing format.
  # https://epsg.io/27700
  SOURCE_CRS = "EPSG:27700".freeze

  # We need coordinates in WGS84 format, the standard Latitude/Longitude
  # coordinate system used in GPS and online mapping tools.
  # https://epsg.io/4326
  TARGET_CRS = "EPSG:4326".freeze

  def initialize
    # cs2cs is provided as part of the PROJ library
    # Mac: brew install proj
    # Linux (Debian): apt-get install proj-bin
    # Usage: https://proj.org/en/9.4/apps/cs2cs.html
    @cs2cs = IO.popen(["cs2cs", "-d", "10", SOURCE_CRS, TARGET_CRS], "r+")
  end

  def transform(easting:, northing:)
    @cs2cs.write "#{easting} #{northing}\n"
    output = @cs2cs.gets.split(" ")
    { latitude: output[0], longitude: output[1] }
  end

  def close
    @cs2cs.close
  end
end

As I said, I'm quite impressed by the performance of this approach. I assume it's because it's cheap to write to the existing process rather than continually spawning new ones.

require "benchmark"

Benchmark.measure do
  ct = CoordinateTransformer.new
  20_000.times do |i|
    ct.transform(easting: i, northing: i)
  end
  ct.close
end

# =>   0.048043   0.050908   0.207849 (  0.200632)

It's frustrating I couldn't get the proj4rb gem to work 😢 thank you for your help digging into this with me though @cfis.

@cfis
Copy link
Owner

cfis commented Apr 9, 2024

Glad you have a solution.

The proj4 gem uses FFI, so does not include any native code (so nothing to compile) and thus is architecture independent. Of course the FFI gem is compiled against the ffi library. But I assume those are well tested and work fine on M3. But I don't have any other explanation for what you seeing.

Unfortunately I don't have a M3 (or M1/M2) computer so can't test for myself.

@ollietreend
Copy link
Contributor Author

@cfis interestingly, I've just tried running the same Docker image under x86-64/AMD64 and it produces different results compared to ARM 👀

❌ ARM64 (e.g. Apple Silicon) produces bad results

$ docker run --rm --platform linux/arm64 $(docker build -q --platform linux/arm64 .)
v0: 533498.0, v1: 181201.0, v2: 0.0, v3: 0.0
v0: 4.4748e-320, v1: 0.0, v2: 1.0e-323, v3: 1.39067016643389e-309

✅ AMD64 (e.g. Intel) produces good results

docker run --rm --platform linux/amd64 $(docker build -q --platform linux/amd64 .)
v0: 533498.0, v1: 181201.0, v2: 0.0, v3: 0.0
v0: 51.513968813644965, v1: -0.0775306317158091, v2: 0.0, v3: 0.0

To be clear: I'm running both these commands on the same MacBook with M3 processor. The only difference is that I've changed the --platform option.

By default, Docker runs under --platform linux/arm64 on an M3 (ARM) chip. But when you specify --platform linux/amd64 it'll run as AMD64 (presumably using Rosetta).

I'm not entirely sure what this means for proj4rb. I don't know enough about Ruby FFI to know whether this is a bug in proj4rb or FFI itself.

@ollietreend
Copy link
Contributor Author

In fact, I think this might just be a general macOS compatibility issue, rather than being specific to the newer Apple Silicon chips. I've forked this repo and added macOS GitHub-hosted runners to the CI workflow, and the test suite fails in the same way it does on my MacBook.

I haven't opened a Pull Request, but you can see what I changed here: master...ollietreend:proj4rb:macos-runners

You can see the results of this CI workflow here: https://github.com/ollietreend/proj4rb/actions/runs/8634089078

Notably, the macos-13 runner is Intel whereas the macos-14 runner is Apple Silicon (ARM). But both test suites fail with the same errors.

@cfis
Copy link
Owner

cfis commented Apr 10, 2024

Thanks for keeping at this. I do have an intel Mac, I will give that a try.

That looks like a great update to the test matrix, can you create an MR?

ollietreend added a commit to DFE-Digital/itt-mentor-services that referenced this issue Apr 11, 2024
I've updated the `Gias::CsvTransformer` to convert the Easting/Northing
values provided by GIAS into Latitude/Longitude coordinates that can be
used for location-based search.

Previously we were using a geocoding service to get the
Latitude/Longitude coordinates for schools based on their address.
However this requires N+1 external API calls (one for each school) which
is both slow and costly.

Turns out it's much more efficient to just make use of the geolocation
data GIAS already provides. It's a bit fiddly to perform the translation
because it requires an external library called [PROJ][1].  However once
that's installed it's pretty straightforward.

I would have liked to have used the gem [proj4rb][2] which is a Ruby
binding for PROJ. However, for [unexplained reasons][3], I couldn't
get that gem to work correctly. So instead I'm calling out to a command
line app that PROJ provides called [cs2cs][4], which converts coordinates
from one 'coordinate reference system' to another.

In my testing, this performs surprisingly well. Even though we need to
call out to an external process, it will easily convert the coordinates
of 20,000 schools in ~200ms. Not too shabby!

[1]: https://proj.org
[2]: https://github.com/cfis/proj4rb
[3]: cfis/proj4rb#23
[4]: https://proj.org/en/9.4/apps/cs2cs.html
@ollietreend
Copy link
Contributor Author

@cfis No problem – I've just opened a pull request. See #24.

Although I'm not sure it'll be possible to merge it, since the tests fail on macOS runners 🙈

ollietreend added a commit to DFE-Digital/itt-mentor-services that referenced this issue Apr 12, 2024
I've updated the `Gias::CsvTransformer` to convert the Easting/Northing
values provided by GIAS into Latitude/Longitude coordinates that can be
used for location-based search.

Previously we were using a geocoding service to get the
Latitude/Longitude coordinates for schools based on their address.
However this requires N+1 external API calls (one for each school) which
is both slow and costly.

Turns out it's much more efficient to just make use of the geolocation
data GIAS already provides. It's a bit fiddly to perform the translation
because it requires an external library called [PROJ][1].  However once
that's installed it's pretty straightforward.

I would have liked to have used the gem [proj4rb][2] which is a Ruby
binding for PROJ. However, for [unexplained reasons][3], I couldn't
get that gem to work correctly. So instead I'm calling out to a command
line app that PROJ provides called [cs2cs][4], which converts coordinates
from one 'coordinate reference system' to another.

In my testing, this performs surprisingly well. Even though we need to
call out to an external process, it will easily convert the coordinates
of 20,000 schools in ~200ms. Not too shabby!

[1]: https://proj.org
[2]: https://github.com/cfis/proj4rb
[3]: cfis/proj4rb#23
[4]: https://proj.org/en/9.4/apps/cs2cs.html
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants