texts: ArcPy to select by location not showing selected features in attribute table

texts:
I have a button where it first calculate no of features in a specific layer and then runs the "are_identical_to" and update "D". After that it checks if the no of counts in two layers is same? If not, then runs the "are_identical_to" on two different layers and should show up the selected features in my data frame.

However it never shows the selected features in data frame or attribute table?

How can I achieve this.

code is shown below:

import arcpy
import pythonaddins

class Check1(object):
    """Implementation for Preproduction_Check_addin.Check1 (Button)"""
    def __init__(self):
        self.enabled = True
        self.checked = False
    def onClick(self):
        lyr1 = ""
        lyr2 = ""
        lyr3 = ""
        result2 = ""
        result3 = ""
        count2 = ""
        count3 = ""
        mxd = arcpy.mapping.MapDocument("CURRENT")
        layers = arcpy.mapping.ListLayers(mxd)
        for layer in layers:
            if layer.isFeatureLayer and layer.workspacePath != "":
                if "edits" in layer.name:
                    lyr1 = layer
                elif "fresh" in layer.name:
                    lyr2 = layer
                    result2 = arcpy.GetCount_management(lyr2.dataSource)
                    count2 = int(result2.getOutput(0))
                    print("Total Features in fresh data:")
                    print(count2)
                elif "D" in layer.name:
                    lyr3 = layer

        Selection = arcpy.SelectLayerByLocation_management(lyr1, "ARE_IDENTICAL_TO", lyr2)
        if arcpy.Describe(lyr1.name).FIDSet:
            with arcpy.da.UpdateCursor(lyr1.name, "Update") as curs:
                for row in curs:
                    row = ("D",)
                    curs.updateRow(row)
                del curs
                del row
                print("Remark 'D' Updated")

        result3 = arcpy.GetCount_management(lyr3.dataSource)
        count3 = int(result3.getOutput(0))

        print("Total features in 'D' Layer:")
        print(count3)

        if count3 == count2:
            print("Everything OK!, Please continue.")
        else:
            Selection = arcpy.SelectLayerByLocation_management(lyr1, "ARE_IDENTICAL_TO", lyr3,0,"SWITCH_SELECTION")
            print("Some features not matched, Please review.")


texts: Information about algorithm used by Terrain models

texts: I have calculated slope and aspect using the DEM (Terrain models) GUI in QGIS Essen. I was also able to calculate the profile and tangential curvature using the toolbox r.slope.aspect. I was wondering if anyone can point me in the direction of documentation that describes the algorithm used to calculate these terrain parameters. I would like to know HOW they are calculated not how to use the plugins. For example, I would like to now how large of an averaging window the program uses, how it choses direction, what it does on the edges of a raster, etc.


texts: How to Determine to world map size in pixel given this data?

texts:
How to Determine the world map size in pixels given this data?

Map size of a particular area say, for example:

Makati city, Philippines with bounds:

Top-Right: 121.0281, 14.56497

Bottom-Left: 121.0144, 14.54637

Image size: width="4095px" height="4095px"

Scale: at zoom level 17 | 1:4,000

Now my question is... how can I know the size of the world map in pixels using these available data?

I'm still new to  geographic information systems.... so please if possible make it understandable for a novice like me... :)


texts: Copying selected (targeted) features to new shapefile using ArcPy?

texts:
I am trying to copy specific attribute table rows from one shapefile to a new shapefile that will contain only the specified rows in its attribute table, that is, if a field GDP_TOP10 = 1.

I am using this at the moment:

#Grab rows using Search Cursor with WHERE clause
rows = arcpy.SearchCursor("C:/CGDP.shp", "GDP_TOP10 = 1")
for row in rows:
    arcpy.CopyFeatures_management("CGDP.shp", "New.shp")
del row
del rows

However, this just copies the entire shapefile, not the rows I want to select with the cursor. If I execute a print statement under my for loop instead, it works and prints only the relevant fields, but the copy to new shapefile command does not seem to copy specifically.

Is there a different ArcPy command to only copy rows that match a query?


texts: Shapefile doesn't fit in CartoDB basemap

texts: When I import some shapefiles in CartoDB they don't fit properly at street level, as you can see in the screenshot. In this case, the shp was generated in QGis with CRS EPSG:29193 - SAD69. I have already tried converting to others CRS with no results. Does anyone have any idea on how to fix this?


texts: How can you identify the geoprocessing tool that was used to create a layer in ArcMap?

texts: In my Table of Contents, I have a layer that is named "symmetrical differences" but I am not certain it was created using that tool.


texts:
I am trying to edit a not-spatial CSV table in PyQGIS 2.18 (i.e. deleting existing fields and adding new one ) for the first time (I can do it with polygonal vector layers but I guess something is different). I just wrote the script below but it doesn't work (i.e. it does not delete fields and returns false when starting editing). I also tried to see the files of editablegeocsv plugin, as suggested in another post by underdark, but I am a beginner and I did not understand it. PS: I cannot try with last Python / QGIS releases.
    from PyQt4.QtCore import *
    from PyQt4.QtGui import *
    import processing, math, numpy, csv

    uri= "file:///D:/.../PROVA.csv?delimiter=,"
    layer = QgsVectorLayer(uri, 'name', 'delimitedtext')
    if layer.isValid() is True:
        print "OK!!"
    else:
        print "NOT VALID!"#it returns "OK!!"

    caps = layer.dataProvider().capabilities()
    if caps & QgsVectorDataProvider.DeleteAttributes:
        res = layer.dataProvider().deleteAttributes([0,1,7])
    layer.updateFields() #it neither does return any message nor delete any column

    #and then,
    layer.startEditing() #it returns "False"
    layer.dataProvider().addAttributes([QgsField("NewField",QVariant.Int)])
    layer.commitChanges()


texts:
I'm trying to do an operation on a polygon layer to ask a different layer and if a condition is met then do an operation. It would be something like that:

Case
when ("DAP-OSM-BuildingsH-0m",'BuildCat') =1
THEN intersecting_geom_count('DAP-OSM-BuildingsH-0mClass1')

where DAP-OSM-BuildingsH-0m is the other layer.
"BuildCat" is a column in the layer DAP-OSM-BuildingsH-0m
and the condition is If in the layer:

DAP-OSM-BuildingsH-0m and the value of the column "BuildCat" = 1 then the operation is:
intersecting_geom_count('DAP-OSM-BuildingsH-0mClass1')


texts: Do Esri developer certifications expire?

texts: Does the Esri Web Application Developer certification expire after some period of time?


texts: Can I have my facilities load just once using Closest Facility when using ModelBuilder ArcGIS?

texts:
I have tried to use Network Analyst as an extension and came across a memory issue when trying to run a Closest Facility query.

The query I want to run is routing zone origin and destinations, as facilities and incidents.  I have 2650 zones, so herein lies the memory issue as the query is routing to the network (as I'm also routing to the network, not straight line), meaning it's 2650 x 2650 routes.

I thought the best to do is to use ModelBuilder, so I could get the model to export the routes after each iteration or after a batch.  Help save on the memory.

My problem is the amount of time this is taking, I have set up a feature selection iterator so one facility is loaded at a time.  However, each time I run the model it loads in the 2650 incidents.

So my question is, can I set up ModelBuilder to not load in the locations each time?  As the loading of the locations takes about 6 hours to run.


texts: Creating unique field in shapefile for editing in ArcMap

texts: How do I make a unique and not null field in my attribute table on ArcMap such that when creating new features(through vectorizing) and populating the fields, it gives and error if the field value already exist or is null and doesn't save that feature?


texts: Rotating shape around defined point instead of centroid in QGIS?

texts:
The rotate feature in QGIS 3 rotates shapes around the centroid of that shape (with our without a defined angle).

When I start to rotate, I see that a cross appears at the centroid of that shape.

I would like to move that cross to rotate the shape around that cross.

I tried to move the cross with both left mouse clicks (no solution) and right mouse clicks (QGIS crashes after this action).

Are there any known ways to rotate shapes (just as in ArcMap) around a given point?


texts:
I have installed the contour plugin on my kubuntu 13 laptop, but when I am opening qgis (2.0.1) I get the following message:

Couldn't load plugin contour due an error when calling its classFactory() method

Traceback (most recent call last):
  File "/usr/lib/python2.7/dist-packages/qgis/utils.py", line 204, in startPlugin
    plugins[packageName] = package.classFactory(iface)
  File "/home/maria/.qgis2/python/plugins/contour/__init__.py", line 47, in classFactory
    from contour import Contour
  File "/usr/lib/python2.7/dist-packages/qgis/utils.py", line 453, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)
  File "/home/maria/.qgis2/python/plugins/contour/contour.py", line 41, in
    from shapely.geometry import MultiLineString, MultiPolygon
  File "/usr/lib/python2.7/dist-packages/qgis/utils.py", line 453, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)
ImportError: No module named shapely.geometry

How can I solve this problem?


texts: Unknown encoding format when reading data from USB COM port using Trimble Pathfinder ProXRT receiver


texts:
Given a projection and coordinates of a single point, how do I calculate the meridian convergence at that point?

Is it even possible? I'm unable to find anything useful using Google.

If it is not possible, is there a generic way this can be done in JavaScript? I need to support different projections, so a single hard-coded formula will not do.


texts: Calculating distance from points to nearest "eastern" polygon using ArcGIS Desktop?

texts:
I am using ArcGIS 10, and trying to do the following:

I have a layer file of points (industrial facilities), and I would like to calculate the distance from each point to the nearest "eastern" polygon (U.S. state).  I know there are tools for calculating distances from points to polygons, but I don't know how (or, if you can) specify direction of the polygon.

Any advice?


texts: How to print a label layer in the ArcGIS Server JS API

texts:
I'm trying to use a PrintTask to print a map, but am finding that a LabelLayer does not correctly print.

I modified the Label Layer sample to include a Print Task, and found that the label layer is simply ignored - only the basemap and polygon layer are included on the printout. See this JS Fiddle to reproduce.

Is a label layer supported in the print task?


texts: OpenLayers WFS GetFeature large amount of data

texts:
I am using WFS GetFeature request to GeoServer from OpenLayers application. And highlight different color that return features on map. This is working fine for 990 records. But not returning for 1800 record.

I am using AJAX request:

        var url = "http://localhost:8080/geoserver/cits/wms?"
            + "REQUEST=GetFeature&"
            + "SERVICE=WFS&"
            + "version=1.1.0&"
            + "typename=pipes&"
            + "outputFormat=text/javascript&"
            + "format_options=&"
            + "callback: JSON_CALLBACK";

No error and no record. I thing it is timed out. How can I solve this issue?


texts: Convert from latitude, longitude to geoTIFF screen pixel x,y coordinates in Python

texts:
I'm trying to use rasterio or pyproj to convert latitude and longitudes into the SCREEN xy coordinates of a geotiff, so that I can superimpose 2D objects on it in a raster plotting library.
Specifically: I have a geotiff that is 512x512 in epsg:3857 representing a small area of satellite imagery.
I have a known lat long that is in that expanse somewhere. If it were in dead center, the transformation I want would yield 256,256 (50%,50% of 512,512).
How to do it? I tried:
import rasterio as rio
latitude, longitude = 39.60378691188893, -104.9429961240613
dataset = rio.open('satellite_downloads/geotiff/home.tif', crs='epsg:3857')
pprint(dataset.crs)
pprint(dataset.transform)
py, px = dataset.index(x= longitude, y=latitude)
print('Pixel Y, X coords: {}, {}'.format(py, px))

But this returns

CRS({'init': 'epsg:3857'}) Affine(0.14929107090938487, 0.0,
-11682239.10936522,
0.0, -0.1492910709366697, 4808569.838743366)
Pixel Y, X coords:  32209094, 782507220

How do I convert that into a value between 0 - 511 (my raster xy range)?


texts: Rubber sheeting polygon collection into another using PostGIS?

texts:
How do I create a 4 x 4 grid over the orange grid like the pink(ish) one?

The pink grid are geometries that have already been provided and are wonderful.  I'm attempting to create something similar over the other orange grids.

I could do this manually, but it would take a long time.

My first attempt was a fishnet function using Creating regular polygon grid in PostGIS?.  The result is the purple in the bottom right of the picture.  I could fiddle with the function for a while to get the grid to be rectangular instead of square but the orange grid varies slightly in places and is completely not perfect in others.  I do not want an overlapping smaller grid.

So my next idea is to take the 'perfect' geometries and a align it with the bottom left of the orange grid and 'squish' it into the orange grid area.

I don't know if this is A) possible or B) the best way to go about solving my problem.

Any suggestions?


texts:
I am wondering how you can save multiple data for one point feature in ArcMap?

For example, I have a shape-file with point features, the data that needs to be included for each point feature consists of different years (2007, 2008, 2009 and 2010). Each year includes data such as total number of students, enrollments, est. How can I save this data in arcMap for each specific point feature?

I need to be able to save all data for individual point features, from which I can load into Geoserver, and when a point feature is selected in OpenLayers (Geoserver), all the (yearly) data will be displayed.

So far I have tried joining the data in arcmap, but when served up in GeoServer, only the year 2007 is displayed when a point feature is selected, and the rest of the years are not displayed (happens for all point features).


texts: Annotation Published to my ArcGIS Server is being displayed as a point when I build my application with ArcGIS Viewer for Silverlight. How do I get the annotation to display as text and not a point feature?


texts:
I have a requirement, upon user action I need to zoomin to the polygon/feature on layer. To achieve this

Getting the feature object based on unique Id by using
    wfs?request=GetFeature
Getting required feature as expected

{"type":"Feature","id":"SEARCH.fid--2a37c556_1505cda2d3b_1593","geometry":{"type":"Polygon","coordinates":[[[3400.95051588509,4083.02548123558],[3502.22770354563,3981.74829357504],[3598.40160615068,4077.87129986026],[3498.61223022689,4180.68719557739],[3400.95051588509,4083.02548123558]]]},"geometry_name":"GEOM","properties":{"LAYER":"CAGE_INFO","TAG":"7234","NUMBER":"000108"}}

Now, I used the below to options and unsuccessful with the both

map.getView().fit(feature.geometry.coordinates,map.getSize());
map.getView().fit(feature.geometry,map.getSize());

Please suggest


texts: How to accurately place geographic object into photo?

texts:
My task is to create before/after photo visualization of a wind turbine which is planned to build.

A photo is taken from a point of known coordinates, and coordinates of a turbine are also known. How to accurately place this turbine into the photo?

What I already do is to identify some object on the photo and on an aerial photo. Then I measure distance on the aerial photo map to this object and estimate its heiht. Next I measure a height of this object in the photo in pixels and by proportion calculate height of the turbine.
It's very rough calculation because real object heiht estimation is biased.

Next problem is with a direction where to place the turbine. I tried some angular calculations, but they were rather misleading and I ended up with drawing line from point to turbine and further, and looking on what background I should place the turbine.

If you have any ideas how do this smarter and more accurately, please share.


texts:
While updating some ArcObjects projects to ArcGIS 10, I realized that I don't have a comprehensive logging strategy.  Sometimes it's log4net.  Sometimes it's System Trace/Debug.  Sometimes custom logging. Sometimes all three within in the same project.

log4net certainly seems to be popular.  I use it and I have had no trouble with it, but I don't feel like I'm experienced enough with log4net to anticipate possible issues with its use.

This question is for ArcGIS Desktop only: AddIn or COM extensibility.  I think the answers to the same question for ArcGIS Server might be different and deserve its own question and its own set of tags.

What issues exist with using log4net as a logging framework for ArcGIS Desktop?

Some related Stackoverflow.com questions:
https://stackoverflow.com/questions/576456/log4net-versus-tracesource

https://stackoverflow.com/questions/98080/what-is-the-best-logging-solution-for-a-c-net-3-5-project

https://stackoverflow.com/questions/2410438/log4net-performance

https://stackoverflow.com/questions/126540/what-is-your-net-logging-framework-of-choice


texts: Changing style of north arrow “ESRI North 4” to “ESRI North 1” using ArcPy?

texts:
I have a Python assignment where I have to write a script to change the style of the north arrow from “ESRI North 4” to “ESRI North 1”. I wrote the following code:

import arcpy
mxd = arcpy.mapping.MapDocument(r"F:/Lab9/Georgia.mxd")

for elm in arcpy.mapping.ListLayoutElements(mxd, "MAPSURROUND_ELEMENT"):

    if elm.name == "North Arrow":
        elm.style = "ESRI North 1"

mxd.save()

del mxd

It gave me no error, but It didn't change the existing style to the one I required.


texts: Convert XYZ to stars curvilinear grid and plot

texts:
I have data in XYZ format (long, lat, elevation) that does not fit on a rectangular grid, but I would like to plot it as grid cells of different sizes.

To illustrate what I mean, the last two plots in this stars vignette illustrate exactly what I would like to do:
https://r-spatial.github.io/stars/articles/stars4.html

In that vignette, they are able to plot the data as points, but then they are able to plot the data as cells of different sizes, which essentially fills in the gaps.

But I am having a hard time figuring out how to convert the data into a stars object and plot.
An example of one of my files can be downloaded here:
https://www.dropbox.com/s/u230o28109hhxee/elevation_8.0.dat?dl=1

Here is what I have so far:

library(stars)
library(sf)

f <- 'elevation_8.0.dat'

temp <- read.table(f)
colnames(temp) <- c('x','y','alt')
sfObj <- st_as_sf(temp, coords=1:2, crs=st_crs(4326))

plot(sfObj, cex=0.15)

dat <- st_as_stars(sfObj)

dat
> dat
stars object with 1 dimensions and 1 attribute
attribute(s):
      alt
 Min.   :-4895.79
 1st Qu.:-1813.55
 Median : 1191.53
 Mean   :   48.76
 3rd Qu.: 1829.30
 Max.   : 3702.25
dimension(s):
         from    to offset delta                       refsys point
geometry    1 43990     NA    NA +proj=longlat +datum=WGS8...  TRUE
                                                        values
geometry POINT (-128.867 42.1795),...,POINT (-102.176 28.3599)

plot(dat, as_points = TRUE)

Error in .image_scale(values, colors, breaks = breaks, key.pos = key.pos,  :
  object 'mar' not found


texts: Geoserver WMS - CORS enabled in web.xml - still error

texts:
I have tried both geoserver 2.12.5 and 2.13.2 with CORS enabled and on both i get an "Invalid CORS request" when accessing WMS from an external webapp.

steps  :

install geoserver, port 8095, i have tried both Linux and Windows
enable CORS in web.xml
restart
access WMS from an angular app that runs on 9000 -> CORS Error
the same call in a separate tab is getting a json response with  the expected data

Accessing it through geoserver "preview layers" is working,
(also the GetFeatureInfo-call as there it is not a CORS-call)

I am using the sample layer nurc:Img_Sample

Using a different IP than 127.0.0.1 does not change anything.

  http://127.0.0.1:8095/geoserver/nurc/wms?SERVICE=WMS&VERSION=1.1.1&REQUEST=GetFeatureInfo&layers=nurc:Img_Sample&time=2018-08-15/2018-09-03&info_format=application/json&BBOX=-523.3886718750001,18.97902595325528,-406.14257812500006,69.3493386397765&FEATURE_COUNT=5&HEIGHT=890&WIDTH=1334&query_layers=nurc:Img_Sample&SRS=EPSG:4326&buffer=15&X=1127&Y=415

This was working in an older version (from 2017 - so maybe 2.11 ?, or did the browsers checks change ?)

From Chrome :

Request Method: OPTIONS
Status Code: 403 Forbidden
Remote Address: 127.0.0.1:8095
Referrer Policy: no-referrer-when-downgrade

Response Headers :

Allow: GET, HEAD, POST, PUT, DELETE, TRACE, OPTIONS, PATCH
Content-Length: 20
Server: Jetty(9.2.13.v20150730)
X-Frame-Options: SAMEORIGIN

!Provisional headers are shown
Access-Control-Request-Headers: authorization
Access-Control-Request-Method: GET
Origin: http://localhost:9000
User-Agent: Mozilla/5.0 (Windows NT 10.0; Win64; x64) AppleWebKit/537.36 (KHTML, like Gecko) Chrome/68.0.3440.106 Safari/537.36

From the console :

`
zone.js:2224 OPTIONS http://127.0.0.1:8095/geoserver/nurc/wms?SERVICE=WMS&VERSION=1.1.1&REQUEST=GetFeatureInfo&layers=nurc:Img_Sample&time=2018-08-15/2018-09-03&info_format=application/json&BBOX=-523.3886718750001,18.97902595325528,-406.14257812500006,69.3493386397765&FEATURE_COUNT=5&HEIGHT=890&WIDTH=1334&query_layers=nurc:Img_Sample&SRS=EPSG:4326&buffer=15&X=1127&Y=415 403 (Forbidden)

zone.js:2224 OPTIONS http://127.0.0.1:8095/geoserver/nurc/wms?SERVICE=WMS&VERSION=1.1.1&REQUEST=GetFeatureInfo&layers=nurc:Img_Sample&time=2018-08-15/2018-09-03&info_format=application/json&BBOX=-523.3886718750001,18.97902595325528,-406.14257812500006,69.3493386397765&FEATURE_COUNT=5&HEIGHT=890&WIDTH=1334&query_layers=nurc:Img_Sample&SRS=EPSG:4326&buffer=15&X=1127&Y=415 403 (Forbidden)

Failed to load http://127.0.0.1:8095/geoserver/nurc/wms?SERVICE=WMS&VERSION=1.1.1&REQUEST=GetFeatureInfo&layers=nurc:Img_Sample&time=2018-08-15/2018-09-03&info_format=application/json&BBOX=-523.3886718750001,18.97902595325528,-406.14257812500006,69.3493386397765&FEATURE_COUNT=5&HEIGHT=890&WIDTH=1334&query_layers=nurc:Img_Sample&SRS=EPSG:4326&buffer=15&X=1127&Y=415: Response to preflight request doesn't pass access control check: No 'Access-Control-Allow-Origin' header is present on the requested resource. Origin 'http://localhost:9000' is therefore not allowed access.

01:00:35.526 vendor.dll.js:8 ERROR HttpErrorResponse {headers: HttpHeaders, status: 0, statusText: "Unknown Error", url: null, ok: false, …}

`

Any ideas what the issue is here ?


texts: How to merge/combine several tables in geodatabase into one table using ModelBuilder?

texts:
I have over 50 csv files that have the same data structure and I need to merge them into one table.

As I need to build a model I have successfully imported all the files into a File Geodatabase.

How to I merge all these files into one table?


texts: Nearest neighbor between point layer and line layer?

texts:
I've asked this question several times on stackoverflow and irc between #qgis and #postgis and I also tried to code it or implement it my self in postgis with no real answer.

Using programming (most preferably python), I'd like to draw a line from a point layer, to its projection on the nearest line of a line or polygon layer.

As of now most of my data is in ESRI's shape and postgis formats; however, I'd rather stay away from a postgis solution as I'm predominantly a shp + qgis user.

An ideal solution would be to implement GDAL/OGR with python or similar libraries

Using the GDAL/OGR libraries where should I start? would it be possible to give a solution plan?
Can I use NetworkX to do the nearest neighbor analysis?
Is this actually possible?

If it's easier, the points could connect to the segment end point instead of a projected point


texts: From areas of several species ranges (polygons) to a map (raster) expressing the count of overlapping polygons

texts:
I am not able to resolve the problem that follows. I have some separate shapefiles files, each one representing a given species range. Some layers have multiple polygons.

I would prefer to solve this problem with R syntax, but as alternative option I could use QGIS.

After import these shapefiles and merging them, I would like to obtain a map representing the count of overlapping polygons. In other terms, my final output would be a sort of “heatmap” in which each cell has the value of the count of overlapping polygons (i.e. the sum of the species falling in a given cell).

Anyone could kind provide me the operative approach to accomplish this?

If required, I can provide the example layers that I'm using during my trials. Despite many many attempts, I have not reached any result.


texts: distance difference in Google Earth with exported GPX or KML compared to source shapefile

texts:
I've converted a shapefile (points) to KML with QGIS.  Points in the shapefile are within 100m distance from each other.  When I measure it in Google Earth, the distance changes to 80m. The same happens with the created GPX file. Why is this?

My original shapefile was created in WGS84/UTM Zone 20S because I needed to work with the metric system. I'm working with data belonging to an area no larger than 6 square kilometers. As to avoid the problem I am having now, I reprojected it to WGS84 (EPSG 4326) by "saving as" and changing its SRC before I then created the KML and GPX. Both KML and GPX show the same surface, the difference comes along with the shapefile.

Sample coordinates are:

    id  lat long
POINT(-61.82795607536099425 -37.83722374457502724)
POINT(-61.82861595101069696 -37.83772772712714527)
POINT(-61.82724009148368083 -37.83767676392982793)
POINT(-61.8278967620074269 -37.83817748597036257)


texts: Efficiently generate Voronoi layers in qgis via Python

texts:
The following python commands both appear to function when called in qgis 2.6 to generate Voronoi diagrams:

processing.runalg("qgis:voronoipolygons", pointLayer, 1, None)
processing.runalg("grass:v.voronoi",pointLayer,False,False, bbox,-1,0.0001,0, None)

but neither of them seems to work well. The first (qgis) is extremely slow, taking 38 seconds with 180 points vs 2 seconds for the grass method. The qgis method also seems to ignore the output variable when it is specified ("C:/tmp/output.shp" in place of None), but the output layer does appear in the Layers pane so I can iterate through the layers and find it.

The second algorithm seems to ignore the output parameter (again, "C:/tmp/output.shp" in place of None), and no result seems to appear in the layer pane. Trying to use memory:name generated an error message. However, the following results in a temp file path being returned:

output = processing.runalg("grass:v.voronoi",PointLayer,False,False, bbox,-1,0.0001,0, None)
print output.output

My datasets will have thousands of points, so the qgis library is much too slow and the grass library requires a workaround to get the result out.

Is there a way to make qgis:voronoipolygons run faster?
Is there a syntax I am missing to get the result layer from the algorithm? It will be an intermediate result used for future calculations so memory storage is preferred.


texts: How to display different SVG icon in the same layer in QGIS based on attribute value

texts:
Based on this question: How to display different SVG icon in the same layer in QGIS
is it possible to add multiple svg files in qgis to the same layer automatically based on attribute value? I have a world dataset and world flags svg files. Can I visualize the flags somehow automatically? (for example adding the path as attribute to each record?)


texts: A random "Could not load cross-domain" message when loading Esri map

texts:
My application displays an Esri map, by accessing the Esri public sampleserver. I use Esri dojo API for this purpose. The Esri dojo I use is of version 1.6 and is locally installed (as a war) on my server. It's important to say that the rest of my dojo code is an older dojo version (1.3.2).
I use Firefox 3.6 and the strange thing is, I sometimes get a "Could not load cross-domain resources" message (at about 10% of my attempts). This problem happens more often on another pc (Firefox 3.6 as well): at about 80% of my attempts. It is usually (but not necessarily) solved after clearing Firefox cookies and cache.
Now, I know the code is ok, as I can see the map in most of my attempts. However, I can't figure out why it sometimes fails on my browser, and why the other browser almost constantly fails.
What do you think? What could be the cause for it?


texts: Getting started with Web GIS development, especially vector editing?

texts:
I would like to start developing a GIS web application for mapping and editing vector data with attributes. I'm familiar with Desktop GIS (ArcGIS, QGIS) and basic SQL. I've also little experience with GeoServer, XML and GML but I've no real programming experience (Java, etc.). I would like to do something similar like in this youtube video. It can be simple, for the start than shown in this video but there should be a possibility to create a feature (point, line, polygon) and enter various attribute data.

So my question is, how to start with the whole thing, what are the basics to read and which software (open source) is recommended?


texts:
I am trying to create a simple script for me to clip the layers of a directory. The moment I create an output directory it no longer creates the interface of the tool. Additionally, I'm not sure I have the clip parameters well created. A greeting and thanks in advance

from qgis.PyQt.QtCore import QCoreApplication, QVariant
from qgis.core import (QgsField, QgsFeature, QgsFeatureSink,
QgsFeatureRequest, QgsProcessing, QgsProcessingAlgorithm,
QgsProcessingParameterFeatureSource, QgsProcessingParameterFeatureSink,
QgsProcessingParameterNumber)
import processing

class clip_dir (QgsProcessingAlgorithm):
    INPUT = 'INPUT'
    INPUT_CLIP = 'CLIP'
    OUTPUT = 'OUTPUT'

    def __init__(self):
        super().__init__()

    def name(self):
        return "Clip de capas"

    def tr(self, text):
        return QCoreApplication.translate("clipeado", text)

    def displayName(self):
        return self.tr("CLIP Geoinnova")

    def group(self):
        return self.tr("clip_capa")

    def groupId(self):
        return "example"

    def shortHelpString(self):
        return self.tr("This process generate un clip with one output and one input")

    def helpUrl(self):
        return "https://qgis.org/api/classQgsProcessingAlgorithm.html"

    def createInstance(self):
        return type(self)()

    def initAlgorithm(self, config=None):
        self.addParameter(QgsProcessingParameterFeatureSource(
            self.INPUT,
            self.tr("Input layer (Required)"),
            [QgsProcessing.TypeVectorAnyGeometry]))
        self.addParameter(QgsProcessingParameterFeatureSource(
            self.INPUT_CLIP,
            self.tr("Clip layer (Required)"),
            [QgsProcessing.TypeVectorAnyGeometry]))
        self.addParameter(QgsProcessingParameterFolderDestination(
            self.OUTPUT,
            self.tr("folder")))

    #def processAlgorithm(self, parameters, context, feedback):
        #layer_input = self.parameterAsVectorLayer(parameters, self.INPUT, context)
        #layer_clip = self.parameterAsVectorLayer(parameters, self.INPUT_CLIP, context)
        #folder = self.parameterAsString(parameters, self.OUTPUT, context)

        #processing.run("qgis:clip", {'INPUT': parameters['INPUT'], 'OVERLAY': parameters['INPUT_CLIP'], 'OUTPUT': parameters['OUTPUT']})
        #return {'OUTPUT': parameters['OUTPUT']}


texts: Why do my kml polygons move during zoom (using leaflet and meteor with OpenStreetMap)?

texts:
I'm using meteor.js and leaflet.js to build a simple application that displays kml polygons on an openstreetmaps layer. I can get the polygons to display fine, but during a zoom operation the polygons move, which is particularly disorienting for users on a mobile pinch-to-zoom device.

Single points do not have this issue. Can anyone point me in the right direction?

A copy of my code is running at newmaps.meteor.com

Please be patient, the dev meteor server is slow to load. You will see a map with a polygon and a marker, note the difference in behaviour between polygon and marker during a zoom operation.

EDIT: replaced site with a non-minified version.  map.js is as follows:

    // on startup run resizing event

Meteor.startup(function () {
    $(window).resize(function () {
        $('#map').css('height', window.innerHeight - 44);
    });
    $(window).resize(); // trigger resize event
});

var numMaps;
var currentMapLayer;
var map;

var info;
Session.set("mapName", "Map01");

Template.map.onRendered(function()
{
    $(window).resize(); // trigger resize event
    map.invalidateSize();
});

function style(feature) {
    return {
        weight: 2,
        opacity: 1,
        color: 'white',
        dashArray: '3',
        fillOpacity: 0.7,
        fillColor: 'blue' //getColor(feature.properties.density)
    };
}

function highlightFeature(e) {
    var layer = e.target;

    layer.setStyle({
        weight: 5,
        color: '#666',
        dashArray: '',
        fillOpacity: 0.7
    });

    if (!L.Browser.ie && !L.Browser.opera) {
        layer.bringToFront();
    }

    info.update(layer.feature.properties);
}

function resetHighlight(e) {
    currentMapLayer.resetStyle(e.target);
    info.update();
}

function zoomToFeature(e) {
    map.fitBounds(e.target.getBounds());
}

function onEachFeature(feature, layer) {
    layer.on({
        mouseover: highlightFeature,
        mouseout: resetHighlight,
        click: function(e){resetHighlight(e);highlightFeature(e);}
        //click: zoomToFeature
    });
}

Template.map.rendered = function () {
    L.Icon.Default.imagePath = 'packages/bevanhunt_leaflet/images';

    //var map = L.map('map', {
    //    doubleClickZoom: false,
    //}).setView([-27.974467340930417, 153.36855959147215], 14);

    map = L.map('map').setView([-27.974467340930417, 153.36855959147215], 14);

    console.log("RENDER");

    L.tileLayer.provider('OpenStreetMap.Mapnik').addTo(map);

    // control that shows state info on hover
    info = L.control();

    info.onAdd = function (map) {
        this._div = L.DomUtil.create('div', 'info');
        this.update();
        return this._div;
    };

    info.update = function (props) {
        this._div.innerHTML = '
Details
' + (props ?
            '' + props.name + '
' + props.description
                : 'Tap on a property');
    };

    info.addTo(map);

    //polygon.getBounds().getCenter();

    /*
     geojson = L.geoJson(statesData, {
     style: style,
     onEachFeature: onEachFeature
     }).addTo(map);
     */

    Meteor.startup(function () {

        if(currentMapLayer) map.removeLayer(currentMapLayer);

            map.setView([51.505, -0.09], 13);

            L.marker([51.5, -0.09]).addTo(map)
                .bindPopup("Hello world!
I am a popup.").openPopup();

            L.circle([51.508, -0.11], 500, {
                color: 'red',
                fillColor: '#f03',
                fillOpacity: 0.5
            }).addTo(map).bindPopup("I am a circle.");

            L.polygon([
                [51.509, -0.08],
                [51.503, -0.06],
                [51.51, -0.047]
            ]).addTo(map).bindPopup("I am a polygon.");

    });

};


texts: Is there a way for QGIS make just the black value (255,255,255) on an .IMG format file transparency. Right now their transparent option makes the whole image transparent.


texts: How to Count point in polygon + on boundary?

texts:
So I'm looking for suggestions in both ArcGIS 10.2 and QGIS on methods on how to count points within a polygon and points that fall on a polygon boundary without adjusting the data/points.

The problem is, if a point is shared by two polygons (on the boundary of both), I'd like this point to be counted to be in each polygon.

I've found that select by location in Arc works as a good visual check, to examine these areas. However, if using a spatial join typically the point is only counted once, or not at all. And if there is a larger dataset visually inspecting isn't always feasible.

A quick googling around shows that you can use the Near Function, but I don't believe this would count the same point more than once in polygons that share a boundary.

Does anyone have any thoughts to improve this methodology in both ArcGIS and QGIS?


texts:
I have a set of polylines that represent fire hydrant laterals and many/all of them are in two peices one from the main to a valve on that line and one from the valve to the hydrant.  Some have small gaps at where the valve is located which causes them to not unsplit when I run that command.

I know functions like intersect have a search radius, is there something similar to this for the unsplit lines tool?  Is there another way to do this?


texts: Dealing with political boundaries in grid cell analysis?

texts:
I have point data for Africa and I want to aggregate these data to grid cells which I create by a fishnet. Some grid cells cover parts of 2 or 3 countries. Normally, such cells are by different methods completely assigned to a single country.

I would like to ask whether there is a method to create grid cells that is more sensitive to national borders?

Can I split such cells with country borders and assign them to different countries or is it possible to duplicate such cells and assign them partially to one country?

I use ArcGIS 10.2 for Desktop.


texts: How to use postgis to reproject 900913 to a different bounding box

texts: The OSM data is in 900913. I need to get 4326, but so that the longitude's range is NOT -180..180, but rather some other, shifted range, based on a parameter. This way, if I need to draw an object near antimeridian (e.g. Alaska), I could request range of -360..0.  The leaflet will be centered at -180, and there won't be any incomplete shapes.


texts: Create a polygon to match non-white (or non-empty) part of a raster

texts:
I have a irregular raster with empty spaces as seen in the attached image, which is the southwest corner of my raster. I would like to create a polygon or a mask that matches the non-white regions of the raster using QGIS.
When I went to "Vector -> Research Tools -> Extract Layer Extent," QGIS created the polygon seen in pink. Note that the left-most edge of the raster is in the image shown. In other words, somehow the layer extent goes beyond the visible part of the raster.
In any event, I would like to create a mask/polygon that covers only the visible part of the raster. In the image shown, it would mean a polygon that does not include the square empty space.


texts: Merging XY coordinate table to shapefile data using ArcGIS Pro

texts:
I am new to ArcGIS Pro. My situation is as follows:
I have a shapefile that contains some data for some variables for different block groups. What I require is the data for those variables for some selected co-ordinate points.
So, what I did is add XY data table to the shapefile as a layer but table just imports to the map but does not extract the relevant variable data for the respective coordinates.
Could you tell me how I can easily do it?


texts:
I'm using GrassGis 7.0.0 for Landsat imagery. In the description of i.Landsat.toar under dos3 and dos4 is:

TAUz = exp[-t/sin(e)]

What is "t"? I cannot find any description about what it is. Can anyone help?


texts: Converting Mapinfo's TAB without knowing its projection?

texts:
I've got a multiple sets of Mapinfo files (*.dat, *.id, *.map, *.tab). But I don't know its projections. When I'm trying ogrinfo -al file.tab I'm getting:

INFO: Open of `file.tab'
      using driver `MapInfo File' successful.

Geometry: Unknown (any)
Extent: (-43996.600000, 28383.340000) - (57101.880000, 167741.110000)
Layer SRS WKT:
LOCAL_CS["Nonearth",
    UNIT["Meter",1.0]]

How can I reproject it to WGS84 or something?


texts: I am struggling on how to calculate tortuosity in QGIS. My data includes numerous lines (the mapping of animal trails) I have tried calculating tortuosity with the field calculator but had no luck is there a plugin that I need to use or have I just missed something?


texts: Adding layer to WMS using GeoServer and OpenLayers?

texts:
I need to add some own layers from geoserver on top of this map with writing functionalities.
The tricky thing is where to embed the code in the present script below.

How should the script look like and where should I embed it in the script below?

The map is projected in EPSG:4326 and my layer in geoserver is in EPSG:3021.
How do I get the layer right?

<!DOCTYPE html>

    <title>Draw Features</title>
    <link rel="stylesheet" href="https://openlayers.org/en/v4.0.1/css/ol.css" type="text/css">

<form class="form-inline">
  Geometry type &nbsp;
  <select id="type">
    <option value="Point">Point</option>
    <option value="LineString">LineString</option>
    <option value="Polygon">Polygon</option>
    <option value="Circle">Circle</option>
    <option value="None">None</option>
  </select>

  DrawEnd - Lon,Lat

</form>


texts: Is it common to save geographic coordinates in minutes instead of degrees?

texts:
Are there any communities or fields where it is common to save geographic coordinates in minutes instead of degrees? (I just received a batch of addresses referenced that way and it made me wonder.)

Are there any advantages in saving lat/lon in 'decimal minutes'?

For example:

2894.511,982.8775

instead of

48.2418,16.3812


texts: IS it possible to use geoserver WFS-T services with ArcGIS android SDK?

texts:
I have to use ARCGIS SDK in my project and I should edit some layers. In fact there is a sample which I can learn how ArcGIS API for android works.In it we the editable layers are from a service like this
here

But I dont want to use ARCGIS for server,I want to use Geoserver instead,Does ArcGIS API for android supports WFS-T or things like that?

this is the link of editing sample using ARCGIS API


texts:
I work with a mechanistic model that uses climate data to simulate some specific output for different time steps and spatial locations. Now I want to investigate the effect of different climate data sources on the model output. To do this, I varied the climate data and ran the model with otherwise constant settings; the groups are:

Climate model that produces the climate data
Emission scenarios with different assumption underlying the climate
data from the same climate model
Spatial resolution of the climate data
Downscaling method to arrive at the specific spatial resolution

I have produced now a bunch of rasters with my model driven by the varying climate data.

Now I would like to extract the effects of the different groups of climate data on the model output and produce raster maps showing at which spatial locations there is the most variation in model output.

I was thinking about applying a principal component analysis to avoid pairwise comparisons of the simulations.

Would a PCA be an appropiate method or is there a better way of comparing the different groups? And what would be the next step after doing the PCA and calculating the loadings?


texts: GeoServer - "tiled=true" doesn't appear to be working for wms layers

texts:
I have a layer which I have seeded in GeoServer using GWC in EPSG:3857. I have made sure that "Enable direct integration with GeoServer WMS" has been selected, however when I run the Openlayers code below, my CPU goes to close to 100% and the GeoServer logs says its tiling. I have the "tiled=true" option in the code, so I am trying to understand why the seeded tiles are not being used.

<!DOCTYPE html>

    <title>Tiled WMS</title>
    <link rel="stylesheet" href="https://openlayers.org/en/v4.6.5/css/ol.css" type="text/css">

I have also tried using the GeoWebCache endpoint URL option, http://localhost:8081/geoserver/gwc/service/wmts?, which works and the seed tiles are used, however there appears to be a loss in quality, see below.

At this point I am happy to use WMTS if the WMS option doesn't work but does anyone know why the quality of the tiles are lower for WMTS than WMS?

UPDATE:

For those having the same quality degradation issue I was having, it was because my resolutions were off. I used the "REQUEST=GetCapabilities" request to solve my problem. Here is the full code:

<!DOCTYPE html>

    <title>WMTS Layer from Capabilities</title>
    <link rel="stylesheet" href="https://openlayers.org/en/v4.6.5/css/ol.css" type="text/css">


texts: Enterprise Architect - domains relationship

texts:
I am building an ArcGIS model in Enterprise Architect.  I have a couple domains inside.

The problem is: I want to make the following relationship between two of it (domains) that after I import the XML file to ArcGIS I will be able to choose one attribute (from the first domain) and the second one (from the connect-relationship domain) will be filled automatically .

Do you know how to make such relationship (connection)?


texts: Using SVG from Wikipedia to fill polygon?

texts:
I'm trying to use the flag of each region in Italy to fill each region's polygon. I have downloaded the shapefile of Administrative areas for Italy from diva-gis and am now attempting to fill each of the 20 polygons with the corresponding flag image taken from Wikipedia at http://en.wikipedia.org/wiki/Regions_of_Italy

When you click on the flag image on the wikipedia page, and try to download a single flag, such as Sardinia:

http://en.wikipedia.org/wiki/Regions_of_Italy#mediaviewer/File:Flag_of_the_Italian_region_Sardinia.svg

you see that the images are already in svg format. However, when I attempt to use them in QGIS, all I see is blank space.

I downloaded Inkscape as I read somewhere this would be the appropriate tool for dealing with SVG, however opening the file in Inkscape and saving as SVG does not result in anything different.

There are 20 regions in Italy and therefore 20 flags, can someone lead me in the right direction?

My plan of attack so far is to change the Style to Categorized and classify the regions of Italy by the NAME_1 Column (which is the region name) and then change each symbol layer type to SVG fill. After that I was thinking I would download each SVG flag individually from wiki and chose it for the SVG fill by clicking on the ellipses (...) button and navigating to the appropriate file.


texts: Clip the result of my QGIS 3.x Processing python script before returning the output?

texts:
I have written a simple QgsProcessingAlgorithm that populates a QgsFeatureSink.
Before outputting the result, I need to call another algorithm to clip my features.
Everything seems to work just fine, the clip algorithm is called and its output shows in the log, except my actual output is not clipped.
It seems the sink is returned no matter what I specify to be actually returned by the function.

What I am doing wrong ?

from qgis.PyQt.QtCore import QCoreApplication, QVariant
from qgis.processing import alg
from qgis.core import (
    QgsField,
    QgsFields,
    QgsFeature,
    QgsFeatureSink,
    QgsGeometry,
    QgsPointXY,
    QgsProcessing,
    QgsProcessingAlgorithm,
    QgsWkbTypes
)
import processing
import random

@alg(name='test_alg', label=alg.tr('Test'), group='test', group_label=alg.tr('Test'))
@alg.input(type=alg.SOURCE, name='INPUT', label='Input')
@alg.input(type=int, name='NB_POINTS', label='Number of points to create', default=1000, minValue=10, maxValue=100000)
@alg.input(type=alg.SINK, name='OUTPUT', label='Output')
def testAlg(instance, parameters, context, feedback, inputs):
    """
    Test script to randomly generate n points unniformly distributed within a polygon extent
    and clip the result with regards to the polygon itself.
    """
    # Define some attributes
    fields = QgsFields()
    fields.append(QgsField('Point_No', QVariant.Double))
    # Initialise variables
    source = instance.parameterAsSource(parameters, 'INPUT', context)
    nb_points = instance.parameterAsInt(parameters, 'NB_POINTS', context)
    # Setup the sink
    (sink, dest_id) = instance.parameterAsSink(parameters, 'OUTPUT', context, fields, QgsWkbTypes.Point, source.sourceCrs())
    # Define the Points generator
    def genPoints(n, extent, attributes):
        """
        Generator to produce the Points yielded as Features with attribute
        """
        feat = QgsFeature()
        xmin, xmax, ymin, ymax = extent.xMinimum(), extent.xMaximum(), extent.yMinimum(), extent.yMaximum()
        total = 100 / n
        for i in range(n):
            feedback.setProgress(i * total)
            attributes[0] = i
            x, y = random.uniform(xmin, xmax), random.uniform(ymin, ymax)
            feat.setGeometry(QgsGeometry.fromPointXY(QgsPointXY(x, y)))
            feat.setAttributes(attributes)
            yield feat

    feedback.pushInfo(f'Generating {nb_points} Points')
    sink.addFeatures(genPoints(nb_points, source.sourceExtent(), [0]),
                      QgsFeatureSink.FastInsert)
    feedback.pushInfo('Clipping to Input Polygon...')
    clipped_points = processing.run('native:clip', {
            'INPUT': dest_id,
            'OVERLAY': parameters['INPUT'],
            'OUTPUT': 'memory:ClippedPoints',
        }, context=context, feedback=feedback)['OUTPUT']

    return {'OUTPUT': clipped_points.id()}


texts: Implementing access control for WFS-T updates?

texts:
I'm implementing a system which has structure like in my question: What's role of WFS-T in webapp?

There's spatial data which is updated through WFS-T and non-spatial data which is updated through traditional (Java) web application.

When user updates non-spatial data, it is easy to implement authorization rules such as

if current_user = sales_region.manager AND sales_region.status = ACTIVE then all spatial fields can be updated

But can I do the same for the WFS-T geometry updates?

if current_user = sales_region.manager AND sales_region.status = ACTIVE then sales_region.area geometry can be updated through WFS-T requests

Isn't the WFS-T interface just an open XML interface? Anyone who has http access to it can modify, add and remove geometries at will? How could I add business rules checks to protect the geometry updates?

I know I could protect the whole WFS-T interface so that only some users can access it altogether, but that's not what I am after. I need to allow or deny some updates based on the data and some defined business rules.

Are there any easy solutions for this?


texts: How can I add dynamic text/properties to a GeoJSON multipoint on each individual feature?

texts:
I'm drawing circles on a map based on lon/lat. These circles are drawn from "subarrays" in an array, e.g:

myArray = [[array(18)]],[array(111)], [array(65)]....]

I map each of these arrays using GeoJSON like this:

var geojsonlatObject = {
        'type': 'FeatureCollection',
            'crs': {
                'type': 'name',
                'properties': {
                'name': 'EPSG:4326'
            }
        },
        'features':
            [{
                'type': 'Feature',
                'geometry': {
                    'type': 'MultiPoint',
                    'coordinates': circlecoordinates().map(function (x) {
                        //console.log(x[0][0]);
                        //console.log(x[0][1]);
                        return ol.proj.transform([x[0][0], x[0][1]], 'EPSG:4326', 'EPSG:3857')
                    })
                }
            }]
    };

Further, I have a style like this:

    var locationandtransitstyle = new ol.style.Style({

        text: new ol.style.Text({
            text: String(resultArray.length),//this.length,//'15',
            scale: 1.7,
            fill: new ol.style.Fill({
                color: '#ffffff'
            })
        }),
        image: new ol.style.Circle({
        radius: 15,
        fill: new ol.style.Fill({
            color: 'rgba(30,144,255,0.7)'
        })
      })
    })

The result of this is that I have the length of myArray in each circle. But I want the length of the subarrays as text inside the circles. Any ideas on how to do this?

My thoughts at the moment is that I can add use the "properties" option in GeoJSON, https://tools.ietf.org/html/rfc7946, but I don't know how to implement it for each feature...

    ......."coordinates": [125.6, 10.1]
  },
  "properties": {
    "name": "***for-circle add array lenght***"
  }

UPDATE

After some more experimenting i came to the following:
I can log length of the subarrays by:

console.log(x.length);

inside my map function. Can i use that to set this value to a properti??:

'type': 'Feature',
                'geometry': {
                    'type': 'MultiPoint',
                    'coordinates': circlecoordinates().map(function (x) {
                        console.log(x[0][0]);
                        console.log(x[0][1]);
                        console.log(x.length);
                        return ol.proj.transform([x[0][0], x[0][1]], 'EPSG:4326', 'EPSG:3857')
                    })
                },
                "properties": {
                    "name": "...value of x.length..."
                }
            }]

Is it possible to do this, or will it not work? Do i need a .map()-function in for the property as well?

Alternatively:

properties": {
                    "name": circlecoordinates().map(function (x) {
                         return    x.length
                            }
                }


texts:
Is it possible to create  a new id serial column inside a view? I need it to be able to add my view in postgis

I have made a full outer join on two tables
I need a unique PK to load the view in QGIS.  Creating a id serial does not help in this case because it will become redundant once the two tables are joined.

This is my view:

CREATE OR REPLACE VIEW myview AS
SELECT
 table1.id
 table1.name
 table2.id
 table2.street
FROM table1
 FULL OUTER JOIN table2
 ON table1.id  = table2.id

I have seen some ROW_NUMBER examples but i dont know how to script it when I am using a join.In PostGIS is it possible to create a View with a unique ID?

a final try before I find another direction to solve my problem. It is a try to use the row_number to create a unique column inside my VIEW. It looks daunting and it does not work.. probably not the best way of solving my problem (Iam not very experienced writing script yet.. but i give it a last try in case it is just a quick fix)

create view testView as (
    select
    row_number() OVER(ORDER BY (SELECT
 table1.id
 table1.name
 table2.id
 table2.street
FROM table1
 FULL OUTER JOIN table2
 ON table1.id  = table2.id) AS OID,
   table1.id
 table1.name
 table2.id
 table2.street
FROM table1
 FULL OUTER JOIN table2
 ON table1.id  = table2.id


texts:
Where can I locate data that can be used to calculate driving distances in the UK using the Network Analyst extension to ArcGIS Desktop 10?

The data that came with my software (I'm in the USA) only covers North America.

I see that ESRI sells StreetMap Premium which covers europe: http://www.esri.com/data/streetmap

How much does this data cost?  I am assuming it is expensive since they don't even list a price and I would prefer some cheaper options.

It looks like you might be able to buy data from Navteq/teleatlas (referenced here http://tereshenkov.wordpress.com/2011/03/08/street-data-for-network-analyst/ ) but that they will require a lot of manual processing.

Is there anything like TIGER for the UK that is provided by the government? The Ordnance Survey website didn't appear to have any streetmap style data.


texts: Refreshing Combobox of Python Add-in of ArcPy?

texts:
I need to refresh the value in the combobox after an event from other tool in the same add-in (ARCGIS 10.1. Python).

I've changed the value of the combobox with the code:

combobox.value = val

and I want to show it in the combobox but I can't refresh this combobox after give a value in other tool.


texts: Disable the "Open Attribute Table" in ArcGIS for some raster layers

texts:
I have made a layer using cell statistics (Spatial Analysis tool). I need to process the attributes before exporting the layer, however the open attribute table is disabled.

How to enable the table or use data?

The layer was converted to other formats including point. Its attribute table is opened while the field is grid code and the original needed field before cell statistics process is not available.


texts: Python loop Polygon to raster not working?

texts:
I want to convert a large number of different shp polygons to raster in ArcGIS 10.1, using an existing field column "class"  for the new raster classes. I am not very fluent in python but I thought this would work, and it doesn´t:

 # Import arcpy module
 import arcpy
... # Local variables:
... project_input = "D:\\RESEARCH_PROJECTS\\NORD-STAR\\AMAZONIA\\1-RAW_DATA\\INPE\\PRODES DIGITAL\\2004\\"
... raster_output = "D:\\RESEARCH_PROJECTS\\NORD-STAR\\AMAZONIA\\1-RAW_DATA\\INPE\\PRODES DIGITAL\\2004\\convert\\"
... # List features
... arcpy.env.workspace = project_input
... fcList = arcpy.ListFeatureClasses()
... # Loop
... for featureClass in fcList:
... # Output
...     raster_outputpath = raster_output + featureClass
...
... # Process: Polygon to Raster
...     arcpy.PolygonToRaster_conversion(featureClass, "CLASS", raster_output, "CELL_CENTER", "NONE", "30")
...     print "finished polygon to raster"
...
Runtime error  Traceback (most recent call last):   File "", line 20, in    File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\conversion.py", line 2436, in PolygonToRaster     raise e ExecuteError: ERROR 999999: Error executing function. Failed to execute (PolygonToRaster).

The script actually calls the first polygon, but that is it:

  Executing: PolygonToRaster
  "D:\RESEARCH_PROJECTS\NORD-STAR\AMAZONIA\1-RAW_DATA\INPE\PRODES
  DIGITAL\2004\PDigital2004_00157_pol.shp" CLASS
  "D:\RESEARCH_PROJECTS\NORD-STAR\AMAZONIA\1-RAW_DATA\INPE\PRODES
  DIGITAL\2004\convert\" CELL_CENTER NONE 30 Start Time: Wed Nov 21
  15:17:11 2012 ERROR 999999: Error executing function. Failed to
  execute (PolygonToRaster). Failed at Wed Nov 21 15:17:11 2012 (Elapsed
  Time: 0,00 seconds)


texts: Unable to edit data in an ArcSDE geodatabase

texts:
I'm dealing for the first time with an ArcSDE geodatabase so I'm facing simple problems to you. ESRI help wasn't of much help (to a begginer at least).

I can't edit my data in ArcMap, most probably because I've already ""registered as versioned" my featuredataset. I did use "Unregister as versioned" from ArcCatalog, but still can't edit my data. And then I got some "schema lock error"" :) So any clues?

Ps: I'm the only user of the geodatabase.

Thank you very much,
Ligia


texts:
I have spent all morning and most of the afternoon today trying to figure out how to work with the geodatabase tool. I need to record 3D information on a 2D plane, and I think this is the way to go about it, however I'm having difficulty.

I am doing I.T. work on network closets at the institution I work at. Each closet has several racks. Each rack has several panels. Each panel has multiple fields of information I need to record; Cable Range, Color, Cable Type, Cable Destination, etc. I could create a shapefile for each panel and then tie that to an excel doc with all of the included data, but that would be tedious and take up a lot of memory space.

Ideally I would have a single point for each rack, which is then tied to data that contains information about the contents of each panel on that rack.

What's the best way to record all of this data?


texts: In ArcMap, can polygon labels always be displayed in a fixed area of the screen?

texts:
If I'm inside an area, I want its properties to be displayed in a fixed area of the screen (top right, let's say).

I am talking about data view, and what I mean is that if I'm zooming into a city, for instance, the country, state, zip code, etc. would show in the fixed area even though the boundaries were not actually visible.


texts: Joining random points to layer using QGIS?

texts: I have one vector layer and one random points layer (generated from random points in extent). How can I join / merge to keep only the random points that align with the map layer?


texts: Is it possible to add a mbtiles file as TileMill layer?


texts:
The projection datum for user defined (Proj4: +proj=longlat +a=6377301.243 +b=6356100.2284 +towgs84=295.00,736.00,257.00,0,0,0,0 +no_defs) works well in Qgis 2.18 but it does not work in Qgis 3.4 version.

in response to :Please provide with some details that prove your case. Some screen shots. – Vadym
Here are some screenshots for reference.
in fig 1, I open a layer in QGIS 3.4.4 which have WGS84/UTM44N projection.

in fig.2, I defined/set custom projection from setting menu named Modified UTM for Central Nepal.

in fig 3, I calculate the value of projection from WGS84 to Custom projection, but the result destination CRS shows the same i.e. WGS84 projection.

in fig 4, I go through project properties menu for CRS, set the project CRS in Modified UTM for Central Nepal (custom projection).

in fig 5, the project projection changed in WGS 84 system instead of desired/selected CRS i.e. Modified UTM for central Nepal.

It's so frustrating. Pls help for this problem and suggest me.

It works well in QGIS 2.18 Version but doesnot in QGIS 3.4.4. I use the MAC, Apple labtop.


texts: Show large information in a map (millions of points)

texts: I have millions of points (census data) about 1.5Gb in a text file. But I only want to show an area based on this points like a heat map, so maybe I have to build polygons based on these points. What is the correct way for do this? How can I start with this?


texts: extract data of polygon shape which intersect with point shape R and create data frame

texts:
I want to extract data of a SpatialPolygonsDataFrame li which intersect with the SpatialPointsDataFrame p R and then create a data frame with the coordinates and values of both objects. I tried

library(sf)

out <- st_intersection(p, li)

I get the error Error in UseMethod("st_intersection") :
  no applicable method for 'st_intersection' applied to an object of class "c('SpatialPointsDataFrame', 'SpatialPoints', 'Spatial', 'SpatialVector')"

Not all points intersect with the polygon shape.

How could I extract the data of a polygon shape based on a point shape and then create a data frame of the result?


texts: Equivalence of protocol.http of openlayers 2 to openlayers 3?

texts:
I am new to openlayers. I have following code in openlayers 2.

for (var i=0;i<wfslayers.length;i++) {
    var object=wfslayers[i];
    var url="/"+object.resourcelocation+"/wfs?token="+object.token+"&layer="+object.name;
    var legendUrl="/"+object.resourcelocation+"/wms?token="+object.token+"&layer="+object.name;

    var layer = new OpenLayers.Layer.Vector(object.name, {
        layerid: object.label_text,
        projection: new OpenLayers.Projection("EPSG:"+object.srid),
        strategies: [new OpenLayers.Strategy.BBOX()],
        protocol: new OpenLayers.Protocol.HTTP({
            url: url,
            format: new OpenLayers.Format.GeoJSON() ,
            geometryName: object.geometry_name
        }),
        legendUrl:legendUrl,
        isBaseLayer: false,
        "visibility": false,
        "chkboxRef":object.chkboxRef
        },{
        "chkboxRef":object.chkboxRef
        });
        wfscatalog.push(layer);
    }

I tried to convert it to openlayers 3. Following is my code. But I need help in the conversion of protocol part and strategies part. Please help.

for (var i=0;i<wfslayers.length;i++) {
    var object=wfslayers[i];
    var url="/"+object.resourcelocation+"/wfs?token="+object.token+"&layer="+object.name;
    var legendUrl="/"+object.resourcelocation+"/wms?token="+object.token+"&layer="+object.name;

    var layer = new ol.layer.Vector({
        projection: "EPSG:"+object.srid,
        isBaseLayer: false,
        visible:false
      });
      layer.layerid=object.label_text;
      layer.legendUrl=legendUrl;
      layer.chkboxRef=object.chkboxRef;

     wfscatalog.push(layer);
    }


texts: Center (and add points too) a map with LonLat

texts:
I'm new to OpenLayers, so I'm sure this will be an easy win for someone :)

I'm trying to center a map on a point (and place markers) using lon/lat but it appears to take everything as meters?

    map = new OpenLayers.Map("map",{
    projection: new OpenLayers.Projection("EPSG:3857"),
    units: "dd",
    numZoomLevels: 18,
    maxResolution: 156543.0339,
    layers: [
        new OpenLayers.Layer.OSM("OpenStreetMap", null, {
            transitionEffect: 'resize'
        })
    ],
    center: new OpenLayers.LonLat(-2.626944, 57.563222),
    zoom: 7
});
var markers = new OpenLayers.Layer.Markers( "Markers" );
map.addLayer(markers);

N.B. this code isn't designed to add markers, I will need to but the problem is clearly with this chunk as I'm being centered (almost) at the 0,0
I've tryed the following code

var size = new OpenLayers.Size(21,25);
var offset = new OpenLayers.Pixel(-(size.w/2), -size.h);
var icon = new OpenLayers.Icon('http://www.openlayers.org/dev/img/marker.png', size, offset);
var i=0;
function Go() {
    x = document.getElementById("x");
    i=i+10;
    x.innerHTML = i;
    markers.addMarker(new OpenLayers.Marker(new OpenLayers.LonLat(i,10),icon));
}

Calling go on a mouse click, and indeed my marker moves a small amount (I'm guessing 10 meters) to the right every click.


texts: How to solve cannot alloc point_buffer error LAStools?

texts:
I have just installed LAStools in my QGIS 3.x version. I have already managed to visualize the point cloud through the lasview tool. Now, I am having problems trying to classify Ground points from the point cloud to extract DEM information.

Can you please help me fix this error (cannot alloc point_buffer) as shown is the image below:

Can you please have a look at the parameters I have as input]1


texts: Set raster properties in arcpy: Multiple no data values for single band raster

texts:
I am using SetRaster properties with arcpy to quickly set the no data value of a raster.

import arcpy
import os

arcpy.SetRasterProperties_management("InputRaster", nodata = "1 0")

I also want to assign value 15 to nodata in my example, so that both values 0 and 15 are assigned no data. In the ArcGIS help for the tool (http://desktop.arcgis.com/en/arcmap/10.4/tools/data-management-toolbox/set-raster-properties.htm) I found the following

  If you want to define multiple NoData values for each band selection,
  use a space delimiter between each NoData value within the
  bands_for_nodata_value parameter.

I haven't found nodata = [[1, 0], [1, 15]] to be working. The reason why I would like to use SetRasterProperties is because it is way faster than for example SetNull. I am using ArcGIS 10.4.1.


texts:
I've been trying to build a graph from a projected (EPSG 27700) clean (no pseudo-nodes) street network. I followed the tutorial here presented employing the QGIS 3.0 functions; now the network analysis tutorial is included in the analysis module of PyQgis.

This is the code:

layer = qgis.utils.iface.mapCanvas().currentLayer()

# undirected graph and distance strategy
graphDirector = QgsVectorLayerDirector(layer, -1, '', '', '', 3)
strategy = QgsNetworkDistanceStrategy()
graphDirector.addStrategy(strategy)

startPoint = QgsPointXY(531389.8,181641.4)
builder = QgsGraphBuilder(lay.sourceCrs(), topologyTolerance = 0.0)
tiedPoints = graphDirector.makeGraph(builder, [startPoint])

graph = builder.graph()
graph.vertexCount()
graph.edgeCount()

The graph is apparently composed by 15358 nodes, even though I was expecting 3799. The edgeCount results in 0 edges instead. I also tried to computed the shortest path and no path was existing between different nodes (as expected, due to the abscence of edges).

Essentially, my question is: what I am doing wrong? Moreover, is the part where I tie the points to the network crucial for creating the graph? I don't understand why that would be relevant for obtaining the graph representation.


texts: Is it possible to georeference raster stored in Postgis?

texts: I'm looking to store a huge raster aerial image (400 Gb) into Postgis, but before storing it i like to know how i could georeference or rectify raster stored in postgis ? is there any Qgis plugins like georeferencer, but could read in-database-raster or any other software ? any hints ?


texts:
What are open source tools for building a web gis application having PAN, Zoom, Query/search, Label etc options?

I have decided to use PostgreSQL for database and GeoServer as middleware but I am a bit confused in selecting the software to design the GUI for the client side.

I am a novice programmer, so will geoexplorer be a good option or are there some other tools which can be more easily customized for my work?


texts: GeoTIFF to KMZ using gdal_translate produce imagery with bad colors

texts:
I have a TIFF file which is georeferenced in UTM coordinates.

So as a first step, I convert the UTM TIFF into a WGS84 TIFF using gdalwarp.
On the next step, I used gdal_translate to convert a TIFF into a KMZ file but I have a problem with the "color" of the image. Here are the commands I used :

gdalwarp test_utm.tif -t_srs EPSG:4326 test_wgs84.tif
gdal_translate -a_nodata 0 -of KMLSUPEROVERLAY -co FORMAT=JPEG test_wgs84.tif test.kml

In fact, the image should be in greyscale (there is only 1 band on the TIFF) but here, the result is fully white with only the shadows which are gray/dark. I think there is a problem with the values of the KMZ ?

Here is a screenshot of the KMZ file :

Here is the TIFF description of the WGS84 converted TIFF :

gdalinfo test_wgs84.tif
Driver: GTiff/GeoTIFF
Files: test_wgs84.tif
Size is 8618, 6296
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 = (-111.186216437919299,38.213784466217568)
Pixel Size = (0.000005191716045,-0.000005191716045)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (-111.1862164,  38.2137845) (111d11'10.38"W, 38d12'49.62"N)
Lower Left  (-111.1862164,  38.1810974) (111d11'10.38"W, 38d10'51.95"N)
Upper Right (-111.1414742,  38.2137845) (111d 8'29.31"W, 38d12'49.62"N)
Lower Right (-111.1414742,  38.1810974) (111d 8'29.31"W, 38d10'51.95"N)
Center      (-111.1638453,  38.1974409) (111d 9'49.84"W, 38d11'50.79"N)
Band 1 Block=8618x1 Type=UInt16, ColorInterp=Gray
  NoData Value=0
  Metadata:
    LAYER_TYPE=athematic
    STATISTICS_MAXIMUM=4744
    STATISTICS_MEAN=836.61375991249
    STATISTICS_MINIMUM=12
    STATISTICS_STDDEV=407.61174893669

Do you have any idea about why I have this problem ? I already convert other images with only 1 band also, and the output KMZ file was great (but I didn't process any conversion with gdal_warp, the tif was already in WGS84). I'm open in any free solution to convert an UTM TIFF into a KMZ file.


texts: ArcIMS equivalent in QGIS for adding images from external client

texts:
I am looking to add some images to the QGIS map canvas from an external server. My client has provided their server URL to access their data (with password credentials). The instructions to access the URL using ArcGIS (via ArcIMS) is also given.

But I want to access the URL and view the image layers from their server in QGIS. Is there any option available to access the URL in QGIS and to add the images in QGIS map canvas?

Or is there any way to develop stand-alone python plugin to achieve this? I had a glance at OpenLayers Plugin python code where they access the Google Physical layer using

OlGoogleMapsLayer.__init__(self, name="Google Physical", html="google_physical.html")

Is there a similar way to access ArcIMS layer using python?


texts: 'Add Field' option is not showing up in the attribute table 'Option' in ArcGIS Online map Viewer

texts: I am trying to add a new field in the attribute table in ArcGIS Online map viewer but for some reason, the option to add field is not showing up. I am using a Cache service and a dynamic service. I can't figure out what could be the problem here.


texts: Importing special/custom csv files to QGIS


texts:
I am working on a project in which I want to take point features( tweets) and compare them to the GDP per capita of countries. I would like to have my end result be such as the following:

GDP Class 1 (55992.21-105044) + high percentage of total tweets tweets =colour 1
GDP Class 1 (55992.21-105044) + low percentage of total tweets =colour 2
...

So on and so forth for each of my 5 GDP classes. I have already tried to work with multiple attributes by that gives me all combinations of single values for each, and I need both in classes to provide my analysis.

Example: USA has A very high GDP per capita, and has tweeted 73% of the total data collected. Canada is in the same GDP/Capita classification, but has contribute only 1% to the total data. Each one shold be coloured in a way that represents this.

EDIT:

More information as requested:

Using : ArcGis 10.2
End Result:a thematic map relating the GDp/Capita of a country to the percentage of tweets it posted during a time frame
What Have I tried?: Using Multiple attributes to relate the 2 attributes together
What have I succeeded in doing: One table contains both the GDP/capita and the number of tweets
How I want it to work: 0-15%,16-55% and 56-75% are my groupings for the percentage of tweets, and GDP/Capita is split into 5 classes. Each GDP/Capita Classificiation should contain sub-classificiations based on the Percents of all tweets from that country.


texts: Find closest point using L.GeometryUtil.closest and open its popup

texts:
How do I use:

L.GeometryUtil.closest(myMap, coords, b)

to find the closest point in a geojson point layer in leaflet. I have tried several things. And the error I receive is:
Uncaught TypeError: Cannot read property '0' of undefined

    let coords=[];

    let closetPt;

    lyrPts = L.geoJSON.ajax('Geojson/pts.geojson',{
        pointToLayer: function (feature, latlng) {
        return L.circleMarker(latlng, markerPts);
      },
        onEachFeature: function (feature) {
                 coords.push(feature.geometry.coordinates);
        }
    }).bindTooltip(function(layer){
      return layer.feature.properties.Layer
    }).addTo(myMap);

    closestPt = L.GeometryUtil.closestLayer(myMap, coords, [-40.5512, 80.9454]);

This part now works when I take out closestPt and put it in a callback.

Is it possible to now automatically open the popup of this closest marker?

I want to use it so that as you are within a distance X of a point, its popup will be displayed.


texts: Sublayers do not exist or are not supported in the geoprocessing tool?

texts:
I'm trying to batch export selected layers that I drag and drop in the tool that I've created by a selecting feature class. It's working when features are input from the source and layers in TOC. However, you get an error when the layers are in group layers. "Parameters are not valid.
ERROR 000732: Input Features: Dataset 'Topography layers\Linear' does not exist or is not supported"

Here's the code:

import arcpy, os, sys

try:
        arcpy.env.overwriteOutput = True

        def selectExport (Layer, LyrSelected, outName):

                for fc in inFeatures:
                        Layer = arcpy.MakeFeatureLayer_management(fc, "outlayer")
                        LyrSelected = arcpy.SelectLayerByLocation_management(Layer,"INTERSECT",SelectLyr,"","NEW_SELECTION")
                        result = int(arcpy.GetCount_management(LyrSelected).getOutput(0))
                        if result > 0:
                                if fc.endswith(".shp"):
                                        bName= fc.split(".")[0]
                                        outName = os.path.join(outPath, os.path.basename(bName))
                                        arcpy.CopyFeatures_management(LyrSelected, outName)
                                else:
                                        outName = os.path.join(outPath, os.path.basename(fc))
                                        arcpy.CopyFeatures_management(LyrSelected, outName)
                                arcpy.AddMessage("*******Geoprocessing done successfully!*******")
                        else:
                                arcpy.AddMessage("*******No features were selected*******")

except:
        tb = sys.exc_info()[2]
        tbinfo = traceback.format_tb(tb)[0]
        pymsg = "PYTHON ERRORS:\nTraceback Info:\n" + tbinfo + "\nError Info:\n     " +        str(sys.exc_type) + ": " + str(sys.exc_value) + "\n"
        msgs = "ARCPY ERRORS:\n" + arcpy.GetMessages(2) + "\n"

        arcpy.AddError(msgs)
        arcpy.AddError(pymsg)

        print msgs
        print pymsg

        arcpy.AddMessage(arcpy.GetMessages(1))

        print arcpy.GetMessages(1)

if __name__ == '__main__':

    outPath = ""
        arcpy.env.workspace = outPath
        inFeatures = arcpy.GetParameterAsText(0).split(";")
        SelectLyr = arcpy.GetParameterAsText(1)
        outPath = arcpy.GetParameterAsText(2)

        selectExport (inFeatures, SelectLyr, outPath)


texts: My selection tool in ArcMap doesn't work

texts: I have been trying to choose some features (rivers) by the selection tool in ArcMap but it seems that it doesn't work. I only have one layer and want to select some features on it. I don't understand why the selection tool can not select the rivers I am putting in its square. I drag the tool over the features and unfortunately, it doesn't highlight them. How can I fix this?


texts: ESRI Local Government Park Locator App Template Deployment Issues

texts:
I am currently trying to implement esri's newest release of the Park Locator | ArcGIS for Local Government and I have been having a few issues that I can't seem to resolve.

First,  If I just change the web map id to the map published on our server in the config.js file on our server and I don't change any of the basemap settings I can get the app to partially work but if I look and the js console on chrome the app gets the following error say a token is required:

Also, when I go to change basemaps the app locks up.

When I search by Activity it throws a subsequent error which looks like it is related to the first:

Second,  if I change my basemap settings to the basemap group and owner in our orgnatation I get a javascript alert saying that I need to set the Geometry Service.  The issue sounds like something similar to this issue.

Here are the links to the two apps I have running to see the issues in real time:

https://gis.linncounty.org/webapps/gis/parklocator/

https://gis.linncounty.org/webapps/gis/parklocator_basemaps/ - This one gets the alert popup saying to provide a geometry service.

I followed esri's tutorial on setting the app up.  I am not sure what settings I have wrong or how to resolve this.


texts: QGIS equivalent to the ArcMap Zonal Statistics as Table

texts:
Is there a QGIS equivalent to the Arc: Spacial Analyst: Zonal: Zonal Statistics as Table?

Given the USDA Cropland Data Layer, Would like to count the pixels by type within a county.


texts: Converting KML with multiple layers to Shapefiles?

texts:
I have a KML-file (available here). It has several polygon paths for Alos/Palsar images.

I am in need to convert all these files from KML to shapefile, saving the layers using the layer "ID".

I am trying to use the piece of code bellow (PyQGIS), but I do not know how to unzip this KML and save the shapefiles using the layer "ID" information

vlayer = QgsVectorLayer("D:/KML-Images/A_Cycle06.kml", "line", "ogr")

writer = QgsVectorFileWriter.writeAsVectorFormat(vlayer, "D:/KML-Images/<layer ID>.shp", "utf-8", None, "ESRI Shapefile")


texts: How to lock the layers into the same location reference

texts: I'm obviously very confused here. So I am using Google Maps as a base layer, as I want to get my new shapes in the right place for a map. Google Maps doesn't seem to allow me to zoom in particularly far, maybe I need to save it as an image rather than use a live map but that is a problem for another day. I select a CRS, enable on the fly, and ensure that the same CRS is used for a new shapeile layer. To my mind that should mean that when I use the Pan Map tool to move about a bit (and I mean only a bit) and my new shapes move position relative to the Google Map. My understanding was that the whoe point was to create maps that retained geographical integrity? Can someone please correct my assumption or my utilisation (preferably the latter)?


texts: Shrinking a line to fit within a raster surface?

texts:
I have a single polyline feature that extends beyond my raster surface. I would like to shrink this line so that the end-points are within the last pixel at the edges of the raster.

I have tried converting the raster to a polygon and clipping but it sometimes puts the endpoint just outside of the raster in no data land.

Is there a simple method that I could use to intersect/shrink the line? ArcObjects is fine.


texts: rename table/feature class in filegdb using ogr or qgis

texts:
According to the OGR SQL docs, it doesn't appear that it is possible to use and ALTER statement to rename an existing table (featureclass) inside of a file geodatabase. RENAME statements can be run on columns, however. Is there a way to rename an existing table using OGR or QGIS?

This works to duplicate an existing table with the new name, and then the original table can be dropped, but it would be nice not to have to do that:

ogr2ogr -f 'FileGDB' -update -append -nln newtblname geodb.gdb/ geodb.gdb/


texts:
Using ArcMap 10.3.

I am trying to change 100+ annotation values from "Placed" to "Unplaced" as quick as possible.

I have tried to select them and use field calculator, but keep getting an error for Empty Values. The only way I can get it to change is to individually click each attribute and use the drop drown box to change it to place/unplaced.


texts: I have a map that has a 12 mile buffer round it already. I have created a buffer inside the original 12 mile line (using a negative value of 9). Is it possible however to stipulate that the new buffer does not overlap the original map polygon? Hope this is clear?


texts: Resolving ArcSDE for Oracle 11g Layer not found (-20) error?

texts:
I am using ArcSDE for Oracle 11g, with Oracle spatial database. I use the sdelayer command ot register my feature class. Surprisingly I get this error of Layer not found(-20) error. The table does exist in Oracle and I had configured feature classes using the same data source.

An Esri Technical Article suggests to increase tablespace, but I do have free spaces in my tablespace.

Please advise me on this?


texts:
Disclaimer: I am basically brand new to VB.NET and ArcObjects.

I am computing some statistics as a component of a larger arcobject script.  I write these components to an ITable.  This ITable is added to the map display as a standalone table.

I find that I am duplicating a fair bit of code to assign a value to each field in the table.  This is not something I would do in Python and wonder if a better methods applies in the VB.Net world.

VB.NET CODE:

   'Populate the stats table with the iteration.
    Dim row As IRow = tTable.CreateRow
    Dim fields As IFields = row.Fields

    Dim fThres As Double = fields.FindField("ThreshDist")
    row.Value(fThres) = lStats(0)
    row.Store()

    Dim fCountR As Double = fields.FindField("CountR")
    row.Value(fCountR) = lStats(1)
    row.Store()

where lStats is a list of doubles, and tTable is an ITable object.

In python I would simply use a dictionary to map the index (as key) to the fields name, as value.  A simple for loop would take care of population using an insert cursor.

Does a similar paradigm exist using VB.Net?  Am I missing something obvious in the documentation?


texts: How to run a simple python script for QGIS from outside (e.g. Sublime Text)?

texts:
I would like to run a few small and simple python scripts for QGIS (Mac) from "outside" of QGIS. With outside I mean in this context, either the normal os command line (terminal.app) or even better, directly out of Sublime Text (text-editor), but definitely not via the inbuilt QGIS python console.

I've read through various tutorial e.g. http://www.qgis.org/pyqgis-cookbook/intro.html#python-applications and I am able to get a reference to the QGIS app, but unfortunately not to qgis.utils.iface or something else deeper. This little code snippet should for instance print out the name of the active layer ... here is what I have:

    import sys
    sys.path.append("/Applications/QGIS.app/Contents/Resources/python")

    from qgis.core import *
    import qgis.utils

    print "helo" # console output: helo

    QgsApplication.setPrefixPath("/Applications/QGIS.app/", True)
    QgsApplication.initQgis()

    print QgsApplication # console output: <class 'qgis.core.QgsApplication'>

    print qgis.utils.iface # = console output: none

    aLayer = qgis.utils.iface.activeLayer()
    print aLayer.name()

    QgsApplication.exitQgis()

Please don't reply me now how to create a plugin or something else "heavy". I am just looking for a quick an easy way to shoot scripts out of a comfortable text-editor to QGIS.


texts: Fixing error in if (any(model$range < 0)) statement when fitting variogram in R?

texts:
I am using R and trying to fit a variogram model that I created:

 mod <- vgm(psill = var(newdfNoNAs$yield_cleaned),range= NA, model = 'Sph',nugget = 0)

And I used

 fit_reml <- fit.variogram.reml(yield_cleaned~1,data = newdfNoNAs,model = mod)

but it turns out

 Error in if (any(model$range < 0)) { : missing value where TRUE/FALSE needed

Can anyone tell me what it means and how to fix this problem?


texts: When does Shapely's "polygonize_full" detect a dangle?

texts:
Shapely's polygonize_full operation, which merges line segments into polygons, returns four geometry collections: polygons, dangles, cuts and invalid ring lines. Dangles are defined in the documentation as "edges which have one or both ends which are not incident on another edge endpoint", while cuts are "connected at both ends but do not form part of polygon".

I assumed that the following code, which tries to polygonize a square with one edge missing, should find two dangles and one cut edge:

import shapely

lines = [
    ((0,0),(1,0)),
    ((1,0),(1,1)),
    ((1,1),(0,1))
]

polygons, dangles, cuts, invalids = shapely.ops.polygonize_full(lines)

However, the polygons, dangles and invalids collections are all empty, while the cuts collection contains the three original edges:

[str(cut) for cut in cuts]
# This returns ['LINESTRING (1 1, 0 1)', 'LINESTRING (1 0, 1 1)', 'LINESTRING (0 0, 1 0)']

When does polygonize_full detect a dangle? Is there a way to get the dangling edges in my test case (0 0, 1 0 and 1 1, 0 1) using Shapely?


texts: Create a .csv of list of number of non null cell from raster in a mapset

texts:
How can I create a .csv file (maybe output of a function) of number of non null cells of raster from a given list in GRASS (e.g. rasters could be: myraster#1... myraster#500)?
Unfortunately g.list cannot print the number of non-null cells, so I was thinking to a loop as:
for Raster in 'g.list rast pattern=myraster* mapset=mymapset'

do eval 'r.univar map="$Raster" '

$n << path/to/my/non_null_cells.csv

done

where n is the number of non-null cells from every raster.
But, even in this case, doing like this I would only overwrite each n with the subsequent in the .csv.
Is there any proper command/function? Should I create a vector and then export that vector to a .csv?


texts: Connecting to Oracle SDE using QGis 2.18 only seeing MDSYS

texts: I am trying to add a layer using the Oracle connection in QGIS. I can connect when I test the connection. But what is returned is only the MDSYS owner (schema) and tables. There is one schema that I am expecting to see but it's not being returned. In other tools I can see the schema, but in QGIS I cannot. I've tried all different variations and checked "also list tables with no geometry." Which I wouldn't want to use because I am trying to grab a spatial table(s) under a schema. I'm wondering if the problem is that I am using Oracle SDE for the spatial tables... Any ideas on how I can return more schemas through this connection tool in QGIS?


texts: Is there an effective limit to the geoprocessing task Geocode Addresses?

texts:
The issue I face has the following background:

We are geocoding an immense directory of records. It is something in excess of 100 million addresses. I have split the addresses up into as small a geographic reason as I feel is feasible, states. Even so, a single state can have in excess of 5 million records to geocode. The arcpy script I composed will loop through each state's addresses and run the proper geocoding process with the appropriate locator.

The problem:

Sometimes it takes 36 hours, sometimes it takes only 20, but the script will stop with an error I have never seen before (this is from the XML log file):
ERROR 001143: Background server threw an excecption.

If the script is run in ArcCatalog or Arcmap's Python window, it will show a red error message Runtime error <class 'arcgisscripting.ExecuteError'>: ERROR 001143: Background server threw an excecption.
If the script is run in IDLE, it will simply cease processing without an error and restart the shell (with the typical =====RESTART===== heading).

I know this arcpy script works with smaller datasets, as I have been using it for months now. What might be the cause of this error? Do I need to split up my address listings into smaller amounts for them to geocode reliably?

EDIT: I tracked down the Event Viewer Logs per @D.E. Wright and this is what I found listed under the most recent failure

Faulting application name: pythonw.exe, version: 0.0.0.0, time stamp: 0x4ba3e4e2
Faulting module name: Geocoding.dll, version: 10.0.1.2800, time stamp: 0x4cbcbb71
Exception code: 0xc0000005
Fault offset: 0x000be1f3
Faulting process id: 0x%9
Faulting application start time: 0x%10
Faulting application path: %11
Faulting module path: %12
Report Id: %13

Not exactly illuminating...


texts:
I have one population raster showing how many 'people' are in a grid cell (i.e. population count per cell)

And I have another raster for water available per cell.

I want an output raster that is classified based on how much water is available for 1 person in a cell. for example:

If 1 person has 1000 cubic meters of water for matching cells in the input rasters, that grid cell in the output raster will have a value of 1.

How do I achieve this in arcmap 10.2. I have the spatial analyst tool activated.

I feel this is a conditional query but I don't know how to do it. It should be something like this in english:

If water raster has cell values of <500, 500-1000, 1000-1700 and >1700 for each 1 person in a cell then output raster will record values of 1, 2, 3 and 4 respectively for cells with these values.

The issue is that population raster shows total number ofpeople in a cell, how do I factor this into the query which is looking at what is available per person.


texts: Calculate distance between points, keeping attributes from input files?

texts:
I'm trying to use Generate Near Table in ArcMap to calculate distance between apartments and public transport.

It works fine since I'm getting the distance in meters, but I would also like to have the information about what the name of the station is.

I.e. from apartment A there is X meters to the nearest station, and here I would also like to have the name of that particular station. Basically another column in the table as you see in the image.

It doesn't seem like this tool allows including any attribute fields from input data.

The only solution I can figure out on my own, would be to break up my shapefile with stations into separate shapefiles, however that would be time consuming and inefficient since there are a lot of them.


texts: How to seralize ArcGIS JSON to .Net core object?

texts: I want to get ArcGIS JSON response (like here) from server and save it to my database. I am using .Net core and entity framework net topologysuite for CRUD operations database. But I coud not convert the Esri JSON object to .Net object. Is there any tool for this?


texts:
I'm working on a model that iterates and calculates a "topo to raster" process for different extents to create a mosaic dataset. The process is very time-consuming even for the small areas and I want to make it so that I can confirm every iteration's execution.
I can't just oversee the model and hit cancel when my desired files are calculated. Once the model executes the process, canceling doesn't stop it. It only makes the process not save the calculated raster which is even worse.
Is there any way to create a pop-up window or something similar that shows a message (for example: "Start calculations for province [XXXX]?") and asks me for a response (for example: Yes or No) and stops the model if my input says so?


texts: Cross validation after kriging, in ArcGIS

texts:
I performed kriging interpolation using the kriging interpolation tool in ArcToolbox, not the Geostatistical Wizard method. However, how do I determine the cross validation without using the geostatistical wizard method. In ArcToolbox, there is a cross-validation tool, but when I try to input the kriging layer I created, it doesn't allow me to.

Can I only use cross validation after creating the krig using the geostatistical wizard method?


texts:
How can i get only attributes names from GeoJSON?

ogrinfo gets attributes with others info:

ogrinfo -al -so rw.geojson

Output is:

    ERROR 4: GeoJSON Driver doesn't support update.
Had to open data source read-only.
INFO: Open of `rw.geojson'
      using driver `GeoJSON' successful.

Layer name: OGRGeoJSON
Geometry: Line String
Feature Count: 31888
Extent: (12.224042, 48.586482) - (18.789360, 51.043735)
Layer SRS WKT:
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        TOWGS84[0,0,0,0,0,0,0],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0,
        AUTHORITY["EPSG","8901"]],
    UNIT["degree",0.0174532925199433,
        AUTHORITY["EPSG","9108"]],
    AUTHORITY["EPSG","4326"]]
FID Column = id
id: Integer (0.0)
railway: String (0.0)
railway_tr: String (0.0)
usage: String (0.0)
i d: Integer (0.0)


texts: What jetty version the geoserver 2.12.1 has?

texts:
An older version of Geoserver (2.1.3) had Jetty, v.6. The problem is that Jetty 6,7 and 8 are End Of Life (here and here).

Do you know what jetty version the Geoserver 2.12.1 has ? Is it a new, stable supported jetty? I am looking at the Geoserver site and cannot find any info.


texts: Image displays in 3 RGB stripes across the surface

texts: I layerstacked the Image in the screenshot above, but it is displaying in 3 stripes of RGB. Anyone who knows what the problem could be?


texts: Orthorectified images have removed relief displacements and tilt displacements. Since the assumption of orthorectification is to keep the light source at infinite distance so that the light rays are parallel to each other and perpendicular to ground. They should not have shadows. Is it true?


texts: Accessing input layer properties from QGIS Modeler?

texts:
So I start my model by dragging an input layer, which in my case it happens to be a vector layer. Let's call it inputVectorLayer. Later on in the model I want to refactor this layer using the Refactor fields algorithm. For each feature I write the correct expression, but for one feature I need an expression that should include something like this:

layer_property('inputVectorLayer', 'name')

The problem is that it doesn't recognize the name that I gave to the representation of the input layer in the Model Builder. How can I achieve this?

If I do something like this:

layer_property('actualLayerNameAsItShowsInTheLayesPanel','name')

it works, but obviously, this name changes. That's the whole point of the module, to use aliases to reference input layers.
I already tried to use double quotes to no avail layer_property("inputVectorLayer", 'name')


texts: Avoid getting power poles when extracting trees in GlobalMapper

texts:
I am using GlobalMapper V21 and trying to extract trees from a municipal dataset.  The classification isn't very good as I see building outline points being considered ground or high vegetation, so I am stepping back to reclassify points.

When I have done that I am trying to extract the medium and high vegetation points but I keep getting all the power poles as well.  I thought about using intensity to differentiate a power pole from tree canopy but not sure if that would work as when I am looking at the intensity of the power pole I am getting similar values to the trees.

Any suggestions on how to get just trees and not power poles?