Fixes opposite pole coordinates

[benjaminmenetrier]

Hi @wdeconinck, I've found a bug in the Rotation implementation, and I think this branch brings a valid bugfix.

Let's take an example: if "south_pole" is specified in the configuration with coordinates (300, -10)

• with develop:

  • l. 129: npole_ = (120, 170)
  • l. 130: npole.lat() > 90 => true
  • l. 131: npole_.lon() = 300

  so npole_ = (300, 170) = (300, 10), which is not correct.

• with this branch:

  • l. 129: npole_ = (120, 10)
  • l. 130: npole.lon() < 0. => false

  so npole_ = (120, 10), which is correct.

[wdeconinck]

Thanks @benjaminmenetrier you're right! I will merge it once all tests have gone through.

[wdeconinck]

This has caused quite a number of tests to fail. Are we saying that the tests themselves are wrong, or that there is something wrong with this change?

[benjaminmenetrier]

I'm looking at it to fix the broken tests.

[benjaminmenetrier]

@wdeconinck: I've fixed the tests with two kinds of modification:

• update of references in many places for tests involving rotations (no problem of auto-consistency, all such tests pass),
• update of the equivalent function in src/tests/projection/test_rotation.cc, where the == about latitude being 90 or -90 was too strict (replaced with is_approximately_equal).

I haven't checked that new references are scientifically correct, but I have done one interesting test: with the bugfix, I manage to reproduce a regional rotated lon/lat grid that we use to download boundary conditions for the AROME model from MARS. So the bugfix makes ATLAS consistent with the MIR interpolation, which was not the case before.

[pmaciel]

Hi, I would like to also review this PR, as it is sensitive to what I'm doing now. I've made a (so far) bullet proof way to normalise point coordinates (including antipodes) and that's actually the different method that worries me, but I'm concerned about the MIR / Atlas comment, and the increase of tolerance

[benjaminmenetrier]

@pmaciel: sure, I would be glad to have your review on this!

[wdeconinck]

I think MIR usually (or always?) uses the South Pole to define the rotation, but still uses the same Rotation class. This is probably why this bug did not show up via MIR.

[benjaminmenetrier]

This is likely. I only noticed the bug when I started to specify the North Pole instead of the South Pole.

[pmaciel]

Hi,

I'm not sure the tolerances change are important here, because I didn't get to test on the CI matrix. The only machine I have fo rthat test is my centre laptoop which is hardly a reference. But I'm surprised the toleracne changed, as we're generally updating the reference results (@benjaminmenetrier)? Is this a multi-platform issue?

I've ported the upcoming eckit::geo PointLonLat::antipode method, which happily generate identical results -- or rather, the tests succeed in the same way, at least on my machine, from the latest version of this branch. So the changes I would propose at this moment aren't strictly necessary, as it wouldn't impose a change on the future adotopn of eckit::geo when the time comes.

If you're curious, the patch I have is attached diff.patch

[pmaciel]

I will add here that PROJ actually has an implementation for rotation, apparently similar to what we try to do, and the reference point is the North pole (and, no option for the South). We could use this as a reference as soon as we validate the rotation this way (which isn't an priority at the moment). For reference: https://proj.org/en/9.2/operations/projections/ob_tran.html