New answers tagged geometry
1
It depends on what you want the heatmap to look like.
If you want areas where the polygons overlap to be brighter depending on how many polygons overlap that area, you can achieve that using transparency and/or blending modes. In the last two images below, notice that the polygons are darker in areas where they overlap. This is not really a heatmap, but it ...
1
You can try to calculate the WKT representation of each geometry and compare them.
To get the WKT representation, you type the following syntax in the field calculator.
With Python, you get the values of the 2 geometries and make a comparison with the operator '==' of the following form:
first_geometry == second_geometry
# Result
True or False
I can ...
2
There should be no line around the outside if layers are identical. Simple test with PostGIS shows a match with a theory.
select ST_AsText(ST_SymDifference(
ST_GeomFromText('POLYGON (( 460 500, 560 420, 460 420, 460 500 ))'),
ST_GeomFromText('POLYGON (( 460 500, 560 420, 460 420, 460 500 ))')));
"GEOMETRYCOLLECTION EMPTY"
It is very unlikely that data ...
0
Here is an example of how you can do this, I'm using the significant urban areas shapefile from the Australian Bureau of Statistics.
Firstly this is what my pnts variable looks like:
x y
<dbl> <dbl>
1 -34.92 138.62
2 -34.93 138.58
3 -34.95 138.52
4 -27.63 152.71
5 -27.57 153.01
6 -33.9 150.73
7 -33.92 150.99
And here is my code:
...
1
To answer yout question about how to set a spatialReference of a polygon layer:
arcpy.CreateFeatureclass_management(path, fc, "POLYGON")
# arcpy.AddField_management(fc, "NAME", "TEXT")
sr = arcpy.SpatialReference("the_name_of_file_of_the_crs")
arcpy.DefineProjection_management(fc, sr)
This is, what I use to define the spatial reference. It takes a name, ...
1
You can use a cross apply. I wouldn't advise running this code on a table with lots of records but in a case where you have a single geometry this works.
select * from yourTable cross apply dbo.STVectorize(YourGeomField)
If you have a table of multiple geometries then consider setting a variable
declare @test geometry
set @test = (select yourGeomField
...
2
It was my understanding that geodesic would (and should) always be larger.
I don't agree.
A polygon feature is defined by its vertices.
The area enclosed by the sides following geodesic lines is not the same as the area enclosed by the sides following straight lines in world sinusoidal projection.
It may be one greater than the other, depending on the ...
3
Geopandas is extremely useful and easy to use for this kind of thing:
http://geopandas.org/
You can load vector data with attributes (e.g. from geopackages, shapefiles, etc.) into ?(geo)pandas dataframes which allow really easy analysis of the data without messing around with ogr/gdal.
import geopandas as gpd
data = gpd.read_file("path.mygeopackage.gpkg")
...
7
You can use the gdal/ogr, fiona (built on gdal/ogr) or geopandas (built on fiona) python libraries.
Below is a fiona example:
import fiona
# No need to pass "layer='etc'" if there's only one layer
with fiona.open('test.gpkg', layer='layer_of_interest') as layer:
for feature in layer:
print(feature['geometry'])
Partial output for one record ...
4
if you have a geopackage file with extension .gpkg you can use ogr2ogr to translate it to a csv file with a geometry column like so:
ogr2ogr -f "CSV" MyLayer.csv MyLayer.gpkg -lco GEOMETRY=AS_WKT
If you use qgis you can do this manually using the convert format tool. You'll just need to input -lco GEOMETRY=AS_WKT in Additional creation options.
3
Run the "Fix Geometries" processing algorithm against your data, and then try once again computing the area on the newly-generated temporary layer.
If it all looks good, save it out to a new shp or gpkg.
1
Here is a code that uses centreline lable (https://github.com/ungarj/label_centerlines) library which can potentially solve your problem. Note this works pretty well with two polygons, but anything over three polygons becomes more difficult. The list of polygons consists of polygons in shapely Polygon format.
import geopandas as gdp
from shapely.geometry ...
2
To resume:
The reprojection seems to be the problem here. If your need is to have an area in meters, the simplest way, as you have 4326 geom, would be to use the geography type:
Add an indexes on your geoms casted as a geography (to speed up the intersection search):
CREATE INDEX idx_table_1_geog ON table_1 USING gist (CAST(geom AS geography));
CREATE ...
1
SELECT ST_AsText('010600...');
returns:
st_astext
-----------------------------------------------------------------------
MULTIPOLYGON ZM (((109.638182917 0.304720466000006 0 -1.79769e+308,...
Top 50 recent answers are included
