Tag Info

New answers tagged

1

For this, you need to create new ogr features and save them to the layer and not simply create geometries ( geom1.AddPoint(linkpoint[0],linkpoint[1]), feature1.SetGeometry(geom1)), with dataSource.SyncToDisk() at the end (SyncToDisk() might be helpful to ensure that a particular feature is flushed to disk but it is not necessary here, look at the script). # ...


0

In your solution you mix geometries and predicates of ogr and Shapely With Shapely only (intersection(), intersects()) from shapely.geometry import Point, LineString, mapping geom1 = LineString([(0, 0), (0, 1)]) geom2 = LineString([(1, 1), (-1,-1)]) intersection = geom1.intersection(geom2) # geo_interface -> GeoJSON mapping(intersection) {'type': ...


2

If you happen to have all the featureclasses you want to convert in a map document, you could just run the Consolidate Map tool on that MXD. Make sure to check the "convert data to file geodatabase option". Inside the output folder will be one or more fgdbs (depending on how many source pgdbs held the original features). If you don't have an MXD with all ...


0

Take a look at the Feature Class To Geodatabase (Conversion) tool. This should hopefully do what you need. You can right click the tool and batch export a number of items.


0

I do this by creating routes from the lines using the Create Route tool to get the lines with common ID attributes (like a river name) all combined together and oriented in a single direction. Then extract the FROM and TO end points of the original lines into separate point feature classes with the Feature Vertices to Points tool. I use the Locate Features ...


0

THe short answer is: I can;t see a way. However, changing attributes on a feature is time consuming at times (I have have a hard time figuring out why it happens). Thus, this should speed up your code a bit: pr.addAttributes( [QgsField(aName, QVariant.Double) ] ) attrs={} #------ # Add Values to Triangle Features in Values ...


0

loaded plugins and their instance are available in the dict qgis.utils.plugins you can manage load and unload with other methods in qgis.utils (do help(qgis.utils) in pythonconsole to show static methods available) this is only the base to solve you problem, you should manage check if running or not using some method in your plugin, or you can check if ...


3

Numbers in any programming language does not have "absolute" precision. As they must be represented as bits in the computer hardware its precision is limited. Anyway review this links, because maybe they can help you: http://comments.gmane.org/gmane.comp.gis.gdal.devel/19331 http://lists.osgeo.org/pipermail/gdal-dev/2011-August/029793.html


1

Your buffer distance is given in the units of your dataset. Your coordinate system is lon/lat geographic. Thus 0.01 = 0.01 decimal degrees. To be able to buffer by 1 km (1000 m) you need to reproject your dataset to a projected coordinate system which has units in metres. There are examples in the Python GDAL/OGR Cookbook of how to reproject layers and ...


0

Move the polygon to the right and down or the raster to the left and up. I would shift the Polygon using QGIS edit tools (or Arc) if you have access as opposed to moving the polygon in GDAL warp.


0

Why don't the polygon and image line up? They both have the same coordinate system but does the image have a world file as well? One way to get an image to be located in a particular coordinate space is to georeference it. For example, if you were using ArcMap you would use the Georeferencing toolbar and then add control points where a location on the image ...


0

Export your image to a format where the coordinates are numeric and text, good choices are Tiff + world file, BIL/BSQ (EHDR) or Esri ASCII Grid (AAIGRID). If you choose Tiff/world then you will need to destroy the embedded georeference using photoshop/gimp. Open the world file/header and edit the coordinates. You should be able to calculate the x & y ...


1

If you don't need to stretch or rotate your raster, you can use the Shift tool (Data Management -> Projections and Transformations). Measure the offset in x and y, then run Shift incorporating these values. edit: this answer assumes ArcGIS, as listed in the original tags


0

Are you using ArcMap? Here's a way to form a polygon to your raster dataset: (1) Multiply your raster by 0 using Raster Calculator. This is to generate a single value raster grid. (2) Use Raster to Polygon to convert the raster dataset to a boundary polygon. You can choose to smooth the polygon or not.


2

Use a self-join of the same table twice with different alias. To find the intersection of one linestring with all others: SELECT b.id FROM trajectory a, trajectory b WHERE a.id = %s AND a.id <> b.id AND ST_Intersects(a.line_geometry, b.line_geometry) And to find all combinations of intersections of linestrings in the table: SELECT a.id || ' ...


1

On the first pass through the loop, the extent environment variable is set to the extent of the projected version of the first shapefile. On the second pass through the loop, the extent environment variable is still associated with the first shapefile. When Project is called, only the features within that extent are processed (presumably no features in the ...


0

It seems to be working now. At the very end of my loop I added the following code. arcpy.env.outputCoordinateSystem = None arcpy.env.extent = None


1

This code has not been tested, but it should get you going in the right direction: import arcpy from arcpy.sa import * arcpy.env.workspace = "Path/To/Workspace.gdb" arcpy.CheckOutExtension("Spatial") #Static Variables field = "VALUE" #Build list of years rasters = arcpy.ListRasters("x_*, "ALL") #Loop through names, get year and use it to select rasters ...


0

This, for future readers. I noted this clue that would be helpful layerSS = QgsMapLayerRegistry.instance().mapLayersByName("TestPoint") print layerSS returns: [< qgis._core.QgsVectorLayer object at 0x000000000FCC6730 >] while: layerSS = QgsMapLayerRegistry.instance().mapLayersByName("TestPoint")[0] print layerSS returns: < ...


0

According to the Help for the Project tool it takes a valid spatial reference object as the output coordinate system, it appears you are providing it a string which is a path to a prj file. Although the help also suggests a prj file can be used try using a spatial reference object?


1

mapLayersByName returns a list-of-QgsMapLayers so you need to use [0] at the end to just get the single one. The defintion for that function is: QList<QgsMapLayer*> QgsMapLayerRegistry::mapLayersByName(QString layerName) Which you can see returns a QList of QgsMapLayer, which is just a list of QgsMapLayers in Pythion. Use ...


2

Works for me in QGIS 2.6 with Processing 2.6 (on Win7) and the test data Spatialite file from the repository. If you are using different data, please share it. Check if network and points are in the same CRS.


0

That script won't work in the latest version of QGIS as the Processing framework has been reworked and some scripts may have been broken. Maybe you will be lucky, you could try to change line 17 (as below): from processing.tools.vector import VectorWriter


1

You can use the osr module (part of the standard GDAL install, so it should come with your Python bindings) and do something like this: import osr import gdal inDS = gdal.open(r'c:\somedirectory\myRaster.tif') inSRS_wkt = inDS.GetProjection() # gives SRS in WKT inSRS_converter = osr.SpatialReference() # makes an empty spatial ref object ...


1

I don't know what could be the difference between console and plugin, I think it should to work either way. Also you try to replace layer.commitChanges() (which makes not sense there) with layer.triggerRepaint() and look if that helps (note it is not necessary to call mapCanvas.refresh()). You can also use custom properties for the layer to enable labeling, ...


0

You can use ogr2ogr to reproject the shapefile. I would probably just write a shell script to call ogr2ogr before calling your python script, but there is a Python interface to ogr2ogr. It's explained here. (That example is how to convert gml to shapefile. Near the bottom of this page you can find an example of how to reproject, though not the python part.)


2

First you will need an iterator to go through your shape files, there are two methods I employ: Method One: a folder full of shape files: import arcpy, sys InF = sys.argv[1] arcpy.env.workspace=InF for fc in arcpy.ListFeatureClasses(): Method Two: a whole tree full of shape files: import sys, os, arcpy InFolder = sys.argv[1] for (path, dirs, files) ...


0

I just discovered a solution, pretty straight forward: I am using QGIS 2.6, and the plugin builder 2.0.3 -- both downloaded within the last couple of weeks. The plugin builder creates a "resources.qrc" file, while creating a main script with this reference: import resources_rc at line 27 Until I removed the trailing "_rc" in the main script, I also did ...


1

Windows users must download an executable installer like it says in the Shapely README.


0

How did you install Shapely ? Shapely is not just a pure Python module, it has some "extras" (extensions written in C,C++, from the GEOS library, that need to be compiled and linked before being used, with the Cython module). This compilation can be done at several stages (pip, easy_install or running python setup.py) but Windows doesn't come with a ...


0

I think user defined functions should do the trick. The concept is well explained here: http://nathanw.net/2012/11/10/user-defined-expression-functions-for-qgis/ from qgis.utils import qgsfunction from qgis.core import QGis @qgsfunction(0, "Python") def custom_function(values, feature, parent): return "hello world"


0

Well is was easier than expected, not ODBC stuff, just the following code: mylayer = QgsVectorLayer(r"D:\Test\test.mdb|layername=points", "mypoints", "ogr") QgsMapLayerRegistry.instance().addMapLayer(mylayer) I hope it will work on the newer versions of QGis.


0

I finally figured it out! Instead of the heatmapplugin I used the kernel density estimation algorithm: import processing processing.runandload("saga:kerneldensityestimation", "path/to/shapefile.shp", "name of the colums with weigths", radius, 0, 0, None, celsize, "path/to/raster.tif")


0

The script is failing to find the Python regular expression library, which is strange because re is part of core python. You could open the Python console in ArcMap, and import re to see if the same error occurs. As @artwork21 suggests, you could browse to your python installation and see if re.py exists. You could start python from the command line and ...


1

The gdalbuildvrt program isn't the right tool for this task. It's for making virtual mosaics and multiband images. I suggest this approach: Define the domain in which you want results (geotransform, rows and columns). We used to call this a "canvas" back in the day. Warp (gdalwarp, nearest neighbor) your source TIFFs to that domain, producing TIFFs with ...


0

I think Nathan W's answer is out of date. I was able to run QGIS (version 2.6) python scripts from the command line (Nathan's option 4) using the following commands. man qgis qgis -nologo --project /path/foo.qgs --code /path/foo.py


0

This is not the answer for your question but may be for your problem. I would populate stream raster, for each point (cell) in streams populate watershad. Within populated watershads, iterate trough all points and calculate HAND something like: cell_value - min(watershad).


0

I don't use QGIS, but with ArcGIS, PDF export has always been tricky, especially with large datasets, at small scales and with poorly designed map documents. In general, I don't think PDF is a good geographic data visualization format because it has no concept of scale-dependent visibility, generalization, etc. But as it is the most popular "print" format, ...


0

I realize this is quite a while after the question, but we just had a big need for a similar capability and ended up building a general, free, open source Python library for the purpose. You can find it at https://github.com/ucd-cws/arcpy_metadata or by running a "pip install metadata" - documentation is still to come, but we tried to keep things relatively ...


0

I have deployed the very same architecture in XenApp + QGIS 2.4. If you digg into QGIS and associated tools, you will find tons of tools that can use all the resources of your virtual servers. I don't think it is possible to disable them all, because users will probably need them, on reasonable datasets. What I do is have big raster data on readonly file ...


2

import arcpy, traceback, os, sys, math from math import radians,sin,cos from arcpy import env env.overwriteOutput = True inFC=arcpy.GetParameterAsText(0) outFolder=arcpy.GetParameterAsText(1) nRows,nCols=4,3 env.workspace = outFolder rectangle=r'in_memory\rectangle' tempf=r'in_memory\many' outFile="fnet.shp" ## ERROR HANDLING def showPyMessage(): ...


1

The "export to Python" option in the modeler has been taken out in the latest version. The modeler has been rewritten and there was not enough time to implement this feature. We have released the new version without it, considering that the improvements that this new modeler have are really worth it. The "export to Python" version will be available again in ...


2

Here is the answer I ended up using: def rotate_point(input_point, center_point, angle=0): angle = math.radians(angle) x = float(input_point.X) - float(center_point.X) y = float(input_point.Y) - float(center_point.Y) xr = (x * math.cos(angle)) - (y * math.sin(angle)) + float(center_point.X) yr = (x * ...


2

Assuming ArcGIS, here is an example Python Addin that calculates areas: class CalculateArea(object): """Implementation for Test_addin.calculateArea (Button)""" def __init__(self): self.enabled = True self.checked = False def onClick(self): # Set mxd to current and data frame to active mxd = ...


2

The best way to print out logging in qgis is to use the QgsMessageLog class. This writes to the debug log window rather then the console. QgsMessageLog.logMessage("message", "name") In the lower right of the main window you will see a icon that you can press to open the debug log.


0

That code are worked: import arcpy,os,sys from arcpy import env env.workspace = r"C:\Project" for mxdname in arcpy.ListFiles("*.mxd"): print mxdname mxd = arcpy.mapping.MapDocument(r"C:\Project\\" + mxdname) df = arcpy.mapping.ListDataFrames(mxd, "Layers")[0] LayerList = ...


2

import arcpy mxd = arcpy.mapping.MapDocument("CURRENT") df = arcpy.mapping.ListDataFrames(mxd, "Transportation")[0] Where [0] refers to the order of the data frame you are looking for. More info: DataFrame (arcpy.mapping)


1

I think the application server will be doing very little so I would just go with whatever web framework you are already familiar with. You can configure GeoServer completely using their REST API and you can manipulate/query the data directly from GeoServer using WMS/WFS (preferably with JSONP) without having to touch PostGIS directly. All the application ...


3

Your expression looks incomplete but I am assuming that the first part of it is trying to find the square of the cube root. My test below on that when R had a value of 8 returned the expected value of 4 to the RESULT field RESULT = [R] ^ ( 2 / 3 ) The above worked when the parser is set to VBScript. With Python, there is hopefully an easier way but I ...


0

I tried solving this by re installing fiona. I was having the same problem. I installed fiona 1.4.8 and it solved the problem.



Top 50 recent answers are included