Skip to content

Instantly share code, notes, and snippets.

@agiudiceandrea

agiudiceandrea/test_geodesic_3.17.py

Last active January 6, 2021 16:25
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save agiudiceandrea/6ce1fa7915e781f355ba624817cfea12 to your computer and use it in GitHub Desktop.
Save agiudiceandrea/6ce1fa7915e781f355ba624817cfea12 to your computer and use it in GitHub Desktop.
Download ZIP
Raw
test_geodesic_3.17.py
QGIS 3.17.0-Master
GEOS 3.8.1-CAPI-1.13.3
PROJ Rel. 8.0.0, March 1st, 2021
SpatiaLite 4.3.0
md = QgsProviderRegistry.instance().providerMetadata('spatialite')
db=md.createConnection(QgsDataSourceUri('dbname=\'C:/TMP/test.sqlite\' table=\'test\'(geometry) sql=').uri(),{})
db.executeSql('select spatialite_version();')
[['4.3.0']]
db.executeSql('select proj4_version();')
[['Rel. 6.3.2, May 1st, 2020']]
db.executeSql('select geos_version();')
[['3.8.1-CAPI-1.13.3']]
db.executeSql('select lwgeom_version();')
[['3']]
db.executeSql('select st_area(st_geomfromtext(\'POLYGON ((20.13293641 59.95688345, 26.94617837 60.47397663, 29.74782155 62.56499443, 27.45254202 68.70650340, 23.75771765 68.24937206, 25.42698984 65.27444593, 21.51545237 63.10353609, 21.40562760 61.12318104, 19.41123592 60.40477513, 20.13293641 59.95688345))\') ,true);')
[[251199344354.4308]]
db.executeSql('select st_area(st_geomfromtext(\'POLYGON ((20.13293641 59.95688345, 26.94617837 60.47397663, 29.74782155 62.56499443, 27.45254202 68.70650340, 23.75771765 68.24937206, 25.42698984 65.27444593, 21.51545237 63.10353609, 21.40562760 61.12318104, 19.41123592 60.40477513, 20.13293641 59.95688345))\') ,false);')
[[249631004700.70947]]
db.executeSql('select st_area(st_transform(st_geomfromtext(\'POLYGON ((20.13293641 59.95688345, 26.94617837 60.47397663, 29.74782155 62.56499443, 27.45254202 68.70650340, 23.75771765 68.24937206, 25.42698984 65.27444593, 21.51545237 63.10353609, 21.40562760 61.12318104, 19.41123592 60.40477513, 20.13293641 59.95688345))\', 4326), 32635));')
[[251038084443.76132]]
db.executeSql('select st_area(st_geomfromtext(\'POLYGON ((20.13293641 59.95688345, 26.94617837 60.47397663, 29.74782155 62.56499443, 27.45254202 68.70650340, 23.75771765 68.24937206, 25.42698984 65.27444593, 21.51545237 63.10353609, 21.40562760 61.12318104, 19.41123592 60.40477513, 20.13293641 59.95688345))\'));')
[[45.25847152075945]]
geom=QgsGeometry.fromWkt("POLYGON ((20.13293641 59.95688345, 26.94617837 60.47397663, 29.74782155 62.56499443, 27.45254202 68.70650340, 23.75771765 68.24937206, 25.42698984 65.27444593, 21.51545237 63.10353609, 21.40562760 61.12318104, 19.41123592 60.40477513, 20.13293641 59.95688345))")
d = QgsDistanceArea()
d.measureArea(geom)
45.258471520759485
d = QgsDistanceArea()
d.setEllipsoid('WGS84')
True
d.measureArea(geom)
249566957499.7546
from shapely import wkt
from pyproj import Geod
pyproj info:
pyproj: 2.6.1.post1
PROJ: 7.0.1
poly = wkt.loads('POLYGON ((20.13293641 59.95688345, 26.94617837 60.47397663, 29.74782155 62.56499443, 27.45254202 68.70650340, 23.75771765 68.24937206, 25.42698984 65.27444593, 21.51545237 63.10353609, 21.40562760 61.12318104, 19.41123592 60.40477513, 20.13293641 59.95688345))')
poly.area
45.25847152075945
geod = Geod(ellps="WGS84")
abs(geod.geometry_area_perimeter(poly)[0])
251199344354.4309
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment