Polygon area calculation discrepancies using PyQGIS

Question

I use the following code to calculate the areas of a multi-polygon shapefile. However, to validate these results, I calculate the polygon’ areas directly in table attribute by creating a new field and applying the formula: $area. I found discrepancies in calculated areas by both methods.

Does the QgsDistanceArea class consider the layer CRS (in my case I use EPSG:32632) and $area does not, that means the geographic WSG84 EPSG:4326 is considered by default?

def area_poly(layer):
  L = []
  features = layer.getFeatures()
  d = QgsDistanceArea()
  d.setSourceCrs(layer.crs(), QgsProject.instance().transformContext())    
  d.setEllipsoid(QgsProject.instance().ellipsoid())
  for i, element in enumerate(features):
    for poly in element.geometry().asMultiPolygon():
        m = d.measurePolygon(poly[0])
        L.append((i, round(m,1))
  return L

Answer 1

It's the opposite - $area will consider both the layer CRS and project ellipsoid setting. But your code does not, as creating a default constructed QgsDistanceArea is not sufficient to get ellipsoid based measurements, and will currently be returning euclidean areas only.

You should revise your code to (for QGIS 3.x):

d = QgsDistanceArea()
d.setSourceCrs(layer.crs(), QgsProject.instance().transformContext())
d.setEllipsoid(QgsProject.instance().ellipsoid())