Opened 4 years ago

Closed 4 years ago

#3369 closed defect (fixed)

Followup to area calculation fix in changeset 71169

Reported by: ndawson Owned by: grass-dev@…
Priority: normal Milestone: 7.2.3
Component: Default Version: unspecified
Keywords: Cc:
CPU: Unspecified Platform: Unspecified

Description

While investigating the QGIS ticket https://issues.qgis.org/issues/16820 I came across changeset 71169 as a fix for a similar issue in grass.

When porting this fix across to QGIS, I had to modify the thres constant value to 0.7e-7, outside of the 1e-4 to 1e-7 range noted in the comment in area_poly1.c. I found this to be the maximum possible value which allowed an existing QGIS regression test to pass:

https://github.com/qgis/QGIS/blob/master/tests/src/core/testqgsdistancearea.cpp#L371 (as a result of https://issues.qgis.org/issues/14675)

So this ticket is just a heads-up that the grass threshold value of 1e-6 may not be suitable, and this may need lowering to fix a real-world area calculation issue.

Change History (9)

comment:1 in reply to:  description Changed 4 years ago by mmetz

Replying to ndawson:

While investigating the QGIS ticket https://issues.qgis.org/issues/16820 I came across changeset 71169 as a fix for a similar issue in grass.

When porting this fix across to QGIS, I had to modify the thres constant value to 0.7e-7, outside of the 1e-4 to 1e-7 range noted in the comment in area_poly1.c. I found this to be the maximum possible value which allowed an existing QGIS regression test to pass:

https://github.com/qgis/QGIS/blob/master/tests/src/core/testqgsdistancearea.cpp#L371 (as a result of https://issues.qgis.org/issues/14675)

About the first example (regression14675): what is SRSid 145L in EPSG or proj4 terms?

About the second example (regression16820): GRASS reports an area of 43.3280133198665 sqm in the original CRS, and an area of 43.2035658006178 sqm reprojected to latlon. The QGIS test has a reference of 43.183369 sqm, i.e. the projected GRASS result is a bit closer to the original than the QGIS reference. The threshold as used in GRASS seems to work fine here.

So this ticket is just a heads-up that the grass threshold value of 1e-6 may not be suitable, and this may need lowering to fix a real-world area calculation issue.

Waiting for SRS info on regression14675 to test in GRASS.

comment:2 Changed 4 years ago by ndawson

About the first example (regression14675): what is SRSid 145L in EPSG or proj4 terms?

Argh, missed that sorry. It's EPSG:2154.

About the second example (regression16820): GRASS reports an area of 43.3280133198665 sqm in the original CRS, and an area of 43.2035658006178 sqm reprojected to latlon. The QGIS test has a reference of 43.183369 sqm, i.e. the projected GRASS result is a bit closer to the original than the QGIS reference. The threshold as used in GRASS seems to work fine here.

Hmm - good point. To eliminate any discrepancies due to projection handling/some other factor I tried using the original threshold of 1e-6 in QGIS and got an area of 43.3280029296875.

I can't find a reliable threshold which satisfies both these requirements. Ubuntu 16.04 with a threshold of 0.8e-7 does, but the same threshold on Fedora 25 fails both tests.

comment:3 in reply to:  2 ; Changed 4 years ago by mmetz

Replying to ndawson:

About the first example (regression14675): what is SRSid 145L in EPSG or proj4 terms?

Argh, missed that sorry. It's EPSG:2154.

Thanks! Here are the results

  • QGIS reference: 0.83301
  • GRASS in EPSG:2154: 0.86121968362156
  • GRASS in latlon: 3.3238796585328

I had to lower the threshold to 1e-8 in order to get

  • GRASS in latlon: 0.86134805688084

About the second example (regression16820): GRASS reports an area of 43.3280133198665 sqm in the original CRS, and an area of 43.2035658006178 sqm reprojected to latlon. The QGIS test has a reference of 43.183369 sqm, i.e. the projected GRASS result is a bit closer to the original than the QGIS reference. The threshold as used in GRASS seems to work fine here.

Hmm - good point. To eliminate any discrepancies due to projection handling/some other factor I tried using the original threshold of 1e-6 in QGIS and got an area of 43.3280029296875.

I can't find a reliable threshold which satisfies both these requirements. Ubuntu 16.04 with a threshold of 0.8e-7 does, but the same threshold on Fedora 25 fails both tests.

This threshold method is not perfect. I found that the threshold must be larger, the closer it gets to the equator, and smaller, the closer it gets to the poles.

comment:4 in reply to:  3 Changed 4 years ago by mmetz

Replying to mmetz:

Replying to ndawson:

About the first example (regression14675): what is SRSid 145L in EPSG or proj4 terms?

Argh, missed that sorry. It's EPSG:2154.

Thanks! Here are the results

  • QGIS reference: 0.83301
  • GRASS in EPSG:2154: 0.86121968362156
  • GRASS in latlon: 3.3238796585328

I had to lower the threshold to 1e-8 in order to get

  • GRASS in latlon: 0.86134805688084

After I applied the advice in my own comment, a threshold of 1e-6 works fine

  • GRASS in latlon: 0.861746563835888

If different latitudes are regarded as nearly identical (fabs(dy) <= thresh), then the average of the two latitudes must be used, not one of them.

This correction is now also in relbr72 (r71260).

comment:6 Changed 4 years ago by martinl

Milestone: 7.2.2

@mmetz Are you planning backport also to relbr70 or the ticket can be simply closed? Thanks.

comment:7 in reply to:  6 Changed 4 years ago by mmetz

Replying to ndawson:

Thanks - I've ported r71260 to QGIS in https://github.com/qgis/QGIS/commit/773b2e7f9e3744a0420f55a03d85fd2adebe14ee

Nice! I'm not sure if a constant threshold is a good idea, but at least for the test cases and synthetic data from equator to pole the results look reasonable.

Replying to martinl:

@mmetz Are you planning backport also to relbr70 or the ticket can be simply closed? Thanks.

Backported to relbr70 in r71261.

comment:8 Changed 4 years ago by neteler

Milestone: 7.2.27.2.3

Ticket retargeted after milestone closed

comment:9 Changed 4 years ago by martinl

Resolution: fixed
Status: newclosed
Note: See TracTickets for help on using tickets.