I have a shapefile with about 2000 polygons. I am using QGIS and I need to write a standalone Python script which opens the shapefile and calculates the area of all polygons separately and outputs the calculated area in a column "Area" in square kilometers. I tried to do it using the example in Calculate area in square meters from degrees?
eligibleLayer = QgsVectorLayer("/home/usr/Desktop/eligible_areas.shp", "eligible_areas", "ogr") if not eligibleLayer.isValid(): print "Eligible areas layer failed to load!" lProvider = eligibleLayer.dataProvider() lProvider.addAttributes( [ QgsField("Area",QVariant.Double) ] ) area = 0 for gFeat in eligibleLayer.getFeatures(): calculator = QgsDistanceArea() calculator.setEllipsoid('WGS84') calculator.setEllipsoidalMode(True) geom = gFeat.geometry() landArea = gFeat.attributes()[gProvider.fieldNameIndex('Area')].toDouble() if geom.isMultipart() is False: # if only simple polygon, calculate only for this polyg = geom.asPolygon() #transform to list of points if len(polyg) > 0: area = calculator.measurePolygon(polyg) else: #is Multipart multi = geom.asMultiPolygon() for polyg in multi: area = area + calculator.measurePolygon(polyg)
The script calculates the variable area but I am unsure if it is doing the right thing since I checked with the geometry tools and I get something completely different. Furthermore, I am unsure about how to update the above code such that I save the information in a newly created attribute column called "Area".