texts: SRID is messed up when google-maps layer is added

texts:
I've some problems with using a PostGIS layer within QGIS.
First of all, I had this same problem on both QGIS 2.4.0 and 2.6.1
I cannot display the Google Maps layer AND my PostGIS layer together. First of all, my PostGIS layer has the same SRID of GoogleMaps (i.e. 3857). I checked it with

SELECT Find_SRID('public', 'l001', 'the_geom')

then I read here and in other blogs that it can be a proxy problem so I checked the "Use a Proxy for web access

but when I do this and try to load again the Google Maps layer, I get the following error:

Traceback (most recent call last):   File "C:/Users/pgcom120/.qgis2/python/plugins\openlayers_plugin\weblayers\weblayer.py", line 97, in addLayer
    self._addLayerCallback(self)   File "C:/Users/pgcom120/.qgis2/python/plugins\openlayers_plugin\openlayers_plugin.py", line 130, in addLayer
    layer = OpenlayersLayer(self.iface, self._olLayerTypeRegistry)   File "C:/Users/pgcom120/.qgis2/python/plugins\openlayers_plugin\openlayers_layer.py", line 330, in init
    self.olWebPage = OLWebPage(self)   File "C:/Users/pgcom120/.qgis2/python/plugins\openlayers_plugin\openlayers_layer.py", line 48, in init
    proxy = getProxy()   File "C:/Users/pgcom120/.qgis2/python/plugins\openlayers_plugin\tools_network.py", line 42, in getProxy
    proxy.setPort(settings.value("/proxyPort", 10, type=int)) TypeError: unable to convert a QVariant of type 10 to a QMetaType of type 2
Versione Python:
2.7.5 (default, May 15 2013, 22:44:16) [MSC v.1500 64 bit (AMD64)]
Versione di QGIS:
2.4.0-Chugiak Chugiak, 8fdd08a
Percorso Python: ['C:/PROGRA~1/QGISCH~1/apps/qgis/./python/plugins\processing', 'C:/PROGRA~1/QGISCH~1/apps/qgis/./python', u'C:/Users/pgcom120/.qgis2/python', u'C:/Users/pgcom120/.qgis2/python/plugins', 'C:/PROGRA~1/QGISCH~1/apps/qgis/./python/plugins', 'C:\PROGRA~1\QGISCH~1\bin\python27.zip', 'C:\PROGRA~1\QGISCH~1\apps\Python27\DLLs', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\plat-win', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\lib-tk', 'C:\PROGRA~1\QGISCH~1\bin', 'C:\PROGRA~1\QGISCH~1\apps\Python27', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\site-packages', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\site-packages\PIL', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\site-packages\jinja2-2.7.2-py2.7.egg', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\site-packages\markupsafe-0.23-py2.7-win-amd64.egg', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\site-packages\pytz-2012j-py2.7.egg', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\site-packages\win32', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\site-packages\win32\lib', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\site-packages\Pythonwin', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\site-packages\Shapely-1.2.18-py2.7-win-amd64.egg', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\site-packages\six-1.3.0-py2.7.egg', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\site-packages\wx-2.8-msw-unicode', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\site-packages\xlrd-0.9.2-py2.7.egg', 'C:\PROGRA~1\QGISCH~1\apps\Python27\lib\site-packages\xlwt-0.7.5-py2.7.egg', 'C:/Users/pgcom120/.qgis2/python/plugins\DataDrivenInputMask', 'C:\Users\pgcom120\.qgis2\python\plugins\mmqgis/forms', 'C:\Users\pgcom120\.qgis2\python\plugins\QuickMultiAttributeEdit/forms', 'C:\PROGRA~1\QGISCH~1\apps\qgis\python\plugins\fTools\tools', 'C:/ProgettoQGIS2015']

I've also downloaded the Openlayers_plugin and put it in

C:\OSGeo4W64\apps\qgis\python\plugins

but still nothing has changed.
The less worse result I obtained till now is this:

where I didn't flag the proxy option.


texts: Is there an easy way to make feature properties case insensitive?

texts:
I have a system where a user can upload a SHP through a web form. I upload this SHP and import it to GeoServer via REST API. Later the user can interact with the data using a web client with OpenLayers.

Even though the property fields in the SHP are specified by the requirements I'm having lots of issues regarding the case sensitiveness of feature's property names. So, some of the SHP include features with a property called LENGTH, or length, or Length...

I'm wondering if there's an easy way to set OpenLayers in some kind of case-insensitive mode or if should I iterate all features when importing to GeoServer to convert all properties to lowercase (or uppercase, it doesn't really matter)

Edit: Answering notkilroy's comment, the problem is that my features' properties appear as undefined if they I check them lowercased but in the shapefile are in any other case form.

/**
 * Assume the loaded SHP has all properties defined UPPERCASE in this case
 */

var Properties = {
    LENGTH: 'length',
    …
}

var featureLength = feature.getProperties()[Properties.LENGTH]; // undefined

var featureLength = feature.get(Properties.LENGTH); // undefined


texts: Custom QGIS expression to filter atlas pages

texts:
I want to filter my atlas pages based on a list returned by a custom function but it won't work.

This is my function:

from qgis.core import *
from qgis.gui import *
from qgis.utils import *

@qgsfunction(args="auto", group='Custom')
def validZones(layerName, f1, feature, parent):

    lyr = QgsMapLayerRegistry.instance().mapLayersByName(layerName)[0]
    if lyr is None:
        raise Exception("Layer not found: " + layerName)

    expr = QgsExpression("\"f_%s\" IS NOT NULL" % (f1))
    uniqueVal = set([feat['Zone'] for feat in lyr.getFeatures(QgsFeatureRequest(expr))])
    okZone = list(uniqueVal)
    return okZone

And this is the filter expression:

"id_zone" IN validZones('myLayer', '32')

Trying to find the problem I've even tested the simplest function. Something like this:

@qgsfunction(args="auto", group='Custom')
def validZones(feature, parent):
    okZone = ['50','51']
    return okZone

With this filter expression:

"id_zone" IN validZones()

But it doesn't work either.

And as this works:

"id_zone" IN ('50','51')

I also tried using a tuple but to no avail:

@qgsfunction(args="auto", group='Custom')
def validZones(feature, parent):
    okokZone = ('50','51')
    return okZone

What am I doing wrong?


texts: Finding compass direction between two distant GPS points?

texts:
I have a Latitude & Longitude for two locations (me in the West, him across the country), and a web-calculated GPS point for each. I'm trying to figure out which direction -- per my compass -- to look out "toward him" from my rooftop.

Google Maps helped me find the two strings of latitude/longitude (I'm 34.abcdef / -118.stwxyz) and a GPS of that point (I'm at 65 47 31, or however you write it), but I hope you can tell me how/where to point my compass so I aim directly at, say, 40.abcdef/-77.stwxyz.

The compass.google tool was beyond me to even get started.


texts: How Does Directed and has_reverse_cost Arguments when set to [true, false] and [true, true] Affect Routing?

texts:
What will be the difference on the routing behaviour when I try to play with the directed and has_reverse_cost arguments on all the three routing algorithms?
Example: the following SQL queries give me a different result each time.

SELECT rt.gid, ST_AsGeoJSON(rt.the_geom) AS geojson, at_2po_4pgr.gid FROM at_2po_4pgr, (SELECT gid, the_geom FROM dijkstra_sp_delta_directed( 'at_2po_4pgr', 66126, 9013, 0.1, false, false) ) as rt WHERE at_2po_4pgr.gid=rt.gid;
SELECT rt.gid, ST_AsGeoJSON(rt.the_geom) AS geojson, at_2po_4pgr.gid FROM at_2po_4pgr, (SELECT gid, the_geom FROM dijkstra_sp_delta_directed( 'at_2po_4pgr', 66126, 9013, 0.1, true, false) ) as rt WHERE at_2po_4pgr.gid=rt.gid;
SELECT rt.gid, ST_AsGeoJSON(rt.the_geom) AS geojson, at_2po_4pgr.gid FROM at_2po_4pgr, (SELECT gid, the_geom FROM dijkstra_sp_delta_directed( 'at_2po_4pgr', 66126, 9013, 0.1, true, true) ) as rt WHERE at_2po_4pgr.gid=rt.gid;

I can understand if both are set to false but what if one is true and the other is false?
Any explanation please?
Thanks.


texts: Spatial resolution of Google Earth Imagery

texts: What is the spatial resolution of imagery used in Google Earth? At various places on the web, I have been reading it varies from place to place - 60 cm, 2.5 m, 15 m and so on. I was wondering if there is any way/tool to find out what is the spatial resolution at a place of interest to me?


texts: Publishing XY Event Layer with ArcGIS Server

texts:
How can I publish an XY event layer and have changes which are made to the table made to the XY  event layer published service as well?  I've made changes to the table, but they were not reflected in the published service. What is the solution so that I am able to see any changes made to the table, that they will be made to the XY event layer as well?

I am able to view the data via XY Event layer, which I initially thought was an in-memory layer, but it is, however,  an instance of the view/or data set that I am using.  The ideal scenario would be:

Make event layer of XY from SDE table
Publish service definition
All changes will be reflected in published event layer.

I am using ArcGIS Server 10.4.1


texts: How to install HouseNumbering Plugin for QGIS 3.2?

texts: I am trying to install the HouseNumbering Plugin for QGIS 3.2 but I am getting an error that this plugin is only valid for QGIS 2-2.9. Is there a specific plugin available for QGIS 3.2 or can I make changes to the current plugin and use it somehow?


texts: Merging raster data in QGIS using wildcards

texts:
I have a problem merging raster data sets in QGIS, and I have searched various forums, but have not found a solution yet.
An "ordinary" merge, for example merging the raster files

height_tile1_average.tif,
height_tile2_average.tif,
height_tile3_average.tif

that are all within one folder works fine.

But, in the folder there are different tifs, such like

a) height_tile1_average.tif,
height_tile2_average.tif,

b) height_tile1_maximum.tif "height_tile2_maximum.tif

So is there an option to merge only the a) *average.tif"-files, and after that only the b) *maximum.tif"-files using a wildcard?


texts: Difference between a *.lyr and layer directly added from *.gdb/sde to mxd?

texts:
What is the difference between a  *.lyr  and a layer directly added from *.gdb/sde to mxd?

Doing it difference in mxd?

Which is better?


texts: MapWindow - Load basemap image given extent?

texts: I'm new to MapWindow and am attempting to display an image given the geographic extents of an area. Tiling would be ideal, but I imagine that retrieving a simple image from a WMS would be easiest. Are there any simple examples available for this? I'm using MapWinGIS and .NET.


texts: QGIS Export to PDF Using Python Console - version 2.14

texts:
I'm trying to automate a process where I'm joining attribute data with a shapefile. I'm able to create the desired map I'm looking for but I am having issues exporting it from the QGIS interface via PDF. When I run the code it only produces my legend in the PDF format and not my map. My code is below. I am new to coding and GIS.

from qgis.core import *
from PyQt4 import QtGui
from qgis.core import *
from qgis.gui import *
from qgis.utils import iface
from qgis.core import *
from PyQt4.QtXml import QDomDocument
from PyQt4.QtCore import *
from PyQt4.QtGui import *
import pyodbc
import os
import sys
import socket
import csv
import processing

county_uri = 'O:/Broker Support/Cat Modeling/QGIS/Shapefiles/US Counties/UScounties.shp'
shp = QgsVectorLayer(county_uri, 'UScounties', 'ogr')
QgsMapLayerRegistry.instance().addMapLayer(shp)
csv_uri = "file:///O:/Path/File.csv?delimeter=,"
csv = QgsVectorLayer(csv_uri, "File", "delimitedtext")
QgsMapLayerRegistry.instance().addMapLayer(csv)
shpField='state_coun'
csvField='state_coun'
joinObject = QgsVectorJoinInfo()
joinObject.joinLayerId = csv.id()
joinObject.joinFieldName = csvField
joinObject.targetFieldName = shpField
joinObject.memoryCache = True
shp.addJoin(joinObject)
True
colorRamp = QgsVectorGradientColorRampV2.create({'color1':'188,255,255,255', 'color2':'0,0,128,255','stops':'0.25;176,196,222,255:0.50;70,130,180,255:0.75;65,105,225,255'})

renderer = QgsGraduatedSymbolRendererV2.createRenderer(shp, 'File', 6, QgsGraduatedSymbolRendererV2.Jenks, QgsSymbolV2.defaultSymbol( shp.geometryType() ), colorRamp )
shp.setRendererV2( renderer )
QgsMapLayerRegistry.instance().addMapLayer(shp)
qgis.utils.iface.activeLayer()
qgis.utils.iface.actionZoomToLayer().trigger()
def savePDF():
mapRenderer = iface.mapCanvas().mapRenderer()
c = QgsComposition(mapRenderer)
c.setPlotStyle(QgsComposition.Print)
x, y = 0, 0
w, h = c.paperWidth(), c.paperHeight()
composerMap = QgsComposerMap(c, x ,y, w, h)
rect = iface.mapCanvas().currentLayer().extent()
composerMap.setNewExtent(rect)
c.addItem(composerMap)

legend = QgsComposerLegend(c)
legend.model().setLayerSet(mapRenderer.layerSet())
c.addItem(legend)

outpath = "\\\\path\File.pdf"
printer = QPrinter()
printer.setOutputFormat(QPrinter.PdfFormat)
printer.setOutputFileName(outpath)
printer.setPaperSize(QSizeF(c.paperWidth(), c.paperHeight()), QPrinter.Millimeter)
printer.setFullPage(True)
printer.setColorMode(QPrinter.Color)
printer.setResolution(c.printResolution())

pdfPainter = QPainter(printer)
paperRectMM = printer.pageRect(QPrinter.Millimeter)
paperRectPixel = printer.pageRect(QPrinter.DevicePixel)
c.render(pdfPainter, paperRectPixel, paperRectMM)
pdfPainter.end()

QTimer.singleShot(3000, savePDF)


texts:
In QGIS v2.6, I was able to save a vector layer as ESRI FileGDB. However, in v2.8 and v2.10 this no longer seems available.

Am I missing something?

Or why is this no longer available in the standard installation for Windows?


texts:
I have a point layer in QGIS. I have created a
graduated color scheme, but it doesn't look good because
the border of all the symbols is showing and are overlapping (so you hardly see the color of the marker itself anymore).

Is there a way to remove the border of a standard symbol (symbol layer type: simple marker)?

The only thing I found with google was http://lists.osgeo.org/pipermail/qgis-user/2012-March/015778.html, but that was not very helpfull for me (I don't find the Outline option to set to None)

I am using QGIS 1.7.4

UPDATE: summary possible answers (thanks to @Marty and @SS_Rebelious):

Change the "Symbol layer type" from "Simple marker" to "Font marker" (in symbol selector window), select the '.' symbol and change the color.
You can set the border color to the same RGB value as the fill color (but this is a lot of manual work if you have a graduated color scheme => you have to do it manually for each class).
Another option is to use Old Symbology where you have the option to set the 'outline' to 'None'. However, the graduated color scheme is not as advanced as in the New Symbology (no color ramp to choose) and you have to set the outline for each class manually.

I used the first option, as this is the most automated way.


texts:
I was searching for some sort of a guide to develop a plugin on GeoServer, but couldn't find it. Is there any good example I could study instead?

Thanks
FD


texts:
I am trying to create a web map using GeoServer, and doing styling via YSLD.

According to YSLD documentation, I should be able to include CQL in order to do filtering.

For a test, I am using OpenStreetMap data from a PostGIS data store.  Here is some YSLD that works:

- name: name
  rules:
  - name: motorwaytrunkcasing
    filter: ${highway='motorway'}
    symbolizers:
    - line:
        stroke-color: '#440054'
        stroke-width: 5

No problem there.  However, when I add the following CQL:

filter: ${highway='motorway' AND other_tags NOT LIKE '%"tunnel"=>"yes"%' AND (other_tags LIKE '%"bridge"=>"yes"%' OR other_tags LIKE '%"bridge"=>"viaduct"%') AND other_tags LIKE '%"layer"=>"2"%'}

It fails with the following error:

${highway='motorway' AND other_tags NOT LIKE '%"tunnel"=>"yes"%' AND (other_tags LIKE '%"bridge"=>"yes"%' OR other_tags LIKE '%"bridge"=>"viaduct"%') AND other_tags LIKE '%"layer"=>"2"%'} symbolizers ^ Lexical error at line 1, column 1. Encountered: "$" (36), after : ""

I've copy-pasted this CQL to test it using OpenLayers on the same GeoServer, and it is valid; it correctly filters the data.

Might this be a bug?


texts: Draw SQL server geometry records as graphics on map using ArcGIS API for WPF 2.4 -can it be done and how?

texts:
I have records (of public festivals changing daily) in the in SQL server 2012 that have a geometry column (polygons where they occur). I need to show the   records as graphics on a map using ArcGIS API for WPF 2.4. I tried researching it in the internet, but haven't got results. Does anybody know if can it be done and if so,  how?

P.S. : I know there is the QueryLayer\QueryTask, but aside the fact that it's a bit too complex(in my opinion), I'm looking for something else


texts: QGIS 2.12 symbology does not work correctly

texts:
I use QGIS for a long times on mac and windows and I have a problem with the 2.12 version on my OSX 10.8.5 when I want to edit style of my vector layer.

When I click on properties of a vector layer and I want to edit the style, QGIS tells me that This layer doesn't have any editable properties. I attach a file to be as clear as possible. Is anyone having the same problem ? Do you find a solution?

I had the same problem with the 2.10 version.


texts:
I need to create a "total" polygon by merging together all the geometries from  about 10+ tables into one big polygon. For a test with just 2 tables I tried ST_UNION and ST_COLLECT, the result was 500,000 rows and took 15 minute to process. To process all 10 will probably take several days..

Any ideas how I can tackle this problem?
Can I get the outermost line of all the polygons and convert that back to a polygon perhaps? Each table contain several polygons.

EDIT:
The system I am running on is i7-4800MQ CPU 2.7GHz with 8GB ram, 64-bit, win7

The number of vertex points on the datasets vary between 25335 to 585723, average vertex points is estimated to about 100000.

I am not very familiar with spacial queries but what I tried so far is this simple query

SELECT
row_number () over () as x_id,
ST_COLLECT(layer1.geom, layer2.geom) as geom
FROM
sensus.layer1,
sensus.layer2

Did the same with ST_UNION.
I have not tried the to run with all the layers yet as I am just trying out with a few to begin with. My initial idea was to run

SELECT
row_number () over () as x_id,
ST_COLLECT(ARRAY[list of all 10 layers here])
FROM all layers


texts:
I am trying to use the selection expression box to highlight counties by state FIPS. When I use the OR conditional all counties across the US get highlighted. So for example, STATEFIPS = 1 will successfully highlight just counties in Alabama, as I want.

But when I do STATEFIPS = 1 OR 12 then I get all counties across the US highlighted, not just the ones in Florida and Alabama. Is this a bug or am I doing the expression wrong?


texts: What to do if PostGIS gives error that "geometry requires more points"?

texts:
I have a .TAB file.

I convert it to ESRI(shp)-->sql and insert it into database.

Using PostGIS when I want to convert the data into ESRI file, it shows an error:

  geometry requires more points

and it didn't generate .prj file...

So how should I handle this problem i.e. what should I do to make the geometry valid?


texts:
I am currently working the QGIS software and I want to extract roads of a certain region from OSM (Open Street Maps) in order to create a black and white .tiff mask as follows:

It is also preferable that the mask keeps the geoData encapsulated.
I do have 2 layers, the vector roads and the black background that is behind them.
My problem is with exporting, as the background is big (33806 * 21462 px) I was not able to export the above picture with this exact desired size (which is crucial for my project).
I do however have the vector roads in a .dxf file and a black background as .tif with the desired height and width, each one by itself.
I am looking for a way to merge them into one file without losing the original size.
I tried using ArcGIS but I am not familiar with it.


texts: How to Buffer data in a Geographic coordinate system?

texts:
I am having an issue with buffer. My data is in Geographic coordinate system. I am trying to do a buffer of 600m and i am finding different results for each case, have tried thse so far:

Map in web maercator
map in UTM zone 12
map in geogaphic coordinate

I am missing some data points due to this and wanted to check how I can resolve this issue. I know that UTM is the most accurate projection for distance but I wanted to know if by using WebMercator map if there is a way to get correct buffer results , any ideas? Thanks jay


texts: Mapbox Studio Classic layers to TileMill

texts: Is it possible to copy a layer from Mapbox Studio Classic and use it in TileMill? Or at least how can I open the layer in say, a text editor and see its content so I can replicate in TileMill. This is the only way I can see that I can get the mapbox.mapbox-terrain-v2 in TileMill and export it as an MBTile.


texts: Why is Sentinel-1 using geographic coordinate system and Sentinel-2 projected one?

texts: I was wondering if there's any particular reason for Sentinel-1 to be georeferenced to WGS '84 (a geographic coordinate system), while Sentinel-2 is using UTM (a projected coordinate system). I couldn't find any justification online, so my initial thought was that this might be the case in order to avoid distortion of the visualization when viewing within a GIS software.


texts:
I have 2 table: 1. parcel 2. building

How to get list of empty parcels (without building)?

I can get parcel with building by use: ST_DWithin(parcel.geom, building.geom, 0)

...so i made my own solution - i call this sql after insert new parcels

--sql start
--first set all free to 1

update parcel  set free=1;

--then set free to 0 for parcel with building

update parcel  set free=0 where id in (
  select distinct on(p.id) p.id
  from   parcel as p,   building as b
  where  ST_DWithin(p.geom, b.geom, 0)  );

-- so next i can call

select * from parcel where free=0;

-- it is super fast ;)


texts: Highlighting terrain only between defined elevation range

texts:
I am trying to highlight terrain between two elevations in order to highlight where certain plants could be found. For example if a plant in found in a certain valley between 600-1200m, I'd like to be able to make a polygon around the valley and then highlight the terrain in a false color that is between the two elevations to show where it could theoretically be found. Above and below these elevations would transparent. The map so far consists of a satellite image layer and a DEM layer.
How can I accomplish the goal?


texts: Difference between Albedo and Surface Reflectance

texts: Can anyone tell me what is the difference between Albedo and surface reflectance? I know for sure that they are not the same. Any links / comments are most welcome.


texts:
i am going through the qgis 2.18 training manual, Forestry Application (georeferencing) I am using Windows OS8.
 When i load rautjarvi_aerial.tif or forest_stands.shp layers, there is no visual.

When i click on Zoom to layer, they will show but the other layers will now show.

I have followed the manual instructions to the letter but no luck. Is there a bug i am unaware of?

Any assistance will be nice.


texts: How to assign a specific IP to Apache Tomcat server

texts: I have installed Apache tomcat 6 on windows server 2008, i want to assign a specific IP to apache tomcat server, where i do that so that Apache should listen on that specific IP ?


texts: Auto-update a feature class with Python + .csv?

texts:
I am hoping to find out how to automatically update a feature class. We have semi-regular updates to an ArcGIS Viewer for Flex website. Looks like a Python script may be the way to go about it (and I know zilch about Python, yikes!).

I am able to get .csv off the website and have been able to import the xy data and create the initial class with that.

Issues: The data from the CSV is in WGS84 and needs to be projected to NAD 1983 StatePlane California III FIPS 0403. Can a script do this as well?

I may be making a mountain out of a mole hill, so in the mean time I can manually edit the information.

The computer has IDLE (Python GUI).


texts: Create a drop-down menu in a QGIS Processing script?

texts:
I'm trying to write a QGIS processing script using the new API (3.4), but I hardly find any documentation. Based on this excellent resource by Anita Graser (https://anitagraser.com/2018/03/25/processing-script-template-for-qgis3/), I could manage to set up a basic script that mostly fits my needs.

But I have some fixed-value (eg: ['walking', 'car']) options for which I'd like to include drop-down menus. I assume that I have to use a QgsProcessingParameterString but using the available documentation, I could'nt figure out how to do it.


texts:
We have a current architecture of FlexViewer and part ArcGIS Flex API based custom tool that's sitting on a IIS server to draw zones.

ESRI is migrating to ESRI Javascript API and we have two options of either using the web app builder or use plain vanilla arcgis js api.

My TA is suggesting that we use plain vanilla ArcGIS js api, and to load the esri js libraries and host our own version of dojo libraries using the dojo module loader. Also, he is suggesting not to point to any ArcGIS URL directly from the client side, but repoint all the ArcGIS layer url's to a Java layer web services url. That would mirror the ArcGIS rest services.

He also wants all the coded domains in a separate table and all the client side configuration (urls and properties file) to be stored in an Oracle table.

He thinks it would give him better control in the Java side and it would align with his architecture (non ArcGIS-based)

I'm looking for community developer inputs, on the benefits/downside of this proposed architecture.


texts: How to getValues for raster with many values?

texts:
I have a raster, which is not very big, only a few hundred megabytes compressed, but it has about 10 billion values, each stored in 1 byte (0..255). Uncompressed it would take about 10gb (lots of empty cells). When I try getValues(myRaster) I get an error that it cannot allocate vector of size 40Gb. Judging by the number, I get it that it tries to put everything in a vector with 4 bytes values. This means that in just one line, when I try to process a 300mb file, I get an error that I need 40gb of ram. Here's some example code:

library(raster)
myRaster <- raster("myRasterFile.tiff") # <--this file is not more than a few hundred mb
getValues(myRaster)
#Error: cannot allocate vector of size 40 Gb

I do understand that R is not made with programmers in mind. But is there any way I can force getValues(myRaster) to return a vector coerced as 1 byte values?

I know I can just split the raster into whatever fits my RAM, and maybe even run things in parallel as many times as it fits my CPU, but that code would not showcase at all the simplicity and elegance of R for statistics.

I also know I can just increase the swap temporarily to more than 40Gb, thus the virtual RAM, and then just wait until it's done, but this code will literally work only on my computer.

Does anyone have any idea of how would this can be solved?


texts: Using Python script to stack layers from multiple directories and produce multiple outputs in ERDAS?

texts:
I am using ERDAS Imagine 2013 and am trying to create a Python script that takes lists created in another script containing the path to each scene's raw data folder and base image ID's (for example LC80230282013183LGN00).  I then use a for-loop to go to each data folder, selecting the layers to stack to produce my stacked images.  Also, I'm working with ETM+ and OLI data so I have a conditional statement to determine which bands to stack (1-7 for OLI and 1-5 & 7 for ETM+).  What I have so far will stack the layers for the first data folder in the list, but then I get an error when it attempts to stack layers in the next folder in the list.  I'm fairly inexperienced with Python and coding in general and I feel like either my logic is wrong or I'm just way off with my code.

What I have:

from imagine import modeler

from getpath import imageName, imagePath, outputPathStack

m = modeler.Model()

numImages = len(imageName)

for i in range (0,numImages):

    #check for landsat 8 image
    if "LC8" in str(imageName[i]):
        ri1 = m.RasterInput(imagePath[i] + "/" + imageName[i] + "_B1.tif")
        ri2 = m.RasterInput(imagePath[i] + "/" + imageName[i] + "_B2.tif")
        ri3 = m.RasterInput(imagePath[i] + "/" + imageName[i] + "_B3.tif")
        ri4 = m.RasterInput(imagePath[i] + "/" + imageName[i] + "_B4.tif")
        ri5 = m.RasterInput(imagePath[i] + "/" + imageName[i] + "_B5.tif")
        ri6 = m.RasterInput(imagePath[i] + "/" + imageName[i] + "_B6.tif")
        ri7 = m.RasterInput(imagePath[i] + "/" + imageName[i] + "_B7.tif")

        StackBands = m.StackLayers(ri1, ri2, ri3, ri4, ri5, ri6, ri7)
        ro = m.RasterOutput(StackBands, outputPathStack + \
                            imageName[i] + "_stack.img")
    #check for landsat 7 image
    elif "LE7" in str(imageName[i]):
        ri1 = m.RasterInput(imagePath[i] + "/" + imageName[i] + "_B1.tif")
        ri2 = m.RasterInput(imagePath[i] + "/" + imageName[i] + "_B2.tif")
        ri3 = m.RasterInput(imagePath[i] + "/" + imageName[i] + "_B3.tif")
        ri4 = m.RasterInput(imagePath[i] + "/" + imageName[i] + "_B4.tif")
        ri5 = m.RasterInput(imagePath[i] + "/" + imageName[i] + "_B5.tif")
        ri6 = m.RasterInput(imagePath[i] + "/" + imageName[i] + "_B7.tif")

        StackBands = m.StackLayers(ri1, ri2, ri3, ri4, ri5, ri6)
        ro = m.RasterOutput(StackBands, outputPathStack + \
                            imageName[i] + "_stack.img")
    m.Execute()

The error I get is:

File "C:\Workspace\stacklayers.py", line 35, in
    m.Execute()
RuntimeError: erdas::sblib::Operator::SetErrorMessage failed
Spatial Model failed in Raster Output 2.  Error: erdas::sbsupportlib::LayerOperator::OnExecute failed
erdas::sbsupportlib::LayerOpOnExecute::WriteBlock failed
Bad response from exchange.


texts: 3D model of rooftops based on LiDAR and cadastral data?

texts:
I´m trying to model building rooftops based on a LiDAR dataset. I've already interpolated using Delaunay's triangulation a DEM with 1 meter resolution. Due to LiDAR data is not accurate in building edges, I managed to get cadastre data so I have building boundaries.

Does anyone know how to build a 3D model of rooftops based on the previously mentioned inputs (DEM and cadastre data).

Any help would be highly appreciated.


texts: Leaflet Omnivore - KML fill and stroke colors are always the default blue

texts:
I'm trying to import kml file and style them.

In my KML, we have different styles for example

<Style id="YellowPolyStyle">
        <LabelStyle>
            <color>00000000</color>
            <scale>0</scale>
        </LabelStyle>
        <BalloonStyle>
        </BalloonStyle>
        <LineStyle>
            <color>ff6e6e6e</color>
            <width>0.4</width>
        </LineStyle>
        <PolyStyle>
            <color>ff00d9ff</color>
        </PolyStyle>
    </Style>
    <Style id="GreenPolyStyle">
        <LabelStyle>
            <color>00000000</color>
            <scale>0</scale>
        </LabelStyle>
        <BalloonStyle>
        </BalloonStyle>
        <LineStyle>
            <color>ff6e6e6e</color>
            <width>0.4</width>
        </LineStyle>
        <PolyStyle>
            <color>ff009961</color>
        </PolyStyle>
    </Style>

And I calling those style for filling my KML layer. Here is my sample

<Placemark id="ID_00000">
  <name>Grid ID: 14</name>
  <Snippet maxLines="0"></Snippet>
  <description></description>
  <styleUrl>#GreenPolyStyle</styleUrl>
  .......
</Placemark>

But it is not working. It is filling default Blue colour in the KML layer.

So how can I fix this issue. I am using Leaflet Omnivore library for loading the KML.

I manually tried to attach style using setStyle function but it is also not working. Here is my sample code

var runLayer = omnivore.kml(kml)
  .on('ready', function() {
      runLayer.eachLayer(function(layer) {
      if(layer.feature.properties.styleUrl == "#YellowPolyStyle") {
              layer.setStyle(
                  {
                    "LabelStyle": {
                      "color": "00000000",
                      "scale": "0"
                    },
                    "BalloonStyle": "",
                    "LineStyle": {
                      "color": "ff6e6e6e",
                      "width": "0.4"
                    },
                    "PolyStyle": { "color": "ff00d9ff" }
                  }
              );
        }

How can this be done?


texts: Repointing and rebuilding address locators?

texts:
I have a few ArcGIS address locators that originally used data from a file geodatabase.  I would like to repoint these locators to the same layers in an enterprise geodatabase instead.  Additionally, as these layers are sometimes edited, I would like to create a script to automatically rebuild the locators to reflect the changes.  I was able to manually repoint the locator data sources to the enterprise geodatabase.  However, every time I attempt to rebuild, I get the error "The address locator needs to be repaired in order to be rebuilt.  Do you want to proceed with repair?"  When I say ok, it opens the locator's properties with a red X next to the data source (which is in the format Database Connections\MyDatabase.sde\fileName).  I will manually repoint it to that exact same data source.  It will be fine with that, allow the rebuild, and the updated features will appear in the locator.  However, the next time I try to rebuild, same issue.

If I try to use the Rebuild Address Locator tool in ArcMap, it says Error 000005: Could not open the address locator.  Then if I try to rebuild or open the address locator, it will say that it can't open the address locator because it was made in ArcMap 9.2 or earlier (it wasn't, it was made in 10.3).  This appears to permanently break it and all I can do is delete it and recreate it.  This happens even when I stop the services that use the locator and make sure all relevant files are closed.

I am able to successfully rebuild the address locators manually, and it will work.  However, I'd like to be able to automate this process so that the locators will be rebuilt whenever the source data is updated.  Everything works fine when using the data from the file geodabase but the problems start when I switch to the enterprise geodatabase.  The address locators are located on the C drive of the server where the geocoding services are published.  The file geodatabase is also on the C drive in the same folder (but the locators aren't in the geodatabase).  The enterprise database is in the Database Connections on the same server.

Do I need to move data around or use different file paths?


texts: Hard wood tree and Soft wood(shrubs) Vegetation Differentiaton

texts:
How I Differentiate Hard wood tree (like Mango, Sal) and Soft wood(shrubs like reed) Vegetation using Multispectral or Hyperspectral Image? I have TM, ETM, OLI, Hyperion, Aster Image and ERDAS 2011, ENVI 5.1, Arcmap 10.1.
Is it possible to differentiate?


texts: categorized style in 2.14 returns faulty values

texts:
i have just tested latest release of qgis with a simple dataset. it is polygon data - type geodatabase with some categorial columns. the one i've tested characterises land-use types in 15 categories. colum-type is integer - value range is from 1 to 15.
layer properties / style / categorized / classify returns the following classification:

these values originate from the first numeric column within the corresponding attribute table - instead of the selected column.

are there any known workarounds yet ?


texts: Converting .LAS to .PCD file using CloudCompare

texts: I want to convert a LAS file to PCD format using CloudCompare software. How I can do that?


texts: Displaying multiple similar overlapping polygons in QGIS?

texts:
This is more of a cartographic style question, using QGIS.

I have a dataset that consists of land parcels (previous archaeological projects) as polygons. Multiple projects from different years share the exact same spatial shape.

In composer, I would like to export a map that somehow can show that there are multiple polygons on top of each other, and also include all of the labels.

Currently, I am unable to show more than one shape except by manually decreasing the offset for the top polygon and making it slightly smaller than the one it is on top of.

Is there a better way to do this or to automate this process in the dataset? I think perhaps by writing a script to identify similar parcels, and then add an offset value, so that the map would look like nested differently-colored shapes,  but this seems like an awkward work-around.


texts: Start postgis (Opengeo suite azure linux vm)

texts:
I have a fresh Azure Linux VM based on Opengeo suite.  The dashboard and geoserver are working.  But I can't connect or configure postgis.  I follow this guide.

Via putty, logged in with the only user azure gives you (wgadmin in my case), I enter: sudo -u postgres psql postgres

and get:

could not change directory to "/home/wgadmin": Permission denied
psql: could not connect to server: No such file or directory
        Is the server running locally and accepting
        connections on Unix domain socket "/tmp/.s.PGSQL.5432"?

PSQL - V returns 9.3.5.

I have changed the postgres users password another way.  I then go on to try and modify: /etc/postgresql/9.3/main/pg_hba.conf but it doesn't seem to be present.

Am I doing something wrong?

A fresh install should be easy to configure.


texts:
I'm iterating over a layer with a python-script in QGIS 2.18. I am going to pan to every feature and ask the user a question regarding this feature an do something depending on his answer.
How is it possible to enable user-pan/zoom etc. while QInputDialog waits for an answer ?
Do I have to use other dialogs to be able to pan/zoom while looking for an answer?

canvas = qgis.utils.iface.mapCanvas()
pLayer=QgsVectorLayer(r'pathToSHP', "polygonLayer","ogr")
QgsMapLayerRegistry.instance().addMapLayers([pLayer])
for po in pLayer.getFeatures():
    pLayer.select(po.id())
    qgis.utils.iface.actionPanToSelected().trigger()
    pLayer.deselect(po.id())
    checked = False
    while not checked:
        answer=QInputDialog.getText(None,'...','....')
        if answer[0] in ('Y','y'):
          #doSomething


texts: ArcGIS Server Object Extension (SOE) RestOperation configuration

texts:
How would you accept a call to a RestOperation with string and reply back with a json array of the column names in the table holding this information?

So, if string is 'owners' then get names of col1,col4,col9,col10, if string is 'zipcode' get name of col9. etc...

New to the SOE, and excited for some guidance, thanks.


texts: How do I overlay points/polygons in a 3D map using QGIS2threejs?

texts:
I am trying to produce a topographical map of an area, but I need my points to sit on top of it. I have created the actual 3d rendering in QGIS2threejs, and I can see the option to include points, but they appear underneath the 3d rendering, even though they do have height data!
What am I doing wrong?


texts: Using vector field for buffer model in QGIS Modeler

texts:
I am creating a model that will create a viewshed and a buffer around an observer layer /point. I want to use a field in the observer layer for the buffer range.
I tried the option 'Field type'. But it won't let me select anything. Probably because the system does not know what to select.
Furthermore, I tried a 'Vector Field' as input. And pre-calculated value. The model then generates no matter what a 10m buffer.
How can I create a model where I'm able to select the column name of the buffer range?
EDIT:
Added the graph with the solution.


texts: Shapefiles are not being added to MXD using ArcPy?

texts:
I want to add a shapefile to an mxd document using arcpy but it generates the following error when used by adding script in arctoolbox

import arcpy
mxd = arcpy.mapping.MapDocument("CURRENT")
dataFrame = arcpy.mapping.ListDataFrames(mxd, "*")[0]
addLayer = arcpy.mapping.Layer(r"D:/Data/input/Line.shp")
arcpy.mapping.AddLayer(dataFrame, addLayer)

ERROR--------------

<type 'exceptions.ValueError'>: Object: CreateObject Layer invalid data source
Failed to execute (Script).


texts: QGIS >2.10: Adding layers to legend groups programatically fails ocassionally with valid layers

texts:
I'm working in a project in which, for some reasons, the beast approach to handle a certain amount of point/vector/raster layers is to load them all into memory, and keep them in the QGIS legend.
For sorting and user sanity purposes, I wanted to keep them neatly ordered in groups in the QGIS legend interface.

Unfortunately, a problem has shown up in this (somewhat simple) piece of code:

    for layer in layerList:
        print("Adding layer: "+layer.name()+" which isValid: "+str(layer.isValid()))
        with self.groupAssignmentLock:
            group.addLayer(layer)
        iface.legendInterface().setLayerVisible(layer, False)

I added that lock to guarantee only one thread would try to make any operation in a group at the same time, to no avail.

I have registered listeners to the signals QGIS emits when a layer will be added and has been added to the QgsLayerTreeNode object (signals willAddChildren and addedChildren respectively) and the result is this:

>>> Node will be added with : 114 from 114 to 114
Currently grouped layers: 475 out of 484
Number of elements in this node: 115
Adding layer: 2016-02-17T14:00:00Z_temperature which isValid: True
 >>> Node will be added with : 115 from 115 to 115
Currently grouped layers: 476 out of 484
Number of elements in this node: 116
Adding layer: 2016-02-19T10:00:00Z_temperature which isValid: True
 >>> Node will be added with : 116 from 116 to 116
Currently grouped layers: 477 out of 484
Number of elements in this node: 117
Adding layer: 2016-02-19T17:00:00Z_temperature which isValid: True
 >>> Node will be added with : 117 from 117 to 117
Currently grouped layers: 478 out of 484
Number of elements in this node: 118
Adding layer: 2016-02-18T13:00:00Z_temperature which isValid: True
 >>> Node will be added with : 118 from 118 to 118
Currently grouped layers: 479 out of 484
Number of elements in this node: 119
Adding layer: 2016-02-15T05:00:00Z_temperature which isValid: True
 >>> Node will be added with : 119 from 119 to 119
Currently grouped layers: 480 out of 484
Number of elements in this node: 120

This is the end block of the log file output. For some reason or another, 4 layers in this example were not added to groups (which is a problem, as I use these signals as a way of finding out whenever QGIS has finished preparing everything for the plug-in or not).

An example of what happens when a layer fails to be added is this:

Adding layer: 2016-02-15T16:00:00Z_temperature which isValid: True
 >>> Node will be added with : 48 from 48 to 48
Currently grouped layers: 409 out of 484
Number of elements in this node: 49
Adding layer: 2016-02-14T23:00:00Z_temperature which isValid: True
Adding layer: 2016-02-15T23:00:00Z_temperature which isValid: True
 >>> Node will be added with : 49 from 49 to 49
Currently grouped layers: 410 out of 484
Number of elements in this node: 50

As you can see, two messages of "Adding layer:.." are consecutive, showing how neither willAddChildren nor addedChildren are being called.

I'm not sure it's a problem of my own, thus I decided to ask you if this problem has been experienced before. Otherwise, I'm going to assume it is a QGIS bug.

As a note, I'll append the code of the listeners:

def onChildrenAddedToNodeGroup(self, node, indexFrom, indexTo):
    """
    Will be called every time QGIS loads a new layer into a
    layer group in the legend interface.
    """
    self.groupedQgisLayers = self.groupedQgisLayers + 1
    print("Currently grouped layers: "+str(self.groupedQgisLayers)+" out of "+str(self.layersToBeGrouped))
    print("Number of elements in this node: "+str(len(node.children())))
    if self.groupedQgisLayers == self.layersToBeGrouped:
        iface.mapCanvas().refresh()
        self.animatorReady.emit(self.animatedLayerCount)

def onChildrenWillBeAddedToGroup(self, node, indexFrom, indexTo):
    print(" >>> Node will be added with : "+str(len(node.children()))
          + " from "+str(indexFrom) + " to "+str(indexTo))
    pass


texts: Interpretation of PostGIS geocode function ratings

texts:
Edit: I am using Postgres 9.2 and recently installed PostGIS along with the Tiger Geocoder from the extras folder and downloaded all the latest Tiger data.

I have been testing the Tiger Geocoder function (http://postgis.net/docs/Extras.html) to process a big database of addresses in order to normalize/validate them and get lat-long coordinates.

I know that the function returns a rating between 0-100 (0 being the best), but I have not been able to find any documentation stating how to interpret the score.

How can I interpret the scores between 0 and 100?


texts:
I have an Image Pyramid in GeoServer 2.2 RC3 (not the latests - I know, but I have other pyramids working here).

The total size is around 1TB and file size is around 1-2GB for each GeoTIFF.

When looking at the results using WMS, it all looks blurred. Going into the GeoServer log, I found that when I request an image, it selects the level which I think is causing it to upscale the result to fit the screen:

Request:
Size in pixels: 1632x960
Size in meters: 1406x827

I have a GSD = 1 meter/pixel at the level GeoServer selects (which is 3). Shouldn't it select level 2 with a GSD = 0.5 meter/pixel to ensure GeoServer downsamples the result - thus loosing information instead of creating it?


texts: Align and georeference GIS layers exported from CAD

texts:
I have two separate CAD files of the same locale. The first one has surface data like street curbs, building footprints, etc. It is not georeferenced, but when I add it to other items in a projected (NAD 1983 FIPS US feet Connecticut) data frame, it aligns perfectly.

The other CAD file contains subsurface utilities and a surface basemap as a reference. When I try adding them, they show up very very small not too far from the properly projected data.

When I try to project subsurface layers, they don't show up at all.

Trying to georeference manually gives me an error after I manually attempt to enlarge it using georef toolbar tools, even if I delete the spatial index.

I have ArcGIS 10.1 SP1 and a trial of CAD 2014.

Screenshot displays properties boxes, the surface data on the right, subsurface on the left.


texts: Level DEM from ArcHydro Tools 0x80040202 error

texts:
When I tried to run Level DEM tool from the ArcHydro Tools, it returns me an error message:
Level DEM System.Runtime.InteropServices.COMException (0x80040202): The index passed was not within the valid range.

App versions:
ArcGIS Desktop: 10.5
ArcHydro Tools: 10.5.0.221

DEM is in the ESRI grid format (raster is inside a folder); the DEM's values are in meters (float) but I tried change it to cm (int). In the DEM_cm raster I have changed in the prj.adf file Zunits value to 100 - without any progress.
lakes polygon feature class (inside geodatabase)

DEM in Projected Coordinate System (WKID: 2180).
Can anyone help me solve this problem?


texts: How should I interpret R terrain/slope output in degrees?

texts:
I'm calculating the average land gradient (in degrees) in Brazil and I do not know how to interpret the results.

As I'm calculating the slope in degrees, 0.25 is 0.25 degrees, 25 degrees or should I use a conversion factor?


texts:
I have some GeoTIFF files that I want to render as an image mosaic layer in GeoServer. The files were created with gdal_grid, and they look fine in QGIS:

And they look fine in GeoServer when loaded as individual GeoTIFF stores. Here's the larger of the two rendering in the GeoServer OpenLayers preview:

But if I attempt to use an Image Mosaic store (using the exact same files, in the exact same directory), GeoServer renders them upside down. Each file is flipped within its own (correct) bounding box:

Is there a way to get the Image Mosaic to render them in the proper orientation? I'd prefer not having to regenerate the images, but if there's a tweak to the gdal_grid command line that would make it work, that'd be fine too.
Both the images are in the "EPSG:4326" coordinate system.

gdalinfo output from one of the images:
$ gdalinfo bath-SCI-3sec.tif
Driver: GTiff/GeoTIFF
Files: bath-SCI-3sec.tif
Size is 1801, 601
Coordinate System is:
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433],
    AUTHORITY["EPSG","4326"]]
Origin = (-118.500416666666993,32.499583333333298)
Pixel Size = (0.000833333333334,0.000833333333333)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (-118.5004167,  32.4995833) (118d30' 1.50"W, 32d29'58.50"N)
Lower Left  (-118.5004167,  33.0004167) (118d30' 1.50"W, 33d 0' 1.50"N)
Upper Right (-116.9995833,  32.4995833) (116d59'58.50"W, 32d29'58.50"N)
Lower Right (-116.9995833,  33.0004167) (116d59'58.50"W, 33d 0' 1.50"N)
Center      (-117.7500000,  32.7500000) (117d45' 0.00"W, 32d45' 0.00"N)
Band 1 Block=1801x1 Type=Float32, ColorInterp=Gray

Command line used to make the image:
$ gdal_grid \
    -ot Float32 \
    -outsize 1801 601 \
    -txe -118.500416666667 -116.999583333333 \
    -tye 32.4995833333333 33.0004166666667 \
    -a nearest \
    -a_srs EPSG:4326 \
    bath-SCI-3sec.vrt \
    bath-SCI-3sec.tif

Contents of bath-SCI-3sec.vrt:
<OGRVRTDataSource>
    <OGRVRTLayer name="bath-SCI-3sec">
        <SrcDataSource relativeToVRT="1">bath-SCI-3sec.csv</SrcDataSource>
        <GeometryType>wkbPoint</GeometryType>
        <GeometryField encoding="PointFromColumns" x="field_1" y="field_2" z="field_3" />
    </OGRVRTLayer>
</OGRVRTDataSource>

Here is the (GeoServer-generated) properties file for the image mosaic:
#-Automagically created from GeoTools-
#Wed Jun 17 10:00:07 CDT 2020
MosaicCRS=EPSG\:4326
Levels=8.333333333351735E-4,8.333333333335202E-4
Heterogeneous=false
Name=bathymetry
HeterogeneousCRS=false
TypeName=bathymetry
Caching=false
ExpandToRGB=false
LocationAttribute=location
PathType=RELATIVE
SuggestedSPI=it.geosolutions.imageioimpl.plugins.tiff.TIFFImageReaderSpi
CheckAuxiliaryMetadata=false
LevelsNum=1

Here's what GeoServer says about the Image Mosaic layer's CRS:


texts: Calling a function for two geojson files in Leaflet

texts:
I have two geojson files :

1- world_countries.geojson : each feature has properties "name" and "relation" as in the following example:

{
  "type": "FeatureCollection",
  "features": [
    {
      "type": "Feature",
      "properties": {
        "name": "France",
        "relation":{
          "Algeria": "Colony",
          "Spain" : "EU",
        }
      },
      "geometry": {
        "type": "Polygon", "coordinates": []
      }
    }
  ]
}

2- African_countries.geojson  :

{
    "type": "FeatureCollection",
    "features": [
        {
           "type": "Feature",
            "properties": {
                 "name": "Algeria",
                 "Population" : "30M",
                 "Capital" : "Alger",
            }
        },
   "geometry": {
    "type": "Polygon", "coordinates": []
    }
  }
]
}

Now when clicking on an african country (algeria) and other country (France) i need windows to show up with a message "Relation between Algeria and France is : Colony".


texts: How to select and move below vertex - Vertex Tool (QGIS 3.6)

texts:
I have a path (tracks) generated by points, which I snapped before to a grid (so there is some vertex at the way back that are overlayed), that at some point overlay each other. When you select in the Vertex Editor the first vertex, and this means a vertex that is below, it does not get highlighted.

When you have a single vertex selected it gets blue, but when you have one above another it is not possible to see where it is, neither move or edit it.
Is there a way to the selected vertex to be highlighted even if it is under another vertex?


texts:
I am trying to set up an animation showing lines which represent erosion contours moving back. Effectively the animation will show 1 whole line then the next line then the next line etc. for about 20 lines.

I've researched the time series and MMQGIS tools but both seem to focus on moving points along a line or a line which extends to its full length.

Is it possible to make a whole line change position again and again in animation.

In the animation I want L1 to appear then L2 then L3 etc for about 20 lines


texts: Adding custom layer to OSM data in switch2osm?

texts:
I'm new to switch2osm and I've created an OSM Tile Server with switch2osm on Ubuntu 16.04 LTS. It uses Apache, Mod_tile, renderd and Mapnik to render tiles. It uses PostgreSQL as DBMS and osm2pgsql to import Open Street Map shape files to database (I've downloaded my shapefiles from Geofabrik).

I want to add my own layers (like my GPS data of my organization) to PostgreSQL and style them through Mapnik to be shown on the tiles. How can I do this?


texts:
I want to implement a raster layer on top of a Leaflet map that includes the following functionalities:

Digital Elevation Model (DEM) terrain height included in each raster cell
GPS-Track (I use this plugin for GPS tracks. Is there another one I should use?)
Assign a color to each raster cell that was visited by a GPS point.

What is the easiest way to go about this? Are there any better plugins I should use?


texts:
We've recently had one of our GIS Techs create some QA/QC check via data reviewer. However no one seems to be able to use the checks he has created.

We placed the single file .rbj on our shared network
Our editing is versioned
Using 10.2
Error reads:

  Warning! There was an exception loading the batch job(s). Loading
  aborted

Any thoughts?


texts: Error 999999 when using .isalnum in ArcGIS field calculator?

texts:
ArcGIS keeps throwing an error when I use the following script in the field calculator, using Python:

pre-logic script code:

y= ()
if !fieldx!.isalnum():
    y= 1
else:
    y= 0

Fieldy=

y

The error in the results window reads:

  'ERROR 999999: Error executing function, Syntax Error, Failed to execture (CalculateField)'

I've checked the classes of the fields: fieldx is a string and fieldy a short integer, so that can't be the problem. I've also enabled the shapefile for editing. Note that I've added fieldy using the addfield function.

Other fieldcalculations, for example:

fieldy = 5

do work, so the problem should be in the syntax.

What am I doing wrong?


texts: QGIS relation_aggregate, how does this work?

texts:
I've been trying to use this feature but I'm not being able to make it work. I'm using this syntax:

relation_aggregate(relation := 'tubo_et_face_et_inicial', aggregate := 'sum', expression := "lenght_linear")

to try to sum up all the values in the "lenght_linear field of the child table, to a virtual field in the parent table. Is this syntax correct, or am I doing something wrong?


texts: Finding all possible adjacent/closest geometrical points from a given point

texts:
So, I have a dataframe like this,
import numpy as np
import pandas as pd
import descartes
from shapely.geometry import Point, Polygon
import geopandas as gpd
import matplotlib.pyplot as plt

df = pd.DataFrame({'Address':['280 Broadway','1 Liberty Island','141 John Street'],
                   'Latitude':[ 40.71,40.69,40.71],
                   'Longitude':[-74.01,-74.05,-74.00]
                    })

%matplotlib inline

geometry = [Point(xy) for xy in zip( df["Longitude"],df["Latitude"])]
crs = {'init':'epsg:4326'}
df = gpd.GeoDataFrame(df,
                     crs=crs,
                    geometry=geometry)
df.head()

I converted the lat and lon to geometry points and I am trying to find all possible closest points for each address using the geometrical points. For example, all possible closest points adjacent to 280 Broadway which lies next to each other for one block.There could be more than one point if the points are adjacent to each other containing in a polygon shape.
This was my approach but didn't really get what I wanted,
df.insert(4, 'nearest_geometry', None)
from shapely.geometry import Point, MultiPoint
from shapely.ops import nearest_points

for index, row in df.iterrows():
    point = row.geometry
    multipoint = df.drop(index, axis=0).geometry.unary_union
    queried_geom, nearest_geom = nearest_points(point, multipoint)
    df.loc[index, 'nearest_geometry'] = nearest_geom

Desired Output:
Address       Lat     Lon    geometry                    nearest_points
280 Broadway  40.71  -74.01  POINT (-74.01000 40.71000)  POINT(NEAREST GEOMETRIC POINT)


texts: Aggregating field on lines where they overlap in QGIS?

texts:
I have a network of lines representing a pedestrian footpath network. Each line intersects and follows another to an endpoint. Each line has a number representing population.

Is there an aggregate tool in QGIS that breaks the lines where they overlap and aggregate the population along that section, adding on each number at each intersection? For instance two lines with 12 and 15 meet, they follow a new segment with 27 until a line with 24 meets, the next section having 51.

I have this working with a SQL stored procedure but it is a bit cumbersome.


texts: Inserting lat/long data into geography column of SQL Server 2012?

texts:
is it possible to insert lat/long data into geography column just by typing the coordinates lat/long as a double or string:

INSERT INTO LOT_EXPLOITATION_LOT VALUES (3621,'oui','Forêt Boréale',560000,'Interdiction d''exploiter le bois, aux alentours du Lac Achinuseskasich, dans un rayon de 15 km.',(**'-75.45','45.36'**))


texts:
Trying to get this to work on a mobile device. Works fine on a computer. Anybody have any suggestions as to why it wont show up on my phone?

    <title>Leaflet Mobile App</title>
    <link rel="stylesheet" href="http://cdn.leafletjs.com/leaflet-0.7.2/leaflet.css">

    <meta name="viewport" content="width=device-width, initial-scale=1.0, maximum-scale=1.0, user-scalable=no">


texts: How to load geojson features outside view with Mapbox GL?

texts:
I am trying to load all features from an imported geoJSON. I want to do so because I need some propeties in order to create a list of all elements.

But when trying to access them it looks like I can only query the features inside the current view.

Is there a way to acces them even if they are not in the current view?


texts: Why is ArcGIS Server 10.0 getFeatureInfo not being buffered?

texts:
Using AGS10.0 I have a point featureclass published as wms. In my Openlayers client, I am not able to click on the point and recieve the GetFeatureInfo response unless I hit the target very accurately. Strangely, when I zoom in to the point it is harder to hit than when the map is zoomed out in a smaller scale.

Does anyone have a trick to be able to buffer the GetFeatureInfo request so that a getFeatureInfo request becomes more user friendly?

For more context I have sent this post to the esri forum (http://forums.arcgis.com/threads/108479-Set-tolerance-for-ArcGIS-Sever-10.0-GetFeatureInfo-query?p=384651&ESRISessionID=__kafq_4ATGCzLXerP68fkOpCCGe1ttHoAcLV9s54Wb-36WXv_w_8fA7JLQQUYXtsZCVwQ==)


texts: I have a kml file which displays correctly in Google Earth but only as a bounding box in QGIS . What can I do to see the associated raster image?


texts: Would a few thousand unmatched/tied NYC addresses be a source data issue?

texts:
Background: Am geocoding thousands of addresses, a process I am quite a novice too, using the NY State address composite locator. Using ArcMap 10.3.

14% remain tied and 12% unmatched. I'm inclined to think this is a source data issue but because of my inexperience I can't be sure.

A few notes:

The unmatched tied addresses are almost exclusively hyphenated addresses. But many of the addresses that were successfully matched are hyphenated as well.
It seems a handful of streets are the worst offenders for unmatched/tied addresses.
None of the matched addresses have a score below 100.


texts: Nominatim where does it get children info

texts:
I am trying to build a reverse geocoder that returns all the areas that a geo location point is located in. Neighbourhoods being the lowest level i am interested in.

I am using the Europe OSM file provided by geofabrik, and import it to PostGIS using Imposm3.

Unfortunately OSM stores neighbourhoods by point mostly, so i cannot really use them in my reverse geocoding result.

So i have started looking around what other services do, and i can see that nominatim has children relationships, but i cannot figure out how it aquires it, as i cannot find any such relationships in the raw OSM file.

An example is here: http://nominatim.openstreetmap.org/details.php?place_id=59349063
Centro Historico is a neighbourhood, and is just stored as a node with id: 4718852002.

If i had the children of that node as nominatim has it, i would be able to calculate the neighbourhood polygon.

So the question is how does nominatim get a hold of the children for a node like this?

Alternatively to get neighbourhood polygons i have looked at google, as they also have these lovely pink areas if you search for neighbourhoods - but as far as i understand they dont provide those in any API.

It is crucial for my project that i get the polygons, as i want to be able to find all locations inside these polygons, restaurants for example.


texts: Geoprocessor: using arc toolbox in vb.net

texts:
i am doing the following in order to use the predefined tool points to line in arctool box

        Dim ptoline As ESRI.ArcGIS.DataManagementTools.PointsToLine = New ESRI.ArcGIS.DataManagementTools.PointsToLine
        Dim pgeo As ESRI.ArcGIS.Geoprocessor.Geoprocessor = New ESRI.ArcGIS.Geoprocessor.Geoprocessor
        ptoline.Input_Features = " D:\Mayur\shapefiles\new_pooint_3.shp"
        ptoline.Output_Feature_Class = " D:\Mayur\shapefiles\new_line123.shp"

        pgeo.OverwriteOutput = True

        Dim geoProcessorResult As IGeoProcessorResult
        geoProcessorResult = pgeo.Execute(ptoline, Nothing)

i dont know the exact syntax of it. after running it shows error on

geoProcessorResult = pgeo.Execute(ptoline, Nothing)

Error HRESULT E_FAIL has been returned from a call to a COM component.

the references geoprocessor and datamanagement tools are included.

the path to the input feature is correct.
the o/p file is created automatically.

is there any prerequisite to run geoprocessor?
do I have to do that workspace thing also ?

all suggestions are welcome


texts: Why ground model from point cloud in Global Mapper is inaccurate?

texts:
I have a LAS file that I assumed was based on LiDAR, but it is actually based on ADS (airborne digital sensor) imagery, according to the body that captured the data.

I load the file in to Global Mapper and generate a point cloud which is totally unclassified. Using the auto-generate ground points feature, about 60% of the points become ground. However, when I generate an elevation grid, the areas that should be flat (in this case, fairways on a golf course) remain quite bumpy.

Is my point cloud inaccurate because it's not LiDAR?

Below are images of the area in question.


texts:
I have a set of 150 sites globally, and I want to map them, and colour by a given variable. Unfortunately they are very non-uniformly distributed, with many of them in small clusters of only a few hundred meters. This means that with any reasonable sized markers on a map, the markers overlap a lot, and many of the markers can not be seen. I would like a way to spread out the markers (with lines connecting them to their original position), so that each marker can be seen without overlap. The exact position is not super important here, but the approximate global position is.

I'm using python with basemap/matplotlib, and I've tried a kind of force-directed approach that kind of works, but not very well. Here is an example:

As you can see, the European points are not very well un-overlapped.
The code for that spread is:

def spread_points(lat_lon, d=5, max_iter=300):

    new_lat_lon = lat_lon.copy(deep=True)

    dists = scipy.spatial.distance_matrix(new_lat_lon, new_lat_lon)
    too_close = (0 < dists) & (dists < d)

    itr = 0
    while too_close.any() and itr < max_iter:
        itr += 1

        neighbour_means = []

        for i in range(len(new_lat_lon)):
            # for each point, find all "too close" points, and take their average.
            group_ind = too_close[i]
            if not group_ind.any():
                # no too close neighbours, use self.
                neighbour_means.append(new_lat_lon.iloc[i])
                continue
            group_ind[i] = True
            neighbour_means.append(new_lat_lon[group_ind].mean())

        max_diff = 0
        for i in range(len(new_lat_lon)):
            dir_grp_mean = (new_lat_lon.iloc[i] - neighbour_means[i])
            grp_dist = norm(dir_grp_mean)
            if grp_dist == 0:
                dir_grp_unit = 0
                grp_force = 0
            else:
                dir_grp_unit = dir_grp_mean / grp_dist
                grp_force = np.sqrt(d / grp_dist)

            # move each point away from the mean and toward home and jitter a tiny bit
            new_loc = new_lat_lon.iloc[i] \
                + grp_force * dir_grp_unit \
                - (new_lat_lon.iloc[i] - lat_lon.iloc[i]) / d \
                + np.random.randn(2) * 0.01 * d
            if not np.isfinite(new_loc).all():
                import ipdb
                ipdb.set_trace()
            max_diff = max(max_diff, norm(new_loc - new_lat_lon.iloc[i]))
            new_lat_lon.iloc[i] = new_loc

        if max_diff < d / 100:
            logger.info("Quitting after {i} iterations - differences are tiny")
            break

        dists = scipy.spatial.distance_matrix(new_lat_lon, new_lat_lon)
        too_close = (0 < dists) & (dists < d)

    if itr == max_iter:
        logger.info("hit max_iter")

    return new_lat_lon

There is another algorithm in Inkscape called removeoverlaps, which seems like it might work, but I'm having a hard time deciphering how it works. Edit: The results really aren't ideal anyway, because it moves points horizontally first then vertically, resulting in an odd layout:

It would be much better if points were spread radially...

Are there any other algorithms for this kind of operation that already exist in or would be easy to implement in python?


texts: What is the distinction and unit of the number written in the slope field

texts:
I added a line shapefile to a surface using the Arc map 10.3
In order to calculate the 3D length and slope  of each feature line
The result was as in the picture.

What is the distinction or the unit of the number written in the slope field?

How can I read it or convert it to a slope degree?


texts: Recovering layer (*.lyr) file that is missing shapefile?

texts: Somebody in the office deleted a point shapefile that has a lyr file as well. How can I recover it? Is there a Script that allows me to recover it?


texts: Does ArcGIS 10 have a "Match Properties" tool like AutoCAD

texts:
My question is short-and-sweet..

I have text, pointers & lines to be added to a series of maps.  ** Creating spatial data via an edit session is not an option **    This project is HUGE, and there will be hundreds of maps for each facet of the project.  Consistency is important.

My question .... is there a tool in ArcView with the same functionality as [Match Properties][1]?

It would save me SO MUCH TIME!!  This "one touch" functionality would be far superior to changing properties via the multi-stage process in ArcView (especially for rotation angles).  I appreciate an feedback.
Take care all


texts: Can't pan around my OpenLayers map .. unless there are no Layers under the mouse

texts:
So I've got my map displayed and some Features added to a Layer. All good, etc.

When I try and pan around my map, I can only pan when I left-click the mouse .. hold .. and drag the mouse ... when the mouse cursor is NOT over a feature/shape on the vector layer.

Otherwise, the panning works great.

What's the trick to say .. can pan the map at any time / anywhere the mouse cursor is?


texts: Find the next vertices at the intersect line

texts: I'm Trying to calculate the angle between two intersected lines at geometric network (simple edge)in order to check if the lines are close to be perpendicular. I already have an Array of the intersection vertices and a function that calculate the angle using Pythagorean theorem. How I can get the next vertices at the intersect line, in order to calculate the angle between this vertices and the intersect point using the function ?


texts: Wrong result while extracting NDVI in Google Earth Engine

texts:
I have here a FeatureCollection containing some productive forests of the company I work for.
The age in month for each forest site is important to me, since my study only looks for forests that are between 1 and 5 month of age.
With this age I calculated, at a desktop GIS software, the maximum date and the minimum date to filter the Landsat 8 images inside a period of the year (September to May).
So, with all these informations together as a FeatureCollection in my Assets, I have to  define a function to be mapped to the whole FeatureCollection containing company's forests.
The function needs to get values from columns MIN_DATE and MAX_DATE to filter the Landsat images, apply cloud masking, calculate NDVI and return the feature to the collection with the mean of NDVI computed as a column.
I've already developed something here, but it keeps returning me some weird result.
The code:
/ Script to get historical NDVI for each productive area

// Import FeatureCollection
var IPTalhoes = ee.FeatureCollection('users/ee_rootDir/toiFinal');

// Define function to map over collection
function calcHistoricalNDVI(feature) {

  // Get initial and end date
  var dateIni = feature.get('MIN_DATE');
  var dateEnd = feature.get('MAX_DATE');

  // Import Landsat 8 ImageCollection
  var imgCollection = ee.ImageCollection("LANDSAT/LC08/C01/T1_TOA");

  // Filter imgCollection
  var imgFiltered = imgCollection
                    .filterDate(dateIni, dateEnd)
                    .filterBounds(feature)
                    .map(function updateCloudMask(image) {
                      var mask = ee.Algorithms.Landsat.simpleCloudScore(image);
                      var imgUpdated = image.updateMask(mask.select('cloud').lt(20));
                      return(imgUpdated);
                    });

  // Calculate NDVI
  var imgNDVI = imgFiltered
                .map(function calcNDVI(image) {
                  return(image.expression(
                  '(b(4) - b(3)) / (b(4) + b(3))'
                  ).select(['B5'], ['NDVI']));
                });

  // NDVI mean inside a period
  var NDVIMean = imgNDVI.mean();

  // Reducer statistics - mean by talhao
  var NDVIMean2 = NDVIMean.reduceRegion({
    reducer: ee.Reducer.mean(),
    geometry: feature.geometry(),
    scale: 30
    });

  // Turn dictionary with NDVI values to FeatureCollection
  var NDVIMean2Feat = feature.set('NDVI', NDVIMean2);

  // Return Feature as function result
  return(NDVIMean2Feat);
}

// Map main function over collection
var NDVITalhoes = IPTalhoes.map(calcHistoricalNDVI);

// Print results
print(NDVITalhoes);

The issue:
When I print the Feature Collection, what I see is that there is a column named 'NDVI' in FeatureCollection with an "Object" within, but when I open each feature in console, I can't find any property with this name, and indeed there are only 27 properties while there are 29 columns for the whole FeatureCollection.
Can anyone help me out?


texts:
I have a raster and need to only identify cells that are adjacent within a 3x3 window.

I tried various tools but without any success. My current cells are 1x1 and currently using ArcGIS.


texts: Local Moran's I on spatially heterogeneously sampled dataset

texts:
I have a dataset where points are sampled in space. Values are attached to the points, and I'd like to check if these values are correlated in space with Geoda's local Moran's I.
However, the sampled data leads to a naturally very clustered spatial distribution of points. As such:

Does this forbid an analysis such as local Moran's I? If so, maybe choosing a large enough distance threshold to compute the weighting matrix could make it work?


texts: Setting composer layout units to inches in QGIS

texts:
I open a new composer and select ANSI A for my paper. The notation after the name is "Letter; 8.5inx11in", but the width and height are displayed in mm, and the units selection is grayed out, so I can't change it. The rulers at the top and left of the layout are in mm, and the parameters for "Snap to grid" are also in mm. If I then select "Custom", I can select inches for units, and the Width and Height display as 8.5 and 11. However the rulers in the layout are still in mm, as are the units for the grid.

Can the entire composer project be converted to inches?

qgis 2.01 64 bit
windows 7, 64 bit


texts:
In the Esri blog Print your Story Map Journal the instructions for doing so start:

  Step 1. Open the Story Map Journal.

  Step 2. Click the share link found in the app header.

  Step 3. Click Print in the upper right corner.

  This will open the Map Journal in print preview mode.

When I try to do this, on a Story Map Journal which I have published, using the steps below, and Internet Explorer 11, I get stuck.

Open William & James Slight - Map Engravers & Test Cricketer which is a Story Map Journal
Click the share link found in the app header (it has tip that says "Get a shot link")

Click Print in the upper right corner.

This should open the Map Journal in print preview mode but what I observe instead is that it opens a web page that says:

  To access the printable version of this Story Maps [sic]: open the story and
  use the print button located in the share dialog.

I think that is precisely what I did to get the same point so I seem to be stuck in a loop.

Does anyone know how I can print this particular Story Map Journal?


texts: Finding QSettings keys from source code?

texts:
I'm looking to change a number of QGIS settings with a python script. Whilst I can find references to some specific keys through various questions/answers on StackExchange, I can't find a method to actually look these up myself.

The question  How do I access QGIS program settings programmatically? indicates that I need to go into the source code to look these up. However, I'm only a novice with this scripting stuff, and don't really know where to start or what to look for.

Can you provide a generalised set of steps for finding these keys in the source code?


texts: How to find lat/lon values for every pixel in a GeoTIFF file?

texts: I want to a store GeoTIFF file pixel by pixel in the database. I need to store pixel values along with the lat/lon of that pixel on the GeoTIFF image. Is it possible to do that?


texts:
Can you build a login window before loading QGIS interface?

I use QGIS 2.8

I just built message windows when loading QGIS interface, but no more.


texts: Calculating areas of polygons that match with overlapping raster threshold value?

texts:
I have a polygon of districts (districts.shp) and a raster of rainfall (rainfall.tiff).

I want to calculate the area of each district that saw rainfall greater than various numbers.

I'd like to produce a table as follows

District    >50mm    >100mm    >150mm
District1    100km2    50km2    20km2
District2    200km2    80km2    10km2

I did use raster calculator to give me rasters of 0s and 1s when this was true. But I think this is cumbersome, and I'm not sure of the next step. I think zonal statistics is one way to go.


texts: Record time for each vertex of line in QField

texts: I use QField 1.9.5, and QGIS 3.16, and I use GeoPackages to record vector data in QField. I use a line layer to record my track. I also use QField to record photographs taken with my tablet as points. I would like to record the time for each vertex of my track, so this can be extracted and used to geotag high-resolution photographs taken on a DSLR, where I have synchronized the camera's time with my tablet. Is it possible to record the time for each vertex in a line?


texts: Trouble shooting ogr2ogr PostGres conversion to shapefile

texts:
I have a postgres database that I am trying to do a bulk export of to individual shapefiles. The database consists of several hundred tables.

If I run this command:

ogr2ogr -f "ESRI Shapefile" dbdata PG:"host=someHostname user=someUser dbname=someDB password=somePassword"

It exports a small portion of the total table's to a shapefiles with no errors.  When I specifically target one of the tables that failed to convert it results in only the *.dbf file being created.  The ogrinfo command returns only the tables that are successfully created.  When I connect ArcMap to the original databse and add the tables that are not being converted to a map they all render, so I know they are "good".

Is there a log or command I can use to output why the tables are failing to be converted?

Version Number: GDAL 2.1.0, released 2016/04/25 postgres 9.3


texts: How to improve performance when using ArcGIS cursors in Python with big tables?

texts:
I have a pretty big point feature class in a file geodatabase (~4 000 000 records). This is a regular grid of points with a 100m resolution.

I need to perform a kind of generalization on this layer. For this, I create a new grid where each point lies in the middle of 4 "old" points:

 *     *     *     *
    o     o     o
 *     *     *     *
    o     o     o
 *     *     *     *

[*] = point of the original grid - [o] = point of the new grid

The attribute value of each new point is calculated based on the weighted values of its 4 neighbors in the old grid. I thus loop on all the points of my new grid and, for each of them, I loop on all the points of my old grid, in order to find the neighbors (by comparing the values of X and Y in the attribute table). Once 4 neighbors have been found, we get out of the loop.

There is no methodological complexity here but my problem is that, based on my first tests, this script will last for weeks to complete...

Do you see any possibility to make it more efficient? A few ideas on the top of my head:

Index the fields X and Y => I did that but didn't notice any significant performance change
Do a spatial query to find the neighbors rather than an attribute-based one. Would that actually help? What spatial function in ArcGIS should do the job? I doubt that, e.g., buffering each new point will prove more efficient
Transform the feature class into a NumPy Array. Would that help? I haven't worked a lot with NumPy so far and I wouldn't like to dive into it unless someone tells me it might really help reducing the processing time
Anything else?


texts: QGIS 2.6 summarizing area (m²) for current atributes

texts: If I have a code for each polygon and I have counted the area in m² for each of them (all this information is in single columns). How can I summarize these areas by codes (for example summarize all areas with code 100 etc.) in QGIS 2.6.1.?


texts: How do I start to debug COM error "HRESULT E_FAIL" for TableToGeodatabase geoprocessing tool?

texts:
I've got some code that creates a table using OLEDB connections to an access table.  I want to convert this table to a Geodatabase table, so I ask a question How do I convert an access table to a proper geodatabase table?  I tried to use the code there, which threw HRESULT E_FAIL, and then after reading Blah's link about using the geoprocessor I came up with this code:

Geoprocessor GP = new Geoprocessor();

ESRI.ArcGIS.ConversionTools.TableToGeodatabase tableToGeodatabaseTool =
    new ESRI.ArcGIS.ConversionTools.TableToGeodatabase();

tableToGeodatabaseTool.Input_Table = _cSM.Geodatabase_FilePath + "\\" + TableName;
tableToGeodatabaseTool.Output_Geodatabase = _cSM.Geodatabase_FilePath;

GP.Execute(tableToGeodatabaseTool, null);

... Which also throws HRESULT E_FAIL.  I'm not sure why, and I'm not sure how to go about debugging it.

This code is being run from a custom toolbar in ArcGIS.  It seems like ArcGIS is locking the access database, but I assume because I'm using ArcObjects that's not the problem. Any ideas?

Edit: Here are the details to the exception that is being thrown.  You will note that "inner exception" is empty; it is returned as null:

System.Runtime.InteropServices.COMException was unhandled by user code

  Message="Error HRESULT E_FAIL has been returned from a call to a COM component."

  Source="ESRI.ArcGIS.Geoprocessing"

  ErrorCode=-2147467259

  StackTrace:

   at ESRI.ArcGIS.Geoprocessing.GeoProcessorClass.Execute(String Name, IVariantArray ipValues, ITrackCancel pTrackCancel)

   at ESRI.ArcGIS.Geoprocessor.Geoprocessor.ExecuteInner(IGPProcess process, ITrackCancel trackCancel, IGeoProcessor igp, IVariantArray iva)

   at ESRI.ArcGIS.Geoprocessor.Geoprocessor.Execute(IGPProcess process, ITrackCancel trackCancel)

   at Reports.Validation.validate_StumpageValues() in C:\Users\Public\Documents\SilvAssist\Reports\Reports\Validation.cs:line 131

   at Reports.frm_Options.run_Validations() in C:\Users\Public\Documents\SilvAssist\Reports\Reports\frm_Options.cs:line 160

   at Reports.frm_Options.bgWorker_Validations_DoWork(Object sender, DoWorkEventArgs e) in C:\Users\Public\Documents\SilvAssist\Reports\Reports\frm_Options.cs:line 190

   at System.ComponentModel.BackgroundWorker.OnDoWork(DoWorkEventArgs e)

   at System.ComponentModel.BackgroundWorker.WorkerThreadStart(Object argument)

InnerException:


texts: Putting parentheses around label value in ArcGIS for Desktop?

texts:
How can I create a label expression to put parentheses around an attribute value I need labeled on a map?

I'm using ArcGIS version 10.3.


texts: I have used the Topology Checker plugin in QGIS to create rules to identify errors in my layer. For example, I have used it to search for all overlapping polygons. This plugin is great for finding such errors, but how can I select them in my original layer. I have a very large data set and doing them by hand is not an option.


texts: Data View/Layout View Changes - ArcMap 10

texts:
I have code that will hide a certain set of points from the view using an ExclusionSet, and then switch to layout view for printing.

I noticed in my 9.1 version, when after I am done printing, and I switch back to data view, and the points I had before running the code are shown. This was my intention.

In 10 on the other hand, the code runs the same, but when I return to data view, the feature exclusions are still in-tact. I have to manually remove them via the display properties of the feature layer.

Is there anyway to have separate properties for both views? Whether it be code or manual?

Thanks!


texts: Do I have to reproject the coordinate system when two layers overlap already?

texts: I have two layers in QGIS, one the coordinate system of WGS84 and another uses GRS80. The two layers overlap, which looks good. My question is do I have to reproject them to the same coordinate system or just leave them as they are? Someone told me that reprojection could introduce inaccuracies to the data. So does it matter if I don't reproject? Under what situations is the reprojection needed?


texts: leaflet - change routing vehicle when clicked on button

texts:
I want to change the routing vehicle when clicked on a button. The default routing vehicle is "car" ("auto") when clicked a button routing vehicle should change to "foot" ("pedestrian").

I am using leaflet.easy.button and leaflet.routing.machine plugins, as routeservice I use mapzen.

Links:

leaflet.routing.machine:
https://github.com/perliedman/leaflet-routing-machine
leaflet.easy.button:
https://github.com/CliffCloud/Leaflet.EasyButton
mapzen:
https://github.com/mapzen/lrm-mapzen

My code:

Routing.control :

var myRouter = L.Routing.control({
    router: L.Routing.mapzen('valhalla-ABCDEFG', {
        costing:"auto"}),
    routeWhileDragging: true,
    geocoder: L.Control.Geocoder.nominatim(),
    formatter: new L.Routing.mapzenFormatter()
}).addTo(map);

Easy.Button :

L.easyButton('<img src="pic.png">', function(btn,map){
    myRouter.route({costing: "pedestrian"});
}).addTo(map);

Any ideas what could be wrong? Or maybe a different way to approach this problem.


texts: Choose what zoomlevel to read from for mbtiles file in GDAL

texts:
GDAL is supposed to read MBTiles files and to do so, it selects the highest zoomlevel available. However, is there a way to choose another zoomlevel? i.e. override the default and go for a lower one?

I'm currently accessing the file like this (and it defaults to the highest zoomlevel):

gdal_translate tiles.mbtiles export.tiff


texts: How to add additional weight to shortest path?

texts: I am calculating shortest bath between two points. I need to add an additional feature eg. air toxicity in my calculation, so that the shortest path calculated takes into account the air toxicity and gives the least toxic path. How can I  do this?


texts: How can I made requests (QgsFeatureRequest) or selections (selectByExpression) from a list of QgsFeature?

texts:
I've got a processing algorithm that does a request to select some features that correspond to the incorrect filling of the attributes. It reads a JSON file containing a description of expression (ruleName), layer's name, and a list of expressions (allRules), just like below.
{
    "ruleName": {
        "layer": {
            "allRules": ["rule 1", "...", "rule n"]
        }
    }
}

I must give to the users an option to run the algorithm only on selected features. In this case, the type of object received change from QgsVectorLayer to a list of QgsFeature, therefore, and I'm stuck here because I can't make any requests through QgsFeatureRequest() or selection using the selectByExpression() method, follow the code below.
    def checkedFeatures(self, rules, layerList, onlySelected, returnIterator=True):
        """
        Select features by conditional rules stored in the rules JSON file.
        :param (dict) rules: dictionary from conditional rules;
        :param (QgsVectorLayer) layerList: list from all loaded layers;
        :param (bool) onlySelected: boolean.
        """

        failedDict = {}
        for ruleName in rules:
            failedList = []
            for lyr in layerList:
                loadedLyrName = lyr.name()
                if loadedLyrName in rules[ruleName]:
                    allRules = rules[ruleName][loadedLyrName]['allRules']
                    for rule in allRules:
                        request = QgsFeatureRequest().setFilterExpression(rule)
                        if onlySelected:
                            selected = lyr.selectedFeatures()
                            # this FOR block below receives a QgsFeature list from lyr.selectedFeatures()
                            # and so I can't make any requests through QgsFeatureRequest()
                            # or selection using the selectByExpression() method
                            for feat in selected: # I cant do the selection here
                                feat.selectByExpression(rule)
                                failedList+=[feat for feat in lyr.selectedFeatures()]
                        else:
                            for feat in lyr.getFeatures(request):
                                failedList.append(feat)
            failedDict[ruleName] = failedList
        return failedDict

The question is: how can I made requests or selections from a list of QgsFeature? Am I missing something?


texts: Getting lat/lng of an object X inches from the GPS?

texts:
I'm working on a project where I have a GPS and two receivers to create a triangle with a transmitter. The problem I'm trying to figure out right now is how can I get the position of the two receivers?

I know it can be done, but I'm blanking on the math and Google hasn't been too helpful, maybe I'm not searching it correctly.

Here's my circumstance:
GPS is reading lat: 42.331429 and lng: -83.045753. My device is sitting perfectly level and pointing due north. My two receivers are 30 inches to the left and 30 inches to the right. Can anyone help me figure out the equation to get their positions? I know these are going to be very small differences, but that's as far as I've gotten so far.

Also, if things were more complex, how would that be solved? Such as pointing North/West, and 5 degrees rotated on each axis, x, y, and z? I realize that's kind of a crazy scenario, but I like to potentially prepare for every scenario. These are based on being in fixed positions from the GPS and the plane that the GPS is on.


texts:
In one folder I can have one or more xml files. I would like (need) to use them one by one as input dataset (reader). Actually each xml file will be used within according fme workbench regardless some values within xml file.
This question is linked with my old question Variable XML dataset as reader in FME wher I got answer for one file in folder.

If I have two xml files,  for instance, fme read them both at the "same time" and use values from both files for one fme wb.


texts: OSM Reverse Geocoding - Find closest node for Coordinates

texts:
currently I am experimenting with Neo4j as Database for OSM Data. In the end I'd like to do routing stuff. I had a look at the Neo4j/spatial, but this is to blown up for me (not optimized for routing). Maybe I am just missing the information on how to use the spatial plugin for my already loaded OSM Neo4j database.

I need to find the closest node to given coordinates. I tried nominatim and photon, but both do not provide the ability to select the node_type. Is there another reverse geocoding provider that provides this feature?

Or would it be the better approach to create something like a locationIndex that allows me to query my own dataset?


texts: Working with ambiguous AddGeometryColumn functions

texts:
I'm using PostGIS to add geometry columns to tables, but keep receiving the error:

ERROR: function addgeometrycolumn(character varying, character varying, integer, character varying, integer) is not unique
SQL state: 42725
Hint: Could not choose a best candidate function. You might need to add explicit type casts.

The statement I am trying to call is:

SELECT addGeometryColumn('schema', 'table', 'geom', 3857, 'POINT', 2);

I have tried type-casting each argument but I received the same message. I have also removed and re-installed the postgis extension, which made no difference.

I have used this function countless times before and am somewhat puzzled as to why it isn't working all of a sudden. I have also tried a few other postgis functions which gave me similar errors.

I am considering reinstalling postgresql/postgis but I will have to spend time reconfiguring the set-up so I am hoping there is an easier way to fix/get around this.

Does anybody have any advice as to what may be causing this?


texts: Randomly subsetting % of polygons by class/attributes using ArcPy

texts:
Referring to this topic How to randomly subset X% of selected points?. I tried in ArcPy and I can't figure out how to change this code to randomly select 50% of the polygons based on the attribute (CLASS) So that each set will have equal attributes (ex. 50% of forests, water, meadow in one set and 50% of forests, water, meadow in the other set).

def SelectRandomByPercent (layer, percent):
#layer variable is the layer name in TOC
#percent is percent as whole number  (0-100)
if percent > 100:
    print "percent is greater than 100"
    return
if percent < 0:
    print "percent is less than zero"
    return
import random
fc = arcpy.Describe (layer).catalogPath
featureCount = float (arcpy.GetCount_management (fc).getOutput (0))
count = int (featureCount * float (percent) / float (100))
if not count:
    arcpy.SelectLayerByAttribute_management (layer, "CLEAR_SELECTION")
    return
oids = [oid for oid, in arcpy.da.SearchCursor (fc, "OID@")]
oidFldName = arcpy.Describe (layer).OIDFieldName
path = arcpy.Describe (layer).path
delimOidFld = arcpy.AddFieldDelimiters (path, oidFldName)
randOids = random.sample (oids, count)
oidsStr = ", ".join (map (str, randOids))
sql = "{0} IN ({1})".format (delimOidFld, oidsStr)


texts:
I'm using rgdal to process rasters and shapefiles that contain CRS information [CH1903+ LV95 (EPSG2056) or CH1903 LV03  (EPSG:21781)]. For some reason part of the CRS information is always discarded (the datum).
Why is it discarding the datum and how can I preserve that information?
> class_ch = raster("ras_ch_25m_class.tif")

Warning: In showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj = prefer_proj):
  Discarded datum Unknown based on Bessel 1841 ellipsoid in CRS definition

> crs(class_ch)
CRS arguments:
 +proj=somerc +lat_0=46.9524055555556 +lon_0=7.43958333333333 +k_0=1 +x_0=2600000 +y_0=1200000 +ellps=bessel +units=m +no_defs

The same when loading a shapefile:
> lfi = shapefile("LFI.shp")

In OGRSpatialRef(dsn, layer, morphFromESRI = morphFromESRI,  ... :
  Discarded datum CH1903 in CRS definition:
+proj=somerc +lat_0=46.9524055555556 +lon_0=7.43958333333333 +k_0=1 +x_0=600000 +y_0=200000 +ellps=bessel +towgs84=674.374,15.056,405.346,0,0,0,0 +units=m +no_defs

> crs(lfi)
CRS arguments:
 +proj=somerc +lat_0=46.9524055555556 +lon_0=7.43958333333333 +k_0=1 +x_0=600000 +y_0=200000 +ellps=bessel +towgs84=674.374,15.056,405.346,0,0,0,0
+units=m +no_defs

I did some research and found some literature about PROJ4 strings being replaced by PROJ6 strings with GDAL>=3.0.0 (I'm using GDAL 3.0.4) and notations being deprecated and so forth. Since I'm not very firm in the field, I couldn't make much sense of it. Apparently the used versions can create problems, however I'm using a fresh installation of R, with GDAL 3.0.4, rdal 1.5.18, sp 1.4.4:
rgdal: version: 1.5-18, (SVN revision 1082)
Geospatial Data Abstraction Library extensions to R successfully loaded
Loaded GDAL runtime: GDAL 3.0.4, released 2020/01/28
Path to GDAL shared files: .../win-library/4.0/rgdal/gdal
GDAL binary built with GEOS: TRUE
Loaded PROJ runtime: Rel. 6.3.1, February 10th, 2020, [PJ_VERSION: 631]
Path to PROJ shared files: ../win-library/4.0/rgdal/proj
Linking to sp version:1.4-4

When I continue to work with the data, I run into errors since areas presumingly don't overlap, however this is most certainly due to the discarded CRS-info (since some data is in LV03 and the other in LV95, however that information got discarded).
I thought maybe I can reassign the CRS-information, however I get the same warning when creating the CRS object with a string.
crs_LV95 = CRS(SRS_string="EPSG:2056")

Warning: In showSRID(SRS_string, format = "PROJ", multiline = "NO", prefer_proj = prefer_proj) :
  Discarded datum CH1903+ in CRS definition

> crs_LV95
CRS arguments:
 +proj=somerc +lat_0=46.9524055555556 +lon_0=7.43958333333333 +k_0=1 +x_0=2600000 +y_0=1200000 +ellps=bessel +units=m +no_defs

I'm lost here... Is there a different way in which I should work with the data?
I haven't worked with rgdal for a long time and not much at all, I'm uncertain whether I'm overseeing something simple here. But something as basic as the CRS shouldn't be a problem/cause warnings on such basic operations such as creating a CRS object or loading data.


texts: Check if Folder exists on ArcGGIS Server

texts:
I am writing a python script to deploy a .sd file with the arcpy.UploadServiceDefinition_server method (http://desktop.arcgis.com/en/arcmap/10.3/tools/server-toolbox/upload-service-definition.htm). The 5th parameter for this method, in_folder_type, request that I specify that if the destination folder is a new folder or an existing folder. I may not know this at the time that the script is executed, so I would like a way to check if the folder I want to deploy my script to exists already and create it if it doesn't. I tried deploying using the "Existing" flag and a folder name that doesn't exist yet in case that would create the folder for me but the script just put my service in the root of the folder.

Does anyone know a good way in Python or using arcpy to resolve this?


texts: Reproject global latlong into UTM slices

texts: I'm working with a global shapefile in R that is latlong projected. I'd like to reproject it into the appropriate UTM projection for individual zones. But, I'd like to chop it up into those zones before reprojection, instead of having multiple global reprojections into UTM most of which have meaningless coordinates. Also, computational efficiency. So, how would I split a latlong shapefile into areas corresponding to UTM zones before the reprojection (and a little overhang is fine)?


texts: Using layer loadstart & loadend events in OpenLayers 3?

texts:
OpenLayers 2 has these layer events "loadstart & loadend."

What is equivalent to them in OpenLayers 3?

While a vector layer loads and is rendered, I need to show a loading icon.


texts: Exporting Imagery in ArcGIS - problems with lines on exported image?

texts:
I am using the pan-sharpened Landsat imagery that is available on ArcGIS Desktop (ie data from ArcGIS Online).

When I export a map with the pan-sharpened Landsat imagery in the background, lines show up on the image. However, this does not show when I am looking at the image in ArcGIS 10. I am enclosing an excerpt of both the exported map (top) and the screenshot of how it looks in ArcGIS (bottom) so you can see the problem. Can you please help me – how can I export the map without the lines?

Someone from Esri suggested that I turn off the DRA however I think that I am unable to do this as it is a file that I brought in from ArcGIS Online. When I access the Image Analysis window, it seems like the DRA box is not clickable. Also, he suggested that I use calculate statistics or build pyramids and statistics tool on the source data. Again, I don’t think I can do this since it is a file from ArcGIS Online.

If someone can help me do any of the above with the ArcGIS Online file or how I might be able to fix this problem in a different way, please let me know what I might be doing wrong.