# Calculate ellipsoidal area for a projected layer in PyQGIS 3

I'd like to loop though each feature in a layer and calculate its ellipsoidal area using QgsDistanceArea() in PyQGIS 3.

Everything works fine on unprojected layers, but I get strange results when input layer is in projected CRS.

So far I have written this:

``````from qgis.core import *

def calculate_area(in_lyr_name, ellipsoid, units):
input_layer = QgsProject.instance().mapLayersByName(input_lyr_name)

area = QgsDistanceArea()
area.setEllipsoid(ellipsoid)

features = input_layer.getFeatures()
feature_area = 0
error_features_counter = 0

for i, feat in enumerate(features):
geom = feat.geometry()
polygon_area = 0
try:
if geom.isMultipart():
polygons = geom.asMultiPolygon()
for polygon in polygons:
polygon_area += area.measurePolygon(polygon)
else:
polygon = geom.asPolygon()
polygon_area = area.measurePolygon(polygon)

feature_area += polygon_area

except AttributeError:  # catches NoneType geometries (broken etc.)
error_features_counter += 1
pass

# calculated area is in sq. metres (see the "else" case)
if units == 'km²':
final_area = feature_area / 1e6
elif units == 'Ha':
final_area = feature_area / 10000
else:
final_area = feature_area
return final area, error_features_counter
``````

How do I get a constant result that does not depend on the presence of projection?

If I´m not mistaken, you will need to provide a transformation context and set the layers source CRS to the `QgsDistanceArea()` object accordingly. Try:

``````...

def calculate_area(in_lyr_name, ellipsoid, units):
input_layer = QgsProject.instance().mapLayersByName(input_lyr_name)
lyr_crs = input_layer.crs()

elps_crs = QgsCoordinateReferenceSystem()
elps_crs.createFromUserInput(ellipsoid)

trans_context = QgsCoordinateTransformContext()
trans_context.calculateDatumTransforms(lyr_crs, elps_crs)

area = QgsDistanceArea()
area.setEllipsoid(ellipsoid)
area.setSourceCrs(lyr_crs, trans_context)

...
``````

Seems to work for me (uses ellipsoid, returns sqm), but pyQGIS is not my expertise.

• Thank you! Exact method name in QGIS3 turned to be .createFromUserInput() – Basile Jun 7 '18 at 19:45