texts: Splitting string in Python parser of ArcGIS Field Calculator?

texts:
I am trying to use Python to split a string field into three new fields at the separator " - ". The field is in the format "a - b - c". So I initially had this in the field calculator for newFieldA:

!myField!.partition(" - ")[0]

This gave me the value a in newFieldA. Similarly, I calculated newFieldB and newFieldC. After a more experienced Python user took a look at it, this was the final result:

' ' if !myField! == 'NULL' else !myField!.split(' - ')[0]

The test for NULL was put in because there were some NULL values throwing errors. Each of the three fields were calculated using this code, just the index was changed (0, 1 or 2).

My problem now is, after changing spelling errors in myField, I tried to split the field again, using the exact same code which worked previously. It now fails to run, giving the error of "wrong field name or unbalanced quotation marks". I am definitely using the correct field, and the quotation marks are not unbalanced. Also, I literally copied and pasted this from the code I used previously, so I don't see why it would be failing now. I feel like I'm missing something really simple. Does anyone have any ideas?

**edit** I've placed the following in the codeblock

define splitField(value):
' ' if value Is None
else:
  value.split(' - ')[0]

and in the expression window

splitField(!myField!)

I am getting a syntax error on line 1, I have not been using Python long enough to see what it could possibly be, as that looks right to me. Does anyone have a suggestion to make the code better?

edit I'm working with crime data, so in the field SUBACTION I have "Policing - Theft - Fraud". I then created three fields MANDATE, CATEGORY, and OFFENCE. The end result should be "Policing" in MANDATE, "Theft" in CATEGORY and "Fraud" in OFFENCE. I changed the codeblock to

def splitField(value):
  ' ' if !SUBACTION! == 'NULL'
      else:
        value.split(' - ')[0]

It is now giving a syntax error in line 2.


texts: Joining layers with respective elements with PostGIS Topology

texts:
I'm currently using PostGIS topology extension, but I have some difficulties in understanding how the structure works:

One of the key points is the use of "layer"s: as of what I understand, features attributes should be stored in a table out of the topology's schema (the one named topo_actualname) and registered as a layer of that topology with AddTopoGeometryColumn.

However, is there a simple way to join the attributes (stored in the layer table) with the respective features (elements in the node, face, or edge_data)?

Now, what I do is:

SELECT whatever
FROM layer_tb l
     JOIN topo_topologyname.edge_data e ON (l.topo).id=edge_id;

But I guess the whole layer concept is rather useless if I have to know both the topology schema name and the layer name to get the information I want.

In fact, I think I understood that the topo column on the layer has enough information to know where the respective topology is, and moreover the topology schema stores a reference to each layer table for each topology.

Is there a short/simple/proper way to join information together? I was looking for something in the topology extension functions, but couldn't find anything useful.


texts:
On the right is a single row feature class polylineZ (with elev data) displayed on a new 2D map with no base map.

On the left is that same feature copied (using drag from contents pane) to a new 3D scene with no base map.

Why does the 3D look radically different? What factors would vary the image after copying from 2D to 3D? or is it something else I'm missing?

All I want is the feature to match. I am not groot. I am noob.

ArcGIS Pro v2.1.3


texts:
I'd like to use dot density to display some data in Arc, but the areas (polygons) I'm displaying it for are so small in some cases that I suspect the pattern is being misrepresented.

If several small areas fall within the size of one 'dot' for example, will ArcGIS suitably aggregate the data from multiple areas to display a dot?  If there isn't enough data (for example because the number of households is below the dot quantity), will those be aggregated into neighbouring areas?

If I am concerned about misrepresenting the data in this way, is there a way to work around this?


texts:
I need to install pgRouting on a Debian 7.5 64-bit server.

I have PostgreSQL 9.1 and PostGIS 1.5 running. But somehow I can not install pgRouting.

1 - To use add-apt-repository ran the command:

apt-get install python-software-properties

2 - I entered the site georepublic https://launchpad.net/ ~ / + archive / ubuntu / pgrouting-unstable ppa to get the source to add to the list but only has to Ubuntu.

Can anyone help me???


texts:
I am trying to create a polyline from the points that I query from a polygon.   I have a polygon, and I am getting its edges, casting them to segments, and then use QueryNormal to get the mid point for that edge. I get 2 different points from 2 edges, I want to create a curve from those 2 points, not a straight line, but follow the polygon path from one point to another. the pictures below explain things better. I tried using the topological operator and cut but thats not doing what I want, I tried using the polygon.Split() but also that does not give me what I want. Is it possible to split a polygon based on points?

New Edits:-
After being able to get the subcurve using the 2 input point distances, it is always giving me the longest part of the polyCurve, here is apicture of what its doing, where I want the short path as per the picture above.


texts:
I'm fairly new to Java and GeoTools and I'm trying to plot this grib of meteorological data; I have the code for shapefiles, but I don't know where to go from here:
File file = JFileDataStoreChooser.showOpenFile("shp", null);
        if (file == null) {
            return;
        }

        FileDataStore store = FileDataStoreFinder.getDataStore(file);
        SimpleFeatureSource featureSource = store.getFeatureSource();
        CachingFeatureSource cache = new CachingFeatureSource(featureSource);
        // Create a map content and add our shapefile to it
        MapContent map = new MapContent();
        map.setTitle("Using cached features");
        Style style = SLD.createSimpleStyle(featureSource.getSchema());
        FeatureLayer layer = new FeatureLayer(cache, style);
        map.addLayer(layer);
        // Now display the map
        JMapFrame.showMap(map);


texts: Changing Annotation Layer styles in a geodatabase

texts:
I have a workflow, using FME, that outputs an annotation layer in a file geodatabase.
At first I create some pre-defined Annotation Classes in an empty geodatabase that will be used to style Water and Road text features etc. I note down the AnnotationClassID of each class.
I then run the workflow which assigns the AnnotationClassID noted earlier to each feature based on a feature class attribute from the source data. This workflow also sets the FontSize, Rotation and Vertical and Horizontal Alignments of the Annotation feature.
This runs fine and I can see the results in ArcMap as expected. However, I was wondering if it is possible to tweak the styles of this output Annotation layer?
When I try editing one feature and changing its AnnotationClassID, called "SymbolID" in the attribute table, from a Water to a Road text style of the pre-defined Annotation Classes, the FontSize automatically reduces from a large number, 93543 (set from the workflow above), to 8, drastically reducing it appearance on the map to almost a dot.
I guess this happens because I left the FontSize in the defined Annotation Class as the default "8" in the first place.
What I want to know is if it is possible to preserve the FontSize when changing its AnnotationClassID?


texts:
After importing OSM-Data to my qgis-project (as SpatiaLite database and/or shapefile-extract from planet.osm), some features are missing or have no attribute values. For example that lake: http://www.openstreetmap.org/relation/552029
Its just missing in the shapefile-import and the SpatiaLite import have the gemoetries of the lake, but no values!

The most other features are fine.


texts: How to show a selected polygon feature, using a single marker in its center

texts:
I've been asked to replicate a map where the selected polygon is shown with a red tick/check mark:

I need to tell ArcMap to use this symbol when a polygon is selected, so that un-selecting the polygon removes the check mark, and selecting a new polygon shows the check mark for that polygon. I don't want to edit a point layer and move the point around, which would be a manual workaround.

I tried setting the layer's selection symbol...

... but the only options here assume that the marker/image should be distributed multiple times throughout the selected polygon:

Is it possible to show a single marker within a selected polygon? I suspect this will require ArcObjects coding (which is probably out of the picture for this) but someone may know of a workaround.


texts: How to convert bathymetry sounding data to DEM (digital elevation model)?

texts:
https://www.ngdc.noaa.gov/docucomp/page?url=http://surveys.ngdc.noaa.gov/mgg/NOS/hsmdb/H04001-H06000/H05455_hsmdb.xml&view=hydro/survey&header=none

I downloaded the bathymetry sounding data in xyz format but noticed that the depth is positive. Why is there positive depth when the gps location is clearly on water?

With ArcGis 10.4 I was able to extract xyz data from the DEM file and plot the image below. The DEM file was used since it covered most of the area needed except for H as shown below, that's why I was trying to use bathymetry sounding data. Orange dots= depth 0, green dots is <0.

If I use bathymetry sounding data there is many orange dots with depth 0.


texts:
I want to do some work with projections in my c++ application with QGIS API v2.14. But canvas and layers didn't return their CRS. Even when I'm setting it to them with setDestinationCrs() for canvas or setCrs() for the layer they returns nothing anyway. QGIS desktop application returns canvas and layers CRS from console at the same time.
OS is CentOS 7. QGIS API v 1.8 on CentOS 6 has the same problem too.


texts: How to reduce all value above certain value using Raster Calculator?

texts:
I have a raster with elevations above 2000m that I want to "cut" to 1800m.

I cannot quite figure out the correct string and most simple way to do this in the ArcGIS 9.3 raster calculator.

I tried to ceiling function without much luck; would an expression similar to IF([my_raster]>2000,1800,[my_raster]) do the trick?


texts: pgr_drivingdistance seems to miss a few lines

texts:
For multiple stations I am running a pgr_drivingdistance function. Its nearly correct but it has a strange behaviour I can not place.

The function seems to leave out some of the road segments around a station. Leaving it with the agg_cost value of null. Then when the station next to it runs through the model it updates the costs that are near another station with a value that is based on the distance towards that station.

I prepare the roadnet with the following code:

    select pgr_createTopology('bikenet', 0.000001, 'geom', 'gid')

I loop the stations through the following code:

    SELECT
         *
    FROM
         pgr_drivingDistance('
         SELECT
            gid as id, source, target, time as cost
         FROM
            bikenet', LOOPEDFILLINSOURCE,
         10000::float, FALSE)

Any idea what is causing this behaviour and how to fix it? In underneath picture you can see the result with road segments labelled with agg_cost.

In the second image the red lines show the roads that are not processed at all. So the green dot is the starting point. The black lines are processed and the red ones are not.


texts: How to Convert Interpolated Map to Excel using ArcGIS for Desktop?

texts:
I have interpolated some SO2 data using IDW and have calculated means on county level using Zonal Statistics. I would like to convert this data now into an excel file. I already got the hint, that I should use the Conversion Tool and Table to Excel. However, I can only convert the two tables I have imported from shape files. That is the map with the counties and the data with SO2 values and measuring stations. Do I have to convert the values I obtained from Zonal Statistics?

I hope I was able to describe my problem sufficiently.

I am using ArcMap.


texts: Applying Attribute Filter and Spatial Filter simultaneously (GDAL - Python)

texts:
quick question on something that has to be simple, but I'm getting stuck on. I want to use an Attribute Filter and a Spatial Filter simultaneously.

I have two vector files (an overlaying grid and a stream). I want to iteratively select road segments that fall within each grid square and that have a specific attribute value for length.

I've been using -

roads.SetSpatialFilter(grid_geometry)

roads.SetAttributeFilter("Length = 'xxx'")

However, regardless of which filter I use first, applying the second filter resets the process (so I either get the attribute subset or the spatial subset, not both).


texts: Understanding percentage drop-raster (flow direction-tool)

texts:
I'm creating a Topographic Wetness Index using a python script (http://arcscripts.esri.com/details.asp?dbid=16750)

As I don't want this just to be a blackbox, I'm trying to understand the different rasters this script produces. One of them is outdrop-raster, which is a by-product of flow direction raster. According to ESRI:

  "The output drop raster is calculated as the difference in z-value divided by the path length between the cell centers, expressed in
  percentages. For adjacent cells, this is analogous to the percent
  slope between cells. Across a flat area, the distance becomes the
  distance to the nearest cell of lower elevation. The result is a map
  of percent rise in the path of steepest descent from each cell."

..which seems reasonable. The problem is, that there are several numbers which pop out. For example 16,666.. % and 25 % are by far the most common numbers, why? I'm sure there is an explanation for this, but I just can't figure it out. It doesn't seem logical to have continuous elevation (thus slope) data and results that are not (that) continuous. Original data is low-density LiDAR DEM (2 m pix), values are 32 bit float (197.132, 198.013 m etc.).

5 by 5 example, darker gray equals flow direcetion to lower left and percentage 16,666  % drop while lighter gray here has 25 % drop and flow direction to the left


texts: Assign colour to my shapefile with matplotlib

texts:
import shapefile as shp
import matplotlib.pyplot as pltsf = shp.Reader("D:/New folder1/21_8_2016.shp")

plt.figure()
for shape in sf.shapeRecords():
    x = [i[0] for i in shape.shape.points[:]]
    y = [i[1] for i in shape.shape.points[:]]
    plt.plot(x, y)
plt.show()

this shows my hsapefile in alot colours how do i assign it a single colour ?


texts: Empty downloaded files from SOILGRIDS 250m

texts:
I am trying to download pH in water mean values for an area long = 48.2 to 48.9 degrees and lat = -4 to -3 degrees. Using the downloading tool of the portal https://soilgrids.org/, an empty file was returned. It was the same with clay or organic carbon content, and with different areas.
It is not the first time I download data for SOILGRIDS, but the first time it failed.


texts:
Is there a way to automate the placement and the text based on a field using callout labels?

One thing I've tried, which is quite cumbersome, is to convert the graphical objects after they have been drawn to an annotation class and then spatially joining the new features to the points they are labeling, then calculating the text field with attributes from the joined table.

I'm really hoping there is a way to generate box/callout labels in a more automatic fashion. I don't want to have to draw all several hundred callouts.

Software Version: ArcGis 10.3


texts: Create buffer around one category of polylines in feature

texts:
I have a feature that contains polylines of fluvial and tidal streams in ArcMap.

I want to create a buffer around the tidal streams only. I tried selecting only the tidal streams using the attribute table and then using the buffer tool in the editor toolbar, but got a notification saying it's not possible to buffer a network of features.

What's the best way to go about this?


texts: Including rivers in least cost analysis in ArcGIS 10

texts:
I am relatively new to ArcGIS.

I am looking to find the least cost path between a  number of points in a landscape (eastern Sicily) that avoids rivers (obviously a Roman cart would not plough through a river). I have DEMs for the elevation and ran the analysis for the slope, and I had a shapefile for the rivers which I used to also create a feature; I've used conversion tools with this to create a raster for the rivers too. I'm just not sure where to go from here!

I'm using ArcGIS 10 on Windows 7.

I do know how to use the cost distance, back link, and cost past tools to find the least cost path; my problem is getting the rivers incorporated into this.


texts: Need to create an XAML GUI to reference a local map (not on the MapServer site) that I'm using in ArcMap

texts:
I know how to create an XAML GUI for a MapServer map . . .

<esri:Map x:Name="MyMap" Extent="-120,20,-100,40" >     <esri:Map.Layers>          <esri:ArcGISTiledMapServiceLayer ID="StreetMapLayer"                Url="http://server.arcgisonline.com/ArcGIS/rest/services/ESRI_StreetMap_World_2D/MapServer" />
     </esri:Map.Layers></esri:Map>

However, I need to know how to set up the syntax for using a map that I created from a Shapefile Feature Class that I have brought into ArcMap.


texts: Forcing ArcGIS Desktop to re-read file geodatabase dataset extent?

texts:
I have feature class sitting in FGDB. I can add it to map, right click on it in the table of content and pick "Zoom To Layer". This is what I see:

I start editing session, remove most of the features in this feature class, save edits, stop editing. Right click on the layer, "Zoom To Layer". This is what I see:

The only option I have is to close ArcGIS session and start a new one.

Is there a better way to tell ArcGIS to use the true dataset extent?


texts: Enabling fullscreen option from cartodb.createVis

texts:
Is there an option to enable the fullscreen toggle through the API?

Seems the docs haven't been updated to reflect this update: http://developers.cartodb.com/documentation/cartodb-js.html#sec-3-1


texts:
I am trying to load a OBJ file in GIScene.js. Unfortunately, the OBJ file is not getting visualized inside scene of GIScene.js, even though the scene of GIScene.js is getting visualized. There are no errors and I don't think that there is any issue in the code.

<meta charset="UTF-8">
<title>GIScene.js Scene</title>

GIScene.js - Scene

The examples (https://giscience.github.io/GIScene.js/examples/index.html) are using either CTM or JSON files for loading inside the scene of GIScene.js, none of the examples are showing how to visualize OBJ files. How can I solve this issue ?


texts: Looking for GeoJSON describing French departments?

texts: Is there any place where I could find a GeoJSON describing the borders of the French departments, or do I have to convert a shapefile?


texts: Zoom to selected feature, without using 'current' in ArcPy?

texts:
I am wondering if it is possible to build a zoom to geoprocessing service without using 'current' for the MXD variable in Python?  Using 'current' my tool works in ArcMap, but not in my geoprocessing service.

Here is my code that work fine locally in ArcMap:

import arcpy

MXD=arcpy.mapping.MapDocument('current')
dataframe = arcpy.mapping.ListDataFrames(MXD, "Web Mercator")[0]
dataframe.zoomToSelectedFeatures()
arcpy.RefreshActiveView()

When I publish my service, the zoom feature no longer works.  I have read online that you can't use 'current' in a geoprocessing service, that's fine.  When I enter the path to the MXD, the zoom function no longer works.


texts: Missing QGIS raster menu tools with error'"osgeo [python-gdal]" module is missing'?

texts:
I recently had to reinstall QGIS (version 2.2.0-9) and the raster menu is missing its tools. I understand that the GdalTools plugin must be installed and enabled. I have already installed the GDal complete framework package. When I try to load the Gdal Tools plguin, the following error message occurs:

A similar error message occurs when I go to 'manage and install plugins' and try to reinstall the openlayers plugin. The following message occurs:

Does anyone know how to install the "osgeo [python-gdal]" module mentioned?

I am using MacOsx 10.8.5.


texts:
I am trying to convert multiple raster layers with hundreds of files in each into web tiles.

qgis2web start processing, at "exporting layer1 to PNG" it shoots the CPU to 100%, hours later it crashes with python error.

here is the complete result:

FileNotFoundError: [Errno 2] No such file or directory: '/var/folders/p6/qvt9dldn13sckjk_2jf2hxf80000gn/T/ttn10rp01545831356_piped.tif'

Traceback (most recent call last):
  File "/Users/igzebedze/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins/qgis2web/utils.py", line 399, in exportRaster
    processing.run("gdal:warpreproject", warpArgs)
  File "/Applications/QGIS3.app/Contents/MacOS/../Resources/python/plugins/processing/tools/general.py", line 96, in run
    return Processing.runAlgorithm(algOrName, parameters, onFinish, feedback, context)
  File "/Applications/QGIS3.app/Contents/MacOS/../Resources/python/plugins/processing/core/Processing.py", line 139, in runAlgorithm
    raise QgsProcessingException(msg)
_core.QgsProcessingException: Unable to execute algorithm
Could not load source layer for INPUT: /var/folders/p6/qvt9dldn13sckjk_2jf2hxf80000gn/T/ttn10rp01545831356_piped.tif not found

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/Users/igzebedze/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins/qgis2web/maindialog.py", line 335, in saveMap
    feedback=self.feedback)
  File "/Users/igzebedze/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins/qgis2web/olwriter.py", line 91, in write
    folder=dest_folder)
  File "/Users/igzebedze/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins/qgis2web/olwriter.py", line 130, in writeOL
    popup, json, restrictToExtent, extent, feedback, matchCRS)
  File "/Users/igzebedze/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins/qgis2web/utils.py", line 225, in exportLayers
    exportRaster(layer, count, layersFolder, feedback, iface, matchCRS)
  File "/Users/igzebedze/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins/qgis2web/utils.py", line 401, in exportRaster
    shutil.copyfile(piped_file, piped_3857)
  File "/Library/Frameworks/Python.framework/Versions/3.6/lib/python3.6/shutil.py", line 120, in copyfile
    with open(src, 'rb') as fsrc:
FileNotFoundError: [Errno 2] No such file or directory: '/var/folders/p6/qvt9dldn13sckjk_2jf2hxf80000gn/T/ttn10rp01545831356_piped.tif'

Python version: 3.6.7 (v3.6.7:6ec5cf24b7, Oct 20 2018, 03:02:14) [GCC 4.2.1 Compatible Apple LLVM 6.0 (clang-600.0.57)]
QGIS version: 3.4.1-Madeira Madeira, exported
Python Path:
* /Applications/QGIS3.app/Contents/MacOS/../Resources/python
* /Users/igzebedze/Library/Application Support/QGIS/QGIS3/profiles/default/python
* /Users/igzebedze/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins
* /Applications/QGIS3.app/Contents/MacOS/../Resources/python/plugins
* /Library/Frameworks/SQLite3.framework/Versions/D/Python/3.6
* /Library/Frameworks/Python.framework/Versions/3.6/lib/python36.zip
* /Library/Frameworks/Python.framework/Versions/3.6/lib/python3.6
* /Library/Frameworks/Python.framework/Versions/3.6/lib/python3.6/lib-dynload
* /Library/Frameworks/Python.framework/Versions/3.6/lib/python3.6/site-packages
* /Users/igzebedze/Library/Application Support/QGIS/QGIS3/profiles/default/python
* /Users/igzebedze/Library/Application Support/QGIS/QGIS3/profiles/default/python/plugins
* /Users/igzebedze/Dropbox (Personal)/devel/zemljevidi/tiles


texts: How to display street names from OS StreetAny layer?

texts: Using QGIS version 2.10.1-Pisa.  I have imported layers from an Ordnance Survey data set. I have a layer called 'Street any' which gives me the boundaries of the streets which I can see.  I can see from the attribute table for the layer that the street names are there. How do I get the street names visible please. Long time since I did this and completely forgotten most of the little I once knew :(


texts: can't find search statistics in geonetwork

texts: I'm using Geonetwork and I need to show search statistics. In admin console, I went to "statistics & statut" but the problem is, I can't find "search statistics" inside.


texts: LiDAR: Seeking information about "Return Values"

texts:
I do not fully understand LiDAR Return Values, as opposed to Class Codes.

Both Arc and LAStools have options for selecting returns based these "values".  For example, here is the Arc 10.7.1 Make LAS Dataset Layer tool, which has a "Return Values" section (LAStools has similar Return Values):

Furthermore, ESRI documentation describes Class Codes, but I can't find anything in their documentation about Return Values.

Importantly, these "values" are not mentioned in the ASPRS LAS Specifications, which I assume (perhaps incorrectly) is the official repository of all things LAS.

I would like to learn more about these Return Values: Are there any published standards?  How are they derived?  Are their values consistent between different software packages?


texts: Calculating range of values for a class for landscape analysis

texts:
I am working on a large scale landscape analysis in Qgis: on several raster layers i made (using SRTM dem as source data) i am doing a reclassification procedure. Reclassification is based on a set of point samples: i use these points on each raster layer (slope, elevation, aspect etc), sample values and use general trends in results as guidance for range of values for reclassification process. Example: range of values is 0 - 66, and major number of sample points fall between 0 and 1 = class 1. However, rest of values are not represented so clearly - i.e. the change is gradual (1/6 of sample points are in the 1 - 2 range, 1/12 of points are in 2-3 range and so on and so forth). I would like to use the sourced values from point sample set to calculate ranges for classes in a way that is not arbitrary but more adapted to the total range of data and sample points distribution throughout the range. Can someone advise how to do this?
I am providing a link to sampled point set in excell that should illustrate data range and distribution: sampled points raster 1


texts: According to my watch that claims to measure altitude, as far as I can see it says that in Stockholm Sweden the altitude is 170 meters but can that be right? How exactly can we know that 0 altitude in Bangladeh is the same as 0 altitude in Stockholm?


texts:
I currently have a 87Mb Compressed GeoTiff with a resolution of 50m. I need to interpolate this to 2.5m

I have used gdalwarp to achieve this using a cubicspline, but it took 12hrs and resulted in a 44Gb GeoTiff. I did use -co COMPRESS=LZW but wondering if gdal_warp uses these creation options?

I then tried -r average which only took a few hours but still resulted in a 22Gb GeoTiff.

If I look at the results the image still looks very course (granular) even though gdalinfo does say it now has a pixel size of 2.5

Am i doing something wrong or expecting too much from gdal

thanks


texts: Using Arcpy to search geodatabase for feature with smallest number of rows and export as .SHP

texts: I have created a bunch of grids for a city boundary area. I would like to know how to use the Arcpy module to look through a geodatabase and find the grid with the smallest number of rows (assuming that these grids are all the same size and the fewest is the most efficient grid) and then export only that feature out as a .SHP.


texts: Any alternative for the plugin Grids for Atlas in Qgis 3?. I use it a lot because I always create atlas in my works but currently I have to create the coverage layer in qgis 2.18 and the continue my work in Qgis3.


texts: import qgis geoprocessing tools in python os x 10.9.5

texts:
I'd like to use the python bindings from qgis, but in a standalone project.
I need the geoprocessing tools for vectors.

I updated QGIS (Wien - 2.8), python (3.4.3), etc.
However, it always ends with the following error :

      import qgis.core
      Traceback (most recent call last):
        File "", line 1, in
        File "/Applications/QGIS.app/Contents/Resources/python/qgis/init.py", line 26, in
          import sip
      ImportError: dlopen(/Applications/QGIS.app/Contents/Resources/python/sip.so, 2): Symbol not found: _PyCObject_Type
        Referenced from: /Applications/QGIS.app/Contents/Resources/python/sip.so
        Expected in: flat namespace
       in /Applications/QGIS.app/Contents/Resources/python/sip.so

I set the following path :
export PYTHONPATH=/Applications/QGIS.app/Contents/Resources/python/
export DYLD_LIBRARY_PATH=/Applications/QGIS.app/Contents/MacOS/lib/:/Applications/QGIS.app/Contents/Frameworks/

The issue is obviously about SIP, but I truly don't know what to do.
The python console inside QGIS is not insulting me when I'm giving it "import qgis.core".
So I 'd bet there is something about the versions of SIP, Python "command line" and Python "QGIS command line".

I can't say what... and I'mm really stuck.

I am using miniconda3 for my libraries (rsgislib, gdal, ogr, ...).


texts: Creating many origin-destination routes with pgRouting

texts:
I'm using a PostGIS / pgRouting database setup with a road network. I'm trying to create an origin-destination table with the complete route information, so not just the sum of the length, like in this post.

I'm trying to get something like this as output table:

Origin | Destination | vertex_id | cost
1        2             415         14.2
1        2             508         4.3
1        3             919         2.4
1        3             1024        6.8

The following code for calculating a single route works fine, but how do I automate this for a lot of origin-destination pairs?

SELECT * FROM shortest_path('
SELECT gid AS id,
start_id::int4 AS source,
end_id::int4 AS target,
length::float8 AS cost
FROM network',
1, -- start_id
2, -- end_id
false, false);

I've tried something with subqueries, like this where I'm getting the start_ids from table nodes:

SELECT id,
(SELECT * FROM shortest_path('
SELECT gid AS id,
start_id::int4 AS source,
end_id::int4 AS target,
length::float8 AS cost
FROM network',
id::int4, -- start_id
450969, -- end_id
false, false))
FROM nodes;

But this doesn't work because the subquery returns multiple rows (as expected): ERROR: subquery must return only one column. Any help?


texts:
I need to accomplish two things in bulk:

Convert SHP to KML
Modify the fill and outline color and opacity either of SHP or KML files

What tools can you recommend for these two tasks?

I've tried Google Earth Pro, which allows me to import in bulk and set the fill/outline color/opacity in bulk.  But, I'm unable to export to KML in BULK.  For the collection of data I have right now I have 193 *.shp files.  In Google Earth Pro, I would have to export each one manually.  In the future I will have 1000's of *.shp files and one-at-a-time won't cut it then.

Can you bulk export to KML in Google Earth Pro?


texts: Calculate attribute height of point data

texts: I have an x (easting),  y (northing) and z (height) values associated with my point shapefile.  I would like to calculate the z values from the DTM I have - how do I do this in QGIS? I am fairly new to QGIS. I can't do it in ArcGIS as don't have the extensions.


texts: How to interpret Lidar scan line problems?

texts:
I am calculating various gridmetrics using the Lidar processing tools in FUSION.  The dataset I am working with contains 154 flight lines of high density data, at a nominal pulse spacing (NPS) of 1.5 meters (full metadata available here).  Specifically, I am using the grid metrics command in FUSION to calculate all returns above 2m (column 53).  The attached figure shows the problem I have been encountering.  You can see there are detectable grid values at the scan line overlap areas.  There are even detectable values in the non-overlap scan line areas over lakes (observed as white blob in screenshot).

How should I be interpreting these apparently faulty results--is it a result of the Lidar dataset, or is it some type of processing issue with FUSION?


texts: How to add a background map to a java project?

texts:
I'm new at GIS. I am trying to make a java application that can display a map with certain points. I have achieved this using the geotools library. However, I want to put the map from OpenStreetMap or Google Streets in the background also. I already know that PostGIS can be used with GeoTools. However I don't know where to start. How is the best way to start once I have installed PostGIS?

Hello, I have also come across KML. I'm sorry I'm a beginner. What KML parsing mean, and how can I use this to show a map in my Java Application?


texts: Using a grid code to visualize habitat values

texts:
I am working on wildlife and agriculture data in a region where I am attempting to layer the different types of habitat to show the total value of a parcel. For example, I have layers for endangered species habitat, migration corridors, and winter range and have scored each of these differently depending on the value they represent. If a parcel has endangered species habitat it scores as 11, migration corridors score an 8, and winter range scores a 5. I have taken these habitat layers and overlaid a grid layer for each of them but this process is new to me so I am wondering how I should display the values to then result in a map that shows the total value of a parcel. If a parcel has all three habitat types it should score a 24.
I am using ArcMap 10.8.1 and working with polygon data.
How do I get a map with the sums of the habitat types using a grid? And do these values go in the grid code in the attribute table?


texts: Iterating Search and Update cursor keeps inserting the same results for all rows

texts:
Ive created a script to find the start XY and end XY of a polyline feature. I know I could do it manually using "Calculate Geometry", but wheres the fun in that.

At the moment, for each row, its inserting the same start and end XY and i dont know why.

import arcpy
fc = 'UKT Unknown Owners with Pipelines'

with arcpy.da.SearchCursor(fc, ["OID@", "SHAPE@"]) as scursor:
    with arcpy.da.UpdateCursor(fc, ["Start", "End"]) as ucursor:
        for srow in scursor:
            for urow in ucursor:
                startpt = srow[1].firstPoint
                startx = startpt.X
                starty = startpt.Y
                endpt = srow[1].lastPoint
                endx = endpt.X
                endy = endpt.Y
                startxy = "{0}, {1}".format(startx, starty)
                endxy = "{0}, {1}".format(endx, endy)
                urow[0] = startxy
                urow[1] = endxy
                ucursor.updateRow(urow)


texts:
If anyone encountered the same problem and come up with the solution please give me an idea if what's wrong with this code.

The InfoWindow is not working when i click the map polygon.

All parts are working except for InfoWindow Pop up.

function initialize() {

var mapOptions = {
    zoom: 5,
    center: new google.maps.LatLng(12.63276, 122.27588),
    mapTypeId: google.maps.MapTypeId.ROADMAP
  };

  var map = new google.maps.Map(document.getElementById('map-canvas'), mapOptions);
  var arr = new Array();
  var polygons = [];
  var bounds = new google.maps.LatLngBounds();

  $.ajax({
    url: baseUrl + "/polygon/xml",
    async: !0,
    type: "POST",
    data: {},
    dataType: "xml",
    success: function(t) {
      var xml = xmlParse(getXmlString(t));
      var subdivision = xml.getElementsByTagName("Subdivision");
      for (var i = 0; i < subdivision.length; i++) {
        arr = [];
        var coordinates = xml.documentElement.getElementsByTagName("Subdivision")[i].getElementsByTagName("coord");
        for (var j = 0; j < coordinates.length; j++) {
          arr.push(new google.maps.LatLng(
            parseFloat(coordinates[j].getAttribute("lat")),
            parseFloat(coordinates[j].getAttribute("lng"))
          ));
          bounds.extend(arr[arr.length - 1])
        }
        polygons.push(new google.maps.Polygon({
          paths: arr,
          strokeColor: '#FF0000',
          strokeOpacity: 0.8,
          strokeWeight: 1,
          fillColor: '#FFF000',
          fillOpacity: 0.8
        }));
        polygons[polygons.length - 1].setMap(map);

        // xxx
        google.maps.event.addListener(polygons[polygons.length - 1], 'click', function(e) {
          var content = "
";
          content += "Latitude: " + e.latLng.lat() + "
";
          content += "Longitude: " + e.latLng.lng() + "
";
          showInfoWindow(e.latLng, content);
        });
        infoWindow = new google.maps.InfoWindow({});
        // zzz

      }
      map.fitBounds(bounds);
    }
  });
}

function showInfoWindow(latLng, content) {
  // alert(latLng);
  // alert(content);
  infoWindow.setContent('Rollie');
  infoWindow.setPosition(latLng);
  infoWindow.open(map);
}

function getXmlString(xml) {
  if (window.ActiveXObject) {
    return xml.xml;
  }
  return new XMLSerializer().serializeToString(xml);
}

function xmlParse(str) {
  if (typeof ActiveXObject != 'undefined' && typeof GetObject != 'undefined') {
    var doc = new ActiveXObject('Microsoft.XMLDOM');
    doc.loadXML(str);
    return doc;
  }

  if (typeof DOMParser != 'undefined') {
    return (new DOMParser()).parseFromString(str, 'text/xml');
  }

  return createElement('div', null);
}

var baseUrl = "http://" + document.domain + ":8080";
google.maps.event.addDomListener(window, 'load', initialize);


texts: Locate points that are within a certain distance of other points in QGIS

texts:
I need a way of locating all points in my dataset that are within a certain distance of any other point. For example, in the following data set in which EASTING AND NORTHING are British National Grid coordinates (to the nearest metre); I need to identify points which are within 5 meters of other points (ID #s 01 and 03 in this example). Can someone please suggest a QGIS technique for accomplishing this?

ID   EASTING    NORTHING
01   297736 703033
02   297991 702934
03   297735 703034
04   299488 702670
05   299476 762627
06   299443 702669
07   299164 702448
08   299210 702444


texts: Where is the server.dat equivalent at 10.1?

texts:
Prior to 10.1, I could locate the server.dat file, and programmatically find out the arcgisoutput folder URL location (virtual directory when using IIS)

e.g.

<Directory>
 <Path>c:\arcgisserver\arcgisoutput</Path>
 <URL>http://lap-303689/arcgisoutput</URL>
 <Description></Description>
 <Type>output</Type>
 <Cleaning>sliding</Cleaning>
 <MaxFileAge>600</MaxFileAge>
</Directory>

At 10.1, the arcgisserver folder default location is now:
C:\arcgisserver\directories\arcgisoutput

This is fine, but how do I find out what URL was set to in the post install?

I understand that there is now an embedded Tomcat within 10.1 (similar to Java edition of AGS), but I am unsure on which file contains the URL info (including the port number).

-EDIT-
Have reopened this Q, as have agreed with FredP that there should be a better approach to get the arcgisserver physical AND virtual paths from ArcObjects.


texts: Copy Features slowing down Model after many iterations

texts:
I would like to make catchments for many thousands of points using Network Analyst. The number of points is too big for one "Solve," so I am iterating through each point, solving and exporting the shapefile

I used Model Builder to start:

It works, but after a few thousand catchments it begins to slow down. Watching the geoprocessing window, I see "Copy Features" is the culprit: after each iteration it seems to hang on it longer and longer. The add location/solve seem to operate at the same speed. I tried adding a delete "in_memory" to the model but it didn't help. I would like to know what could be causing the "Copy Features" to drag, and if it is possible to fix it.

With model builder it is hard to figure out what is going on. I would like to switch to python but am stymied there as well (see this question)


texts:
This is in ArcMap 10.0.

When using data frames to label UTM coordinates, I'm sometimes asked to add leading zeros (e.g., 590 would become 00590) to the tick or line labels of the data frame.  I'd like to automate this.  I'd be happy using an appropriate style, exporting and editing an attribute table, or editing a shapefile--but I can't figure out how to even access this data outside of ArcMap.


texts:
I have a raster which I want to do a clip with a vector layer (polygon) using its atttribute values.
What I have done so far is the folliwing:

from qgis.core import *
from qgis.gui import *
import processing
from processing.core.Processing import Processing

VECTOR = QgsVectorLayer('MyShapefile.shp','MyShapefile','ogr')
RASTER = QgsRasterLayer('MyRaster','MyRaster')
row_info = VECTOR.getFeatures()
for row in row_info:
   row_name = row['NAME_FIELD']
   clip = processing.runalg("gdalogr:cliprasterbymasklayer", RASTER, VECTOR, '0', False, True, True, 2, 0, 1, 1, 1, False, 2, False, "", 'OUTPUT\Clip_{}.format(row_name)')

I have to use a 'select features' using NAME_FIELD value but I can't find how to iterate throw values, select each value , clip according to that value and go to the next value to clip again.

I am using Python Editor in QGIS 2.18


texts: Creating specific sub-study area in ArcGIS Desktop and classifying symboloy?

texts:
I am using ArcGIS Desktop 10.1.

My study area is quite large (sometimes spanning across a few grid squares) and when I am running a viewshed analysis or calculating slope or euclidean distance it calculates for the entire region. But at times I just want to concentrate on a very small area (e.g. within a 10 km radius).

In these cases is there anyway to make sure that these analyses just focus on the required area rather than the entire region?

When it calculates for a really large area I am having difficulties with symbology - I can't classify results into anything meaningful!


texts:
Currently I am working on a project where I have to draw boundaries of specific region of UK. Exactly like this

I want to know where to retrieve this kind of GIS data. I have tried openstreetmap data for uk, in that some of the data are unavailable for UK cities.


texts: Is there demand for higher resolution alternative to GPS tracking in GIS market?

texts:
If there was an alternative to GPS that provided higher resolution tracking, what applications would there be for it in GIS?

Note that this means higher resolution, but not necessarily higher accuracy. For instance, the absolute position (lat/long) would not be more accurate than GPS, but the relative movement during tracking could be within a few inches. This means you could get a very detailed track recording of a device's movement, allowing you to see small changes in position or elevation. If you were to carry this device while hiking a trail, you could see elevation changes with 1-5 inch resolution. I imagine this might be useful for GIS, but I am not familiar with the market.


texts: Display Table of contents of the map in ArcMap in the online published Map/Web map application

texts:
The following is the sample layer information of the map in ArcMap 10.1:

How can I have the same Table of contents in the published map online or in web map application, so that user can select/unselect layer's content?

Is there any feature I need to enable to have all contents appear for online user?
Do I need to implement JavaScript API?Is there any sample available?

Thanks


texts: Replacing nodata values by a constant using GDAL CLI tools

texts:
I have a raster DEM file containing 1 band: f='DEM.tif'
It has some NODATA values which I want to replace by a true, actual constant value (i.e. no more nodata stuff), let's say, zero.
Here are all the things I've tested:

gdal_translate; not working:

gdal_translate -of GTiff -a_nodata 0 ${f} ${f%.*}_cleaned.tif

well... it actually replaces the nodata value of 3.399999^38 to a nodata value of 0, but it doesn't replace it by a true 0 (which would be an actual zero value, no a nodata value).

I also tried the accepted answer here using gdalbuildvrt: Redefining nodata value into zero in QGIS? but with no success. Nodata values still show as ... well, nodata values.

Also, there is the gdal_fillnodata tool, which seemed promising: https://gdal.org/programs/gdal_fillnodata.html but it actually interpolates values from the closest borders, which doesn't make sense in my case.

Finally, gdalwarp was my last candidate, but it is not working better:

gdalwarp -srcnodata 3.39999995214436425e+38 -dstnodata 0 ${f} ${f%.*}_cleaned.tif

because, here again, it replaces the nodata value by another nodata value, which is not what I would like to achieve.
I've done that successfully in no time with QGIS using the "Fill NoData cells" on my local machine as a test, but I must work on my full dataset on a headless server, without QGIS installed.
So, do you know how to fill nodata values by a constant using any of the GDAL CLI (bash) based tools?  I hope I simply missed something obvious...
Edits

Removed as it was a malformed version of @user2856's working answer.

I forgot to mention this when I wrote the post, but I also tried @Mike T's first part of his answer:

gdal_edit.py -unsetnodata ${f}

(Disclaimer: this is a non reversible operation as it's actually working 'inplace', so backup your raster before running this command!)
Here, on the contrary, gdalinfo doesn't report NoData anymore. But do not rejoice too quickly, as a loading into QGIS shows:

The second part of Mike T's answer:
gdal_calc.py -A ${f}.tif --outfile=result.tif --calc="numpy.where(A>=3e34, 0, A)"

is actually giving a result.tif raster with:
gdalinfo result.tif | grep No
  NoData Value=3.40282346600000016e+38

That's funny? Probably some weird numerical issue under the hood.
But the most funny part is that if I first run gdal_edit.py -unsetnodata
and then the latter gdal_calc.py command, the resulting raster has so nodata again (heck, where do these actually come from?)!

General info on my setup just in case:
# gdalinfo --version
GDAL 3.3.0dev-33cf0e31a992be112b3091f012368d15605ed51d, released 2021/03/17

from osgeo/gdal:ubuntu-small-latest Docker image available at https://hub.docker.com/r/osgeo/gdal


texts: Smallest GIS library for embedded devices?

texts:
I need to use very simple geometry operations (read SHP or other standard format and simple searches) within an embedded device running a simple linux based os.

What is the smallest and simplest API I could use?

It could be written in C/C++ or Python in order to allow the integration with my system.

Thanks,

Samuel


texts:
I have the following code that works well to print out the individual values of pixels on a image (first band), using C-GDAL.

#include "gdal/gdal.h"
#include "gdal/cpl_conv.h"
#include <stdio.h>
#include <time.h>
int main()
{

    // to store execution time of code
    double time_spent = 0.0;
    clock_t begin = clock();

    GDALDatasetH raster;
    GDALAllRegister();
    raster = GDALOpen( "/home/roger/Documents/98_GDAL_C++_snippets/data/selected_files_l2a/S2A_MSIL2A_20180408T155531_N0206_R011_T17PLK_20180408T204658.SAFE/GRANULE/L2A_T17PLK_A014594_20180408T155557/IMG_DATA/R10m/T17PLK_20180408T155531_B02_10m.tif", GA_ReadOnly);
    if(raster==NULL)
    {
    printf("Invalid Raster Dataset");
    }
    else
    {
    GDALRasterBandH raster_band;
    int XSize, YSize;

    raster_band = GDALGetRasterBand( raster, 1);
    GDALGetBlockSize( raster_band, &XSize, &YSize);
    //int *pafScan;
    int nXSize = GDALGetRasterBandXSize (raster_band);
    int nYSize = GDALGetRasterBandYSize (raster_band);
    //int line = CPLMalloc(sizeof(int) * nXSize * 1);
    float *line;
    line = (float *) CPLMalloc(sizeof(float)*nXSize);
    int row, col;
    int value;

    // looping through the pixels of this line;
    for (row = 0; row < nYSize; row++) {
        GDALRasterIO(raster_band, GF_Read,0, row, nXSize,1, line, nXSize, 1, GDT_Float32, 0, 0);
        // this is where you'd do your processing on the pixel
        for(col=0;col<nXSize; col++){
            value = line[col];
            printf("%d\n", value);

        }

        }

        CPLFree(line);

    }
    clock_t end = clock();
    time_spent += (double)(end - begin) / CLOCKS_PER_SEC;
    printf("Time elapsed is %f minutes\n", time_spent/60);

    return 0;
}

Is there a way to only get the first pixel value(or a specific value), specifying the row and col coordinates in C?


texts:
My ArcMap 10.2 colour palette seems to be glitching out. What I think are the 'base' colours (those which are always available even when you don't have any style references) keep getting added. I think they might get added again each time I open arc, until the palette reaches capacity. I can add in style references such as the ESRI-optimized if there's room in the palette, but they eventually get overridden by the base colours. I've found I can wipe the palette by deleting my personal styles file but then the same thing happens again.

Does anyone have any idea how to fix this?

Edit: I think the colours that keep getting added might be the ESRI-CAD style? Always seems to be checked in the Style References window.


texts: Strange Polygon Shape Rendering in ArcGIS JavaScript API 3.22

texts:
When I upload a shapefile using the existing WAB Geoprocessing widget Shapefile uploader. Below is the screenshot. What could be the cause of the problem?
The image below was supposed to show the shape of Kenya. The issue happens in other polygons too.


texts: Using Checkboxes in Layer Tree of GeoExt?

texts:
I am trying to generate a customized layer tree using GeoExt 2.0.2
My code is as follows:

var store = Ext.create('Ext.data.TreeStore', {
            model: 'GeoExt.data.LayerTreeModel',
            root: {
                expanded: true,
                children: [
                    {
                        plugins: ['gx_baselayercontainer'],
                        expanded: true,
                        text: "Mappe di base"
                    },{
                        expanded: true,
                        text: 'Dati amministrativi',
                        leaf: false,
                        children: [
                        {
                         text: "Area progetto",
                         layer: "area_progetto",
                         nodeType: "gx_layer",
                         iconCls: 'area_progetto',
                         checked: true
                        },
                        {
                         text: "Comuni italiani",
                         layer: "comuni_italia",
                         nodeType: "gx_layer",
                         iconCls: 'comuni_italia',
                         checked: true
                        }
                       ]
                     },{
                        expanded: true,
                        text: 'Percorso Via Regina',
                        leaf: false,
                        children: [
                        {
                         text: "Tratte",
                         layer: "via_regina_tratte",
                         nodeType: "gx_layer",
                         iconCls: 'via_regina_tratte',
                         checked: true
                        }
                       ]
                  }
               ]
             }
        });

tree = Ext.create('GeoExt.tree.Panel', {
        border: true,
        region: "west",
        title: "Layers",
        width: 250,
        split: true,
        collapsible: true,
        collapseMode: "mini",
        autoScroll: true,
        rootVisible: false,
        lines: false,
        store: store
    });

My purpose here is to create multiple folders with the layers that I select and I want to assign different icons for each layer. This code created subnodes as I want, and I can assign a different icon for each layer in the tree, but the checkboxes are not working. Do you have any ideas or suggestions about why they are not working and how I can make them work?

NOTE: radio buttons for base layers are working

NOTE2: Setting the "store" as follows actually work for overlay layers, but in this case one big folder is created containing every layer except the base maps and I cannot create different folders and assign different icons for each layer.

   var store = Ext.create('Ext.data.TreeStore', {
        model: 'GeoExt.data.LayerTreeModel',
        root: {
            expanded: true,
            children: [
                {
                    plugins: ['gx_baselayercontainer'],
                    expanded: true,
                    text: "Base Maps"
                },{
                    plugins: [{
                        ptype: 'gx_overlaylayercontainer'
                    }],
                    expanded: true,
                    text: 'Overlay Layers'
                }
            ]
        }
    });


texts: Getting the middle of a way over Overpass API

texts:
I´m trying to import several Ways and Nodes from the OSM. But i dont want to get the Ways as a way with different nodes, only the middle of the way.

In my case, all ways are an area (Parks, Parking etc.) and i only need one point - the middle of the area. I could calculate this after the request, but maybe theres a more elegant way?


texts: Empty geometry error in GeoPandas dataframe

texts:
I want to reproject my ground truth data (land reference map) to UTM.

# Load Reference Data (reference map = ground truth data): labeled vector data in a shapefile (crops in Jendouba)
LULC_filename = "Jendouba.shp"
land_use_ref= gpd.read_file(str(LULC_filename))

# Convert CRS (WGS 84) to UTM
# CRS is the Coordinates Reference System
land_use_ref = land_use_ref.to_crs(Jendouba_crs.pyproj_crs())

#plot country
land_use_ref.plot()
plt.axis('off');
land_use_ref.crs

which raised the following error:

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
<ipython-input-87-fe96bab05714> in <module>()
      6 # Convert CRS (WGS 84) to UTM
      7 # CRS is the Coordinates Reference System
----> 8 land_use_ref = land_use_ref.to_crs(Jendouba_crs.pyproj_crs())
      9
     10 #plot country

4 frames
pandas/_libs/lib.pyx in pandas._libs.lib.map_infer()

/usr/local/lib/python3.6/dist-packages/shapely/ops.py in transform(func, geom)
    222     also satisfy the requirements for `func`.
    223     """
--> 224     if geom.is_empty:
    225         return geom
    226     if geom.type in ('Point', 'LineString', 'LinearRing', 'Polygon'):

AttributeError: 'NoneType' object has no attribute 'is_empty'

I checked the content of the dataframe:

print(land_use_ref.info())

And geometry had a missing value:

class 'geopandas.geodataframe.GeoDataFrame'>
RangeIndex: 249 entries, 0 to 248
Data columns (total 11 columns):
 #   Column      Non-Null Count  Dtype
---  ------      --------------  -----
 0   Type        249 non-null    object
 1   Observatio  100 non-null    object
 2   Crop_heigh  249 non-null    int64
 3   Previous_y  7 non-null      object
 4   Irrigation  28 non-null     object
 5   Sewing_dat  249 non-null    object
 6   Soil_type   99 non-null     object
 7   Soil_col_1  95 non-null     object
 8   FID222      249 non-null    int64
 9   Label       249 non-null    int64
 10  geometry    248 non-null    geometry
dtypes: geometry(1), int64(3), object(7)
memory usage: 21.5+ KB
None

So I tried to drop the row that has this missing value:

polygons_clean = Polygon[Polygon.geometry.type == 'POLYGON']

But I got the following error:

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
<ipython-input-97-062acdbf0eeb> in <module>()
----> 1 polygons_clean = Polygon[Polygon.geometry.type == 'POLYGON']

AttributeError: type object 'Polygon' has no attribute 'geometry'

This is

Edit: I used this:
land_use_ref['geometry']=land_use_ref['geometry'].dropna
and the problem persists and the weird thing is that the number of non-null samples become 249 (previously 248).


texts: Conflicting scale in layouts between AutoCAD and QGIS?

texts:
I have a recurrent problem.

Whenever I import data from AutoCAD (through dxf, properly referenced) into QGIS and then compose a layout at the same scale as the layout I had in AutoCAD, they appear to be at different sizes.

Not a huge difference, but definitely noticeable.

I checked, for instance, the size of a random line, and it's the same in both programs, the problem is only in the layout.


texts: How to use SAGA ShapesToGrid via Processing in a Python plugin?

texts:
I run these lines in my code:

vectorlayer_vector =  ftools_utils.getVectorLayerByName(inputLayer)
Elevation = self.lineAttrib.currentText()
cellSize = int(self.linePix.value())
outPath = self.inputLayerCombo3.text()
output = ftools_utils.getRasterLayerByName(outPath)

Processing.runAlgorithm("saga:shapestogrid", vectorlayer_vector, Elevation, 0, 0, 4, None, cellSize, output)

When I run, I get an tiff image but with nothing. The .aux file is not created and when I try to open the tiff image that I created appears this message:

c:...teste3.tif is not a supported raster data source

So, I guess that the output is wrong. What I have to do to get a correct output? What I missing?

I update my code but I have no results yet.

Please see if you can help me.

My code are:

inputLayer = self.inputLayerCombo.currentText()
Elevation = str(self.lineAttrib.currentText())
cellSize = int(self.linePix.value())
outPath = self.inputLayerCombo3.text()

Processing.initialize()
Processing.runAlgorithm("saga:shapestogrid", inputLayer, Elevation, 0, 0, 4, "10000.0,14000.0,50000.0,54000.0" , cellSize,outPath)

Is a simple code but I don't have any results. I don't know what else to do.
It doesn't create a raster.


texts: How do I rotate a coordinate system (ArcGIS/GDAL driver)?

texts:
I am trying to get oblique tiles into a mosaic dataset. The tiling scheme works perfectly well for vertical and north imagery (using google tiling scheme) through the GDAL driver, however I am struggling to get 45 degree tiles in.

The tiles are indexed very similarly to google tiles (starting at top left) as per example below (East orientation):

The crux of the problem is the CRS as defined by the GDAL XML. This is what it looks like for vertical imagery:

<GDAL_WMS>
<Service name="TMS">
<ServerUrl>http://example.com?x=${x}&amp;y=${y}&amp;z=${z}</ServerUrl>
</Service>
<DataWindow>
<UpperLeftX>12523442.7125</UpperLeftX>
<UpperLeftY>-1878516.40687</UpperLeftY>
<LowerRightX>13149614.8481</LowerRightX>
<LowerRightY>-2504688.5425</LowerRightY>
<TileLevel>21</TileLevel>
<SizeX>8388608</SizeX>
<SizeY>8388608</SizeY>
<TileX>1703936</TileX>
<TileY>1146880</TileY>
<YOrigin>top</YOrigin>
</DataWindow>
<Projection>PROJCS["WGS_1984_Web_Mercator_Auxiliary_Sphere",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433],AXIS["Lat",SOUTH],AXIS["Long",EAST]],PROJECTION["Mercator_Auxiliary_Sphere"],PARAMETER["False_Easting",0.0],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",0.0],PARAMETER["Standard_Parallel_1",0.0],PARAMETER["Auxiliary_Sphere_Type",0.0],UNIT["Meter",1.0],AUTHORITY["EPSG",3857]]</Projection>
<BlockSizeX>256</BlockSizeX>
<BlockSizeY>256</BlockSizeY>
<BandsCount>3</BandsCount>
<UnsafeSSL>true</UnsafeSSL>
<MaxConnections>10</MaxConnections>
</GDAL_WMS>

My theory is that I need to rotate the coordinate system so that, in the example of East orientation, Y axis gots towards East, and X axis goes towards South. After trying a few things I could not get it to work. The best candidate was an affine transform as described here https://gis.stackexchange.com/a/159134/52713, however that just resulted in an invalid dataset (using identity transform or no-op transform), showing up in ArcMap like so:


texts:
I've recently began trying out basic geoprocessing tools by using arcpy (clipping layers, buffering, etc).

However, this seems slower than using the ArcMap GUI.

Are there any advantages to using ArcPy for simple tasks that can otherwise be done directly through the GUI?


texts: network distance from subjects to park within buffer

texts:
I have multiple subjects with buffers of 0.5 miles and parks in the city. I want to calculate a network distance between subjects and parks from each subjects with its buffer area. In some buffers there are 3 or more parks, so I want to know the distance of all the parks within the buffer. Is there any tool in ArcGIS 10 or higher and or already build model that I can use?
In figure below, circles are 0.5 miles buffer of a subject (blue dots), red lines are the street, and colored polygons are the parks. I want to know the distances of each park within the polygon from the subjects (dots). I have more than 100 subjects and their buffer. Any Idea to do easily?


texts: Missing a feature after digitalization of a geopackage vector in QGIS

texts:
I am trying for the first time the geopackage format and I am facing a weird behavior, for me at least.

I hade a file with two polygon features and I have added a third one by digitatizing. After editing I can see the third feature in the attribute table but not in the canvas. Performing a check validity returns a layer called valid_output with the three features.

What could it be wrong?

Here the link of the file: http://s000.tinyupload.com/?file_id=09465121551351065440

I also see a "[3]" in the layer name that puzzels me.

What is it?


texts: Removing one to many links when performing spatial join in GeoPandas

texts:
I am spatially joining census data to building data using GeoPandas spatial join gpd.sjoin(left_df=gdf_building,right_df=gdf_census_data, how="left", op="intersects"). In the results, I observed that if a building is intersecting two census blocks, GeoPandas generates two copies of the census properties in the building shape. How can I remove one of the duplicates and only keep the data corresponding to the census block with the largest area of intersection?
Red polygons shows buildings and the yellow ones shows census blocks:

This builidng intersects two polygons (census blocks)

Two separate links in the building attribute


texts: How can I turn a python list into a table?

texts:
I am trying to create a table of coded values so each textual description has a numerical attribute. The order of numbers is not important - they just need to be unique.

What I have so far is a python list of all unique values, so all I need to do now is add this to a table and create numbers. Is there a method in arcpy that will allow me to save my list as a table with an objectID?


texts: Concatenating export path in Google Earth Engine

texts:
My goal is to store different exports from Earth Engine projects in separate subfolders of a main GEE_export folder on my Google Drive. I'm using  Export.image.toDrive:

Export.image.toDrive({
    image: export_img,
    folder: export_path,
    fileNamePrefix: filename,
});

However, if I pass a path like GEE_export/subdir or GEE_export\subdir or even GEE_export\\subdir it just creates a folder with a corresponding name in my Drive root folder.

What is wrong?


texts:
I want to use Time Manager plugin with QGIS 2.18.9 and unable to correctly format the date and assign it to a field within a table.

I have an imported table with a field called "Year Built" containing a year value as a long integer (looks like "2015"). I want to use that value in Time Manager. Time Manager needs Gregorian format YYYY-MM-DD HH:MM:SS.ssssss.

I have tried a number of ways to convert it in the expression dialog and place it in a new field called TM_Date. While in the Expression Dialog the "Output Preview" shows . But after clicking on OK and UPDATE ALL, the result in the field is a NULL value.

My last attempt was the following:

to_datetime( format_date(concat(to_string( "Year Built") ||'-01-01 01:01:01.001'),'yyyy.MM.dd HH.mm.ss.zzz'))

Not sure what I'm doing wrong here.


texts: ArcGIS points to raster conversion with multiple fields in the attribute table

texts: I have 2,000 points each containing fields like day, month, year, time, brightness, name, serial, etc. But when I try to convert them to raster, only 1 field will be considered. Is there a way I can convert it to raster containing at least 4 fields? I use ArcGIS in the conversion


texts: OpenStreetMap Style Editor besides Mapbox TileMill and CloudMade

texts:
I'm looking for an online OSM style editor that lets me style OSM ways and nodes by tags. I know CloudMade's editor is no longer available, and Mapbox's does this but has zoom level rendering limits (eg, won't render buildings below zoom level 13).

My goal is to get a map of the USA (or world if I can) that shows where every OSM building has been mapped (ways with building=).  And a second map that shows where every OSM building + street address is (ways with addr:housenumber=).

Here's an example map I did using a state OSM extract in Kentucky, processed with Osmosis, and QGIS to render:

https://twitter.com/CivicDataAlly/status/644963898463023108

But unfortunately the North America extract is too large for me to work with, and QGIS can't even open the large Osmosis processed file on my machine.  I also set up osm2pgsql and PostgreSQL/PostGIS databases, but again the file sizes and processing times are too huge.

Is there an online tool where I can style buildings at a low zoom level?

SOLVED

I finally got the image I needed of all the buildings and buildings with street addresses for all of North America, by using the info in the answer below.

I basically setup postgres and used pgAdmin3 to manage it and then osm2pgsql in the command line with the correct options and styles to load the raw .pbf export into the database, then connect to it with QGIS to visualize it.

Once the data was in postgres it was easy and fast to deal with, even displaying it in QGIS on my modest computer.

Buildings

Here are the stats and what I did for the building footprint layer:

osm2pgsql -U mschnuerle -d mschnuerle --create --slim --flat-nodes nabuildings.bin --cache 18000 --number-processes 3  --style buildings.style  north-america-latest.osm.pbf;

Processing: Node(794783k 1292.3k/s) Way(52669k 32.51k/s) Relation(471970 89.52/s)  parse time: 7507s (2 hours)

node cache: stored: 794783721(100.00%), storage efficiency: 66.31% (dense blocks: 695842, sparse nodes: 243029485), hit rate: 100.00%

Osm2pgsql took 88452s overall (24 hours)

Backup DB size: 4GB

with only this in the buildings.style file (everything else deleted):

way         building           text    polygon
way         way_area           real    linear

Addresses

And here is what I did for only the buildings that had street addresses attached to them:

osm2pgsql -U mschnuerle -d mschnuerle --create --slim --drop --flat-nodes naaddress.bin --cache 6000 --number-processes 4 --verbose --exclude-invalid-polygon --prefix address --style address.style  north-america-latest.osm.pbf;

Processing: Node(794783k 1271.7k/s) Way(52669k 31.67k/s) Relation(471970 37.13/s)  parse time: 14998s (4h)

Osm2pgsql took 17455s overall (5h)

with this in the address.style file (everything else deleted):

way    addr:housenumber       text    polygon
way    way_area               real    linear
way    building               text    delete

And my machine specs are this, for those asking:

Mac Specs

OS X El Capitan 10.11 Beta
iMac 27" Mid 2011
3.1 GHz Intel Core i5 (4 cores)
25 GB 1333 MHz DDR3
1TB SATA HD

This is the only map online that shows building footprints at this level, as far as I know.  If you find any other ones, please let me know here.

I also plan to create all the Mapbox tiles needed using Tippecanoe and post the map online over the next few weeks.


texts: Converting from EPSG:4978 to EPSG:4326 using pyproj always produces 0.0 latitude

texts:
I am trying to convert between EPSG 4326 and 4978. I have wrote the following codes:
import pyproj
lon = -74.6604194
lat = 40.3450248
transformer = pyproj.Transformer.from_crs("EPSG:4326", "EPSG:4978", always_xy = True)
print(transformer)
xy = transformer.transform(lon, lat)
print('Projected XY', xy)
lonlat = transformer.transform(xy[0], xy[1], direction='INVERSE')
print('Lon Lat', lonlat)

This is the result it produces. As you can see it converts the lon/lat to the right xy. However, when I try to convert xy back to lon/lat, I always get a 0.0 lat value. I have tried with other xy values and I always get a 0.0 lat. What seems to be wrong with my code?
proj=pipeline step proj=unitconvert xy_in=deg xy_out=rad step proj=cart ellps=WGS84
Projected XY (1287775.7028415564, -4694570.895022585)
Lon Lat (-74.6604194, 0.0)


texts: Querying/getting EID ( Edge id's of network ) in ArcGIS API for JavaScript?

texts:
I did a server-side(c#) network analysis and the output is the EID's(NetEid's) which needs a INETelements to query for features and Object Id's.
It takes long time especially when you found many features.
Is there anyway to send the Eid's to Java-script( since I have to mark them and ...) and query through the Eid's on client side or any quicker way to do the query in server-side( maybe just for one layer? )   ?

My code :

 INetElements INET2   = ipNetwork as INetElements;

     ipTraceFlowSolver.FindFlowElements(oFlowMethod, oFlowElements, out m_ipEnumNetEID_Junctions, out m_ipEnumNetEID_Edges);
           _ResultEdges = m_ipEnumNetEID_Edges;
       for (int i = 0; i < Ipcount; i++)
                    {
                         var Eid= m_ipEnumNetEID_Edges.Next() ;

                         INET2.QueryIDs(Eid,esriElementType.esriETEdge,out UserClassIDs,out UserID,out userSubId);
    .
    .
    .


texts: Android application development by using OpenStreetMap

texts:
I have managed to install my own tile server for one country and we can see the map, it's working.

The END Goal is, we will build a vehicle tracking software by using our tiles.
Currently we are on learning stage as we never worked with OSM or its API. Since we are going to host everything in-house, I just wanted to get a basic concept,

From an administrator point of view, what tools I need to configure to allow our developer to work with OpenStreetMap for vehicle tracking?

I refer to this site https://wiki.openstreetmap.org/wiki/Develop, so do just Need to install Rubyport? That will allow developer to start working? I already have Tiles working.


texts: Lower distance limit in Regular Points tool

texts:
I need to define a regular grid of points in QGIS, which later I'll edit further based on its interaction with some data. Up until now, I was using the "Vector Grid" of the same Vector Research Tools, because by using "Regular Points" I had some complications. Then, I simply took the data of the polygons and computed their center later in Excel.

A few days ago, I managed to get the desired interaction with Regular Points, which I'd prefer to work with. The only problem is that, for some reason, QGIS won't let me select a distance low enough for my needs. I'm attaching a screenshot of the tool:

The CRS I'm using is WGS84 (EPSG: 4326, WGS 84). I would need a distance in degrees of 0,000010, which is what I was using before with the polygons, but the value in the screenshot is the lowest it will go. The points should be closer than 1.5 meters, but exact distance is not required (i.e. both 0.5 meters and 1 meter are fine).

Any ideas as to how to lower the value? I don't get why this tool has a lower limit, while the polygon one doesn't. Sorry if this is trivial, I searched for a similar question but found none.


texts:
I'm working with a large dataset and I prefer to handle it within a postgres database with postgis. The data is looped in python into a long string (or multiple if one string gets too big). I would like to add a lot of data if the data (in my case column "pos") is within a polygon.

The dataset is much bigger, hence this is an example of what I would like to do:

INSERT INTO I_test_test (pos, row_id)
VALUES(ST_PointFromText('POINT(13.55259726 55.39571705)', 4326), 1,
       ST_PointFromText('POINT(13.55387947 55.39514755)', 4326), 2)
WHERE st_intersects(pos, ST_GeomFromText('POLYGON((13.5542 55.3956,13.5543 55.3952,13.5586 55.3949,13.5587 55.3954,13.5542 55.3956))',4326))

I could write an sql question and check every "pos" before I add them into the INSERT INTO, hence since it is 10 000 - 1 000 000 rows that I would like to add it would take forever...

Any suggestions?


texts: Multiple Json files conversion to feature class using ArcPy

texts:
I am a beginner in Python.
I am migrating from visual basic.
I have multiple .json files that I want to convert them all in one script using python to a file Geodatabase:
What is the issue in the following code
    # Esri start of added imports
    import sys, os, arcpy
    import arcpy
    import os.path
    import os

   try:
for x in range(100):

    gdb=arcpy.GetParameterAsText(0)
    arcpy.env.workspace = "C:\buildings_cabinets"
    arcpy.JSONToFeatures_conversion("C:\buildings_cab\ x + ".json",
    os.path.join("C:\buildings_cab\yu.gdb", x + ".json))

 except Exception as e:
    arcpy.AddMessage(e.message)


texts:
I'm using mongo to store some geolocalized data, with the goal of retriving them using $geointersect. Specifically I have this document stored in my db:

  {
    "loc": {
      "geometry": {
        "type": "Polygon",
        "coordinates": [
          [
            [
              -179.875,
              -89.875
            ],
            [
              179.875,
              -89.875
            ],
            [
              179.875,
              89.875
            ],
            [
              -179.875,
              89.875
            ],
            [
              -179.875,
              -89.875
            ]
          ]
        ]
      },
      "crs": {
        "type": "name",
        "properties": {
          "name": "urn:x-mongodb:crs:strictwinding:EPSG:4326"
        }
      },
      "type": "Feature",
      "properties": {}
    },
    "_id": "576af8e31d41c87fa1f1d04f"
  }

Which as you may notice covers almost the entire World. Now it doesn't matter what coordinates I input for my $geointersect query, it will never be returned... Does anyone know why?

An example of the query I'm using could be:

[
  {
    "loc.geometry": {
      "$geoIntersects": {
        "$geometry": {
          "type": "Polygon",
          "coordinates": [
            [
              [
                13.4307861328125,
                41.599013054830216
              ],
              [
                13.9801025390625,
                41.599013054830216
              ],
              [
                13.9801025390625,
                41.80407814427234
              ],
              [
                13.4307861328125,
                41.80407814427234
              ],
              [
                13.4307861328125,
                41.599013054830216
              ]
            ]
          ]
        }
      }
    }
  },
  {
    "loc": 1
  }
]

Of course the coordinates of the query are contained inside the polygon in the saved document (I mean... it's as big as the the world, duh!) but for some reason it doesn't find any match... I'm kinda lost.


texts: Convert classified point matrix to shapefile

texts:
I've use R to generate a matrix of points tagged to separate categories. Each point has a lat, long. and category. When I plot them in QGIS (using colored points of the appropriate size) it looks like a territory map - which is the application. The points cluster together in discernible market areas.

I would like to use R to generate shapefiles which follow the color boundaries and which contain the appropriate classification label. I would use those market boundary shapefiles in a variety of other applications.

The size of the matrix is 10000 x 10000.


texts: File format for specifying GPS motion as a function of time?

texts:
What file format can specify GPS lat/lon/elev time history?
If I want to specify vehicle or pedestrian motion along a waypoint path, what is a compact way of specifying this?
.kml has waypoints but no time.
.gpx might work, but how?
This seems to be a glaring omission in file formats.


texts: I ma begining with vector tiles and geo data in general. I am running a T-rex server to serv data post a postGIS database. i have an issue with my browser part because i can't display anything ! Something should be wrong with my .toml config file.


texts:
I'm new to QGIS, and am trying to figure out if there is a way to view metadata within QGIS 3.10.1.

I'm adding layers downloaded from various state and federal websites, and the shapefiles include a .xml file (ISO-8859-1) - when opened with Notepad, it includes all the info I am looking for (namely, what the attribute names mean). It seems like there should be a way to view it from within QGIS - I found the Metadata tab under Layer Properties but nothing is there. At the bottom is an option to "Load Metadata", but when I browse to the .xml file associated with the shapefile it is greyed out.

I can copy and paste text into the Metadata tab from the file I have open in Notepad, but that would take way too long!

Any suggestions?


texts: Compare Two Tables and Delete matching Records based on a Field in Original Table

texts:
I am having an issue with this and think there is a simple solution.

I have a large table, upwards of 500,000 records and need to add new records to it daily. Occasionally a new record will just be an update to a field or two and I want the record from the new table to replace the old record with the same Unique ID.

I have tried appending the new (smaller) table to my original (large) table and then running a delete identical, however sometimes the old row is kept and the new record is deleted.

I would like to be able to compare the tables based on Unique ID, delete all identical ID's in my original table and then append the new table.

Either in ModelBuilder or python as this needs to be done daily.


texts:
I am working in QGIS and have been wanting to do a simple Least Cost Path Analysis using the SAGA tool. I have created a Slope raster for it, as well as an Accumulated Cost Surface. I have two points which I have imported from Google Earth and transformed into two separate point layers. I have tried numerous times to run the plug in, but I keep getting the same error:
ERROR: The start-point layer contains no legal point.
I tried to re import the points, change the CRS, researching on the internet and watching tutorials, but nothing seem to work or to give a reason for this issue. Is there a way to fix this?


texts: OpenLayers 2.12 + KML from String + Change markers'icons

texts:
Reading this post How to add KML data but from variable - not from url?
I changed the way of reading KML files (from HTTP to a string loaded using jQuery)

/**
 * The OpenLayers object
 * @type {OpenLayers.Projection}
 * @var
 */
var geographic = new OpenLayers.Projection('EPSG:4326');

/**
 * The OpenLayers object
 * @type {OpenLayers.Projection}
 * @var
 */
var mercator = new OpenLayers.Projection('EPSG:900913');

****
success: function(data) {
    if (data != '') {
        var format = new OpenLayers.Format.KML({
            'internalProjection': mercator, //EPSG:3785/900913
            'externalProjection': geographic, //EPSG:4326
            extractStyles: true,
            extractAttributes: true
        });
        sFeatures = format.read(data);
    }
}

I would like to use the icons that come from the KML string, but when I click on a feature I would like to increase the size of the icon or even change it.

My Layer definition:

var hydrographsLayerVector = new OpenLayers.Layer.Vector('hydrographs', {
    //renderers: renderer,
    styleMap: new OpenLayers.StyleMap({
        "default": new OpenLayers.Style(OpenLayers.Util.applyDefaults({
            externalGraphic: "http://openlayers.org/dev/img/marker-green.png",
            graphicOpacity: 1,
            rotation: -45,
            pointRadius: 10
        }, OpenLayers.Feature.Vector.style["default"])),
        "select": new OpenLayers.Style({
            externalGraphic: "http://openlayers.org/dev/img/marker-blue.png"
        })
    })
});

If I remove the lines extractStyles and extractAttributes = true I see the effect using the styleMap defined in the new OpenLayers.Layer.Vector('hydrographs'*.

How can I use the styleMap:"select": if I use extractStyles and extractAttributes = true?
Ideally, I would like to increase the size of the feature onmouseover and if the user clicks on it I want to keep the bigger size until the "unFeatureUnSelected".

Thanks!


texts:
I have features with multiple themeID, for example, themeID=49 and themeID=53, and when I load them into qgis it shows attribute themeID=(2:49,53). The problem is that I use themeID in sld to choose what svg sign to add to it, and now I can't find a way to use this. It would nice to solve this with use of sld but to leave possibility if there is only one themeID (than it won't be in form of(1:49))

feature in gml

<app:poMember>
<app:po gml:id="point1">
    <app:idElement>1001914</app:idElement>
    <app:desc>null</app:desc>
    <app:status>1</app:status>
    <app:markPKO>602346</app:markPKO>
    <app:themeID>49</app:themeID>
    <app:themeID>53</app:themeID>
    <app:xPoint>4963127.35</app:xPoint>
    <app:yPoint>6436882.64</app:yPoint>
    <app:zPoint>0.0</app:zPoint>
    <app:geoTopoMember>
        <gml:Node gml:id="n2">
            <gml:pointProperty>
                <gml:Point gml:id="p1" srsDimension="2" srsName="urn:ogc:def:crs:EPSG::31276">
                    <gml:pos>4963127.35 6436882.64</gml:pos>
                </gml:Point>
            </gml:pointProperty>
        </gml:Node>
    </app:geoTopoMember>
</app:po>

sld that I use

<Rule><Name>DPoint</Name><Title>DPoint</Title><ogc:Filter> <ogc:PropertyIsEqualTo><ogc:PropertyName>themeID</ogc:PropertyName><ogc:Literal>49</ogc:Literal> </ogc:PropertyIsEqualTo> </ogc:Filter> <PointSymbolizer><Graphic><ExternalGraphic><OnlineResource xlink:type="simple"   xlink:href="/DPoint.svg"/><Format>image/svg+xml</Format></ExternalGraphic><Size>16</Size> <Rotation>0</Rotation></Graphic></PointSymbolizer></Rule><Rule><Name>BPoint</Name><Title>BPoint</Title><ogc:Filter> <ogc:PropertyIsEqualTo><ogc:PropertyName>themeID</ogc:PropertyName><ogc:Literal>53</ogc:Literal> </ogc:PropertyIsEqualTo> </ogc:Filter> <PointSymbolizer><Graphic><ExternalGraphic><OnlineResource xlink:type="simple"    xlink:href="/BPoint.svg"/><Format>image/svg+xml</Format></ExternalGraphic><Size>16</Size> <Rotation>0</Rotation></Graphic></PointSymbolizer></Rule>


texts: Load vector layer directly from mongodb with PyQGIS

texts:
I want to load my vector layer directly from a MongoDB database using Python in QGIS. I found code here but it dumps the layer before loading it. Here is the code I found:
client = MongoClient()
client = MongoClient('localhost', 27019)
db = client['pb']
cursor=db.section.find()
feature = QgsFeature()# create layer

os.chdir( "C:\Users\Public" )
feature =[]

for data in cursor:
    geom=(data['geometry'])
    id  = (data['properties'])

    feature.append({"geometry":geom, 'properties':id})

json_str = json.dumps(feature)
final={"type": "FeatureCollection","features":feature}

with open('data_3.json', 'w') as f:
     json.dump(final, f)

vlayer = QgsVectorLayer('C:\Users\Public\data_3.json|geometrytype=LineString',"mygeojson","ogr")
QgsMapLayerRegistry.instance().addMapLayer(vlayer)

But I want to know if it is possible to load directly the layer from database without saving it on a file? And how to do it?


texts: Geometry column isn't visible in the attribute table in QGIS while it is in Python?

texts:
This code turns a dissolved shapefile to singlepart but the geometry column although it can be shown in Python when I exported as a shp and imported in QGIS the geometry column is not there. What should I do?

import geopandas as gp
import pandas as pd

def multi2single(dissolved):
    gpdf_singlepoly = gpdf[gpdf.geometry.type == 'Polygon']
    gpdf_multipoly = gpdf[gpdf.geometry.type == 'MultiPolygon']

    for i, row in gpdf_multipoly.iterrows():
        Series_geometries = pd.Series(row.geometry)
        df = pd.concat([gp.GeoDataFrame(row, crs=gpdf_multipoly.crs).T]*len(Series_geometries), ignore_index=True)
        df['geometry']  = Series_geometries
        gpdf_singlepoly = pd.concat([gpdf_singlepoly, df])

    gpdf_singlepoly.reset_index(inplace=True, drop=True)
    return gpdf_singlepoly

here is the call:

multi2single(dissolved)


texts: how to remove small water bodies in binary mask array that we get from ndwi values?

texts:
I've detected water bodies in satellite images (modis, sentinel, landsat) in special area polygons, by ndwi values and based on special thresholds. area polygons and image tiles are coordinate restricted, so there is only one main water body in image (as you can see in the picture below). i need to identify water body borders so i defined a binary mask that represents water body area (1s) and land area (0s). the problem is, there are some small water bodies (as noises) in the image. i want to get access to them and turn their values in binary mask array from 1 to 0, without effecting the main water body. i found some codes can count such like binary areas base on adjacency, but i dont have access to those areas to change the values!i need to access those areas seperately as an distinct objects (or something like this).

what should i do?


texts: Deleting holes of contour polygon in QGIS

texts:
I'm trying to create hypsometric areas from a DEM. My current workflow in QGIS is:

Extracting contours with SAGA, as I need contours from specific elevations (250, 500, 1000...) and not just every x meters.
Converting the contours from multipart to singlepart, and then to polygon after fixing geometries.

The problem is, the resulting polygon has other polygons inside it that should be holes of the first polygon (image). I can't seem to find a way to "substract" them, and I've been searching all over Google for days. Do you have any ideas?


texts: Refresh of a customized WMS-Layer (with cql-filter) generates a new map instances on website

texts:
I've two simple problems, but I can't fix it by my own :o(

My OL-map features a basemap and a cql-filtered WMS from
my local geoserver (topp:states).

The map is just shown, when i click the filterform
-> how can i load the map befor filtering?
the filter wors fine but when i change the parameter and click
the filterform again a new map-instance comes up in my website!
-> how can i reload (just) the map with the filtered wms?

Any Help? I'm so thankfull for a "code-supplement" of my code:
this is my code:

document.getElementById("submitBtn").addEventListener("click", getFilteredMap);

function getFilteredMap(){

    var a = document.getElementById("persons").value;
    var b = document.getElementById("cars").value;
    var c = document.getElementById("penners").value;

var filterstring = 'PERSONS>' + a + 'AND CARPOOL>' + b + 'AND UNEMPLOY>' + c;

    var layers = [
        new ol.layer.Tile({
            source: new ol.source.MapQuest({layer: 'osm'})
        }),
        new ol.layer.Image({
            source: new ol.source.ImageWMS({url: 'http://localhost:8082/geoserver/wms',
            params: {'LAYERS':'topp:states', 'cql_filter': filterstring},
            serverType: 'geoserver'})}
        )
    ];

var map = new ol.Map({
    layers: layers,
    target: 'map',
    view: new ol.View({
        center: [0, 0],
        zoom: 2
    })
});}

I think i can't make it the layer variable a global variable beacuae
i neet to compute the new cql-statement wich is a part of my layer-variable.

Maybee i'm wrong ?!

can you correct my code snipplet?


texts: Converting grid_code values to lat and longs or UTMs using ArcGIS Desktop?

texts:
I have made a text file from raster to points via a DEM in ArcGIS Desktop 10.1 where I have inserted values in an another column that I want to display. Next to point ID is grid_code, which I believe substitutes for lat and long.

How do I convert these grid_code values(which are z-scores of the study area) to lat and longs or UTMs?


texts: ArcGIS Tool Script Validation...browser window default location for multivalue option box

texts:
I am trying to do some Tool Validation for a script I'm writing in ArcGIS 10.2. I have a Multivalue File field which the user can browse to a directory to select the files they require to run. When opening the browser I would like it to automatically browse to a specific folder where the files reside so it's quicker/easier for the use to select the files required.  I've looked though the http://pro.arcgis.com/en/pro-app/arcpy/geoprocessing_and_python/programming-a-toolvalidator-class.htm documentation but can't see anything related to achiving this.

Any Thoughts on if this is possible?


texts:
I have coordinates in *GPS (WGS84) (deg) format and I need convert them to South African Cape/Lo27 format. I know I can do this on this website (https://twcc.fr/en/#), but unfortunately I need to convert a few hundred coordinates, are there any conversion formulas that I can use?

Here is an example of my data:

What I have:

Longitude: -32.9220

Latitude: 27.7586

What I need:

X: -70991.542

Y: 3644059.052

Convergence: 0.4123


texts: Method to find "if map data is available or not"

texts:
I am using Arcgis javascript api 3.2. I use the mapservice,

 https://server.arcgisonline.com/ArcGIS/rest/services/World_Street_Map/MapServer.

I realized that map data is not available when the scale(slider) is zoomed to the fullest at certain places on the map. Is there a method to find out if map data is available or not? Could you please share the steps in doing so.


texts: Is it possible to reproject KML using ESRI's Javascript API?

texts:
Specifically into NAD83.

Google tells me that KML only supports WGS84, but I can't seem to find much info on being able to reproject KML using ESRI's Javascript API.

Basically, I've got this base map in Nad83 and I need to overlay some kml.


texts:
I am trying to compare the contents of two columns named C1 and C2, each of them having some 250 values. I want to create another field, named C3, with TRUE if a value from C1 is found anywhere in C2, or FALSE, if not.

I found something which is very close, but it works only for the same entry in each column.

I am using QGIS 2.12 and I am curious if there is any SQL expression that I can use.

edit: C1 is called NAME, C2 is DIST_ORIG


texts: How to rename fields from postgis when using ogr2ogr to create a GPX file?

texts:
I want to export a GPX file from my postgis database. Luckily I have a name field in postgis, so the 'name' field automatically maps directly over. But I want to export another field into the GPX comment field. How do I do that?

I tried using a -sql statement with 'select name, change as cmt', but when I did that, ogr2ogr lost the geometry field: ERROR 6: Geometry type of `None' not supported in GPX.

I tried adding the geometry field to the sql statement and got ERROR 6: Cannot create GPX layer sql_statement with unknown geometry type. Adding an explicit -geomfield command didn't help - I got the same error.


texts: I am recently learning SWAT modeling in QGIS and I have a problem in delineating the watershed as the DEM file is not working. I downloaded SRTM from USGS Earth Explorer, created mosaic and clipped raster by mask layer to create DEM of Sample area. Then I started a project in SWAT and after selecting the DEM warning appears saying 'DEM does not seem to be projected: Its units are unknown'.


texts: Cropping GeoTIFF to square at specific lat/long for use as Unity 3d terrain heightmap?

texts:
I have some GeoTIFF files that contain elevation data. I have an area defined by latitude and longitude, and I would like to extract the elevation as a heightmap from the GeoTIFF files to use in Unity. Unity requires the heightmap to be both square and the dimensions to be 2^n + 1. The terrain is represented in Unity as a literal square mesh in x[0-2159.9], z[0-2159.9] with elevation in y dimension.

I've been able to import my files into QGIS 3.6 without any trouble and they all line up as expected (NOAA source and SRTM). I've even been able to use gdal_translate to convert to a raw heightmap (ENVI). What I'm struggling with is how to define a square region in QGIS (I can't even figure out how to add a layer that will let me add arbitrary shapes, let alone how to figure the dimensions of that shape). I can't just do delta latitude = delta longitude because that's really a parallelogram, or at least I assume it is.

UPDATE: to simplify things, I have updated the region of interest to be a square, and the lat/long is now:

NE: -82.67761,28.00410
SW: -82.69958,27.98470

Ultimately, what I'm trying to do is create a terrain in Unity that matches up exactly to the transformations that Mapbox does (the source of my coordinates is the the coordinates of the Mapbox tiles) so that I can map points of interest from Mapbox onto my terrain in the correct locations.


texts: Why did my Print Widget stop working in Web AppBuilder?

texts: I have a Web Mapping Application designed in Web AppBuilder Developer Edition for ArcGIS Enterprise.  The print widget uses a custom print service that, up until a couple weeks ago, was working just fine.  But now, when attempting to print from one of our templates, it either throws an error with no description as to what went wrong, or sits and sits like it's processing, but without generating the print.  Our print templates are on a network location which has a Data Store pointing to it.  That Data Store validates with no issue.  What can I do to fix this?


texts: Where can I find historic country border data?

texts: I am looking for border polygons from World War I (from 1914-1918). Is there a resource out there that maintains this kind of information?


texts: How illustrate cells in attribute table?

texts:
I have an attribute table with 2883 rows which are the OA of Leicestershire with their Easting and Northing column. I have four more columns in the attribute table (four land use classes) that carry the population density of each land use class in each OA (i.e. 4*2883=11532 cells).

I want to create a polygon shapefile with 11532 polygons which every polygon's attribute is the population density mentioned in the original attribute table (it should be a table with 11532 rows and three columns of Easting, Northing and population density).

Does any one know how to sort this problem out?


texts: QGIS Modeler - drag and drop data from ToC

texts:
Is it possible to drag and drop data from the ToC to a graphical modeler and then draw connections between data and algorithms?

Also can we do the same with items in the history panel?

This would speed up the process of creating a model.


texts: Travelling salesperson in pgRouting without starting point

texts:
I'm a beginner with PostGIS and it's the first time I try pgRouting.

I work for a regional transit agency and I'm trying to develop a mean to evaluate the easiness or difficulty to offer a bus service for the population of a certain area.

"Traditionally", classic dispersion statistics where used, like the average nearest neighbor.  In theory, it would be easier to offer a good service to a more dense population. But in reality, the road network have a great deal of importance. You can have a quite dense population, but on each side of a river or train track and if the nearest crossing is miles away, you will need to travel a lot longer to be able to pick up the maximum of persons.

We call it the "Pacman index" :) On the left, the theoric bus would pick up 6.7 persons/km. On the right, 5.3 p/km. The two populations have a average nearest neighbor of 1 km.

Once we have the results, we'll compare them to our knowledge of the area and the service already offered (and see if a relatively "easy" area have many buses or if someone, somewhere is doing wonders to serve a difficult area). That's our theory and now we'd like to test it.

I have a set of locations (every intersection within a city that could be a bus stop) and I am trying to find the shortest path between them. Looks like a classic TSP problem.

My problem is that I really don't care about the starting point and it can pass twice on the same path too, as long as the result is the shortest path to "pass by" all locations.

Network Analyst can do it, but it takes forever (my test dataset only have 2 dozen points, but my final one will have 1000+).

With pgRouting, I have to assign a starting point and an end point (if I don't specify an end point, it's doing a loop and I don't want that.

From the documentation, I've found this

CREATE TABLE tsp_route AS
SELECT seq, id1, id2, round(cost::numeric, 2) AS cost
   FROM pgr_tsp('SELECT id, x, y, geom FROM vertex_table ORDER BY id', 6, 5);

It there a way to do this with pgRouting?

Any help would be appreciated.


texts: qgis-bin.exe has stopped working after changing symbol selector?

texts:
i am having troubles with QGIS 2.6.1.  It crashes when i try to change the symbol selector within each vector layer.  If it doesnt crash it seems that the changes i have made in symbol selector do not display.  see screen shot below for more info.

crash error message: qgis-bin.exe has stopped working.  a problem caused the program to stop working correctly. Windows will close the program and notify you if a solution is available.


texts: ArcPy/Python method to ascertain VBA requirements in mxd

texts: I need to be able to look at an MXD and ascertain if it has any VBA code embedded within it.  It's ArcGIS 10.2, and the VBA extension tool is installed on all machines on the estate, but not all mxd files will have VBA code embedded in them.  I've been looking all over but cannot find a method in ArcPy to check the mxds