texts: How to get Database name for layer in ArcObjects?

texts:
I am writing an ArcGIS ArcMap Add-In and need to determine what database a given layer from.  How can I determine this information? I attempted to use the FeatureLayer to get a FeatureClass, and then get the FeatureDataset from that.  From FeatureDataset I can get the name and parse out the first component that seems to be the DB name.

However this approach seems really hacky and does not work if the layer is not part of a FeatureDataset.  When I look at the layer's properties, it clearly shows the Database (see the highlighted line in the screenshot) so I think I can get the DB from a layer but I don't know what interface to use.


texts: How to quick-fix drive time in residential area with Osm2PgRouting?

texts:
I'm currently testing pgRouting using OpenStreempMap XML dataset loaded via Osm2PgRouting. My goal is to solve shortest path problems using time driven costs instead of path length.

I noticed that since V 2.1 osm2pgrouting include two nice columns namely cost_s and reverse_cost_s. They aim to express cost in seconds. Theses columns are computed using maxspeed_forward maxspeed_backward properties.

When the OSM dataset include maxspeed attribute you can expect reasonable results. However when theses tags are not filled a default value is applied (controlled by the config.xml file).

This lead to very high maxspeed values in residential area when the road type is set to primary/secondary/tertiary. But 90Km/H inside a town really doesn't sound realistic.

Has anybody already have a fix for that kind of behavior?

I was thinking about setting a lower maxspeed for road that intersect landuse=residential in OSM. Am I reinventing the wheel or is it something to look at?


texts: GeoServer: Cannot create PoolableConnectionFactory

texts:
I am trying to add a PostGIS database a data store in my GeoServer installation but I keep getting this error:
Error creating data store, check the parameters.
Error message: Unable to obtain connection: Cannot create PoolableConnectionFactory
 ( The authentication type 10 is not supported. Check that you have configured
   the pg_hba.conf file to include the client's IP address or subnet,
   and that it is using an authentication scheme supported by the driver.)

Here is the code in my pg_Hba.conf file.
# TYPE  DATABASE        USER            ADDRESS                 METHOD

# "local" is for Unix domain socket connections only
local   all             all                                     md5
# IPv4 local connections:
host    all             all             127.0.0.1/32            md5
# IPv6 local connections:
host    all             all             ::1/128                 md5
# Allow replication connections from localhost, by a user with the
# replication privilege.
local   replication     all                                     md5
host    replication     all             127.0.0.1/32            md5
host    replication     all             ::1/128                 md5
host    all     all             30.19.120.49/31            md5

Has anyone had this issue before? Was it solvable by adding an entry to pg_hba.conf?
Since my GeoServer resides on http://localhost:8080/ I tried adding:


texts: Need help getting started with QgsComposition

texts:
I am using QgsComposition for the first time, but I fail at the first step :-(

mapRenderer = iface.mapCanvas().mapRenderer()

  line 31, in mapRenderer = iface.mapCanvas().mapRenderer()

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

PS: The Script

from PyQt4 import *
from PyQt4.QtCore import *
from PyQt4.QtGui import *
from glob import glob
from qgis.gui import *
from qgis.core import *
from qgis.utils import*

QgsApplication.setPrefixPath("C:/OSGeo4W/apps/qgis/", True)
QgsApplication.initQgis()

name="test.tiff"
fileName = name
fileInfo = QFileInfo(fileName)
baseName = fileInfo.baseName()
rlayer = QgsRasterLayer(fileName, baseName)

QgsMapLayerRegistry.instance().addMapLayer(rlayer)

mapRenderer = iface.mapCanvas().mapRenderer()
c = QgsComposition(mapRenderer)


texts: Stereo Analyst & ArcGIS 10.1 - volume calculation

texts:
I use to perform analysis and imagery exploitation with the software above.

I have created a Stereo Analyst report with some 3d polygon by a stereoscopic block.
How do I calculate the volume of some soil piles with irregular 3d polygon, with different height for each point?


texts: Defining Winkel Tripel projection in QGIS

texts:
I'm new to QGIS and used to ArcMap.

I'm trying to project to Winkel III and it's not included in
my version 2.18.2

I came across a thread outlining how to define a North American Conic, which included a link to http://spatialreference.org/ where the user apparently got the proj4 line that they plugged into the custom CRS in QGIS with success.

However, there is no proj4 line on spatialreference.org for Winkel III.

Does anybody have it/willing to share it?

Then again, another thread mentions how Winkel cannot be custom defined in this manner:


texts: Precision of point coordinates in polygon geometry is different from its array

texts:
I was using ArcPy to create new polygons from coordinates but when I create the polygon geometry from my ArcPy Array, the coordinates of the polygon corners are slightly different from that of the points in the array used to create it.
See the code below
>>> anc = 0,0
>>> arr = arcpy.Array()
>>> arr.append(arcpy.Point(*anc))
>>> arr
<Array [<Point (0.0, 0.0, #, #)>]>
>>> arr.append(arcpy.Point(45,20))
>>> arr.append(arcpy.Point(25,10))
>>> pg = arcpy.Polygon(arr)
>>> pg
<Polygon object at 0x168f0410[0x168f0b20]>
>>> p = pg.getPart(0)
>>> p
<Array [<Point (4.76837158203e-07, 4.76837158203e-07, #, #)>, <Point (45.0001004934, 20.0001004934, #, #)>, <Point (25.0001004934, 10.0001004934, #, #)>, <Point (4.76837158203e-07, 4.76837158203e-07, #, #)>]>
>>> arr
<Array [<Point (0.0, 0.0, #, #)>, <Point (45.0, 20.0, #, #)>, <Point (25.0, 10.0, #, #)>]>
>>>

Anyone has any idea why it's like this or how I can correct this?


texts: Troubleshooting pythonw.exe has stopped working with Trace Geometric Network?

texts:
I have not yet pulled this down to a code snippet that has just enough in it to reproduce the error below that I am getting when I run my code in IDLE.

pythonw.exe has stopped working
Windows is checking for a solution to the problem

However, I am seeing it when my script repeatedly runs Trace Geometric Network which I am doing to try and identify isolation valves for almost 1,000 locations on a water network.

These are a few lines from the code:

        print_message("Trying to trace")
        arcpy.TraceGeometricNetwork_management(
            r"C:\Users\101791\test.gdb\WaterNetwork\WaterNetwork",
            r"in_memory\WaterDistribution_Net",
            r"C:\Users\101791\test.gdb\Flag",
            "FIND_CONNECTED", "", "", "", "", "wSystemValve",
            "TRACE_ENDS", "", "", "", "AS_IS", "", "", "", "AS_IS")
        print_message("Traced OK")

and this is the ouput I get.  The error in the image above occurs after it has run Trace Geometric Network eight times successfully immediately after printing "Trying to trace" the ninth time.

1
2074
Trying to trace
Traced OK
Isolation Valves for KCA 2074: [u'RV71818']
2
2074
Trying to trace
Traced OK
Isolation Valves for KCA 2074: [u'RV71818']
3
2074
Trying to trace
Traced OK
Isolation Valves for KCA 2074: [u'RV130266']
4
2074
Trying to trace
Traced OK
Isolation Valves for KCA 2074: [u'RV12053', u'RV415348']
5
2074
Trying to trace
Traced OK
Isolation Valves for KCA 2074: [u'RV12052', u'RV12051']
6
2074
Trying to trace
Traced OK
Isolation Valves for KCA 2074: [u'RV12060', u'RV12059']
7
2074
Trying to trace
Traced OK
Isolation Valves for KCA 2074: [u'RV12060', u'RV12061']
8
2074
Trying to trace
Traced OK
Isolation Valves for KCA 2074: [u'WSF405337']
9
2086
Trying to trace

If I try to start the process on the ninth flag pythonw.exe stops immediately.  If I set my code to skip trying the 9-10th and 12-14th flags then all the others work until it gets to the 28th.  The first failure occurs when the KCA value changes for the first time, but the 11th flag (which works) has the same KCA value of 2091 as the 12th-14th (which work) so particular KCA values do not appear to be the cause.

Does anyone have any thoughts on how I may be able to troubleshoot this?

I have a suspicion there may be some interplay between how much memory is needed to perform a geometric network trace and the way the Trace Geometric Network tool writes its result as an in_memory group layer.

I am using ArcGIS Desktop 10.2.1.01 on Windows 7 SP1.


texts: ArcGIS equivalent of "Update Selection Set"

texts: In MapInfo's MapBasic window I can perform a mass update of a field in a selection by entering UPDATE SELECTION SET TYPE = "APARTMENT". How do I do that in ArcGIS now that I have a selection?


texts:
I am searching for default location for QGIS 3 Plugin folder

I checked here but could not find :

C:\Program Files\QGIS 3.0\apps\qgis\python\plugins
C:\Users\Abhijit.qgis2 (cannot locate .qgis3,  should it be there ?
)


texts: Is it possible to export a PDF/SVG as raster with vector labels?

texts:
Is it possible (or will it be possible in the future) to export a PDF/SVG as raster  with vector labels?

I'm using effects (blending, transparency etc.) which I can not export as vector. So I have to export my maps as raster. But since the map will be printed, the raster labels will never gets sharp as they get if they are printed as vectors.

My solution for now is to export once as raster and once as vector and merge them manually via Inkscape. Thats a lot of work to do with 50 different maps!


texts: Moving Versioned Data between different geodatabases

texts:
Can you please let me know how to migrate feature class which is registered as version and archiving is enabled on it

From: Database Environment: ArcSDE 9.3, RHEL 5, Oracle 11.2.0.2
 To Database Environment: ArcSDE 10.1, RHEL 6, Oracle 11.2.0.3.0

I used import feature class tool and copy tool to migrate these feature classes, but they succeded only in feature class data and no versioned/archiving stuff was found in destination database.

Please let me know the best practice for this....


texts: Method to buffer points using euclidean distance and change in elevation

texts:
I am looking for a method to buffer a set of thousands of points with known elevations by 150m Euclidean distance, and then limit the area in each buffer to 50m higher or lower in elevation than each point using a DEM.

What methods and tools (proprietary or open-source) would you recommend with the goal of developing an automatic workflow, tool or script to create a Euclidean buffer with an elevation condition?


texts: How can I add an action from a QGIS toolbar to my custom toolbar programatically?

texts:
I have a toolbar and I want it to have an already existing action, but checking the documentation I only found methods to add newly created actions, not already existing ones. Is there a way to do this?

The action I want to add is the Touch, Zoom and Pan tool from Map Navigation toolbar. I know I can reference it this way:

qgis.utils.iface.mapNavToolToolBar().actions()[0]

But there's nothing I can do with it.


texts: GeoServer - Conversion formula between scale and resolution

texts:
I am trying to figure out what is the conversion formula between scale and resolution used by GeoServer, specifically when creating a new gridset for Tile Caching. One would expect that the conversion formula for a given resolution expressed in m/px would be:

ScaleDenominator = Resolution * (PPI / 0.0254)

However, I have discovered that the following formula gives exactly the correlation between scale and resolution found in GeoServer:

ScaleDenominator = Resolution * (90 / 0.0252)

The value of 90 for PPI is understandable because you can find in GeoServer's documentation (here) that

  The OGC standard output resolution is 90 DPI.

However, I can't understand the value 0.0252. From my point
of view there are 2 possibilities:

Either there's a bug in GeoServer because 0.0252 is used for conversion between meters and inches, instead of 0.0254.
Or GeoServer is correctly using a conversion ratio of 0.0254 between meters and inches, but it is using a PPI of 90.71428571428571, which doesn't make much sense to me either.

Can someone help with this? I am using GeoServer 2.13.1, and similar results have consistently been obtained when using different coordinate systems.


texts: Add Local .osm xml File as a Map Source in MOBAC

texts: I am generating custom OSM XML files using JOSM. I now need to add these as a Custom Map Source in MOBAC. Can I do this from the local file or do I need to create a local web server to load them?


texts: Merging a geodataframe and pandas dataframe based on a column

texts:
I have a simple geopandas dataframe that includes geometry and a column called 'MUKEY':

merged_spatial_df.head()

I want to merge it to a tabular (.csv) pandas dataframe (which also has a column called 'MUKEY') based on 'MUKEY'.

merged_tab_df.head()

There are 31,000 rows in merged_spatial_df and about 391 in merged_tab_df, but each unique MUKEY value in merged_tab_df corresponds to one in merged_spatial_df. I tried the following but can't seem to merge them together and .sjoin requires 2 geodataframes. I tried:

merged_master = gpd.GeoDataFrame(merged_tab_df.merge(merged_spatial_df, how='right'))

but I got an error - ValueError: You are trying to merge on int64 and object columns. If you wish to proceed you should use pd.concat, so I tried:

merged_master = gpd.GeoDataFrame(pd.concat([merged_tab_df, merged_spatial_df], join ='inner'))

None of the concats I tried produced viable results. What am I doing wrong here?


texts: QGIS Semi Automatic Classification SCP Plugin Runtime Error

texts:
I am using Landsat 8 satellite imagery for doing Supervised Classification using Semi Automatic Classification Plugin (SCP) in QGIS.

I have created ROI's and done all the procedures for doing the Supervised Classification. But, when I am finally executing the SCP plugin for creating the output file of classification, I am getting the given below error.

An error has occurred while executing Python code:
RuntimeError: Deleting C:/Users/Sreeraj/AppData/Local/Temp/semiautomaticclassification/20180817_125852974000/20180817_125944159000/sig_1_1_20180817_130833289000.tif failed: Permission denied

Traceback (most recent call last):

File "C:/Users/Sreeraj/.qgis2/python/plugins\SemiAutomaticClassificationPlugin\classificationdock\classificationdock.py", line 520, in runClassificationAction
self.runClassification()

File "C:/Users/Sreeraj/.qgis2/python/plugins\SemiAutomaticClassificationPlugin\classificationdock\classificationdock.py", line 599, in runClassification
ok, opOut, mOut = self.runAlgorithm(cfg.algName, img, sL, cfg.clssPth, cfg.macroclassCheck, None, None, None, compress)

File "C:/Users/Sreeraj/.qgis2/python/plugins\SemiAutomaticClassificationPlugin\classificationdock\classificationdock.py", line 432, in runAlgorithm
oRL, opOut = cfg.utls.createSignatureClassRaster(signatureList, rD, cfg.tmpDir, cfg.NoDataVal, None, previewSize, previewPoint, "No")

File "C:/Users/Sreeraj/.qgis2/python/plugins\SemiAutomaticClassificationPlugin\core\utils.py", line 846, in createSignatureClassRaster
oRL = self.createRasterFromReference(gdalRasterRef, 1, outputRasterList, nodataValue, "GTiff", cfg.rasterDataType, previewSize, previewPoint, compress, compressFormat)

File "C:/Users/Sreeraj/.qgis2/python/plugins\SemiAutomaticClassificationPlugin\core\utils.py", line 2744, in createRasterFromReference
oR = tD.Create(o, c, r, bandNumber, format)

File "C:\OSGeo4W\apps\Python27\lib\site-packages\osgeo\gdal.py", line 1502, in Create
return _gdal.Driver_Create(self, *args, **kwargs)

RuntimeError: Deleting C:/Users/Sreeraj/AppData/Local/Temp/semiautomaticclassification/20180817_125852974000/20180817_125944159000/sig_1_1_20180817_130833289000.tif failed:

Permission denied

Eventhough I tried to solve this by many different ways, all my attempts were unsuccessful.

I can see many files generated in the folder C:/Users/Sreeraj/AppData/Local/Temp/semiautomaticclassification/ when I execute this plugin and I can also see a tif file named sig_1_1_20180817_130833289000.tif (this tif file is having 0 KB size).

How can I solve this issue and what is the exact reason for this error ?

Note: When I am executing the plugin for the continuous 2nd attempt (even after getting above error during 1st attempt), I am getting the given below error in the Python console :-

Exception RuntimeError: 'TIFFResetField:C:\\Users\\Sreeraj\\AppData\\Local\\Temp\\semiautomaticclassification/20180817_135719707000/sig_2_2_20180817_142402476000.tif: Can not read TIFF directory entry.' in <built-in function delete_Dataset> ignored


texts: Individual Feature Class that jumps TO individual Feature in Webmap

texts:
I am using the WAB(3) (Web App Builder, Beta 3) doing development of ESRI-based maps using FeatureClasses stored in our own City Government REST Services.

We currently have a program called 'Adopt-A-Spot' where corporations/individuals can adopt a section of town for various maintenance and cleanliness provisions. What I would like is to offer a webmap with two different kinds of dots.

Dot #1 : SPOT - These are locations on the map that are either ADOPTED or AWAITING ADOPTION.
Dot #2 : ADOPTER - These are locations where the Adopter is based out of (home for individual, corporate HQ for companies)

I would like for someone clicking on an a SPOT to have a link that will jump to the ADOPTER and vice versa.

Using the Javascript API for ESRI / Dojo / WAB(3) is there a best way of doing this?


texts: REST query to determine which polygons geometry is inside and outside of

texts:
I am experimenting with feature services (ArcGIS Server) using REST.  I have a specific problem now I'd like to solve but not sure of the best way or syntax.

I have an input geometry that I query to find out which polygons in a feature services it intersects.  For this I send in a polygon as my input geometry and use 'esriSpatialRelIntersects' as the value for 'spatialRel' in the REST query.  Very nice, it works well.

However, what I want to also identify is that some of my input polygon is outside of any other polygon.
In the example below, I have three polygons. POLY C is my input geometry.  It intersects POLY A and POLY B.  However some of the polygon is over 'empty space'.

How can I determine that?

I need to know if some of the inputGeometry does not intersect the target features.

Perhaps there is a function I can use that will calculate the area that the polygon intesects others and I can compare this to the total area.


texts: Producing MapServer mapscript map inset?

texts:
How do I produce a map similar to the one in http://www.nexusenergy.com.au/6468050/nexus-energy-oil-gas-exploration-and-productio.htm using geodjango to drive mapserver?

I would like to dynamically create the insets, based on user supplied input.


texts: I would like to display Google maps (Raster and Vector) in Java using GeoTools. I am stuck because I am new to this library. Are there any tutorials out there that can help me ?


texts: Determining future location of Fire Station using ArcGIS Desktop?

texts:
I have been tasked with creating a map and doing some analysis that shows the future growth (10 years) of our county. They would like to see where the future growth is likely to occur and also show where fire station(s) should be located to support the new growth.

My first thought is to look at the census. They give me the 2010 numbers and I can also get a current estimate, however, finding anything beyond that is tricky.

In short, I'm not exactly sure the best methodology to use for this analysis. Can anybody give me some pointers?


texts: Scatterplot derived from lat/lon data and census block group

texts:
The Guardian  plotted homicide locations on top of city maps here:

https://www.theguardian.com/us-news/ng-interactive/2017/jan/09/special-report-fixing-gun-violence-in-america

Given lat/lon data for homicides, I'd like to create a census-block-group scatterplot of population density versus homicide prevalence.

I'm familiar with R, python, and linux commandline tools, and wouldn't mind learning some QGIS.

Having never done anything like this, there is probably a better way (?), but my thinking is to convert homicide prevalence to a heatmap, then calculate the average heatmap value over the area of each census block group.

Are there GIS utilities/libraries that would simplify the task?


texts:
Is it possible to read raw data into a geopandas GeoDataFrame, a la a pandas DataFrame?
For example, the following works:
import io
import pandas as pd
import requests
data = requests.get("https://data.cityofnewyork.us/api/geospatial/arq3-7z49?method=export&format=GeoJSON")
pd.read_json(io.BytesIO(data.content))

The following does not:
import geopandas as gpd
import io
import requests
data = requests.get("https://data.cityofnewyork.us/api/geospatial/arq3-7z49?method=export&format=GeoJSON")
gpd.read_file(io.BytesIO(data.content))

In other words, is it possible to read geospatial data that's in memory without saving that data to disk first?


texts:
I'm working on a GIS application that requires to work with Big Data technologies on a cluster (Cloudera CDH 6.2.0). One node will host a GeoServer instance to provide maps to a web application (built with OpenLayers).
My problem is how to provide rasters (GeoTIFF format) to GeoServer. I have two main solutions for the storage:

Keep raw raster files on hdfs,
Yse a tool (HBase?Accumulo?) to save raster data in some way.

Which one is the best?

Once this decision is made, the second point is how to connect GeoServer to data. I need to build layers to provide maps to the webapp but i don't want to store rasters on the local filesystem of the cluster node where GeoServer runs. I want to keep raw data only on HDFS or another store.

Is that solution possible? Is it a good architectural choice? What is a good way to connect GeoServer to HDFS/HBase/other datastore?


texts: QGIS Nautical Chart Symbols for sector lights

texts:
This is I know going to be of limited interest to many of you but I teach navigation and I am working on producing test questions and charts for assessing learners knowledge. My problem is solvable with post editing but as I am finding QGIS very versatile as I am learning it I think there may be another way to achieve this without needing post editing.

The effect I am looking for is this:

And in QGIS I have a layer with length of the vector lines and angle, and using the Vector Field marker symbol I get this:

This I am really impressed with and think I am probably asking to much to be able to either add colour bands as in image 1 or arrows between the lines labeled with the colour abbreviation.

One thought that I have is that perhaps there is some way to draw a radius and clip it to the sectors.

I am using Qgis Lyon and have a Mac so GIS alternatives are not easily possible.


texts: ArcGIS WebApp Builder not Displaying on Mobile

texts:
We create an MXD in desktop and publish the data to ArcGIS Server. Then upload to ArcGIS Online as a Map Image Layer, create a Web Map, and create a Web Mapping Application. If we are on a desktop computer and using Chrome, Firefox, or IE, the Web Mapping Application works just like it should and we have the look and feel of the original MXD. However if we use on a large tablet like Samsung Galaxy Note 4 (Android based) and use Chrome or Firefox, we can only see the ESRI base map and the widgets; we cannot see the Map Image Layer that should have been there from our Web Map.

We were successfully able to upload a zipped shapefile to AGOL, create a Feature Layer, create a Web Map, and then create a Web Mapping Application. The Web App works just like it should, on both desktop and the tablet; the basemap, widgets and map service data all display correctly.

The service, map and app are all shared to Everyone (public).

We cannot figure out why the tablet is having issues with the Web Application that uses a Map Service.


texts: Do OpenLayers maps always show/flash an empty screen while zooming?

texts:
When double-clicking on an OpenLayers map, the screen is momentarily blank while the next zoom level draws.

Is this avoidable? It's not ideal, given that many other mapping APIs show a smooth zooming effect.


texts: How many users can I create in GeoServer?

texts:
I need to know how many users can I create in GeoServer (2.6x)?

Is it infinite or finite?


texts:
I use PostGIS 2.0 and Python API (GeoAlchemy) to work with it.

I want to convert topological informations (points, lines, faces) from PostGIS tables to topological network (points, lines) to study it with network analyses.
Initial data is a planar graph (shapefile, POLYGON).In PostGIS, my first steps were :
 - to create a topological environment (with 'CreateTopology')
 - to convert a shapefile (geometry) into a topological table (topology, with 'AddTopoGeometryColumn', 'toTopoGeom')

Then, I've got :
 - a topological extension ('postgis_topology')
 - a topological schema (called 'my_topo')
 - a table with topological lines ('my_topo.edge_data')
 - a table with topological nodes ('my_topo.node')
 - a table with topological surface ('public.es', topological information is in 'topo' column)

I want to study topological surfaces coming from the shapefile as with network analysis, so, surfaces must be considered as nodes of the network.
I know that 'NetworkX', a Python module, permit to read/write network data from/to PostGIS (module 'pg_net').

How can I convert my tables ('my_topo.edge_data', 'my_topo.node', 'public.es') to a network ? How can I convert polygons ('public.es.topo') to nodes ?


texts: GRASS : r.resamp.rst error : "Not enough disk space--cannot write files"

texts:
Today, I am confronted to a difficulty I cannot manage.

When I try to use the tool "r.resamp.rst", it gives me as answer:

*r.resamp.rst input=map@PERMANENT ew_res=11 ns_res=11 pcurv=curv
Processing all selected output files will require
261320 bytes of disk space for temp files.
Temporarily changing the region to desired resolution ...
Changing back to the original region ...
Percent complete:
Not enough disk space--cannot write files
Not enough disk space--cannot write files
Not enough disk space--cannot write files
Not enough disk space--cannot write files
.....
interpolate() failed
dnorm in mainc after grid before out1= 12.961481
ERROR: split_and_interpolate() failed*

I have 14GB left on my hard disk so I don't understand why they are saying me that :/

I tried to find solution, and I found some answers (e.g. here : http://lists.osgeo.org/pipermail/grass-dev/2005-April/018152.html)
but I don't understand ^^
They speak about a "ulimit -f" function, but where should I write it ? It doesn't work in GRASS.
I am on Windows, maybe is that the problem.

SO in conclusion, I am more than lost so if someone could help me it would be perfect!
Thank you,
b.raoul


texts: Making the selection transparent in QGIS 3.0

texts:
Despite setting the layer or the symbology of the features to semi-transparent, when performing a selection, the features become solid. I want to edit those features based on properties of the features/background image below them.
How can I set a transparency when features are selected?


texts: QMetaTiles plugin in QGIS doesn't work on Linux Mint

texts:
I am trying to use the QMetaTiles plugin in QGIS on Linux Mint 18.3 OS, but it doesn't work.

The error is ImportError: No module named pydevd.

Do somebody know, how to install this module or how to solve this problem?


texts: Field calculate sequential values, but assign same number to duplicates

texts:
I work at a company where I have to constantly update parcel maps and generate tables representing landowner parcels. I have figured out how to sort the landowner name field in the table alphabetically, then generate sequential numbers for map values, but I'd like to have the landowner ID remain the same number if it's the same landowner. This picture shows what I mean.

I have a basic knowledge of model builder and while I don't know python well, I know how to copy and paste code into the field calculator.


texts: Customized legend not showing in geoserver ol3

texts:
Recently we were stuck with updating legend with customized sld but we were not able to get customized colors in legend.

The code was as below

if (layers[i].get('showLegend') === true) {
    try {
        var url = layers[i].getSource().getUrls()[0];
    } catch (err) {
        var url = layers[i].getSource().getUrl();
    }
    var legendImg = document.createElement('img');
    //var str = url + '?REQUEST=GetLegendGraphic&sld_version=1.0.0&layer=' + layers[i].getSource().getParams().layers + '&format=' + format+'&SLD_BODY='+sld+'&legend_options=fontSize:13;fontName:san-sarif;bgColor:0xdddddd;forceLabels:on'+'&WIDTH=12&HEIGHT=12&Scale=1';

    var str = url + '?REQUEST=GetLegendGraphic&sld_version=1.0.0&SLD_BODY:' + sld + '&layer=' + layers[i].getSource().getParams().layers + '&format=' + format + '&legend_options=fontSize:13;fontName:san-sarif;bgColor:0xdddddd;forceLabels:on' + '&WIDTH=12&HEIGHT=12&Scale=1';
    // str=str.replace("<","&lt;").replace(">","&gt;").replace("/","%2f");
    alert(sld1);
    legendImg.src = str;
    legendDiv.appendChild(legendImg);
    var br = document.createElement('br');
    legendDiv.appendChild(br);
}

We are getting square box with text Grey polygon with black outline which is default sld, new sld is in other color.


texts: single feature probability in ERDAS object-oriented classification

texts:
I'm using ERDAS for my UAV image classification (consists 6 land cover types) and I performed the first step of Object-oriented classification, Single Feature Probability (SFP). It outputted a probability layer with value ranged from 0 to 1 indicating the similarity of every pixel to the training set.

However, it did not tell me the formula of how they calculated the SFP.

Could anyone give me a clue of it?


texts: GEE: Reduce region function is not producing all properties for all features

texts:
I have a feature collection of sub-districts in sub-Saharan Africa:

var Africa = ee.FeatureCollection
('users/salem043/Africa_Districts').select("ADM3");

I'm also using the Global Forest Change data:

var gfc2014 = ee.Image('UMD/hansen/global_forest_change_2015').clip(Africa);

I am trying to produce annual estimates of square meters of forest cover loss by subdistrict. The following code worked when I did this for countries, not subdistricts:

//Identifies the year in which the pixel transitioned
var lossYear = gfc2014.select(['lossyear']);

var lossByYear = lossAreaImage.addBands(lossYear)
  .reduceRegions({
  reducer: ee.Reducer.sum().group({groupField: 1}),
  collection: Africa,
  scale: 30
});

var addGroups = function(feature) {
  var statsFormatted = ee.List(feature.get('groups'))
  .map(function(el) {
    var d = ee.Dictionary(el);
    return [ee.Number(d.get('group')).format("20%02d"), d.get('sum')];
  });
  var statsDictionary = ee.Dictionary(statsFormatted.flatten());

  return feature.set(statsDictionary);
};

print("Final forest loss df", lossByYear.map(addGroups).limit(10))
lossByYear = lossByYear.map(addGroups)

But when I run it for the subdistricts, there are many missing observations. The code should produce thirteen new variables for each subdistrict: forest loss in 2001, 2002, 2003, ...2014. But the number of new variables produced varies across sub-districts! For some, only a few year variables are produced, while for others, we're a lot closer to having a complete panel.

Is this a missing data issue? Should I assume forest loss is zero in those years for those subdistricts where no variable is produced?

Moreover, when I exported this data, only two year forest loss variables exported: 2013 and 2014.

Code in GEE is here.


texts:
I have a feature collection of conflict events for Iraq in 1999:

var events = ee.FeatureCollection('users/salem043/ACLED19Iraq');

I would like to transform this feature collection into an image where each band is a day and where the pixels are coded as the number of conflict events at that pixel on a given day. Does anyone know the function that will allow me to do this?

My objective is to use this image to obtain tabular data at the grid-day level of the number of conflict events.

My code is here


texts: Rotating geometry and shifting angles in clockwise direction in PostGIS?

texts:
I have point (cp) and line layer (one row each) in my PostgreSQL database. The point (cp) is on line layer and the direction of line is 112 degrees.

For cp, following code generates a buffer of 50 meters and divides the buffer (circle) into 12 sections (pie-like wedges):

WITH circle AS
    (
    SELECT gid, ST_Buffer(geom, 50,25) as geom FROM cp
    ),
dump AS (
    SELECT gid, (ST_DumpPoints(ST_Buffer(geom,60,3))).geom as geom
    FROM cp
    UNION ALL
    SELECT gid, geom FROM cp
    ),
triangles AS    (
    SELECT gid, (ST_Dump(ST_DelaunayTriangles(ST_Collect(geom),0, 0))).geom geom
    FROM dump
    GROUP BY gid
    )
 SELECT a.gid as pt_id, ST_Intersection(a.geom, b.geom) geom
 FROM circle as a, triangles as b

The output of above code (with line and point layer) is:

I need to rotate circle segments geometry i.e., I need to align 1st segment with line direction (112 degrees) such that starting from this direction, all other segments follow a step of 30 degrees till last segment. For example, 112, 142, 172 and so on. The desired output scenario could be like:

I am aware that ST_Rotate() function rotates the geometry but it does in counterclock wise according to this link while I need to rotate circle segments according to line direction (clockwise).

Can anyone suggest me how to use ST_Rotate() for clockwise rotation?


texts: Fiona throws KeyError when attempting to read a shapefile

texts:
I'm attempting to read the AADT_2014.shp file from the New York State Department of Transportation average traffic data (TDV_Shapefile_AADT_2014.zip specifically through fiona (by way of geopandas). Here's the equivalent minimum working example:

import fiona
import io
import zipfile
import requests

r = requests.get("https://www.dot.ny.gov/divisions/engineering/applications/traffic-data-viewer/traffic-data-viewer-repository/TDV_Shapefile_AADT_2014.zip",
                 verify=False)
zipfile.ZipFile(io.BytesIO(r.content), 'r').extractall()
for feature in fiona.open("AADT_2014.shp"):
    print(feature)

geopandas wraps a simple fiona call for reading this data in. It successfully works out localizing the similar NYC Lion dataset, for example.

But with this data fiona raised a KeyError:

---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
<ipython-input-9-d9cf07a334ea> in <module>()
      1 import geopandas as gpd
----> 2 aadt = gpd.read_file("../data/aadt/AADT_2014.shp")

/Users/Honors/anaconda/envs/nyc-traffic-model/lib/python3.5/site-packages/geopandas/io/file.py in read_file(filename, **kwargs)
     26         else:
     27             f_filt = f
---> 28         gdf = GeoDataFrame.from_features(f, crs=crs)
     29
     30     return gdf

/Users/Honors/anaconda/envs/nyc-traffic-model/lib/python3.5/site-packages/geopandas/geodataframe.py in from_features(cls, features, crs)
    185         """
    186         rows = []
--> 187         for f in features:
    188             if hasattr(f, "__geo_interface__"):
    189                 f = f.__geo_interface__

fiona/ogrext.pyx in fiona.ogrext.Iterator.__next__ (fiona/ogrext.c:19524)()

fiona/ogrext.pyx in fiona.ogrext.FeatureBuilder.build (fiona/ogrext.c:4795)()

fiona/_geometry.pyx in fiona._geometry.GeomBuilder.build (fiona/_geometry.c:3406)()

KeyError: 2002

Why is that?


texts: Extracting geometry data of classes in classified image or masking class layers using Google Earth Engine

texts:
I made an unsupervised image and extracted to classes as a new layers. I want to get and use those layers (subset1, subset2, subset3) geometry data. I will use this geometry datas in ".clip(subset1_geometry)" function.
How can I obtain geometry data of those layers?
var bands=['B4','B3','B2'];

var training=clipdata3.sample({
  region:geom,
  numPixels:101,
  scale:10
});

var grupo= ee.Clusterer.wekaKMeans(3).train(training)

var classified= clipdata3.cluster(grupo)

var subset1 = classified.select("cluster").eq(0).selfMask();

var subset2 = classified.select("cluster").eq(1).selfMask();

var subset3 = classified.select("cluster").eq(2).selfMask();

OR
I want to masking my layers with my collection data but my collection has 16 different band. The program tells me only one band can be masking but I want all bands of the collection data on subset layer.
var subset1 = classified.select("cluster").eq(2).mask(clipdata3);


texts: Split 3D line feature at vertices based on slope change

texts:
I've already asked Splitting polylines in ramps? about this subject, but now I think I'm almost getting to the point.

I have to split roads in ramps to calculate LS factor (RUSLE - soil erosion)

These are the steps that I'm trying to do:

Create a DEM using Topo to Raster
Convert the roads (polyline features) to 3D features using Interpolate shape (ArcGIS 10.1)
Split line at vertices.

but now the problem begins...

There are too many vertices in each polyline.

Is there a way to "clean" these vertices leaving only the ones where the slope abruptly changes (from going up to going down) or maybe create points at the Z min and Z max?


texts: Mapping several maps with the same colorscheme in tmap

texts:
How to control the dot colors in tmap? I have the following code:

png(filename = "test.png",
    width = 210, height = 297, units = "mm", pointsize = 1,
    bg = "white", res = 1200)
for (j in 2012:2016){
t <- tm_shape(n2000) +
  tm_borders("black") +
  tm_shape(dk) +
  tm_fill("grey30") +
  tm_borders("grey70") +
  tm_shape(out3[out3$Country==i & out3$Year == j,]) +
  tm_dots("LE_GEAR",
          size = 1, title=paste0("Type , ", i, ", ", j)) +
   tm_legend(text.size=10,
            title.size=12,
            position = c("left","top"),
            bg.color = "white",
            bg.alpha=.2,
            frame="gray50",
            height=.6,
            hist.width=20,
            hist.height=20,
            hist.bg.color="gray60",
            hist.bg.alpha=.5) +
  tm_shape(Labels) +
  tm_text("Name", size = 4)

assign(paste0("t", j), t)
}

tm <- tmap_arrange(t2012, t2013, t2014, t2015, t2016)#, asp = NA)

graphics.off()

Which gives me the map below. The problem is that the colors for for example GNS is different from 2012 to the other years. How can i make sure the colors are the same?


texts:
I wanted a proportional symbols map of a quantitative variable. Following this question Proportional Symbol Formatting I tried with centroid fill and changing the size of the marker. Nevertheless this what I get using the 10m_admin_0_map_units layer from natural earth data:

As I understand each polygon should have just one marker, nevertheless they have more than one marker. I wonder if I am doing something wrong. I tried this solution since you can only work with color for the symbology method. There is no size option for polygons. I am working with Qgis 2.14.


texts: Raster menu does not have any tools in QGIS

texts: In the QGIS raster menu only raster calculator is available, how do you get other tools, like wrap (reproject)? Please see the image below.


texts: CartoDB Geocoding Fails for Connecticut Cities

texts: I uploaded a dataset that includes a bunch of cities in the state of Connecticut. One example is Farmington. I set the state column to the column with the state abbreviations (it uses "CT") and used text "USA" for the country. All the city and state combinations outside Connecticut geocoded fine, but the Connecticut ones, which make up a majority of this dataset, all failed. Is this a technical error or am I missing some kind of magic thing I should be doing?


texts:
I am having trouble trying to merge two shapefiles in R. One of them has 5 polygons and the other one has 1, so I would like to obtain a final shapefile with 6 polygons.

I try to replicate the merge function from ArcGIS, unsuccessfully though. I have tried with several functions (gUnion, unión from raster package, merge), but all of them create more columns instead of adding a new row.

Could anybody help me to solve this issue?

I have tried with the methods proposed in Merging multiple SpatialPolygonDataFrames into 1 SPDF in R? but they do not work in my case.


texts: Controlling Categories in Script Validation Tools - Expanding Groups By Default

texts:
The validation tab in your script tool properties allows you to put parameters in categories. The help states:

  Put parameters in different categories. Categories are expandable
  window shades with a collection of parameters.

I am using the categories to organize my parameters, so I don't care about the expandable/collapsible feature;  it just adds a lot of time to the user.  This is all within a normal toolbox, not a script toolbox.  v10.1

Does anyone know how to have these categories expanded by default?

Here is my only validation code.  I am only utilizing categories.  The other parts are the default.

def initializeParameters(self):
    """Refine the properties of a tool's parameters.  This method is
    called when the tool is opened."""
    self.params[0].category = "General"
    self.params[1].category = "General"
    self.params[2].category = "Inudaciones"
    self.params[3].category = "Cienegas"
    self.params[4].category = "DrenajesPrincipal"
    self.params[5].category = "DrenajesSensillos"
    self.params[6].category = "Inudaciones"
    self.params[7].category = "Cienegas"
    self.params[8].category = "DrenajesPrincipal"
    self.params[9].category = "DrenajesSensillos"
    return


texts: Reduce image collection to get daily sum from hourly precipitation

texts:
I'm using GSMaP image collection in Google Earth Engine, which has hourly precipitation data from 2000 to present. I have a script which helps getting data (as a table) from points from each available image in given period. Could anyone please assist me on how to make the function that it will sum the values of hourly data to daily. Please, find the script which should be improved below:

 //Points
var p1 = /* color: #d63000 */ee.Geometry.Point([69.15, 54.8331]),
    p2 = /* color: #98ff00 */ee.Geometry.Point([70.917, 54.433]),
    p3 = /* color: #0b4a8b */ee.Geometry.Point([66.967, 54.367]);

    // Collection of point
    var pts = ee.FeatureCollection(ee.List([ee.Feature(p1),ee.Feature(p2), ee.Feature(p3)]))

    // Import GSMaP data
    var gsmap = ee.ImageCollection('JAXA/GPM_L3/GSMaP/v6/reanalysis')
    .filterDate(ee.Date.fromYMD(2014,1,1), ee.Date.fromYMD(2015,1,1))
    .select('hourlyPrecipRateGC');

    // Empty Collection to fill
    var ft = ee.FeatureCollection(ee.List([]))

    var fill = function(img, ini) {
      // type cast
      var inift = ee.FeatureCollection(ini)

      // gets the values for the points in the current img
      var ft2 = img.reduceRegions(pts, ee.Reducer.first(),30)

      // gets the date of the img
      var date = img.date().format()

      // writes the date in each feature
      var ft3 = ft2.map(function(f){return f.set("date", date)})

      // merges the FeatureCollections
      return inift.merge(ft3)
    }

    // Iterates over the ImageCollection
    var newft = ee.FeatureCollection(gsmap.iterate(fill, ft))

    // Export
    Export.table.toDrive(newft,
    "KZ_GSMaP_14-15",
    "KZ_GSMaP_point data",
    "KZ_GSMaP_14-15")


texts: ArcGIS Server Error - ClientCommunicatorAdmin restart

texts:
I am seeing some strange behavior each night with our ArcGIS Sever installation where it will throw the following error and become unresponsive.

Sep 01, 2014 4:03:06 AM ClientCommunicatorAdmin restart
WARNING: Failed to restart: java.rmi.ConnectException: Connection refused to host: 192.168.30.14; nested exception is:
  java.net.ConnectException: Connection timed out: connect
Sep 01, 2014 4:03:27 AM RMIConnector RMIClientCommunicatorAdmin-doStop
WARNING: Failed to call the method close():java.rmi.ConnectException: Connection refused to host: 192.168.30.14; nested exception is:
  java.net.ConnectException: Connection timed out: connect
Sep 01, 2014 4:03:27 AM ClientCommunicatorAdmin Checker-run
WARNING: Failed to check connection: java.net.ConnectException: Connection timed out: connect
Sep 01, 2014 4:03:27 AM ClientCommunicatorAdmin Checker-run
WARNING: stopping

I am still able to navigate to the the manager page, however, each service is listed as "stopping" or "stopped". Restarting the service for ArcGIS Server seems to resolve the problem until approximately the same time on the next day

We have tried re-installing but the problem is still reoccurring..

Has anyone run into this problem before and, if so, what did you do to resolve it?


texts: How to convert OSM polyline data to polygons?

texts: I just now downloaded OSM data .I am having the data as polylines. I want that polylines as polygons is it possible in Qgis?below I also posted the screen shot of the polyline also.


texts: ESRI File Geodatabase-Import Table without Geometry

texts: I have an excel spreadsheet that I need to use as the attribute table for a file geodatabase feature class. I know I can import a table to a file GDB but can I import the table as a feature class or to a feature class? If I have it as a feature class I will use the information in the attribute table to digitize my features.


texts: QGIS 2.2 expression string builder removes brackets

texts:
in QGIS 2.0.1 is used to calculate minimum building distances (Abstandsflächen in German) with an ellipse marker with data definded properties. the expression for this is quite complex:

(([Value 1] + [Value 2] / '2' * '1,01' * '0,464' )* '0,25' + sqrt([Value 3] ^ 2 + ([Value 2] / '2' * 1.01) ^ 2))*2

When I do the exact same thing in QGIS 2.2, the expression string builder changes the expression when saving it into this:

Nabenhoehe + Rotor / '2' * '1,01' * '0,464' * '0,25' + sqrt("Nabe Turm" ^ 2 + Rotor / '2' * 1.01 ^ 2) * '2'

As you can see, the brackets are removed, therefore the result of the calculation is wrong.
I have no idea why this works in QGIS 2.0.1, but doesn't in QGIS 2.2 ...

any ideas?


texts:
I have found this post:
QGIS ->ArcGIS Server configuration possible?
My question is more specific >>
Can I publish feature layers to ArcGIS server from QGIS using WMS?


texts: Appending multiple shps from folder to one dataset using ArcPy in standalone script?

texts:
I have users creating data (shp) daily into a specific folder, I already know how to invoke GP tools to append these data sets into one. What I am trying to do is run a task (which I already know) to look into this folder and take shapefiles and append them together (same schema) into one dataset in another location.

If I have a folder (call it date) 042115 with 5 shapes inside 1,2,3,4,5.shp I want python to look into the folder and take these 5 shapefiles and append them together to another location. I know how to do this with the GP tools, but each day the files names change and I have to continually tweak the model to accomodate the changed file names. What I was hoping for was a python script that would take these (at a folder level) and start the operation without worrying about the specific names. I hope this makes sense, at a lower level I am trying to invoke an append.management operation on shp within a folder and export to another dataset.


texts: Recommend topics to be included in a Computer Science for Geospatial Technologies course

texts:
I will be instructing a course at the local university titled Computer Science for Geospatial Technologies. This is an introductory course meant to introduce computer science concepts to geospatial technologies students (GIS & Remote Sensing). In the past I have introduced programming concepts, but I found this went over many of the students' heads.

Currently, I am planning to discuss computer hardware, spatial data types (i.e. shapefiles vs geodatabases), ESRI Geodatabase Model (the university works on an ESRI platform), basic database theory with ArcSDE Personal.

Could anyone recommend some other computer science related topics that practitioners of GIS and Remote Sensing should know before entering the workforce?

UPDATE:  Last years curriculum included:

Google Maps Javascript API/HTML/Google Earth/KML - 5 weeks
Python Scripting - 6 weeks
Database Theory/MS Access - 2 weeks
VBA - 2 weeks

The response I received from students was that not enough time was spent on each topic.  I am speaking to the university to offer a next level course in GIS Programming using Python.


texts: How to remove layers data source M-awareness to publish editable featureServer

texts:
I would like to publish an editable feature layer (featureservice)  on my ArcGIS server. When I want to publish the data it warns however that the layer is m-aware. In the arcgis documentation I read that you can add or delete m-aware data but not edit the existing points. I do need this editability.

Is there a way to delete the m-awareness of a feature dataset?

The dataset is a layer of points, stored on a database (sqlexpress). I am using ArcGIS 10.2


texts:
ogr2ogr -spat_srs EPSG:31370 -spat 179600 167350 179800 167500 -f "ESRI Shapefile" test.shp "WFS:http://ccff02.minfin.fgov.be/geoservices/arcgis/services/WMS/Cadastral_Layers/MapServer/WFSServer?" Cadastral_parcel

There is also a date field in this WFS-web service [LasUpdDTS].
Is it possible to filter on this field: between 2018/05/01 and 2018/08/01
I suppose it is possible with -sql option?


texts: A simple way to convert multipoint to point in sqlite table, using gdal ogr2ogr?

texts:
I have recently spent a good deal of time attempting to navigate my way around sql command line code using ogr2ogr, within the current release of gdal. As per the title to this post I need to convert multipoint to point for a sqlite table using open source command line tools. So qgis or arcgis tools, etc are not possible solutions. I am not a programmer so I am looking for as direct a command line solution as possible, which can be scripted. I have explored saga gis command line options but the appropriate algorithm seems to not accept sqlite input? My best hope so far appears to be ogr2ogr in gdal. I have seen many examples of postgis code online (indeed, on stack exchange) which is supposed to facilitate this relatively simple operation, but I have been unable to rejig the code in sqlite to make it work for my input/desired output (again, sqlite), if indeed this is possible? My current command is:

ogr2ogr -gt 65536 -nlt POINT -f SQLite -spat -0.4 51.42 0.29 52.65 -dsco SPATIALITE=YES output.sqlite input.sqlite -nln outputname -dialect sqlite -sql "SELECT (ST_DumpPoints(geometry)) AS geometry from input"

I get a 'function does not exist' error.

Why do I get this?


texts: Indication arrow or line showing two or more merged polygons

texts:
Is there a plugin or quick way of indicating that two or more polygons have been merged. I often produce maps that might reflect one compartment number and overall one area, but there may be two or three polygons that form that one compartment or overall area. Sometimes its not clear which polygons have been merged if you have many compartments.

I am using QGIS 3.6.3. I make forestry maps for plantation and sometimes a single compartment might have 3 polygons


texts: Why is "calculate geometry" grayed out when data frame properties shows map units are in meters not degrees?

texts: I'm using ArcMap 10.0 and am trying to get the length of a polyline by adding a field to the attributes table and calculating geometry, but that option is grayed out. Someone suggested that it's because my map was in degrees but according to the properties for the data frame, both display and map units are in meters, which means I should be able to get the lengths of these lines, right?


texts:
My point is to use pgrouting with planet data. osm2pgrouting fails at importing whole planet, so I tried osm2po. Looked fine, it generated sql file in around 2,5h. After importing it to the db though I found there are many ways excluded from the db. Looking closer at osm2po logs:

INFO   46,421,095 of 235,414,299 ways extracted - 25,018M

So I commented out every wtr.finalMask entry in osm2po.config file. It helped a bit:

INFO   74,660,756 of 235,414,299 ways extracted - 15,761M

but still there are many ways missing. So I uncommented this one line:

disablePreFilter = true

But then right before way extracting I got this:

INFO   2,362,909,765 of 2,362,909,765 nodes extracted - 29,114M
FATAL  Exception at Thu Jul 02 07:34:50 CEST 2015- 8,932M
       Message: java.lang.RuntimeException: java.lang.ArrayIndexOutOfBoundsException: Array index out of range: -2147483648
       Class: java.lang.RuntimeException
       Stacktrace:
       de.cm.osm2po.converter.OsmXmlParser.parse(OsmXmlParser.java:83)
       de.cm.osm2po.Main.runTiler(Main.java:214)
       de.cm.osm2po.Main.main(Main.java:159)

After reading all the answers here I understand osm2po doesn't extract ways with tags it doesn't know. Do I understand correctly that disabling prefilter should extract all the ways?
If so, is there any way to configure osm2po to build sql file with all the ways from osm?

I'm using osm2po 5.0.0 version. It was executed with:

java -Xmx32G -jar osm2po-core-5.0.0-signed.jar cmd=tjsp tileSize=x,c ./planet-140522.osm

on system having 64GB RAM. Map file is a little bit old, but I need exactly this version.
My osm2po config is there: osm2po.config

Example of missing ways (or sections, not being sure what exactly "way" in osm means) is there: part.of.osm.file (this is part of newer osm file, not the planet, to be able to grep it reasonably fast, but it looked similar when I grepped it on planet file) leaving in the db only these 3 rows:

osm_id;osm_source_id;osm_target_id
111729593;252343523;1370502399
111729593;1370502399;252343594
111729593;252343594;252343601


texts: ArcGIS 10.5 subtype and domains in database

texts:
I can not find ArcGIS domains and subtypes of feature classes in my database table (PostgreSQL 9.5, postgis 2.2).

For example I have a domain named material types:

1  PVC
2  IRON
3  HDPE

I need to access it from web applications.  So I need to send an sql query to databse and get the list of material types from the attribute domain. Is this possible? Or should I not use the domains and subtypes.


texts: Automatic selection of pixels, from a raster series, greater than minimum image of the series

texts:
Given a series of rasters, and the raster image of the minimum values of the series, I am looking for a python script capable of selecting, for each raster, the pixels whose  values are greater than the minimum, i.e. "A/(A>minimum)", being A a multiple raster.

I made the following script, it seems like is calculating something, but at the end I get the message: 'NoneType' object is not iterable. Any idea for solution?

##NAME=name
##Raster=group
##A=multiple raster
##minimum=raster
##X=output raster

import qgis, os, processing
from PyQt4.QtCore import QFileInfo
from qgis.core import QgsMapLayerRegistry, QgsRasterLayer
from qgis.analysis import QgsRasterCalculatorEntry, QgsRasterCalculator

# Split rasters
layers = A.split(';')
output = X + "/"
suffix = "_suffix.tif"
# Get layer object
layer1 = processing.getObject(minimum)

def x(A,minimum,output):
    for ras in layers:
        fileInfo = QFileInfo(ras)
        baseName = fileInfo.baseName()
        lyr = QgsRasterLayer(ras, baseName)

        # Get layer object
        lyr = processing.getObjectFromUri(ras)
        entries = []
        ras = QgsRasterCalculatorEntry()
        ras.ref = 'lyr@1'
        ras.raster = lyr
        ras.bandNumber = 1
        entries.append(ras)

        #define raster 1 ("minimum")
        raster1=QgsRasterCalculatorEntry()
        raster1.ref='minimum@1'
        raster1.raster=minimum
        raster1.bandNumber=1
        entries.append(raster1)

        calc = QgsRasterCalculator( 'lyr@1/(lyr@1>minimum@1)', output + lyr.name() + suffix, 'GTiff', lyr.extent(), lyr.width(), lyr.height(), entries )
        calc.processCalculation()

for results in x(layers, layer1,output + "*" + suffix):
    fileInfo = QFileInfo(results)
    path = fileInfo.filePath()
    baseName = fileInfo.baseName()
    layer = QgsRasterLayer(path, baseName)
    QgsMapLayerRegistry.instance().addMapLayer(layer)

In the case of solution proposed by Ali, I don't understand if the solution is suitable only for a fixed number of input rasters (7). Actually I need a flexible tool which work independently of the number of input, being them 3 or 10 for instance. Anyway I was not able to make a useful script due to my basic knowledge about python.

On the other hand I also tried the solution as proposed by xunilk but I get the following messages:

xmin = layers[0].extent().xMinimum()
   IndexError: list index out of range

In both cases I would need to use the script within the "processing toolbox" as a standalone tool as in the image:

Actually in the case of using just two layers I have a solution which is the following:

##NAME=name
##A=raster
##Minimum=raster
##X=output raster

from qgis.analysis import QgsRasterCalculator, QgsRasterCalculatorEntry

# Get layer object
layer1 = processing.getObject(A)
# Get layer object
layer2 = processing.getObject(Minimum)

def x(A,Minimum,output):
    entries=[]

    #define raster 1 ("A")
    raster1=QgsRasterCalculatorEntry()
    raster1.ref='A@1'
    raster1.raster=A
    raster1.bandNumber=1
    entries.append(raster1)

    #define raster 2 ("Minimum")
    raster2=QgsRasterCalculatorEntry()
    raster2.ref='Minimum@1'
    raster2.raster=Minimum
    raster2.bandNumber=1
    entries.append(raster2)

    #X Processing
    calc=QgsRasterCalculator('A@1/(A@1>Minimum@1)',output,'GTiff',Minimum.extent(),Minimum.width(),Minimum.height(),entries)
    calc.processCalculation()

x(layer1,layer2,X)

The problem here is that I need the factor A as a multiple raster instead of a single raster in order to iterate the same operation, but I don't know how to solve it.

Is there the possibility of using the script for making a tool from the "processing toolbox" by entering the variables directly from the QGIS workspace by using the following way:

##NAME=name
##Raster=group
##A=multiple raster
##minimum=raster
##X=output raster

Instead of referring the script to folder path?


texts:
Where I can find the QML specification which is used for storing symbology information in QGIS?

I tried to Google it - but got no result.


texts: Adding raster to PostgreSQL database using raster2pgsql?

texts:
I added a raster to my PG database using the raster2pgsql and when I try to load it from QGIS it says that I need to specify it's geometry. Is that a common problem, and if so, how can I solve it?

The command that I used was

raster2pgsql -s 2100 -c C:\Users\user\Desktop\data\raster_2100.tif public.raster2 | psql -h localhost -U user1 -p

and SELECT * FROM raster_column prints this:


texts: Arcpy: creating a line on an existing line based on points' coordinates

texts:
I have a point which is some meters away from an existing line:

This point should be snapped to the line (to the closest part) and then, based on the points' new coordinates, a new line Feature Class should be created, which is on the existing line 1,5 meters left and 1,5 meters right (so 3 meters total) of the point. Therefore it should look like the following picture:

Existing Datasets shouldn't be used, the result should be in a new FC.

Here is the code is used so far:

    import arcpy
    arcpy.env.workspace = arcpy.GetParameterAsText(0)
    point_feature = arcpy.GetParameterAsText(1)
    line_feature = arcpy.GetParameterAsText(2)

    # Make backup copy of point feature class, since modification with
    #  the Editing tools below is permanent

    arcpy.CopyFeatures_management(point_feature, "point_feature_backup")

    # Snapping-Tool

    arcpy.Snap_edit(point_feature,[[line_feature,"EDGE","50 METERS"]])

    # Getting the coordinates of the new snapped Point

    with arcpy.da.SearchCursor(point_feature, ['SHAPE@XY']) as cursor:
        print ("Coordinates of the snapped point: {}".format(row[0]))

    # Creating the new line, I used phantasy coordinates here, it was just a test, the new Coordinates should be on the already existing line 1,5m left and 1,5m right of the snapped point
    cursor = arcpy.da.InsertCursor(line_feature, ["SHAPE@"])

    array = arcpy.Array([arcpy.Point(1823123.213, 6139654.819),
                         arcpy.Point(1823123.244, 6139654.869)])

    polyline = arcpy.Polyline(array)

    cursor.insertRow([polyline])

The hardest part for me is to find out how to get the 1,5m on both sides...which tools could be used to do that or is there an alternative way?


texts: How do I query a string in Arcobjects that has a single quote in it?

texts:
I'm trying to query values using a cursor and loop for the process of removing diacritics. The problem I am running in to is that when I try to create a WhereClause like so...

qf.WhereClause = "\""+fieldname+"\" = "+"'"+address+"'";

I get a COMException about my query if the value of address has a single quote in it (i.e. Land o'Lakes). I have tried escaping the ' with address.Replace("'","\'") but the resulting string is unchanged. When I attempted to double-escape with address.Replace("'","\\'") , the string ends up with both backslashes.

What should I do?


texts: "No module named numpy" Python Error in ArcGIS Label Expression?

texts:
I have developed code to label points that pull formatting from a dictionary. I developed this code and it was working fine two weeks ago, today I try to open my map and I am now met with an error which I can't seem to trace in the code.
Before anyone marks my question as a duplicate, I have previously asked unrelated questions about this code
here: Using Python Dictionary to label subscript formatting ArcMap 10.6 [duplicate]

and here: Labeling using chemical formulas in subscripts that contain a mask in ArcMap 10.6
The error I am currently stuck on is:
"Error 0 on line 0
Traceback (most recent call last):
File "< string >", line 1, in < module >
File "c:\programfiles(x86)\arcpy\arcpy_init.py", line 20, in < module >
import numpy
Import Error: No module named numpy

The code I am working with is:
def FindLabel ( [WellShape.GWMW_NAME], [2018ExceedanceFigureData-PEMLF_Mhuber.csv.SPRING_EXC], [2018ExceedanceFigureData-PEMLF_Mhuber.csv.FALL_EXC] ):
    a = [2018ExceedanceFigureData-PEMLF_Mhuber.csv.SPRING_EXC]
    b = [2018ExceedanceFigureData-PEMLF_Mhuber.csv.FALL_EXC]
    formats = {
        'SiO2': 'SiO2',
        'CO32-': 'CO32-',
        'HCO3-': 'HCO3-',
        'CaCO3': 'CaCO3',
        'NH3': 'NH3',
        'NO3-': 'NO3-',
        'NO2-': 'NO2-',
        'NO3-+NO2-': 'NO3-+NO2-',
        'OH-': 'OH-',
        'NO2-': 'NO2-',
        'Br-': 'Br-',
        'SO42-': 'SO42-'
    }
    c = [str(i.strip(',')) for i in a[1:-1].split()]
    d = [str(i.strip(',')) for i in b[1:-1].split()]
    e = ','.join([formats.get(i,i) for i in c])
    f = ','.join([formats.get(i,i) for i in d])
    if a == "(-)" and b != "(-)":
        return "{}\r\n <FNT scale='70'><CLR red='0' green='122' blue='192'>({})</CLR></FNT>".format( [WellShape.GWMW_NAME] ,f )
    elif a != "(-)" and b == "(-)":
        return "{}\r\n <FNT scale='70'><CLR red='255' green='0' blue='255'>({})</CLR></FNT>".format( [WellShape.GWMW_NAME] ,e )
    elif a != "(-)" and b != "(-)":
        return "{}\r\n <FNT scale='70'><CLR red='255' green='0' blue='255'>({})</CLR>\r\n<CLR red='0' green='122' blue='192'>({})</CLR></FNT>".format( [WellShape.GWMW_NAME] ,e ,f )
    elif a == "(-)" and b == "(-)":
        return [WellShape.GWMW_NAME]

Do you have any suggestions for finding/fixing this error?
I've already restarted my computer and application thinking the glitch may not be in the code (as it was working fine just a few weeks ago).
I am using ArcMap 10.6 and havent changed recently.
However my python updated on March 4th from 3.7.2 to 2.7.16... maybe this new version has an issue with numpy?


texts: Unexpected Result After Reprojection Attempt

texts:
I'm a member of a student team that is currently researching how to query raster information using PostGIS v2.0. Our raster's coordinate system is EPSG 2236 (NAD 83 / Florida East) and as we will be using the Google Maps API extensively, we believe that we should use gdalwarp to reproject the raster's coordinate system to EPSG 3857.

The WKT of our raster is:

PROJCS["NAD83 / Florida East",
    GEOGCS["NAD83",
        DATUM["North_American_Datum_1983",
            SPHEROID["GRS 1980",6378137,298.2572221010042,
                AUTHORITY["EPSG","7019"]],
            AUTHORITY["EPSG","6269"]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433],
        AUTHORITY["EPSG","4269"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",24.33333333333333],
    PARAMETER["central_meridian",-81],
    PARAMETER["scale_factor",0.999941177],
    PARAMETER["false_easting",656166.6666666666],
    PARAMETER["false_northing",0],
    UNIT["US survey foot",0.3048006096012192,
        AUTHORITY["EPSG","9003"]]]

When we execute gdalwarp -t_srs 'epsg:3857' input.tif output.tif, the raster's size expands from 25,500 x 30,000 to 25,532 x 30,164. Is this normal? I know I can enforce a size on the output using gdalwarp's -ts switch - is that common practice?

Lastly, the WKT obtained using gdalinfo of the output GeoTIFF has too many unnamed and unknown parameters as shown below:

PROJCS["unnamed",
    GEOGCS["unnamed ellipse",
        DATUM["unknown",
            SPHEROID["unretrievable - using WGS84",6378137,298.257223563]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433]],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AUTHORITY["EPSG","3857"]]

And the coordinates are way off:

Origin = (-9003811.830927321687341,2993729.091127494350076)

Can someone enlighten me in what we are doing wrong?

Thanks in advance!


texts: Does ESRI ArcGIS API for JavaScript support vector tiles with fill-extrusion Style Property?

texts:
ESRI ArcGIS API For JavaScript’s support GL JSON Stylesheets/Vector Tiles that use fill-extrusion property?
https://docs.mapbox.com/mapbox-gl-js/style-spec/layers/#fill-extrusion
{
"id": "building-3d",
"type": "fill-extrusion",
"source": "openmaptiles",
"source-layer": "building",
"paint": {
"fill-extrusion-color": "hsl(39, 41%, 86%)",
"fill-extrusion-height": {
"property": "render_height",
"type": "identity"
},
"fill-extrusion-base": {
"property": "render_min_height",
"type": "identity"
},
"fill-extrusion-opacity": 0.6
},
"filter": [
"all",
[
"!has",
"hide_3d"
Trying to show 3D Buildings from either MapBox or OpenMapTiles/MapTiler


texts: How to get points which are reachable in given amount of time(pgrouting)

texts:
I am trying to get points which are reachable in given amount of drive time from a point. I have a table of hospitals points. And I have OSM map data of an area which is converted to work with pgrouting using osm2pgrouting tool.

Can someone suggest me what is the right way to achieve this?

Approach I am thinking is as bellow.

Find the alphashape polygon for given N minutes.
Use the polygon to search for points inside that polygon.

I also want to show drive time for each individual hospital. I am not getting how to do this.

So far I have done this.
I have updated the ways table by adding one more column as 'time' and updated the values by calculating approximate amount of time to cover the edge.

 UPDATE ways SET time = round(((length_m * 60) / 12000)::NUMERIC, 2);

Calculating alphashapes for some time like bellow.

SELECT 1 As id, ST_SetSRID(pgr_pointsAsPolygon( $$SELECT dd.seq AS id,
 ST_X(v.the_geom) AS x, ST_Y(v.the_geom) As y FROM pgr_drivingDistance($sub$SELECT gid As id,
source, target, time AS cost FROM ways$sub$,(SELECT n.id FROM ways_vertices_pgr AS n
ORDER BY ST_SetSRID(ST_Point(77.54653690,12.97587840),4326) <-> n.the_geom LIMIT 1) ,  4, true ) AS dd
INNER JOIN ways_vertices_pgr AS v ON dd.node = v.id$$), 4326) As geom;

And then using st_within to find points inside the polygon.

Can someone suggest me how and what is the correct way to achieve my requirement?

I am using POSTGIS="2.4.2 r16113" and pgr (2.5.2,v2.5.2,60585f1f7,master,1.58.0)


texts: QGIS Memory Usage Eventually Crashes Program

texts:
I am using the long-term release of QGIS (3.4.12). When I first start the program, it takes up about 625MB of memory. After I have been working with it for several hours, the memory usage often increases up to 15GB, and everything on the computer becomes sluggish. I have 16GB of memory in my computer. When it becomes sluggish, I have to save everything because if I continue past this point, QGIS usually crashes shortly thereafter. Even if I can close the project I am working on, the memory usage doesn't change. When I then try to exist QGIS, it hangs up, and I have to kill it using the Task Manager.

Is there a reason why the memory usage climbs so significantly and stays high? Is there anyway to reclaim it without exiting the project and the program?

In response to comments, yes, the memory was off by an order of magnitude. 625,000KB to start, and 15GB after using for a while. 16GB computer memory. The installed plugins are DB Manager, Geometry Checker, Georeferencer GDAL, GPS Tools, Lat Lon Tools, Meta Search Catalog Client, Processing, and Topology Checker. I haven't really paid a whole lot of attention to what I was specifically working on when it crashed, but the last time I had been working on something and it crashed, I had been using "Offset lines" and "Split Parts" and "Merge Features". I haven't just left the program open without doing any analysis, but I doubt if doing this would result in any increased memory usage.

Operating system is Windows 7.

Update 03 Feb 2020: Still having the same problem. I haven't been able to identify one particular tool or action that causes it. The memory usage just climbs until it is about 15GB. Even if I close the project, the memory usage remains high (see screenshot)

If I try to exit QGIS, it just becomes unresponsive, and I have to kill the program in Task Manager. There is never a graceful exit. I have to be careful to save every so often or I will lose everything I have been working on. Any ideas?


texts: Merging small polygons with biggest adjacent polygon using ArcGIS Desktop?

texts:
I have a polygon shapefile classifieds in different regions. In each region there are many small polygons that I'd like to merge to the biggest polygon adjacent to each one of them, as long as these polygons belong to the same region.

I don't want that polygons neighbors located in different region be united.

I'm looking for an automatic solution. I tried eliminate tools but ArcGIS does not differentiate the union between regions.

Any ideas?


texts:
I have created a WMS layer in geoserver. I have added that WMS layer using openlayers and geoext.  The WMS layer is displaying, but when I zoom in, the layer disappears and when zoomed out the layer is again visible as a dot.

My code is:

 new OpenLayers.Layer.WMS("Layers",
                "http://localhost:8080/geoserver/cite/wms",
                {
                    "LAYERS":'layer1',
                   srs:'EPSG:3857',
                   transparent: true,
                    format: "image/png"
                },

                {
                    isBaseLayer: false,

                    buffer: 0,
                    // exclude this layer from layer container nodes
                    displayInLayerSwitcher: false,
                    visibility: false
                }
            )


texts: Using GDAL/C++ to calculate distance in meters

texts:
I'm using GDAL in C++, with WGS84. I want to calculate distance between two points in meters. The code snippet below tries to find the distance between 25N/75E and 25N/76E. The code, as is, returns 1.0000, which is just a Cartesian number. I need the answer in meters or feet or miles. How do I do that?

OGRSpatialReference  wgs84;
OGRPoint   point1,
           point2;

assert(wgs84.SetWellKnownGeogCS("WGS84")==OGRERR_NONE);
point1.assignSpatialReference(&wgs84);
point2.assignSpatialReference(&wgs84);
point1.setY(25.0000); point1.setX(75.0000);
printf("Distance %f\n",point1.Distance(&point2));

point2.setY(25.0000); point2.setX(76.0000);
printf("Distance %f\n",point1.Distance(&point2));


texts: Method 'to_crs' to reproject with GeoPandas is missing

texts:
I am trying to reproject data using GeoPandas' to_crs method but it seems to be missing. When I type hasattr(gpd, 'to_crs') I get FALSE. Also, when I type dir(geopandas) in the command prompt it doesn't appear in the list of attributes and methods.

>>> dir(geopandas)
['GeoDataFrame', 'GeoSeries', '__builtins__', '__cached__', '__doc__', '__file__', '__loader__', '__name__', '__package
__', '__path__', '__spec__', '__version__', '_version', 'base', 'datasets', 'geodataframe', 'geopandas', 'geoseries', '
gpd', 'io', 'np', 'overlay', 'pd', 'plotting', 'points_from_xy', 'read_file', 'read_postgis', 'sjoin', 'tools']

I have reinstalled GeoPandas several times, to no avail.


texts: Extracting 1st level boundaries from GADM worldmap in QGIS?

texts:
I need a world map, broken down to GADM 1st level. The obvious way would be to download all single country .shp files and merge them. But that is far too tedious of course.

I tried the dissolve function in QGIS on the 500mb GADM world map, but its not dissolving the entities properly. It fails in 1 out of 10 because for example there are some single points remaining which QGIS recognises as borders.

Is there any solution to the dissolve problem?

Or do you know another approach?

Or is there a 1st lvl administrative boundary .shp file out there on the internet to download?


texts: Increasing resolution of 3D image in QGIS GRASS NVIZ?

texts:
I'm trying to display 3D data in NVIZ within QGIS, but the result appears very low resolution in comparison to the same data in ArcScene.

How do I improve the quality?

I've read the question on 'How to increase the resolution of 3D image in GRASS NVIZ' - but this is based on GRASS not in QGIS and the solution didn't appear to work for my data. Both raster and elevation model have been resampled to 0.5 m resolution.

How the data looks in ArcScene:

How the data looks in NVIZ:


texts: Passing the name of the highlighted layer (in the TOC) to a script?

texts:
Users here add a picture to a map using the 'add data'.
Their first task is to right-click the jpg in the Table of Contents, select 'properties', select 'symbology' and then 'Classify' and then choose a new colors for the different ranges (to make the background transparent).

I can create a model to do this and, alternatively, can write a arcpy script to do this but in both cases, I have to hardcode the name of the input layer.

How do I use the selected layer in the TOC (i.e. clicked via the mouse and is now highlighted) as the input (and output) to my script?

Ideally, the user would add the picture, click on it to highlight it (in the TOC) and then hit a key combination to execute the script.


texts:
I recently changed over to a laptop with a 3200x1800 resolution screen, and QGIS no longer displays properly. Many items and font are too small and almost unreadable. Slight improvements can be made by changing style, icon size and font under Settings -> Options -> General, however issues still remain, particularly with pop up menus like attribute tables, etc.

Surely there are many users now with screens >1920x1080 resolution, experiencing similar problems?

I have read this question: Making icons and text larger on high resolution screen for QGIS?

And tried this solution: http://www.danantonielli.com/adobe-app-scaling-on-high-dpi-displays-fix/

It does not work for me. Are there any fixes? I'm finding this a major problem (particularly working from a 13" screen).

Edit: picture uploaded. You can see the inconsistency between the sizing of each element. The drop down plus symbols next to the layers can barely be seen let alone clicked accurately, ditto some of the other icons. Similar issues occur with pop up menus


texts: Calculation of average point distances from one another within a polygon grid

texts:
I am working with processed LiDAR data of forested areas in ArcMap 10, with various area based statistics extracted for each 30x30m polygon grid-cell in a shapefile.  There are approximately 20,000 of these cells.  I also have a layer of points relating to individual tree locations for the whole site.

What I need to do is calculate the average distance between the tree-points and their standard deviations within each grid-cell polygon.  And then write those two values to the corresponding grid-cell.

I have used the GME toolbox (http://www.spatialecology.com/gme/) to derive this data for areas corresponding to my field site locations; however this is simply not feasible over such a large area.

Does anyone have any ideas on how I can solve this issue?


texts: Replacing values in one field with values from another using ArcGIS Pro Calculate Field tool?

texts: I am in the infancy stage with python and have no idea where to start. I simply want to replace values in a field from a feature layer with the values in a field from a table I joined (to the feature layer) using the Calculate Field tool in ArcGIS Pro.


texts: Accessing secure ArcGIS Online map service from custom web mapping application?

texts:
I'm building a custom web mapping application for a client. The application will be hosted externally to ArcGIS Online, but the application is registered with AGOL. The application will access the client's secured AGOL web map services.

I'd like for the user to be able to log in to AGOL, open the application, and not have to log in again. However, each request for a secured resource from AGOL results in a login challenge.

I've set up the web proxy to avoid this problem, but it seems that the web proxy allows anyone access to the secured AGOL web map services.

So, it seems that my options are:

Allow all users access to the client's secured AGOL web map services (hint: not an option).
Require users to provide their login information each time a request is made to a secured AGOL service (hint: again, not an option).

Am I missing something? Is it not possible to do this?


texts:
I am using ArcMap.
I am using Kriging interpolation in order to show the heat-map on my study area. I have total of 6 variables which ranges between 0.119~0.300 when I select kriging from arc toolbox for 2 of them everything goes right but it won't split my region to different colours.it only shows one color. the same setting work fine with other variables.by the way I use all settings for Semivariogram models but no help.


texts: How do I set a CRS when it isn't listed in QGIS? (Mac, not PC)

texts:
I want to work with a shapefile that is in a Coordinate System that I can't locate in QGIS for Mac. When I load this shapefile data, it shows up fine, on-the-fly has generated the CRS data most likely from the .prj file that has the info in it. However, if I add a coastline from WGS-84 coordinates, or any other data from another coordinate, my data ends up in Africa or even farther than that! So I try to set a coordinate system for it so QGIS will have a reference point.  I type in the filter 102631 because that is the number for the Projection data, but that projection isn't listed.

What do I need to do?  Can I manually alter something so QGIS has a reference point for this projection data?

The shapefile was created in a NAD-83 Alaska Coordinate System (NAD 1983 StatePlane Alaska 1 FIPS 5001 Feet ESRI 102631 http://spatialreference.org/ref/esri/102631/).

On QGIS for the Mac this coordinate reference system is not listed.
How can I correct this on the Mac? Is there a way to download more Projections?


texts: My geoserver layer is not updated always correctly. Old features are still displayed in different zoom levels

texts:
It would be great to get an advice on this topic as its something it troubles me regularly in different projects.

I always find it hard to work efficiently with geoserver and geowebcache when it comes to cached layers etc.

I use a stack of technologies, which are difficult to change (in the current project). I know that my geoserver is quite old (2.5).

I also use geoserver api and gsconfig to create new layers, update the bounding box of a layer etc. I also use a set of POST requests (with python) to add new features on my point/polygon and line layers.

It works but not always. Although I send a seed request to make sure that I update my layer based on the new data:

<seedRequest>
 <name>workspace:layer_name.xml</name>
  <srs>
   <number>4326</number>
  </srs>
  <zoomStart>1</zoomStart>
  <zoomStop>20</zoomStop>
  <format>image/png</format>
  <type>truncate</type>
 <threadCount>2</threadCount>
</seedRequest>

I don't always get to see what it actually exists on the layer and is visible through the layer preview functionality of geoserver. Most of the time if I zoom in/out the features appear (in some random zoom level).

I understand that something goes wrong with the seeding process. Tiles are not seeded correctly/ Old features, which were moved or don't exist keep popping up in different zoom levels.

I also suspect that this might be the browser's cache. Of course I cannot ask from the users to hard reset their browser.

Are there specific strategies to follow here? To bypass the cached file etc?

Settings in geoserver

EDIT

Checking the HTTP headers in developers tool I see that my geowebcache-cache-result equals to MISS. There must be something wrong with the WMS request.

This is my RequestURL:

http://geonode.wfp.org/geoserver/wms?SERVICE=WMS&LAYERS=geonode%3Awld_wfp_ipc_phases&FORMAT=image%2Fpng&TRANSPARENT=TRUE&VERSION=1.1.1&REQUEST=GetMap&STYLES=&TILED=true&SRS=EPSG%3A900913&BBOX=-5009377.085,-10018754.17,0,-5009377.085&WIDTH=256&HEIGHT=256

EDIT II

The above image displays the Response Headers for the request of a specific layer. I see now that in other layers there is no issue the geowebcache-cache-result parameter (it equals to HIT).

But if I understand correctly the issue is that while geowebcache in this case works, sometimes it returns outdated tiles. Does someone can understand something more from the above screenshot?


texts: OpenGeo/Geoserver: can I have a store that access more than one PostGIS schema?

texts:
I have an existing PostGIS database that I use with QGIS.  I am now looking at options for getting some of this data up on the web. Hence OpenGeo.

I have Geoserver installed and running and have succeeded in loading a simple layer from my PostGIS database.  My current problem is that in QGIS/PostGIS I have used the base schema for my geometry tables and put all the tables relating to my layers into the layers schema.  I build my layer using sql.

When I create a new 'store' in geoserver I am forced to specify a single schema within the database so this breaks my modus operandi that I use in QGIS.

I am also struggling to figure out how to define "SQL Views" and whether or not these are dependent on 'stores' or whether there is some other mechanism that will allow me to do arbitrarily complex SQL queries on the database.


texts:
I found Why are the Google basemaps no longer appearing in QGIS? which explains a simple method to reach the google tile server and get the Google Physical tiles. Has anyone a similar code with which I can get the Google Streets tiles? I'm unfortunately not really into programming and just a user of QGIS which is why I have to ask this question.

I'm using QGIS 2.14.8 and I'm unfortunately not able to upgrade right now. So using XYZ tiles or QuickMapServices is not possible. I am using OpenLayers plugin but this one is way too often not working the way I'd like it to work.


texts: Will upgrading my RAM and Processor speed up functions on QGIS Desktop 2.12.0?

texts:
I am using QGIS 2.12.0- Lyon (64-bit) at work and have been finding with especially large amounts of data (many vector layers) that all processes on the system have been running slower. I am running on a Windows 7 machine, 64-bit system, with 4GB RAM and an Intel(R) Core i5 processor, at 3.20GHz.
I was wondering if increasing the RAM would make a difference to the speed of loading and manipulating my map on QGIS? Is it worth increasing?

I had also heard that ArcGIS versions before 10 were not able to use multi-core processing, and wondered if this may be the case for QGIS also, therefore there is no need to upgrade processors?


texts: Creating random points from attribute field

texts:
In QGIS 3.6, I wanted to create different layers of random points reading different fields of an attribute table. For example:

City  #cinemas    #malls     #stores
A       15         10          12
B       20         21          10
C       16         20          16

If I have the layer for the map of the cities, how would I plot the distribution of cinemas, malls and stores as random points in each city. reading the values from this attribute table? Without inputting values myself though.


texts: How to add neighborhood names from shp to csv data points QGIS

texts:
I'd like to add the neighborhoods listed in a shp file to the data points in a csv

Link to shp

Link to csv

Here are the data points overlaid on the shapefile map of the neighborhoods.

I'd like to get these neighborhood names

Into this dataset, so that each data point has the correct neighborhood name

Adding what I've tried.

Clicked 'no'

Opened new shapefile, seems empty.


texts: How to create polyline feature from two input points by following/tracing another polyline?

texts:
I have a list of things each have two grid references, start and stop points. I can turn these into points and join them with polylines.

But what I need to do is have the lines that join the points the follow the route of another polyline, in this case it is a river.

I have 38 lists (for various rivers), each with 10-20 items on them, so hopefully I can find a quick soloution.


texts: Make NumPy wrapper for laspy memory view with more than one attribute (i.e., inFile.XYZ)?

texts:
I'm trying to use lidar point X, Y, and Z coordinates (from .las files) in calculations.  I'm using Python 2.7.13 (x64), with laspy 1.5.0.

Reading attributes one at a time with the following is amazingly fast:

inFile = laspy.file.File(las_file, mode="r")

x = inFile.X

But things slow down considerably when concatenating different single-attribute ndarrays into one array with the following (for example, when getting nx3 array of positions):

coords = np.array((inFile.X, inFile.Y, inFile.Z)).transpose()

I see in Inserting LiDAR points (from laspy) in GeoDataFrame without using a numpy array? that getting inFile.X, inFile.Y, and inFile.Z is as fast as it is because they're numpy-wrapped memory views with little copying going on.

Is there a way to create numpy-wrapped memory views containing multiple attributes, so time doesn't need to be spent copying multiple single attributes to one array, on top of doing something with that array?


texts: Linking new feature violates foreign check constraint in PostgreSQL database when working in QGIS3

texts:
I have set up a spatial database in PostgreSQL in which I have some data/tables that need to be linked by M to N relations. For those occasions I have two tables that contains geometry and a pivot table that contains foreign keys of both tables.

I tried to set up QGIS (3.10) in the way the manual insist for working with M to N relations. I have loaded the relations in project properties and I can see it in the attribute form when I identify the feature.  I have already created features in one layer and now I would like to create features in other layer and link it. When I try to set up the link in the attribute form of the new feature I encounter an error:

  ERROR: insert or update on table "park_usek" violates foreign key
  constraint "fk_park"
               DETAIL: Key (fk_park)=(18) is not present in table "park_pol".

I understand that it's because the foreign key in the pivot table can't be directed to the ID of the new feature as it hasn't yet been created. Is there any way I can set QGIS and/or PostgreSQL (via PgAdmin4) so it should fire this only after the creation of the feature or something like this? I understand that I can probably write a script which would handle this, is there is an easier solution?


texts:
When trying to edit a wfs-t dataset served through GeoServer (2.5.2) - I edit my data and save, however everytime I get an empty response. Using Fiddler2 I modify the request and the transaction goes through fine.

Bad part of request
 xsi:schemaLocation="workspaceURI http://webaddress/geoserver/workspace/wfs?SERVICE=WFS&amp;VERSION=1.0.0&amp;REQUEST=DescribeFeatureType&amp;TYPENAME=workspace:prop_bound&amp;SRSNAME=EPSG:3308&amp;username=***&amp;password=***"

Modification of request
xsi:schemaLocation="http://webaddress/geoserver/workspace/wfs?SERVICE=WFS&amp;VERSION=1.0.0&amp;REQUEST=DescribeFeatureType&amp;TYPENAME=workspace:prop_bound&amp;SRSNAME=EPSG:3308&amp;username=***&amp;password=***"

Removing the workspace uri from in front of the schema location address is solving the problem, however is there a way I can modify the request that is made through QGIS or is this a GeoServer issue? I have had this issue using QGIS 2.4,2.6 and 2.8


texts: How to Parse Through a ESRI GDB File/Folder After Opened with Fiona?

texts:
I should point out that I'm a noob.

I'm playing around with making a python gui app and believe I have opened an ESRI GDB file/folder using:  with fiona.open(GDBPath, 'r') as sourceGeoData

What I'd like to do is inspect / parse through the sourceGeoData to verify I have indeed opened up / connected with the GDB.  I don't know if this is still done via Fiona or if I must pass this on to Shapely which is what I believe I need to use for ultimately doing some spatial joins.

What I don't follow is the use of drivers I see sprinkled through Google searches I've done.  How is this related?


texts:
I have some markers bind with popup. I need to give offset for the popup to display on the top of the markers.Now the popup hides the marker. Below is the snapshot of marker and its popup content.

Below is my code and this is not working. Stillth popup hides the marker. :

 var popup = new L.Popup({ offset:[0,6]}).setContent(contentString);
                      markerdata.bindPopup(popup).addTo(m);

could any one help me on this?


texts: How do I organize labels on top or bottom, or right, or left?

texts: How do I organize the labels in the layout, for example, I want the label Year to be on the top of the Acres. I tried variety ways to do this using the Label Manager but can't seem to get it right  Anyone have a better or suggest of a good way to solution this ?


texts: Can anyone point me to any resources regarding developing a custom WMS service in java? I am seeking for a tutorial for an experienced programmer, but including the details of setting up the server and the like as well as some code samples.


texts: Split a large GeoTIFF into smaller regions with Python and GDAL

texts:
I have quite a large GeoTIFF-file (compressed 50MB but unpacked 6GB) and want to split it into several smaller regional files.

I am using Python, and as I am a new GDAL user I do not really have an idea how to start. Is gdalwarp ChunkAndWarpImage the function I am looking for? And how do I apply it in Python? Can I or must I define the region size?


texts: How to make a raster gradient in any direction in R

texts:
How can I create a raster, where the values increase in a certain direction?

I know how to do the "main" directions with full 90 degree steps, as illustrated in the image. But I would need some help to do lets say degree 27, or 154?

Is there a function for that?

My idea was to always create the "North" gradient raster, rotate the input raster to "face Norh", crop/mask the gradient and assign values to the input raster and rotate back. But rotating a raster also doesn't seem to be too trivial.

Anyone some suggesions or magical functions I am not aware of?

Below is some code to produce the plot and play around.

library(raster)
library(sp)
# library(maptools)
Polygon1 <- Polygon(rbind(c(4498482, 2668272), c(4498482, 2669343),
                          c(4499991, 2669343), c(4499991, 2668272)))
Polygon1 <- Polygons(list(Polygon1), 1)
Polygon1 <- SpatialPolygons(list(Polygon1))

## make raster
r <- raster(extent(Polygon1), ncol = 180, nrow = 180)
polyraster <- rasterize(Polygon1, r)
values(polyraster) <- 1

par(mfrow = c(2, 2))

## Make Raster gradient for direction 0
polymat <- as.matrix(polyraster)
polymat <- apply(polymat, 1, function(x) seq(1, ncol(polymat), 1))
polymat <- raster(polymat,
                  xmn = extent(polyraster)[1],
                  xmx = extent(polyraster)[2],
                  ymn = extent(polyraster)[3],
                  ymx = extent(polyraster)[4])
plot(polymat, main="North 0°")

## Make Raster gradient for direction 180
polymat <- as.matrix(polyraster)
polymat <- apply(polymat, 1, function(x) rev(seq(1, ncol(polymat), 1)))
polymat <- raster(polymat,
                  xmn = extent(polyraster)[1],
                  xmx = extent(polyraster)[2],
                  ymn = extent(polyraster)[3],
                  ymx = extent(polyraster)[4])
plot(polymat, main="South 180°")

## Make Raster gradient for direction 270
polymat <- as.matrix(polyraster)
polymat <- apply(polymat, 1, function(x) seq(1, ncol(polymat), 1))
polymat <- raster(t(polymat),
                  xmn = extent(polyraster)[1],
                  xmx = extent(polyraster)[2],
                  ymn = extent(polyraster)[3],
                  ymx = extent(polyraster)[4])
plot(polymat, main="East 270°")

## Make Raster gradient for direction 90
polymat <- as.matrix(polyraster)
polymat <- apply(polymat, 1, function(x) rev(seq(1, ncol(polymat), 1)))
polymat <- raster(t(polymat),
                  xmn = extent(polyraster)[1],
                  xmx = extent(polyraster)[2],
                  ymn = extent(polyraster)[3],
                  ymx = extent(polyraster)[4])
plot(polymat, main="West 90°")

## Rotate Polygon with maptools::elide
poly_rotate <- maptools::elide(Polygon1, rotate = 27)
r <- raster(extent(poly_rotate), ncol = 180, nrow = 180)
polyraster <- rasterize(Polygon1, r)
values(polyraster) <- 1
polymat <- as.matrix(polyraster)
polymat <- apply(polymat, 1, function(x) seq(1, ncol(polymat), 1))
polymat <- raster(polymat,
                  xmn = extent(polyraster)[1],
                  xmx = extent(polyraster)[2],
                  ymn = extent(polyraster)[3],
                  ymx = extent(polyraster)[4])
plot(polymat, main="North 0°")
plot(poly_rotate, add=T)
poly_rot_clip <- crop(polymat, extent(poly_rotate))
poly_rot_clip <- mask(poly_rot_clip, poly_rotate)
plot(poly_rot_clip)
# How to rotate back?


texts: Importing data with SRID 4269 to GEOGRAPHY type of PostGIS?

texts:
In PostGIS, we are importing various US census bureau shapefiles (with SRID = 4269) using shp2pgsql loader.  I have medium-sized database with several tables containing fields of type geography (lat/lon points, polygons, multipolygons, etc) that span locations across the US on which we do various distance calculations, “within” polygon, etc type of geographic functions as well as display data using Google Maps Api.  We currently have no geometry type fields (only SRID used in database is 4326 via geography type).

Because they cover a fairly large area and to keep things simple and consistent (avoid projection details etc), I would like to store the census shapes (ie boundaries for states, metros, counties,block groups..) in geography data type instead of geometry  as well (unless there is a reason NOT to do this or someone can explain why we should store as geometry instead).

If I use shp2pgsql loader with –G option to create a MULTIPOLYGON geography field for these boundaries (and thus use a 4326 SRID) does the loader “convert” the data from a 4269 projection to 4326?

Is a conversion really necessary?

Do I need to first load into temp table with geometry type then perform update to transform/cast into geography type and then store that result?


texts: How to include data from other tables in the info-window in CartoDB

texts:
I have 3 layers(table_1,table_2,table_3) of data. In the info-window, which was created for the second layer.

I would like to display one column of data from table_1 . When I used  {{table_1.name}}to display the data, I am getting the error , "Error in line 40: Cannot read property 'table_1' of undefined

Is there any way to display data from one table in another layer of data?


texts: Activating 'missing' standard GRASS utilities in QGIS Processing toolbox

texts:
I have been using QGIS and Grass for ages and while I've noticed this issue before, I've not been bothered by it until now.  It is (hopefully) a bit of a noob question but when I open GRASS in the QGIS Processing toolbox, I notice that a lot of standard utilities appear to be missing (e.g. r.in.bin and others  - There are only 161 commands listed for Grass 6.4.3 and this figure should be much higher).

I am not talking about add-on plugins but standard parts of the GRASS toolset (according to the Grass documentation).  Looking in the Grass bin folder, however, I can see the are there.  It is not that the bin folder is not recognized as I can see many of the other utilities there in Processing.  I have tried updating to the latest version of QGIS (2.18) with OsGeo4W and re-installing both Grass 6.4.3 and 7.2.  I have checked the paths in Processing-options and all looks fine (in as much as the paths point to the top-level directory in each case as they always have).

How can I make the missing ones appear in the Processing toolbox?