New answers tagged python
0
I think this question is answered by A.A. Changing local saved files is slow with changeAttributeValues but on a layer in memory it's faster. I compared different ways of changing attributes of a local layer.
Methods
changeAttributeValues
change attribute values of Shapefile with dbf package
change attribute values by creating new layer with parts of old ...
0
Returning after some time to share my findings. The reproject method cannot be used in this manner (or it can, but will not result in a raster translation as I was expecting). I've tested this again with a mockup raster. The following code should run as-is in an ipython notebook or as a standalone script:
This first part creates a mock single-band raster at ...
1
To make a map with more than one layers as is your case (the base worl_map layer and the overlay map_bays_gulf) there are at least two posible options as it's explained in the documentation of GeoPandas.
Assign the output of the plot method to a variable and use it as the ax parameter of the second call to plot
Create a Axes object calling plt.subplots(), ...
0
Found the Answer! Basically, when been run through the batch file it can't find the GDAL data file unless it's in the System Variables. Intellij must do something in the background that makes it work.
How to set the GDAL Environemt.
2
The issue is that reading the CIM returns an array of None values in the queryFields property of the data connection. When setting the definition back on the CIM with the array of None values, the data connection on the layer breaks.
A solution is to remove the queryFields attribute from the CIM and then setting the definition.
layerCIM = layer.getDefinition(...
0
I found a work around for this that finally worked to get a token. I believe part of the problem is how outdated and wonky ESRI's sample codes were, which is what I was working off of. The code that ended up working is below.
Basically, I just used a completely different library, the requests library, and was able to get the token without issue.
# Create ...
2
The QgsFontMarkerSymbolLayer must be added as a subSymbol to your QgsCentroidFillSymbolLayer.
I just needed to add 4 lines to your code to obtain something like this :
Here is the lines I added :
symbol_parkcent = QgsCentroidFillSymbolLayer.create()
# Make the centroid transparent
symbol_parkcent.setColor(QColor("transparent"))
#
colorP = QtGui....
1
I think you need to use the UpdateDefinition URL parameter and a JSON object.
Try this link for more information
https://developers.arcgis.com/rest/services-reference/online/update-definition-feature-service-.htm
Also the format of parameters being pass in the post() method need to be changed, such that what is being passed in this method in correct URL ...
0
you can approach the problem in a easy way.
made buffer from the the point (5ft diameter)
trasform buffer in line
intersect line buffer with your conductor line to have a point at 5ft from the end, done.
take care of attributes during the operation.
0
As suggested by Kevin Reid I filed a bug report and it turns out the task was failing internally with the following message:
The member bindings allUsers and allAuthenticatedUsers are not allowed since public access prevention is enforced.
Which means that the project didn't allow public access to the bucket. Changing the project permissions solved the ...
1
I found the mistake. The PrefixPath for QgsApplication was wrong. The PrefixPath that is working for me under Ubuntu is the following:
QgsApplication.setPrefixPath("/usr", True)
4
Get a handle to the driver, then check its metadata:
from osgeo import gdal
raster_drivers = []
vector_drivers = []
for i in range(gdal.GetDriverCount()):
drv = gdal.GetDriver(i)
md = drv.GetMetadata_Dict()
d = [drv.ShortName, drv.LongName, drv.GetMetadataItem(gdal.DMD_EXTENSIONS)]
if 'DCAP_RASTER' in md:
raster_drivers....
4
Here's how to use gdal.VectorTranslate, note that some ogr2ogr options map to slightly different arguments, but they should be easily figured out:
from osgeo import gdal
import numpy as np
import pathlib
gdal.UseExceptions()
## example dummy ids
id1 = np.linspace(start=1, stop=1, num=1, endpoint=True, dtype=int)
id2 = np.linspace(start=1, stop=1, num=1, ...
0
I was seeing the same issue. The environment we are attempting to work in has a Host Intrusion Prevention System(HIPS) rule that is preventing anything from being executed from the TEMP or TMP System directory (C:\windows\Temp). A firewall or antivirus may also cause the same issue.
When we disable HIPS the environment switch succeeds.
The event being ...
1
In addition to the answer of CodeBard:
When Pandas read a csv file, the columns are strings or numbers, never shapely geometries. A Pandas DataFrame with a geometry column (string) is not a GeoPandas GeoDataFrame (geometry column, shapely geometry)
import pandas as pd
df = pd.read_csv("lat_long.csv")
df
lat lon geometry
0 ...
1
You have to include 'OUTPUT' : 'TEMPORARY_OUTPUT' in parameters.
parameters = {
'ALPHA_BAND': False,
'CROP_TO_CUTLINE': True,
'DATA_TYPE': 0,
'INPUT': input,
'KEEP_RESOLUTION': True,
'MASK': mask,
'MULTITHREADING': False,
'OUTPUT': 'memory:',
'SOURCE_CRS': mask.crs(),
'TARGET_CRS': mask.crs(),
'NO_DATA': -9999,
...
1
The close method does not have feature input and your code, as it stands, has no stored knowledge of any feature.
Below is the official example modified to include what you are seeming to try and do (this is most likely assumption based upon your question). It stores all the features in a list and then outputs them all 'en masse' when the transformer is ...
1
You need to convert your pandas DataFrame to a geopandas GeoDataFrame to access the geometry column as shapely object. This post describes how to achieve it.
one solution:
from shapely import wkt
df = pd.read_csv('myFile.csv')
df['geometry'] = df['geometry'].apply(wkt.loads)
gdf = gpd.GeoDataFrame(df, crs='epsg:4326')
0
Line 17 is not needed at all. You do not have to / cannot add an empty row to the layer first. With your line 15 you get the desired vector layer.
The minimal required code could look like this:
from qgis.core import QgsApplication, QgsVectorDataProvider, QgsVectorLayer, QgsVectorFileWriter, QgsFeature, QgsGeometry, QgsPointXY
QgsApplication.setPrefixPath(r&...
1
I tested it myself.
In conclusion, this error is a problem that occurs when osm_id 32947351 is changed to EMPTY after the clip.
And if you check the geometry in the original file, you can see that the value of the linestring is strange. (The start and end points of the line are the same.)
You can check empty rows after operation by using the following ...
2
When you read a csv file with Pandas the result is a Pandas DataFrame, therefore why create a new DataFrame: dataframe = pd.DataFrame() ?
data = pd.read_csv("lat_long.csv")
list(data.columns)
['lat', 'lon']
type(data)
<class 'pandas.core.frame.DataFrame'>
data.head()
lat lon
0 41.389474 2.156421
1 41.383093 2.181116
2 41....
0
It depends on the format of the coordinates you want
import geopandas as gpd
df = gpd.read_file("points.shp")
df.head(3)
FID geometry
0 0 POINT (-0.18411 0.47266)
1 1 POINT (0.16016 0.34106)
2 2 POINT (-0.46957 0.39932)
df['lon'] = df.geometry.x
df['lat'] = df.geometry.x
df['coords'] = df.geometry.apply(lambda p: list(p.coords)...
0
For those with the same issue, I am now using the gdistance package with good results.
Note that for r.cost the input raster values should refer to the time taken to traverse the whole cell.
For gdistance the input raster values should correspond to minutes required to move a metre in that cell.
0
For a reference timing from my test with gdal_grid executable (GDAL 3.3) on Windows. I took the Liechtenstein OSM extract from Geofebrik.de and converted data into SpatiaLite database. Layer "points" had 15286 features
My command
gdal_grid -a invdist:power=15 -of GTiff -outsize 976 1966 -l points liechtenstein.sqlite points.tif
Processing took 6 ...
0
I am using QGIS 3.10 and there i connect to the mapCanvas signals this way:
iface.mapCanvas().renderStarting.connect(self.start_time)
iface.mapCanvas().renderComplete.connect(self.stop_time)
This is how the class looks like. The 2 methods count the time and print a message:
import time
class MyPlugin:
def __init__(self, iface):
self.tic = None
...
1
Windows 10, Python 3.9, QGIS 3.20 (installed with OSGeo4W). This code works fine:
import osgeo_utils.gdal_merge
output_file_path = r'insert path here'
input_files_path = r'insert path here'
parameters = ['', '-o', output_file_path] + input_files_path + ['-separate', '-co', 'COMPRESS=LZW']
osgeo_utils.gdal_merge.main(parameters)
0
If the existing features are all straight lines with no vertices in between the end points, then you could construct python code that would do the math of it for you. It's a little involved, but doable. Figure out the angle of the line by theta = math.atan((y2 - y1) / (x2 - x1)). (You'll have to precede that by testing if x2 == x1 and if so, theta = pi/2 if ...
0
I had the same problem. Its a permissions issue. The way I fixed it was to go to the folder
C:\Users{User Name}\AppData\Local\ESRI\conda and right clicked on the envs folder and turned off Read-Only. Then let that change for all subfolders.
My clone folder is about 1 GB, while my failed ones were about 29mb, so it wasn't pulling all the data for some ...
0
On GIS SE, you should ask one question at a time, so I will answer problem 1.
Without a template raster, geocube determines the area of the raster by the bounding box of the input geometries and resamples the geometries to the raster. That is why it looks how it currently looks
If you want the raster cells to be centered on each point, I recommend doing a ...
0
You might be able to use the complex_split function from here (or here) and test the resulting geometries for which are inside of the rectangle. Like this:
# Split the braid on the rectangle w/o splitting on self-intersections
geom_collection = complex_split(braid, rectangle)
# Filter for the result inside the rectangle
linestrings_within = [ls for ls in ...
0
From ogr2ogr- GDAL documentation
ogr2ogr can be used to convert simple features data between file formats. It can also perform various operations during the process, such as spatial or attribute selection, reducing the set of attributes, setting the output coordinate system or even reprojecting the features during translation.
The problem here is that ...
0
Did you try to remove the quotation marks around the crs - Argument (and only type "crs = wktproj"), as it is a variable and not a string? Looking at your outcome a connection of the problem to the crs seems to be possible to me.
1
The GDAL LIBKML driver relies on the external libkml library which you have installed, however the LIBKML driver is not built in to GDAL by default.
Just installing the libkml library will not automatically add the LIBKML driver to GDAL, you will need to either install a GDAL build that includes the LIBKML driver or compile GDAL from source yourself with the ...
0
GDAL has two KML drivers:
KML https://gdal.org/drivers/vector/kml.html
LIBKML https://gdal.org/drivers/vector/libkml.html
Your orginfo test confirms that you do not have the LIBKML driver. You can use the KML driver that you have
ogr2ogr -f KML out.kml in.shp
but KML driver has some limitations. Alternatively you must find out what is wrong with your LIBKML ...
0
I faced the same issue. The easiest solution, in this case, is not using rasterio. Instead, using gdal will be much easier, without any need to modify the source code. Follow the second step of this video carefully.
Basically, they are creating a virtual raster XML file using gdal.BuildVRT() function of gdal. The virtual raster XML stores all the individual ...
1
This should work out of the box with fiona.__version__ : '1.8.20' (Python 3.6.9):
import fiona
from shapely.geometry import Point, mapping
my_schema = {
'geometry': 'Point',
'properties': {
'id': 'int',
'attribute1': 'float'},
}
my_id = 1
a_shapely_point = Point(9.17044, 45.45340)
with fiona.open('/path/to/shapefile.shp',
...
0
Have you tried adding the path to ArcGIS\Pro\bin\Python\envs to the condarc config file?
You need to add the path before creating the environment and then, when you create a new environment in that location, you should be able to activate it using only the name. Also, its name should appear along with the path when you run "conda env list".
https:/...
0
Since QGIS 3.20 you can use the Help -> About menu to see the Python version.
0
You have probably resolved this, but the solution is to provide a list of identifiers rather than the direct string of the variable you are seeking i.e:
wcs.getCoverage(identifier=[wcs['dop20'].id],
Styles='fc',
bbox=subset_bbox,
format='image/tiff',crs='EPSG:25832', width=200, height=200)
This is a key ...
0
You Can use Intersect tool in the Geoprocessing tab in Arcmap.
0
I recommend to use rasterio if you want to process in python.
I didn't test the following code with multiband raster, but I hope the code will work with your data.
This code stretch raster values 0 to 254 range, set nodata 255, dtype Byte (uint8), and save without compression.
import rasterio
import numpy as np
import numpy.ma as ma
with rasterio.open('data/...
2
Solved! This issue has to-do with Landsat 8 using a different scale then Landsat 7. Resampled so all values go from 0-255
0
Your dims and coords need to have the same name. This should work:
da0 = xr.DataArray(
data = data,
dims = ["x","y"],
coords = dict(
x = (["x","y"], lons),
y = (["x","y"], lats)
)
)
da = da0.rio.write_crs("epsg:4326", inplace=True)
0
If you want to work on rasterio, try as follows.
The following code assumes all rasters have same crs, resolution, extents.
from rasterio.plot import show
import matplotlib.pyplot as plt
import numpy as np
import rasterio
im1 = 'im1.tif'
im2 = 'im2.tif'
im3 = 'im3.tif'
# Read source raster data and plot
fig, axs = plt.subplots(1, 3)
with rasterio.open(im1)...
1
It seems that your code specifies the same CRS for src_crs and dst_crs.
If you don't want to reproject to a different coordinate system, but just want to rasterio.open with the CRS information added, you can write to memory and open instead of writing to a file.
import rasterio
import rasterio.plot
from rasterio import MemoryFile
from rasterio.profiles ...
0
First, you need to see if your layer's CRS is defined. You can check it by:
import geopandas as gpd
gdf = gpd.read_file(<file path>)
print(gdf.crs)
If CRS is not defined, then you need to set the CRS:
gdf = gdf.set_crs(epsg=26918)
Then you can reproject it:
gdf = gdf.to_crs(epsg=4326)
Now you can perform your spatial join:
gdf2 = gpd.read_file(&...
1
If you know the locations you want to route from and to, you can use the built-in method nearest_nodes:
import osmnx as ox
place = "Piedmont, California, USA"
G = ox.graph_from_place(place, network_type="drive")
Gp = ox.project_graph(G) # convert to UTM
orig_point = (567866.7056592009, 4186342.41614376) # UTM coordinates
dest_point = (...
1
To get the distance between two points you could use the geometry_length function from pyproj. It returns the geodesic length in meters of a given shapely geometry.
from pyproj import Geod
from shapely.geometry import Point, LineString
line_string = LineString([Point(-1.154321, 55.124412), Point(-1.2354352, 56.2345235)])
geod = Geod(ellps="WGS84")
...
1
With geopandas 0.6+ & pyproj 2.2+, the axis order shouldn't be an issue: https://github.com/geopandas/geopandas/pull/1122
So, you should be able to do:
gdf4326 = gdf.to_crs('EPSG:4326')
And it should work properly and the warning should go away.
1
Not sure what you were doing but the following do the reprojection without any issue (axis inversion). I've used the following sample https://labs.webgeodatavore.com/partage/nyc_census_blocks.zip (contains nyc_census_blocks.shp and related associated files that use EPSG 26918 like in your case)
import geopandas
# Read input shp
gdf_26918 = geopandas....
Top 50 recent answers are included
