New answers tagged

2

GDAL has no problem with GeoJSON http://geojson.org/geojson-spec.html. The following GeoJSON encodes one multilinestring: { "type": "FeatureCollection", "features": [ { "type": "Feature", "geometry": {"type":"MultiLineString","coordinates":[[[11.20558631,46.48251782],[11.2058444,46.48280049]],[[11.20578705,46.48252192],[11.20596731,46.48275133]]],"crs":{"...


0

If you are looking for pure python with gdal on your cmd you need to do follow steps: Step 1: Install pure python from python.org Feel free to download the latest 2.7x version of python (rather than the 3.x python version). Install python with the default options and directories. Step 2: Next Install the GDAL Binaries Head over to Tamas Szekeres’ Windows ...


1

Um, the Sentinel-1 sensor is C-band SAR (radar) data, not spectral. Given that NDVI is a spectral index of photosynthetically active radiation using the ratio between Red and NIR spectra, what makes you think that you can get NDVI from SAR polarisation? There is absolutely no way to derive NDVI from SAR data! Please, go research some basics on active and ...


0

As per http://python-visualization.github.io/folium/quickstart.html "Folium supports passing any Leaflet.js compatible custom tileset" let lets try that and bypass it's built in support: folium.Map(location=[45.372, -121.6972], zoom_start=12, tiles='"http://{s}.tiles.mapbox.com/v4/wtgeographer.2fb7fc73/{z}/{x}/{y}.png?access_token=pk.xxx', ...


1

Column and row along with the zoom level are needed to identify a specific tile. This section in the osm wiki shows how they are calculated.


0

I'm three years late, but hopefully whoever else might be stuck on an issue like this will benefit from my answer. # adds a LineLayerStyle (Arrow, dotted) # layer being a vector layer object of QgsVectorLayer() def setLineLayerStyle(self, layer, style): if(style == "arrow"): sl = QgsSymbolLayerV2Registry.instance()....


0

I just realized that the problem was in dictoxml or dom parsing. I changed my code a bit everything works nicely: xml = open(path, "r") org_xml = xml.read() dict_xml = xmltodict.parse(org_xml, process_namespaces=True) out = xmltodict.unparse(dict_xml, pretty=True) with open("path/new.qml", 'a') as file: file.write(out.encode('utf-8'))


0

You need to have "write" access to the data source. In this case the data source is a web page resulting from a query. One solution is to create a copy of the data source in the computer's memory, and, then add to that memory: #create an output datasource in memory mem_driver=ogr.GetDriverByName('MEMORY') source = mem_driver....


0

I downloaded the files from the website mentioned by hlinak and then extracted them to the folder nominated by Daniel. After a quick reboot all was well.


1

Thanks to user30184! You helped me solve it. The SrcLayer must specify the same name as the CSV file. Here is the correct, working VRT file. The change is on line 4: <OGRVRTDataSource> <OGRVRTLayer name="example"> <SrcDataSource>D:\Temp\example_short.csv</SrcDataSource> <SrcLayer>example_short</SrcLayer&...


1

This is not a problem of Matplotlib but your script and the module you use for reading shapefiles 1) You know that there are points in the geometries of the Polygon shapefile thus eliminate try... except 2) you load and read the shapefile twice for x and y (memory) for shape in shp.shapeRecords(): xy = [i for i in shape.shape.points[:]] x = [i[0] ...


2

I didn't knew geopandas, it is very easy to add a field, a great library indeed. You just need to read your file add it and save. dataSrc = gpd.read_file('my_shp.shp') dataSrc['new_field'] = 1 dataSrc.to_file('newfile.shp') Selecting is also straightforward: dataSrc = gpd.read_file('my_shp.shp') dataSrc[dataSrc['id']<4].to_file('new_shape.shp')


2

OK, I think I've solved my problem. The issue is noted at the github page for laspy here: https://github.com/grantbrown/laspy/issues/51 I needed to adjust the "start_first_evlr" property to be the full length of the LAS file as shown below. Not sure if that is the best way, but it seems to work. I needed to read the mmap size for the input LAS file (before ...


0

If you installed Arc before Anaconda, your default IDLE, Python Window (in Arc), and Python console will always pull from your Arc-installed modules, i.e. C:\Python27\ArcGIS <version>\Lib (and/or) Lib\site-packages. If you have packages installed with Anaconda that you want to use with Arc's install of Python you'll need to write .pth files in order to ...


-1

This tutorial shows you how to set up conda environment to work with ArcGIS 10.4 and ArcGIS Pro 1.3. Workflow to set up Anaconda with ArcGIS 10.4 Install Anaconda without fouling the Windows environment (paths, registry) to break Esri's python stack Configure Anaconda with the particular add-ons you want, and Configure ArcGIS's Python so that ...


1

This would be extremely simple in PostGIS: SELECT line.unique_name, line.geom, COUNT(circle.geom), SUM(circle.attribute_name) FROM line JOIN circle ON ST_Intersects(line.geom, circle.geom) GROUP BY line.unique_name, line.geom; In QGIS 2.14 there is the option of creating a "Virtual Layer" which might work for this. But the performance is not as good with ...


1

You can use the srsName= string to get the polyline from the WFS in lat/lon coordinates from the start. Alternatively, provided that you've read the coordinate strings into Python as arrays x and y, you can do import pyproj prj = pyproj.Proj(PROJ4_STRING) lon, lat = prj(x, y, inverse=True) where PROJ4_STRING is the libproj SRS string for the original ...


3

In your if line you need to change your or to an and in the last part. if f4 == f5 and (f1 == f2 or f1 == f3) and (f6 != f7 and f6 != f8): This is because with f6 != f7 (6 != 5) it returns 1, even though f6 != f8 fails (6 = 6) in the second part. With the OR there it returns true if one of them is true. With an AND it would return true only if both are ...


0

I'm not an expert regarding python + GDAL/OGR, but looking into something similar I found a possible answer to your question. import ogr, osr import fileinput from datetime import datetime srs = osr.SpatialReference() srs.ImportFromEPSG(4326) driver = ogr.GetDriverByName("GPX") out = driver.CreateDataSource('/tmp/test.gpx') # layer creation: if you use '...


0

To clarify, @AndreJ describes creating a local copy of a plugin, and renaming it. What you cannot do is rename a plugin on the QGIS plugins repository website, even if you are that plugin's author. If you change the name of a plugin you maintain yourself, it will become a new plugin in the repository. Old versions and download stats from the previously named ...


0

There are a few questions I had about your code. Opening the files in 'rb' mode means to read it as a binary. Is this what you meant? I find your indexing a bit confusing. Your lines converting to float are mapping every item from an index on. That means that map(float, row1[16:]) converts every item from 16 to the end, map(float, row1[17:]) converts every ...


0

You can use the KD-tree based method that you propose. I've split the answer into two parts, addressing your questions about projection and how to interpolate with the KDTree. Projection Provided that the SMOS and SMAP datasets are in a cylindrical projected coordinate system, you will need to first convert them to longitude and latitude to use the linked ...


1

The following PyQGIS line should do the trick: self.contourLayer.setCustomProperty("labeling/labelPerPart", "true") As soon as you run such line, you should get the 'Label every part of multi-features' checkbox toggled in the properties dialog. I'm using QGIS 2.14.2.


0

I don't have enough reputation to post this as a comment, but since nobody else has contributed yet, I'll just put this here. Are you sure you are passing your lat/longs as numbers and not strings? Maybe try something like this just to make sure: print haversine(float(floats1[1]),float(floats1[0]),float(floats3[1]),float(floats4[0])) Also, is there a ...


2

Signal can be emited by an object. Ex. your canvas: canvas = iface.mapCanvas() custom_canvas = QgsMapCanvas() # if you want to create new canvas instead of this from main window Canvas can emitt a signal renderComplete after rendering is complete, mapCanvasRefreshed after you refresh this canvas. So you can connect this signal to the proper function and ...


0

There are two modules in QGIS Python included: Scipy csgraph: http://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.sparse.csgraph.dijkstra.html NetworkX: http://networkx.readthedocs.io/en/stable/reference/algorithms.shortest_paths.html Both modules are efficient and written in C.


1

This has been reported as a bug on QGIS (bug number 15279). The reporter says: "the problem is the use of the character "ø" (o-slash) in the directory name for the QGIS program directory: "C:\Program Files\QGIS Nødebo". If you make a name change of the installation directory to ex. "C:\Program Files\QGIS Noedebo" it works without a hitch." I tried this by re-...


0

After a night without sleep, I fixed it myself: selectedLayerIndex = self.dlg.cbox_layer.currentIndex() selectedLayer = layers[selectedLayerIndex] # get unique values fni = selectedLayer.fieldNameIndex('firebrary') unique_values = selectedLayer.dataProvider().uniqueValues(fni) # define categories ...


2

In C:\Users\<username>\.qgis2\python\plugins\, just rename the folder of the plugin. It will appear in the plugin list with that name. Once renamed, it can not be overwritten from the plugin repository. You have to activate it manually in the installed list to get it running. Furthermore, you can change name and version number in the metadata.txt file....


0

This example work perfectly. from PyQt4.QtCore import * lyr = iface.activeLayer() lyr.startEditing() provider = lyr.dataProvider() provider.addAttributes([QgsField("newAttribute", QVariant.Double)]) lyr.updateFields() lyr.commitChanges()


0

It's possible PostGIS is doing something under the covers that you need to do by hand in Shapely. We use prepared geometry automatically in predicate loops, and in Shapely you might need to marshall that on purpose. Depends on what your process actually entails. Advantages of python are the lower install overhead and the library of potential other features ...


0

I was thinking about difficoult soulution and I just overwrote it all. It was so simple. I must just import dialog for second window, add a method to open new dialog and add a line in init to react when button pushed.


1

The quickest way is to keep the geometry in hex-encoded WKB. import psycopg2 from shapely import wkb conn = psycopg2.connect('...') curs = conn.cursor() shps = {} # key: gid, value: Shapely geom curs.execute('select gid, geom as geom from pysal.stl_hom;') for gid, geom in curs: shps[gid] = wkb.loads(geom, hex=True) There's no need to use GeoJSON, ...


0

So, it was a ssl certificate validation problem. The python used by QGIS didn't recognize the organism that had signed my API's ssl certificate. There were a possible work around here (telling urllib2 to ignore certificate validation) but it calls urllib2.HTTPSHandler() which has been implemented in a python version posterior to the one QGIS uses. The ...


1

You can create a new dialog as a class, even in your main .py file. class NewDialog(QtGui.QWidget): def __init__(self, parent): super(NewDialog, self).__init__(parent) Than in your main class you can put function like that: ... def open_new_dialog(self): self.nd = NewDialog(self) self.nd.show() ... This function creates a new dialog and ...


0

I think your error is found in the for-loop syntax. Try the following: ValList = [] for j in lstFields: ValList.append(layer.GetFeature(1).GetField(j)) Note that in your case, if you have already chosen a feature, change the last line to: ValList.append(Feature.GetField(j)) The results in my case are given below. > lstFields = ['data1', 'data3'] # ...


1

I can not reproduce your error, maybe the installer has been fixed by now. Similar to the linked question, which deals with Linux installations, you can run the OSGeo4W installer again, and choose Advanced installation to re-install the python-future package. This might not help with the standalone installer, but copying missing folders from OSGEO4W to ...


1

Personally, I've started using rasterio's windowed read/write However, you can use x and y offsets when writing your array to the output dataset. dataset.GetRasterBand(1).WriteArray(array, xoff, yoff)


2

You are overwriting your data 'matrix'. First you reclassify the negative to 1, then you do other steps... but at the end yo are doing : matrix[np.where((matrix > 1.5)) ] = 6 Which is reclassifying your already changed values to 6. That is why you have only 1 (negative) and 6 the other values (2,3,4,5) -> 6 I would suggest to copy your matrix like ...


0

Part of the question mentions not knowing how to create a schema. You can create a schema with pgAdmin (GUI) or with this SQL statement: CREATE SCHEMA tiger_data; I recommend you install the PostGIS TIGER Geocoder. This extention creates scripts (windows and linux) to download and import US Census shapefiles. You can learn a lot from the process that ...


2

In the FieldPyculator, you will need to: Insert your function as a Global expression and remove the AP parameter as you can choose from the GUI which field you want updated. Slightly change your Field expression by specifying fields using <Mag> instead of !Mag!. This is the global expression I used (not sure what Vector_Mag was so changed these ...


1

You need to iterate through the segments of the line to find the one that touches/intersect the polygons (as in Get the vertices on a LineString either side of a Point) I show here the results with one line and one polygon as example # transform the line to a list of points linegeom = liengeom.asPolyline() # extract the segments seg_start, seg_end of ...


1

This can be done in far fewer lines of code src = gdal.Open(path goes here) ulx, xres, xskew, uly, yskew, yres = src.GetGeotransform() lrx = ulx + (src.RasterXSize * xres) lry = uly + (src.RasterYSize * yres) ulx, uly is the upper left corner, lrx, lry is the lower right corner The osr library (part of gdal) can be used to transform the points to any ...


2

Try out with 'intersects' and 'intersection' methods. Next code produces only two points. mapcanvas = iface.mapCanvas() layers = mapcanvas.layers() feats_line = [ feat for feat in layers[0].getFeatures() ] feats_poly = [ feat for feat in layers[1].getFeatures() ] for feat_l in feats_line: for feat_p in feats_poly: if feat_p.geometry()....


1

If you look at the docs for cursor.fetchall you will see that the return value is a list containing many tuples, each tuple representing a row in the returned query (cursorset). So, in order to test for that a return value [(True,)] is true, you need to use double indexing, [0][0], once to return the first element of the list, and once to select the first ...


0

Sample code snippet is given below uri = QgsDataSourceURI() uri.setConnection("xxx.xxx.xxx.xxx", "port", "dbname", "username", "password") uri.setDataSource("public","tablename","geom") uri.setWkbType(QGis.WKBPolygon) #geometry type has to be mentioned uri.setSrid('4326') vlay=QgsVectorLayer(uri.uri(),"test_poly","postgres") ...


1

The easiest solution is to get a dictionary from query (psycopg2.extras: Real dictionary cursor or Dictionary-like cursor ) and the GeoJson format (ST_AsGeoJSON) import psycopg2 conn = psycopg2.connect("dbname='strati 'host='localhost' user='me'") import psycopg2.extras dict_cur = conn.cursor(cursor_factory=psycopg2.extras.RealDictCursor) dict_cur.execute("...


1

Next code is similar to the cookbook code for testing the OGR 'Transform' method (from EPSG:32612 to EPSG:4326). It also includes 'AutoIdentifyEPSG' and 'GetAuthorityCode' methods for determining EPSG code for projection in my GeoJSON layer (EPSG:32612). from osgeo import ogr, osr import os driver = ogr.GetDriverByName('GeoJSON') # get the input layer ...


0

Can be done fairly easy. First, create the dataset: drv = gdal.GetDriverByName("GTiff") ds = drv.Create("name.tif", width, height, 6, gdal.GDT_Float32) You then need to set the geotransform using ds.SetGeotransform where the argument is a six element tuple: (upper_left_x, x_resolution, x_skew, upper_left_y, y_skew, y_resolution) Set the projection ...


1

There are a lot of different libraries for geoprocessing. Depending on your particular needs could be convenient to use one or another. Even though, it is possible to combine many of them using using a common format or datatype. Here is a small list of the libraries I'm using at the moment. Python/GDAL The so called "Swiss knife of GIS" https://...



Top 50 recent answers are included