Skip to content

Difference between area calculations by QGIS and sf package

4 messages · Rik Ferreira, Roger Bivand

Original
Rik Ferreira
# 
Dear list members,

Why does sf::st_area returns a different area than the output of QGIS $area
for the same geometry?

A reproducible example (state of Minas Gerais, Brazil):

library(sf)
sf::st_area: 588358849314 [m^2]
QGIS $area: 586521037474,860 m?

The layer used in QGIS was exported with st_write() in GeoJSON format.

The error is: 1837811840 m?.

Thank you for the attention.

Att.
Roger Bivand
#
On Mon, 27 Sep 2021, Rik Ferreira wrote:

            

      
      
      
    

        
[1] TRUE

            

      
      
      
    

        
[1] TRUE

            

      
      
      
    

        
Spherical geometry (s2) switched off

            

      
      
      
    

        
586520883040 [m^2]

            

      
      
      
    

        
Spherical geometry (s2) switched on

            

      
      
      
    

        
588358849314 [m^2]

            

      
      
      
    

        
1837966275 [m^2]

See https://r-spatial.org/r/2020/06/17/s2.html

When sf_use_s2() is FALSE, see ?st_area ("ellipsoidal distances are 
computed using st_geod_distance which uses function 'geod_inverse' from 
GeographicLib (part of PROJ)"). You'd have to establish what code QGIS 
uses to know why it and GeographicLib differ.

Hope this clarifies

Roger

            

      
      
      
    

        

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Rik Ferreira
    

    

      
      #
    

  


  

      
Dear Roger, thanks for the explanation.

The error between QGIS and sf drops a lot with sf_use_s2(FALSE) but it's
still bigger than 0. I believe it is due to different parameters passed
because GEOS and PROJ are QGIS dependencies.

I will run more tests to find out what may be causing this issue but the
problem is solved for now.

Thank you!

Em seg., 27 de set. de 2021 ?s 12:32, Roger Bivand <Roger.Bivand at nhh.no>
escreveu:

            

      
      
      
    

        

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Roger Bivand
    

    

      
      #
    

  


  

      
On Mon, 27 Sep 2021, Rik Ferreira wrote:

        

            

      
      
      
    

        
No, the differences have nothing to do with GEOS, which only handles 
planar geometries, or PROJ, which expresses the coordinate reference 
system (here "EPSG:4674", a GEOGCRS in decimal degrees). The differences 
are coming from difference ways of estimating area on the sphere (s2) or a 
chosen ellipsoid (lwgeom, QGIS) all using different approaches. Try 
projecting to UTM zone 23; then GEOS would be in play, because the 
geometry is planar.

Hope this clarifies,

Roger