texts: Finding SDE permissions on feature class using ArcPy?

texts:
I am writing an arcpy script where an SDE feature class is removed then recreated.  I need to set permissions on the recreated feature class that match the permissions before it was removed.

Is there a way to identify existing explicit permissions using arcpy?  I've taken a look through arcpy.Describe() and in the Geodatabase Administration toolbox, but can't see anything that will do this.

I'm going to use arcpy.ChangePrivileges_management() to set the new permissions once the feature class has been recreated.

The only way I've come across to get the permissions is to use use arcpy to pass a SQL query to find permissions:

myDB = "database.sde"
myFC = "FCName"
permissionsQuery = permissionsQuery = """
    SELECT permission_name, state, pr.name
    FROM sys.database_permissions pe
    JOIN sys.database_principals pr ON pe.grantee_principal_id  = pr.principal_id
    WHERE pe.class = 1
    AND pe.major_id = OBJECT_ID('{}')
    AND pe.minor_id = 0
    ORDER BY pe.permission_name""".format(myFC)

dbConn = arcpy.ArcSDESQLExecute(myDB)

permissions = dbConn.execute(permissionsQuery)
print permissions

which will output something like

[[u'DELETE', u'G', u'midavalo'], [u'INSERT', u'G', u'midavalo'], [u'SELECT', u'G', u'midavalo'], [u'SELECT', u'G', u'some_role'], [u'UPDATE', u'G', u'midavalo']]

which I can then use to pass back to arcpy.ChangePrivileges_management() but would like to make use of an arcpy tool to do this properly if available.

ArcGIS (desktop and server) 10.3.1
SQL Server 2012


texts: Adding label to QGIS map composer through Python script?

texts:
I need help in adding a label to QGIS map composer through python script. I could add a label but it is incomplete. The label is like this: 27-02-15 11-25. It only shows the first part: 27-02-15. I don't know why.

Also when I create the composer map, the other elements were put outside of the composer map rectangle. It is possible to overlap elements in map composer?

EDIT
I see that the problem is in the width and height of label box. I can't find a way of change this parameters. The adjustSizeToText should do it, but somehow is not doing that.


texts:
To the best of my understanding, setting a specific day under Analysis Settings for the VRP (vehicle routing problem) uses the Routing_ND historic traffic information.

Is there some step I'm missing to enable traffic? I've routed some orders during PM hours, then changed the times to AM, and the solver returned the exact same times/distance/stops.

I'm using ArcMap 10.2 with Network Analyst and StreetMap 14_1


texts: How to create a zoomed in portion of a map in print composer?

texts: I posted a question yesterday about troubles I was having with labeling (which are still plaguing me), but it hit me today that it would probably be easier for me to just create zoomed in sections of my map to be able to show all of the labels properly (they are very close point locations). The problem is that I have no idea how to do this! I know that I can add another map in a smaller window and zoom from there, but I'm not sure how to get the boxed area perfectly matched up (if that makes sense). Can anyone help?


texts: Column names have moved in attribute table QGIS

texts: I've joined an excel spreadsheet to one of my point shapefile. After performing the join, I can see that the names of the columns have not been set as the name in the attribute table column; instead the names are one row below as shown in the picture. Any insight on how I could fix this?


texts: A way to update CartoDB's info window in real time?

texts:
With mapbox, I'm able to populate a JS array and draw a layer in realtime. I'm wondering if, once I pulled a map with createVis in cartodb, I can populate the infowindow with data per shape.

For example I'm pulling (from the PHP) the latest statistics per shape. And I draw/render the cartodb, and I want, on infowindow open, to show that information (so if I have shape id 60, I want corresponding statistics from my local mysql database). Is this possible? Or do I need to store all the data into cartodb and then redraw the map?


texts: Vector difference between two polygon layers by attribute

texts:
Is there a way to do a 1:1 difference between two polygon layers based on a similar attribute? The layers have the same number of features with a shared attribute. I've tried the native QGIS difference algorithm but it subtracts from the features of the first layer, the dissolved/merged area of the second layer.

Edit: Added photo showing sample layers.


texts: How to convert GOES-R (15) IR image from gray-scale to brightness temprature

texts: How do I convert satellite IR images taken by GOES-R (15) from gray-scale to brightness temperature? We can download the IR satellite images, and read its RGBA values. But how can we convert this RGBA value into cloud temperature?


texts: I would need to use data from this page: http://dosoremi.hu/int_agyag_0_30.html. However, changing windows all the time isn't too effective. Do you think I can get the same images into my QGIS maps?


texts:
I'm running GeoServer on my own server at home. Everything works great usually but when I try to connect to WFS from QGIS over https I get a warning saying SSL Handshake Failed. This is probably because I have a self-signed certificate.

Do I have to pay for an authorized certificate or is there a way around this?

https works fine in the web browser.  I just have to click to accept the certificate. I really don't want to pay $70 a year for a certificate just to use my GIS server for my own personal experimental use.


texts: How to count the number of buses that pass through a stop during a day?

texts:
GTFS spec for reference: https://developers.google.com/transit/gtfs/reference/

I'm trying to create a program that counts the number of buses, trams, rails etc that pass through a certain stop during a day, for a certain route, and do a split-count for those that are end of lines.

My algorithm is as follows:

Given a route R
Given a day D

foreach R.trips as trip
    if trip.calendar_dates != D
         // skip all dates that we are not interested in
         continue
    endif

    endOfLines = 0
    notEndOfLines = 0

    foreach trip.stopTimes as stopTime
        // get all stop times of this trip, ordered by stop_sequence
        stopTimes = stopTime.trip.stopTimes.orderBy(stop_sequence)

        if stopTime.isFirstIn(stopTimes) or stopTime.isLastIn(stopTimes)
            endOfLines++
        else
            notEndOfLines++
        endif
    endforeach

    results[route.route_type][end] = endOfLines
    results[route.route_type][notEnd] = notEndOfLines
endforeach

Explanation of what I'm trying to do:

Take all the trips of the given route. Use calendar_dates to link the trips to the specific day and remove the trips that are not in our day of interest.
Loop all the trips and take their stoptimes. Then, determine if this stop time is the first or last in the trip. If it is, then we count it as an end of line, otherwise we count it as a simple stop.

Since the GTFS spec is somehow vague about how the tables should be used, I'm trying to figure if this logic is correct. Can someone validate and point out if I have any mistakes in my logic?


texts: How to obtain GeoServer layer's shapefile name using REST API

texts:
I have a GeoServer store of type "Directory of Shape Files".  I have layers defined in that store corresponding to specific shapefiles. Some of my layers have been renamed.  I'm looking for how to pull from GeoServer the name of the shapefile a particular layer from this store is referencing.  Seems like it should be easy but I cannot find it anywhere...

Do I have to look in the XML files?

This is an old GeoServer running 2.8.2, I'm planning to move some layers to a newer version, and I want to make sure I get the right data files...


texts: OpenLayers layer-switcher for WMS layers

texts:
I have configured a WMS layer for OpenLayers like:

var wmsSource = new ol.source.ImageWMS({
  url: 'http://<servername>/geoserver/<namespace>/wms',
  params: {'LAYERS': [ 'Layer 001', 'Layer 002', ... 'Layer xxx' ]},
  ratio: 1.2,
  serverType: 'geoserver'
})

var wms = new ol.layer.Image({
  title: 'Basemap',
  type: 'base',
  source: wmsSource
})

The layer-switch I use is this one: https://github.com/walkermatt/ol3-layerswitcher

It is configured like

var layerSwitcher = new ol.control.LayerSwitcher({
    tipLabel: 'Legend'
});
map.addControl(layerSwitcher);

Now the layer-switcher only shows one "Basemap" entry. But I want to have the WMS-Layers listed - best in a group "Basemap".

Any ideas how to do this?


texts: QGIS: Error [24]: The following signature does not match the band set. MC_ID: 1 C_ID: 1

texts:
I have run my classification map in QGIS using the semi-automatic classification plugin or SCP.
I want to use the same ROIs for different sentinel imageries. But I faced an error

Error [24]: The following signature does not match the band set. MC_ID: 1 C_ID: 1

How can I use the same ROIs for different band sets?


texts: Is there Shape Metrics Toolbox or equivalent for ArcGIS 10?

texts:
I need to compute the geometric properties (e.g. contiguity index, dispersion index etc.) of a series of polygons and I would like to use the Shape Metrics Toolbox.

I am using ArcGis 10 and the version of the toolbox that I could find on the website is said to be for ArcGis 9.3. When I try to use the tool I get an error message which might be related to an incompatibility with ArcGis 10:

ERROR 999999: Error executing function.
<type 'exceptions.SystemExit'>: 1
Failed to execute (ShapeMetrics).

I found several papers on line in which the authors mention an ArcGis 10 version of the Shape Metrics Tool but I couldn't find this version anywhere. I have contacted the authors but have received no reply yet.

I was wondering whether any of you have experience with the Shape Metrics toolbox in ArcGis 10,  have encountered the error message above or know of other user-written tools for ArcGis 10 which compute shape metrics.


texts: Is there a vector affine transformation tool in QGIS2?

texts:
In QGIS 1.8, I can transform my point layer with the plugin "Affine transformation". So, this plugin don't exist in QGIS 2.0.

Do you know a method to simply translate (i.e : X -> X+30000) in QGIS 2 ? It's possible in QGIS or I must use ogr2ogr, QGIS 1.8 or openjump ?


texts: Extracting new data with ArcGIS Pro ModelBuilder?

texts:
I'm very, very new to ModelBuilder.

I have a database that's being continuously updated with new addresses. I'm trying to set up ModelBuilder such that, each time I log in to my computer, ArcGIS geocodes any new addresses that have appeared since I last logged in. So I suppose that, with each log in, I'd have the computer scan over the previous log-in's list and geocode any thing that's new.

I've already figured out how to schedule running the model with each new log in, but how would I command ArcGIS to do this comparison?


texts:
I am trying to Add a layer with Youtube in the info window as follows

HTML

  x

Javascript

$('#map').ready(function(e){
    cartodb.createVis(map, 'https://internetadventureclub.cartodb.com/api/v2/viz/4d794d4a-8b43-11e5-ada8-0ecd1bab222de5/viz.json')
        .on('done', function(vis,layers) {
          // get sublayer 0 and set the infowindow template
          var sublayer = layers[1].getSubLayer(0);

          sublayer.infowindow.set('sanitizeTemplate',false);

         //  sublayer.infowindow.set('template', $('#infowindow_template').html());
         }).on('error', function() {
           console.log("some error occurred");
         });
  });

The moment we add other data layers in addition to the video layer. This seems to make the video layer stop working.


texts: Text cutting off when exporting map using ArcMap

texts:
For some reason my labels keep cutting off; the same happens with my legend when converting to pdf or jpeg.
I am using ArcMap 10.8 and have never experienced this issued until my ArcMap was reinstalled on a new laptop.
I have attached a map - the label are suppose to read Welverdien, Utah Camp, camp b and so on, from the picture it is seen that only letters are showing and legend is also scrambled.
How can I fix this?


texts: Importing multiple layers in QGIS: ERROR: Too many connections: max 64

texts: I get this error in the terminal when I open more than 40 layers in a project and QGIS eventually crashes. It doesn't happen when I open them into a completely new project. How can I avoid this?


texts: ArcGIS Composite Bands tool does not recognize temporary layer input

texts:
I am trying to extract monthly layers from a multi-year stack and then re-stack them by month. I am using make raster layer to make temporary layers and append them to a list to use Composite Bands tool for re-stacking. However, it appears that the composite band tool does not recognize the input (ERROR 000732: Input Rasters: Dataset ... does not exist or is not supported). I used a similar method to create temporary layers before, however not from a stack but from a list of individual files, the composite band tool recognized the input and functioned well.

How can I identify the problem?

import arcpy, os
arcpy.env.workspace = 'D:/GIMMSdata'
from arcpy.sa import *
arcpy.CheckOutExtenstion = 'Spatial'

list_NDVI = []

for i in range(1, 205, 12):
    lyrName = 'GimmsNDVI' + str(i) + '.tif' + '_temp'
    print(str(lyrName))
    list_NDVI.append(arcpy.MakeRasterLayer_management('D:/GIMMSdata/GIMMS50N.tif',lyrName, '', '', i))
    i = i + 12
    print(i)

arcpy.CompositeBands_management(list_NDVI, 'd:/GIMMSdata/GimmsNDVIJan.dat')


texts: How to batch process .txt files with X,Y data (in Arc or QGIS)?

texts:
I have a long list of .txt files that I need to convert to feature class. I can bring the table into ArcMap, right-click, Display X,Y Data, and it converts to a feature no problem.

Is anyone aware of a way to convert all of them at once, either in QGIS or ArcGIS? I was hoping to do this in ModelBuilder but there is no "tool" in the search option for Display XY Data.


texts: v.net.steiner tool in QGIS not producing output

texts:
I have a shapefile of a network of points mostly on roads (a county shapefile from TIGER downloads). I want to find a minimal length tree along the road network that hits all points within 100m of the road. I've tired many different ways of using the v.net.steiner GRASS tool but it doesn't produce any output, saying the layer was not correctly generated. I've tried snapping, using random points instead, changing projection, etc, but always the same result.
Is there anything else I can try? Or is there any other way to approximate a Steiner tree out of two shapefiles?


texts: Get a single cluster from cloud of points with specified maximum diameter in postgis

texts:
I have on the order of 2M Geometry Collections with N (median: 2, stddev: 6.6, 75th: 6, 95th: 21) Point objects (repeated points allowed and important) in a PostGIS table. I want to identify clusters within each collection, with a specified maximum diameter. From these I would identify a main cluster (which I assume exists), based on the number of points it contains (and in case of a tie, choose the cluster with the smallest area).

Since I will be using these clusters further in the database, I'd rather use a technique within PostgreSQL/PostGIS.

K-means clustering is not ideal because it requires specifying a priori the number of clusters, which is impossible due to the heterogeneity of the data. And specifying 1 cluster would return the entire collection.


texts: ogr2ogr creates too many columns with LEFT JOIN

texts:
I am creating a GeoPackage using several shapefiles using ogr2ogr (GDAL 2.2.4):
ogr2ogr -f GPKG /path/to/out.gpkg -append -nln out -a_srs EPSG:4326 /path/to/in1.shp
ogr2ogr -f GPKG /path/to/out.gpkg -append -nln out -a_srs EPSG:4326 /path/to/in2.shp
...

So far it does exactly what I want. All files are merged in a single GeoPackage file as long as all Shapefiles share the same column structure.
Now I would like to add columns based on a relation with another table, which can be done using a JOIN from the SQL dialect:
SELECT * FROM in LEFT JOIN '/path/to/other.dbf'.other ON in.a = other.a

I simply use the -sql option to pass the SQL query and voilà!
ogr2ogr -sql "SELECT * FROM in LEFT JOIN '/path/to/other.dbf'.other ON in.a = other.a" -f GPKG /path/to/out.gpkg -append -nln out -a_srs EPSG:4326 /path/to/in.shp

It does what I asked, but it also adds unexpected columns.
Let’s say in.shp contains the following columns: a and b
And other.dbf contains the following columns: a and c
I would except to get columns: a, b and c
However I get columns a, b, c, in.a, in.b, other.a and other.c
(for clarity I did not mention the mandatory columns geom and fid)
I thought that SELECT * was a bit too inclusive, so I tried to select a limited set of columns:
SELECT a FROM in LEFT JOIN '/path/to/other.dbf'.other ON in.a = other.a

And I get the following columns: a, in.a, in.b, other.a and other.c
It looks like the in.* and other.* columns are always inserted as long as I perform a LEFT JOIN.
Is there a way to get rid of them?


texts: gdalwarp's minimum resampling has no effect (possible bug)?

texts:
I'm attempting to use the minimum operator with gdalwarp v2.1.0 but it always mosaics as if I had selected nearest neighbor.

From the docs:

-r resampling_method:
min:  minimum resampling, selects the minimum value from all non-NODATA contributing pixels. (GDAL >= 2.0.0)

The command is staightforward:

gdalwarp.exe -r min -srcnodata -9999 -dstnodata -9999 a.tif b.tif ... z.tif output.tif

But it gives me this output:

When what I expect is this (mosaicked through ArcMap with the minimum operator):

This totally looks like a bug to me but I figured I'd ask here first in case I'm missing something.


texts: ModelBuilder not passing information through Iterate Feature Classes using ArcGIS Pro

texts:
I am using ArcGIS Pro.
I am trying to create a model to iterate through three feature classes in a file geodatabase, all with the same schema. The model performs several joins and field calculations before eventually returning a feature class with several more fields in addition to the original fields. My problem is that the original fields aren't showing up in the model. For example, when specifying a join field, I get:

Warning 000728: "Field ___ does not exist within the table"

As a result, they don't show up in the final feature class either.
What am I doing wrong?
Here's a screenshot of the relevant part of the model:


texts: Apply same coordinate system to raster image and GeoJSON with rasterio

texts:
I have a dataset containing raster images and corresponding features in geojson files. I would like to use rasterio.features.rasterize to create a footprint image of the features on the raster image, but the resulting footprint array stays empty. I think the issue is that the raster image and the features are defined for different coordinate systems.
The raster image uses EPSG-3857:
>>> data = rio.open('image.tif')
>>> data.crs
CRS.from_epsg(3857)

and the underlying transformation seems to be defined in bounding box coordinates:
>>> data.transform
Affine(4.77731426716, 0.0, -1942112.014669758,
       0.0, -4.77731426716, 1658377.7656751834)

The geojson file on the other hand seems to be defined using WGS 84 coordinates:
{
"type": "FeatureCollection",
"crs": { "type": "name", "properties": { "name": "urn:ogc:def:crs:OGC:1.3:CRS84" } },
"features": [
{ "type": "Feature", "properties": { "Id": 0 }, "geometry": { "type": "Polygon", "coordinates": [ [ [ -17.441005042849905, 14.689811850766825 ], [ -17.441172607380693, 14.689811850766825 ], [ -17.441172607380693, 14.689967040094865 ], [ -17.441005042849905, 14.689967040094865 ], [ -17.441005042849905, 14.689811850766825 ] ] ] } },
...

If I apply the affine transformation from the image to pixel coordinates, I end up with coordinates in bounding box coordinates:
>>> data.bounds
BoundingBox(left=-1942112.014669758, bottom=1653485.7958656116, right=-1937220.0448601863, top=1658377.7656751834)
>>> data.transform * (100, 100)
(-1941634.283243042, 1657900.0342484673)

My two questions are:

How do I transform these coordinates into world coordinates? And...
Is there a way to bring the geojson file coordinates on the same system as the raster image to generate the footprint?

As a reference, I wrote the following script to rasterize the features on the image:
>>> import rasterio as rio
>>> from rasterio import features
>>> import geopandas

>>> data = rio.open('image.tif')
>>> df = geopandas.read_file('labels.geojson')

>>> feature_list = list(zip(df['geometry'], [255]*len(df)))
>>> fpt = features.rasterize(shapes=feature_list,
                             out_shape=data.read(1).shape,
                             transform=data.transform)

>>> print(fpt.sum()) # to show that the footprint raster is empty
0


texts: Delete records from table using ModelBuilder

texts:
I'm very new to ModelBuilder so forgive me if there is a relatively simple tool to perform this. I have looked into Delete Field/Table but it does not seem to achieve the goal.

Goal: Create a tool that will update the master parcel feature class on one geoDB with combined data from a feature class and table in another location.

How do I go about deleting the data that is in the master parcel feature class without deleting the headers? I basically want it to be a blank template with headers to later append/join the newest information with. Or would it be simpler to delete the file and just save the joined information in its place?

I'm also not sure if I have permissions to delete files, only modify, so finding a way to go about this without deleting would be great.

edit:

I realize that I may have have explained this poorly. Our master parcel feature class is the result of a manual combination of a geometry feature class and a parcel data table that is provided by our county assessor office. The goal of the tool would be to automate this process.

My current workflow:
BlountParcels - Master Parcel feature class (what needs to be updated/overwritten)
Parcels - Geometry feature class
camagisdata - table parcel data


texts: Labelling random point with it's constraining feature name

texts:
I was wondering if after a random point creation and constraining it to a feature class (with several features). Can I add a field (to the output table) which displays (names) the correlating feature the point was constrained to? I am assuming since the point was constrained to a feature, that information is available somewhere - though I might be completely off on this assumption.

For example, I create 5 random points in a feature class that has 5 features (5 polygons). In the output, am I able to indicate/label in a new field each random's constraining polygon?


texts: I have two layers, one is point layer and other one is a line layer. I want the point features to snap on to the lines, but the points are away from the lines. Is there any way to move the points onto lines?


texts: QGIS Georeferencer add point tool doesn't function

texts:
Using the Georeferener, after selecting the add point tool (to add a GCP) I am unable to place a CGP onto the image. Clicking the mouse pointer does not add a GCP as is has done when I've used previous versions. It just doesn't seem to function.

I'm running QGIS 2.14.3 Essen on a Lenovo thinkpad running Windows 10.

I have successfully georeferenced many images before. However trying with a number of images, JPG and BMP the add point tool does nothing.

Any ideas on a fix?


texts: An example Tissot ellipse for an Equirectangular Projection?

texts:
I am trying to calculate distortion of an equirectangular projection via Tissot indicatrices.  I have tried following the directions on this post (among other things) but it is beyond my understanding as an amateur.

So, I am wondering if someone would be so kind as to calculate a single Tissot ellipse for an example equirectangular lat/long (whichever is your favorite and is distorted on an equirectangular projection).  I don't understand exactly what the variables are and where they come from, so seeing the equations in action would be very useful.

I am trying to understand these equations so I can plug them into a mapping program I am coding.  I asked a bunch of general questions in this thread, but I think a specific example will help me figure the rest out.

Many thanks, as always.

NCashew


texts: Appending to an existing feature class attribute and batching to the rest of the attributes in the feature QGIS 3.10

texts:
I am trying to figure out how to add (before and after the existing data) to a specific feature class and batch it for the rest of the attributes.
For instance, if I have a set of numbers, "12345678", and I need to add a prefix to this to make it "FIB_AER::123456789". I hope explained it in an understandable way.


texts: Using PyScripter with ArcGIS for Desktop?

texts:
Can someone please let me know how I can set the PyScripter 2.5.3.0 x64 works with Python 2.7.3 which comes with ArcGIS
what I did was:

I installed the ArcGIS 10.2
I installed PyScripter 2.5.3.0 x64 which forced my to install Python 3.3
I installed Python 3

Now I can't use the PyScripter for ArcGIS Python. I already add the ArcPy in to PyScripter IDE Options--> Special package as: os, wx, SciPy, ArcPy but still not working.

I also check the "desktop10.2.pth" at c:\python27\ArcGIS10.2\Lib\site-packages which looks like:

C:\Program Files (x86)\ArcGIS\Desktop10.2\bin
C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy
C:\Program Files (x86)\ArcGIS\Desktop10.2\ArcToolbox\Scripts

Can you please let me know how I can fix this issue?


texts: How to create an area or surface from point objects in Smallworld with MAGIK?

texts:
In Smallworld i have a geometry-set with simple-Point objects.

Is there a way to create a polygon, a surface or an area from these simple-Points with Magik? I prefer a pseudo-area, but these is not essential.

Later i would like to Intersect these area with other Point-Objects.


texts: QGIS OpenLayers plugin needs to be reloaded

texts:
Using QGIS 2.12, I successfully add basemaps via Web > OpenLayers plugin.

However, if I save the project (EPSG:3857) and later reopen it, the OpenLayers basemap will not display.  As a workaround, if I delete the "blank" basemap from the table of contents and then add it again via Web > OpenLayers, the basemap will display correctly.  However, this must be redone every time I reopen the project.

This problem is not a show-stopper, but it is aggravating.  Has anyone else experienced this oddity?


texts:
I am trying to use SAGA from within QGIS.

I can run SAGA standalone, but not from within QGIS - nothing appears in the Processing Toolbox.

I tried following Why is SAGA from the Processing Plugin not recognised? and updating the path, however I cannot update the path.  If I change it, it simply reverts to what it used to be.  I've tried closing and reopening QGIS but that does nothing.

I'm also not sure what I should be pointing to - is it the SAGA app? Something else?  I can't seem to point to the App directly, as I cannot navigate to the folder where it is.  If I move the app - it breaks.

I've found a few old instructions online, but no luck following them.  I'd prefer not having to move things back and forth between QGIS and SAGA if I don't have to.

I possibly (probably) have multiple versions of saga on my computer, as I've tried following many different sets of instructions for how to get this working.

QGIS version 2.14
SAGA version 2.31
Mac OS X version 10.12.16


texts: Labelling point elevations of raster (DEM) using QGIS

texts:
I'm trying to create a map with spot levels using an underlying 1m resolution DEM TIFF file. My thought process so far is to create a grid of points with a height value field over the area of interest (say 6m x 6m) and then select the highest in each cell. Then draw a 5 buffer around all these points. If any of these buffers contain more than 1 point keep the highest value and delete the rest. This should result in nicely spaced point that I can label neatly without overlaps or the chance of the label not being the highest.
I've tried making a template for this using the Graphical Modeller, but I'm not great with expressions (I have more of a survey/cartography background and its been a while since I did this at school!)
What is the best way to go about achieving this?
Images attached showing graphically what I'm trying to achieve.
Alternatively is there a way to use labels to give greater weight to points with a higher elevation field so that they are displayed instead of the lower values?

6x6m grid cells

Highest value in each cell

5m radius from high point in each cell. Keep highest point in each circle buffer that contains >1 points

remaining points for DEM elevation labelling

example of what I'm trying to achieve


texts: Finding free building height dataset for major cities around world?

texts: I am searching for building/property data which include building height, building area, lot size, and building type classifications(preferably as complete as possible). I am working on an Earthquake Model for highly urbanized locations, especially those around fault lines.


texts: Automating Multi Data Driven Pages Maps Export to PDF?

texts:
I´m currently working in this project where we are creating a time wise analysis of a certain region with many areas of interest.

We are using the data driven pages tool to address the many AOIs but we are still missing the automation part of exporting each different map (each map is from a different date).


texts: How to resolve a 'Received NULL pointer' ArgumentOutOfRangeException for OSGeo.OGR.Geometry.TransformTo?

texts:
I am attempting to project a geometry using the GDAL c# bindings and I'm getting a ArgumentOutOfRangeException ('Received NULL pointer') exception on the Transform methods of Geometry. I'm also getting the same exception when attempting to project a point using a OSGeo.OSR.CoordinateTransformation (e.g. https://gis.stackexchange.com/a/61574/31).

I am using the DLL in FWTools 2.4.7\csharp.

Has anyone seen this before and know what the fix is?

Here's a test that creates the exception for me:

        [TestMethod()]
        public void TransformPolygonReturnGeoJSON()
        {
            string jsonLambert = string.Empty;
            string projectedJson = string.Empty;

            string wkt = "POLYGON ((1831321.945291 412165.262452,1831324.478486 412152.596478,1831312.445811 412108.265569,1831298.513239 412108.265569,1831287.747162 412108.265569,1831275.714486 412102.565881,1831271.914694 412100.032686,1831266.215006 412098.766089,1831253.549032 412188.694503,1831321.945291 412165.262452))";
            int srid = 102009;

            int outputSRID = 26918;

            string sridWKT  = string.Format("EPSG:{0}", srid);

            string outsridWKT = string.Format("EPSG:{0}", outputSRID);

            var g = Geometry.CreateFromWkt(wkt);

            var sr = new OSGeo.OSR.SpatialReference(sridWKT);

            g.AssignSpatialReference(sr);

            jsonLambert = g.ExportToJson();

            Assert.IsFalse(string.IsNullOrEmpty(jsonLambert));

            var outSR = new OSGeo.OSR.SpatialReference(outsridWKT);

            try
            {
                var result = g.TransformTo(outSR); // <--- Exception occurs here

                Trace.WriteLine(string.Format("result={0}", result));

                projectedJson = g.ExportToJson();

                Assert.IsFalse(string.IsNullOrEmpty(projectedJson));
            }
            catch (Exception ex)
            {
                Trace.WriteLine(ex.Message);
                Trace.WriteLine(ex.StackTrace);

                throw;
            }
        }


texts: Understanding Promise deprecation warnings in ArcGIS API for JavaScript?

texts:
I am trying to update legacy map code to the latest ArcGIS version of 4.6 and noticed this warning in the console when querying a feature layer.

[esri.core.Promise] DEPRECATED: then() -- use .when(callback, errback) instead

So I tried to change then to when, but got when is not a function errors. I also tried running the sample code from the site and got the same errors.

view.whenLayerView(featureLayer).then(function(lyrView) {
          lyrView.watch("updating", function(val) {
            if (!val) { // wait for the layer view to finish updating

              // query all the features available for drawing.
              lyrView.queryFeatures().then(function(results) {

                graphics = results;

                var fragment = document.createDocumentFragment();

                results.forEach(function(result, index) {
                  var attributes = result.attributes;
                  var name = attributes.ZIP + " (" +
                    attributes.PO_NAME + ")"

                  // Create a list zip codes in NY
                  var li = document.createElement("li");
                  li.classList.add("panel-result");
                  li.tabIndex = 0;
                  li.setAttribute("data-result-id", index);
                  li.textContent = name;

                  fragment.appendChild(li);
                });
                // Empty the current list
                listNode.innerHTML = "";
                listNode.appendChild(fragment);
              });
            }
          });
        });

I also tried setting the "esri-promise-compatibility" flag to 1 but that didn't help. Can't see anything in the docs or community forums about this. Anyone come across this before?


texts: Open Source Tools for Mobile GIS app development?

texts:
I want to create a complete Map application with the facilities of edit geometry and attributes of layers. The server and client for the application should be some open source tools. The application need to run with Windows Phone 7, iOS and Android.

Please give me suggestions on the Open Source Tools selection for my development. As a first step i am planning to start with Windows Phone 7.


texts: How can I extend L.mapbox.legendControl?

texts:
I'm trying to override legendControl.onAdd, but I can't seem to do it early enough. I can override the other functions, since they're fired after map init, but I can't seem to change the prototype?  How should I do it?

This doesn't work:

L.mapbox.legendControl.__proto__.onAdd = function(map) {
    console.log('map.legendControl._onAdd');

    this._container = L.DomUtil.create('div', 'map-legends-wrap-container');
    this._subContainer = L.DomUtil.create('div', 'map-legends wax-legends');
    L.DomEvent.disableClickPropagation(this._container);

    this._container.appendChild(this._subContainer);

    this._update();

    return this._container;
};


texts: Splitting lines with PostGIS and view result in QGIS?

texts:
I have two tables, one contains lines (test_line) and the other contains the points at the vertices of the lines (test_line_vertex).

I want to split the line into individual segments based on the vertices.
I tried the following query:

CREATE TABLE test_intersect AS
SELECT ST_Split(a.wkb_geometry,b.wkb_geometry)
FROM test_line as a,
    test_line_vertex as b;

It returns successfully the rows with the line segments, but when I try to display it in QGIS it shows me that something is wrong. The spatial type is not defined and I can't view the lines in QGIS. So I assume something is wrong with my query.


texts: ERDAS mosaicPro from 2D view reference layer is not displayed on top

texts:
I am having an issue using the ERDAS mosaicPro from 2D view tool.  I am using ERDAS 2015.  I have two separate Landsat 8 file bands.  Band 1 from one date and Band 1 from another date.  I have ran the histogram match tool in (Raster \ Radiometric \ Histogram Match).  I want to "stack"/mosaic the two .tif files to fill in missing gaps.

Image one shows Band one from December 20th.  The black blobs used to be cloud/shadow coverage.  All black values have been set to "NoData".  Image two shows Band one from April 15th.  Same characteristics as image one.  I have set image one as the reference layer for the mosaic tool, but image 2 is being placed on top in the output.  I tried switching the reference layer, but I get the same result.  I have also tried using the Mosaic to Raster tool in Arcmap 10.4 and get the same result.  Is there something I need to change in the output properties?  Any suggestions?


texts: Which sensitivity analysis is most suitable on cost factors in route planner?

texts:
I am working on developing a route planner with alternative impedance/cost factor (safety of the route).

I am building the application mainly in pgRouting.

As a final step I would like to perform a sensitivity test. This sensitivity test should identify the influence of the different safety parameters (lighting, speed limit, etc) on the routing. The way to do this is to select random origin destination pairs (OD pairs) and see how varying values influence the routing. I did some research on the different types of SA and I cannot conclude which method best suits this particular problem. I've read about rationalized SA (Monte Carlo Filtering), Variance-based SA (Sobol method), Elementary Elements Test/EET (Morris Method), other One-at-a-time (OAT) methods, regression, etc. I am not sure which one is most applicable in this situation.

Can someone steer me in the right direction?


texts: town, village, prefecture shapefiles for China?

texts:
I'm looking for boundary files for towns, villages and prefectures in China that are free. I've came across http://chinadataonline.org/, but they are limited.

Can anyone tell me where to look?


texts: Can I work directly with spatial (PostGIS/SDE) DB with SQL queries?

texts:
There is some instruments like ArcPy for working with spatial DB. But can I do same work with simple (for example) PostgeSQL driver for Python?

Can I do select data column from DB A and insert data to DB B (only specified tables that include my feature class)? I have heard that PostGIS/SDE add some tables and if I will do it without touching PostGIS/SDE tables (that possible store some of indexes for my data) I will damage my data.


texts: Moving project, settings & plugins from Mac to Windows PC?

texts:
I am a complete novice with QGIS.

I have completed my first project in QGIS 2.18.9 on a Mac. I need to transfer it to another Windows computer.

How can I best do this while retaining all my settings and plugins etc.?


texts:
I'm trying to use the ST_Distance_Sphere function in PostGIS to work, but it's not behaving.

When I run this SQL command:

select id, name, location, ST_Distance_Sphere(location, ST_GeomFromText('POINT(51.51 -0.12)')) as dist_london from stations order by dist_london;

I get this error:

SQL Error [42883]: ERROR: function st_distance_sphere(geometry, geometry) does not exist
  Hint: No function matches the given name and argument types. You might need to add explicit type casts.
  Position: 28

In case it's relevant: location is a Geometry field. It contains only points.

I can use ST_Distance fine. But as I'm looking at spots very far north, it would be better to use something that outputs meters rather than degrees. What's going on here? Is there a PostGIS extension I need to add?

PostGIS_full_version output:

POSTGIS="3.0.0 r17983" [EXTENSION] PGSQL="120" GEOS="3.8.0-CAPI-1.13.1 " PROJ="6.2.1" LIBXML="2.9.10" LIBJSON="0.13.1" LIBPROTOBUF="1.3.2" WAGYU="0.4.3 (Internal)"

(I'm on Postgres 12.)


texts: Converting from SAFE (Sentinel-1) to BEAM-DIMAP?

texts:
I need to convert Sentinel-1 data from the SAFE format (used by the ESA's Sentinel project) to the BEAM-DIMAP format.

I know that GDAL can convert GeoTIFF files (which are used by Sentinel-1) into ENVI, but it adds an XML file containing some metadata the .img and .hdr files. From this link, I understand that the Coordinate Reference System should be in the header file which is not the case using GDAL.

I've also tried to use SNAP, the Sentinel software, to export SAFE to BEAM-DIMAP but there seems to be a loss of data as I can't retrieve the CRS.

Does someone know how to perform the conversion without any loss?

Ideally, I'd like something like SNAP where I give the whole SAFE directory as input and it outputs the .dim file and the .data directory of the BEAM-DIMAP format and doesn't require to do the conversion by hand.


texts:
I have a column (see below) where I want to label only part of the row.

Example:
16-23-1-01-7   -a  -00

I want to label only: 7    -a

Is it possible in QGIS?


texts: Euclidean distance outside a Polygon (ArcMap)

texts:
I got some university coursework to analyse the best potential areas for a residential development project and I'm currently doing multi-criteria evaluation. One of the shapefiles in the dataset given is high voltage power lines, I created a buffer around them of 200m (a recommended safe distance) and want to do a Euclidean Distance around that buffer closer being better because of cheaper land value.

Is there a way I can do Euclidean Distance going away from that buffer and not within it?

Should I standardise it, change it into points delete points within buffer then change it back? Seems a bit complex. Please say there's an easier way.

Please forgive me if this is simple my degree isn't GIS!


texts:
I have a set of point features which have been snapped to a waterways polyline feature class.

I need to select the waterways where they are within 6km upstream or downstream of any given point and create a new feature class from this selection. I do only have a basic licence on ArcMap 10.3.

I've looked a bit into routes and linear referencing but thus far haven't seen a way that work for what I want, and the waterways feature class isn't currently formatted as a route so any way I can avoid converting it would be ideal.


texts: Is there a standard way to show navigability of waterways?

texts:
Is there a standard (or common) way to indicate whether a river is navigable?  Even better, are there standard ways of showing various grades of navigability (accessible to small craft only, accessible to deep water vessels, etc.)?

Roads tend to be well marked to show their usefulness, from dirt trails to paved multilane highways.  I'm doing a project where it would be very handy to indicate the usability of rivers and creeks as well.  I could make up some indicators, but I'd prefer to follow existing models if possible.

I'm looking for a software-agnostic answer.


texts: Calculating proportion of polygon by value in QGIS?

texts:
Using QGIS, I have a single-band raster with with multiple counties.  I have reclassified the pixels into  2 categories, using a singleband pseudocolor, based on their value ( > or < 3800), corresponding to non-forest or forest land.

How do I calculate the proportion of pixels that is in each of the 2 categories, in each county?

I have a lot of images over a lot of years so I would prefer to not take the route the I figured out - converting to shp, exporting as csv, sorting csv and calculating percentages in Excel.


texts: In Quantum GIS, how do I specify the units for a shapefile?

texts:
I'm trying to create a shapefile of fields, drawing polygons around each field, and having the acreage calculated in the attribute table. I believe I have the projections correct. It's not coming up with accurate figures. How do I correct this?

Thank you!


texts: Is there any way to create user login for single project in QGIS? I want to track changes on database when users update features. I know that I can do this on database level but is there way to do this from project level?


texts: Convert Mapbox Vector Tiles to Longitude and Latitude

texts:
I am using Mapbox's QA Tiles and Python to access OSM data. While the tile-based approach is very handy I got stuck with how the geometry is saved. Apparently, coordinates are optimized for Mapbox's rendering needs and are 'vector tile grid coordinates'. To my surprise, I can't find any reference that would help me to understand what they actually did.

In my case, I want to check if a location (as latitude and longitude) are within an OSM area (Polygon or Multipolygon). This is what I got from the QA tile:

{
  "osm": {
    "version": 2,
    "features": [
      {
        "geometry": {
          "type": "Point",
          "coordinates": [
            95810,
            1041260
          ]
        },
        "type": 1,
        "properties": {
          "@uid": 144382,
          "backrest": "yes",
          "@id": 3072121792,
          "@timestamp": 1410453825,
          "@changeset": 25372130,
          "amenity": "bench",
          "@user": "SHuehne",
          "@version": 1,
          "@type": "node"
        },
        "id": 0
      },
      …
    ]
   "extent": 1048576
  }
}

This is the example node in OSM, the QA tile is z=12, x=2163, y=1435.


texts: How to export and import gps points from Garmin 62S and QGIS

texts: I'm pretty new to using QGIS and I'm curious about how to import and export GPS points created in QGIS. I have QGIS 2.0.1 Dufour on Windows XP and a Garmin 62s. Any suggested walk throughs would be great.


texts: Export from WMS Source to ECW file using QGis 2.18

texts:
I create a new WMS connection using QGis 2.18 from this URL:

http://wms.pcn.minambiente.it/ogc?map=/ms_ogc/WMS_v1.3/raster/ortofoto_colore_12.map

Now I want export to ECW file in according to specific coordinates.
I read some article about this but a the moment I'm not able to reproduce .ECW file output.

References:
Using the python console to prepare offline maps from WMS (tile)
Generate Map Tiles with QGIS

Any suggestions?


texts: GeoServer / OpenLayers : Front-End interface

texts:
I am using GeoServer as cartographic server.

I'd like to make a client side interface using OpenLayers, Ext, GeoExt.
That can select by attributes the layers published on GeoServer as WMS, and display the result on a map (OpenStreetMap).

Which IDE to use for this purpose: NetBeans or Eclipse?


texts:
I'm trying to get the QGIS OpenLayers plug-in to work from behind a proxy.

I've set up all my proxy details and this has had the effect of allowing QGIS's Manage & Install Plug-Ins tool to work with no problem. One oddity here is that although my log-in details are stored in Settings->Options->Network->Use Proxy for Web Access, QGIS asks me for my details again when I access Manage & Install Plug-Ins.

The problem is that when I try to add an OpenLayers layer, I just get blank tiles returned or in some cases white tiles with a blue box containing a question mark. I've spoken to our IT department who have monitored what happens when I try to add an OpenLayers layer and they say that the issue appears to be that my credentials aren't being passed to the proxy, it sees it as an anonymous request.

We've tried HttpProxy and DefaultProxy settings with no luck.

I was wondering if anyone had any suggestions how to resolve this?

Below is an example of some of the details that our proxy is returning:

Forefront TMG requires authorization to fulfill the request. Access to the Web Proxy filter is denied.

Client Username: anonymous

URL: http://tile.openstreetmap.org/2/3/2.png

Action: Denied Connection

UPDATE

As advised, I've installed fiddler2. This gives me some unusual results!

When running QGIS and fiddler2 together, with QGIS set up to use the same proxy settings as my browser, Fiddler2 doesn't 'see' any activity coming from QGIS. Therefore I change QGIS proxy settings to point at localhost, which I understand means Fiddler2 acts as a proxy.

Immediately Fiddler2 is able to identify QGIS traffic and OpenLayers begins to work partially. OpenStreetMap layers still fail, Stamen and Google maps load some tiles but not all (Google returns 'Sorry, we have no imagery here' messages for missing tiles).

I'm going to have a dig around in the Fiddler2 logs to see if I can spot what's causing the problem.


texts: Custom format for a grid with degrees, minutes, seconds with suffix using letter from another language

texts:
I am designing a map in French where West is "Ouest", and thus the map should read 45°23′67′′O instead of 45°23′67′′W. I am using a version of QGIS in English (update: in the French version of QGIS, the default option is O). I am currently using QGIS 2.18.11, but intend to switch to 3.0 soon.

1) Is there a way to change the longitude and latitude suffix to another language or just set a different letter?

2) Using the custom option, how can I get the minutes and the seconds?  I found an example here that gets me a decimal degrees version with my own suffix, but I would like to keep degrees, minutes, seconds format. If @grid_number is degrees... what is the equivalent for minutes and seconds?

abs(format_number( @grid_number ,2))
|| ' °'
|| CASE
    WHEN  @grid_axis = 'x'
       THEN
         IF (@grid_number > 0, 'E' , 'O')
       ELSE
           IF (@grid_number > 0, 'N' , 'S')
       END


texts: Using composer manager in QGIS 2.8.1 to export the current canvas as PNG?

texts:
I am a newbie to QGIS Python. I have a shape-file which I want to load to QGIS 2.8.1 and manipulate them by zooming to the layer of the preferred shape-file.

Once I have my canvas ready, I want to use that canvas in the composer manager and export it to a png file.

Below is my script:

import ogr,os
from PyQt4.QtGui import *
from PyQt4.QtCore import *
from qgis.core import *
import qgis.utils
import glob
import time
import math
import processing
from processing.core.Processing import Processing
from PyQt4.QtCore import QTimer
import sys

Processing.initialize()
Processing.updateAlgsList()
OutputFileName = "ABC"    # Temprory global placeholder for filename
canvas = iface.mapCanvas()

#Add List of address below for heatmap shapefile
heatmapfile = ["C:/Users/prashant.kumar/Desktop/QGIS python/Shapefiles/Map_00721.shp"]

# Add array of address below for  company shapefile
cspFile = ["C:/Users/prashant.kumar/Desktop/QGIS python/Shapefiles/HM_13_June_19.shp"]

for lop in range(filesLen):
    flag = 0
    currentShpFile = heatmapfile[lop]
    currentCompShpFile = cspFile[lop]

    #enter the common path for each shapefile here
    currentShpFileName = currentShpFile.strip("C:/Users/prashant.kumar/Desktop/QGIS python/Shapefiles/")
    #enter the common path for each shapefile here
    currentCompShpFileName = currentCompShpFile.strip("C:/Users/prashant.kumar/Desktop/QGIS python/Shapefiles/")

    # enter the output path here
    OutputFileName = "C:/Users/prashant.kumar/Desktop/QGIS python/ImageOut/" + currentShpFileName + ".png"

    # adding the shapefiles, changing the orders and zoom to layer
    wbc = QgsVectorLayer(currentCompShpFile, currentCompShpFileName, 'ogr')
    wb = QgsVectorLayer(currentShpFile, currentShpFileName, 'ogr')
    wbTex = QgsVectorLayer(TexLayer, 'CountyGrid', 'ogr')
    QgsMapLayerRegistry.instance().addMapLayer(wbc)    # Add the shapefile
    QgsMapLayerRegistry.instance().addMapLayer(wb)    # Add the shapefile
    QgsMapLayerRegistry.instance().addMapLayer(wbTex)    # Add the county shapefile
    qgis.utils.iface.setActiveLayer(wbc)   # Makes wb as active shapefile
    qgis.utils.iface.zoomToActiveLayer()    #Zoom to Layer
    print "Ready for composer manager"

    #Starting the canvas manager
    comp = iface.createNewComposer('currentcomposer')

Now I need help so I can use the current state of my canvas to be in the composer manager and export it as a PNG file.


texts: Map service query regarding old version of SQL Server client communications

texts:
I just upgraded to ArcGIS Desktop 10.4.1 on my ArcGIS Server machine. I could not connect via an SDE connection at first but I uninstalled the 2008 SQL Native Client and installed Microsoft® ODBC Driver 13 for SQL Server (we also just upgraded to SQL Server 2014 a few weeks ago) on our ArcGIS Server machine and the SDE connection worked fine after that.

However, I tried a query on our map services and get the error regarding the old client (see screenshot). I did the same thing on the web server machine (uninstalled 2008 and installed ODBC Driver 13), but it's still not working. Do you think I should just uninstall all references to 2008? I'm not sure why all those SQL Server 2008 program files are installed on our web server since we have a separate machine for SQL Server as well.


texts: Polygon Surface Area Calculator in different map projection(s)

texts: I'm simply looking for an online tool or a utility to calculate the surface area of the same polygon (i.e shapefile)  in different projection (projected and geographic projection system) systems. It would be great if you could share your suggestions.


texts: How to create a gradient color expression in data defined properties

texts:
This question is stemming from my other question on data defined attributes: How to create a gradient point on top of a shapefile?

Anyway, can anyone help me create an expression within the data defined properties -> fill color -> expression string builder. See following screenshot

My field is an integer value and I'd prefer light blue to dark blue gradient fill. Can anyone help with the input into the expression builder that would result in gradient filled centroid markers?


texts:
Does anyone know if it is possible to make a range label go over more than one line in MapInfo when doing thematic maps?

My legend has a limited width and I have not been able to find a way to get MapInfo to do this. The customize legend dialog doesn't seem to allow you to add text over two or more lines.

In the example below each item is wider than the frame on the layout that contains the legend. The frame is 7cm wide, and the legend items need 10cm when the text is on the same line.


texts: Is the ArcMap style database format available somewhere?

texts:
The *.style files are MS Access databases where each style type is a table.
The table has one binary column where the actual style is written.
Is this binary format available?
I would like to import some custom color ramps, but I can't see a way to do it... and I hoped in the style database!


texts:
I've created a map and I am using Openlayers Vector Tile Info example https://openlayers.org/en/latest/examples/vector-tile-info.html
to retrieve the point information from GeoServer.
My code is as follows:
// WFS Layer
layer_wfs = new ol.layer.Vector(
{
source: new ol.source.Vector(
{
format: new ol.format.GeoJSON(),
url: 'http://[URL as required]
}
),
title: 'WFS Layer',
visible: false,
style: new ol.style.Style(
    {
        image: new ol.style.Circle(
        {
            radius: 6,
            stroke: new ol.style.Stroke(
            {
                color: [255, 219, 77],
                width: 2
            }
            ),
            fill: new ol.style.Fill(
            {
                color: [255, 255, 0, 0.6]
            }
            )
        }
    )
    }
    )
}
);

// Create a view
myview = new ol.View(
{
center: ol.proj.fromLonLat([-6.67, 54.15], "EPSG:3857"),
zoom: 7,
projection: 'EPSG:3857'
}
);

// Create a map
map = new ol.Map(
{
target: 'map',
interactions: ol.interaction.defaults({ altShiftDragRotate:false, pinchRotate:false }),
layers:[layer_wfs],
view: myview,
controls:[
new ol.control.Zoom(),
new ol.control.ZoomSlider(),
new ol.control.ScaleLine(),
new ol.control.LayerSwitcher()
]
}
);

map.on('pointermove', showInfo);

var info = document.getElementById('info');
function showInfo(event) {
  var features = map.getFeaturesAtPixel(event.pixel);
  if (features.length == 0) {
    info.innerText = '';
    info.style.opacity = 0;
    return;
  }
  var properties = features[0].getProperties();
  info.innerText = JSON.stringify(properties, null, 2);
  info.style.opacity = 1;
}

This brings up the information perfectly, however it seems to be bringing up additional information which I'm not sure how to get rid off (or why it's there as it isn't stored in the attribute table). I only want the Address, Postcode, X and Y coordinate information to appear - how do I get rid of the rest?


texts: Applying categorized symbol to each feature using PyQGIS

texts:
I want to apply unique symbol to each feature. I have used following code for doing that
from qgis.core import *
from PyQt4.QtGui import *
# supply path to where is your qgis installed
QgsApplication.setPrefixPath("/path/to/qgis/installation", True)
# load providers
QgsApplication.initQgis()

# Get the active layer (must be a vector layer)
layer = qgis.utils.iface.activeLayer();
iter = layer.getFeatures();
for feature in iter:
    geom = feature.geometry()
    print "Feature ID %d: " % feature.id()
    qgis.utils.iface.mapCanvas().setSelectionColor( QColor("yellow") );
    layer.setSelectedFeatures([feature.id()])
    qgis.utils.iface.mapCanvas().zoomToSelected( layer )
    qgis.utils.iface.mapCanvas().refresh()
    selected_features = layer.selectedFeatures()
    for i in selected_features:
       attr =i.attributes()
    chaltano= str(attr[layer.fieldNameIndex('test')])
    renderer = QgsCategorizedSymbolRendererV2("test")
    layer.setRendererV2(renderer)
    symbol = QgsSymbolV2.defaultSymbol(layer.geometryType())
    symbol.setColor(QColor("red"))
    cat = QgsRendererCategoryV2(feature.id(), symbol,str(feature.id()))
    renderer.addCategory(cat)

Symbol is applying to only last feature.
Here, "test" is attribute.


texts: Replace Polygon geometry in SDE Database using FME

texts:
I'm still trying to wrap my head around using FME and encountering another issue.

This is my scenario:
I have a master SDE database that i am trying to update using the data from a GDB file. I was able to successfully update the SDE database on the Points data but when i tried to do the same process on Polygon's, the update breaks my SDE data.

FME Translation Setup:

I have 2 readers (GDB file reader and SDE DB reader)
I created 3 writers (INSERT, DELETE, UPDATE) *there could be a better implementation
I am just reusing these readers and writers and created an update process for each feature inside the GDB file.

Process:

Read through the SDE and GDB data.
Run CRC on each dataset.
Merge resulting datasets using FeatureMerger

Type: Attributes and Geometry

Run GeometryCoercer and 2DForcer.
Tester, check if both CRC fields are the same.

If yes, do nothing
If no, update the row.

The translation succeeds but when i try to view the feature in ArcCatalog, it gives me this error:

I might just be doing/setting something stupid that is causing this issue.

UPDATE: Also tested the updated geometries before the insert using the Data Inspector and everything seems to be fine. I think my issue is somewhere in the WRITER.


texts: ArcGIS incredibly slow when opening any toolbox tool

texts:
I've got a fresh install of ArcGIS 10.0, and whenever I try and open a tool from the toolbox (for example, Project Raster) it takes literally 5 minutes to open the dialog, and then at least 30 seconds for any clicks I make in the dialog (for example, to select input files) to register.

Does anyone have any idea why this is? The rest of ArcGIS isn't too slow (at least, no more than normal) but using toolbox tools is just impossible! This doesn't happen on my laptop, and this computer is significantly faster than my laptop (this has 12Gb RAM, quad-core processor etc).

Are there any known problems that might cause this?


texts: ST_Intersection point of two 3D linestrings - how is the Z value calculated?

texts:
I have two 3D linestrings (geometry type LINESTRING Z) in a PostGIS database that intersect in 2D space, but not in 3D space. Given that there is technically no 'shared portion' of geometry in 3D, how is the Z value returned by ST_Intersection calculated?

Here is the query I am using:

SELECT feature, str1, (ST_Dump(geom)).geom.ST_Z as z_value
FROM
(SELECT ST_Intersection(linework.geom, align.geom) as geom, linework.feature, linework.str1
FROM align
INNER JOIN linework ON ST_Intersects(align.geom, linework.geom)
WHERE ST_isvalid(align.geom)='t' AND ST_isvalid(linework.geom)='t') q1;

What I would like is a query that returns the z value of the point of intersection but on one of the linestrings only.

linework contains surveyed features that are crossed by an alignment in align - I would like to 'lift' elevations off the surveyed features at the points of intersection but I'm not sure how to go about it.


texts:
I am working on a project involving roads and traffic signs. I can retrieve Tag:highway=traffic_signals data which maps to individual signals present on roads.

How can I group related traffic signals which are operated in coordination with each other?


texts: How to create random points in overlapping polygons?

texts: I have a shapefile with ~1840 polygons, many of which overlap.  I'm trying to find a way to create 1000 random points per polygon, based on polygon ID.  I don't see a way to this using the GUI Random points or the Processing tool "Random points inside polygon (fixed)." Both tools will create the random points, however they are not linked to a polygon (meaning when I open the random points attribute table there is one column only with values 1:1.840mil).  I've thought of doing a spatial join after creating the random points, but I'm not sure how I would deal with all the points that fall within the overlapping area of the polygons.


texts: Find longitude for point between two coordinates with pyproj?

texts:
I've written a function to get the longitude of a point of known latitude between two known points (might not be working correctly yet), but I feel like there must be a better way of using pyproj to do this.

There isn't a built in method for doing this, but doing this binary search method I'm using seems possibly unnecessary as opposed to just solving for the other half of the coordinate directly.

Any suggestions? Ideally I'd like to generalize this to solve for latitude  when the longitude is known, too.

My function thing (maybe not working, will update once I'm sure):

import pyproj
import random

def unknownlon(lon1,lat1,lon2,lat2,lat_known):
    ''' Binary search by distance to find unknown longitude between points

        latitude and azimuth are known'''

    # Get positive for north movement and negative for south
    vert_change = (lat2 - lat1) / abs(lat2 - lat1)
    horz_change = (lon2 - lon1) / abs(lon2 - lon1)

    # Calculate angle and distance between known points
    fwd_azimuth,back_azimuth,dist = g.inv(lon1,lat1,lon2,lon2)

    # Random distance between known points
    dist_result = random.uniform(0.0, dist)

    # Loop while difference between desired lat and result is high
    lat_diff = 1.0
    while lat_diff > 0.00000001:
        lon_result,lat_result,back_azimuth = g.fwd(lon1,lat1,
                                                   fwd_azimuth,dist_result)

        # Calc difference between desired lat and result
        lat_diff = abs(lat_known - lat_result)

        # Determine whether to use a distance higher or lower than last
        if (lat_result * vert_change) > (lat_known * vert_change):
            dist_result = random.uniform(0.0, dist_result)
        else:
            dist_result = random.uniform(0.0, dist - dist_result)
            lon1 = lon_result
            lat1 = lat_result

    return lon_result, dist_result, lat_diff

if __name__ == '__main__':

    lat_known = 65.0

    lat1 = 68.0
    lon1 = 3.0
    lat2 = 63.0
    lon2 = -12.0

    g = pyproj.Geod(ellps='WGS84') # Use WGS84 ellipsoid

    # Calculate an unknown longitude
    lon_result, dist_result, lat_diff = unknownlon(lon1,lat1,
                                                   lon2,lat2,
                                                   lat_known)


texts: Extracting point values from a raster using .csv files

texts:
My research is about determining forage quality parameters using random forest regression in R. I have the GPS values for 5 different farms(FINCA) as well as the data quality values of Patia. I'm calculating the vegetation indices of a mosaic imagery I have of Patia to input as predictor variables for the grassland forage quantity and quality random forest regression models. All files can be found here:
DATASET
The problem I'm facing is in extracting raster values to point. My dataset is of PlanetScope imagery which has 4 bands: BGR and NIR. When I try to extract the points for the bands or the indices using the datosPatia.csv file, I get NA values for all 65 locations:
###########################################

setwd('D:/uni/project MSc/Patia/mosaic')
library(raster)
library(rgdal)
library(sf)
library(ggplot2)
library(dplyr)

#open your ndvi raster
r<- raster("NDVI.tif")
r

#read the csv file for sample points
pts<- read.csv('datosPatia.csv')
coordinates(pts)= ~ Longitude+ Latitude

#extract point values
ndvi= extract(r, pts)
ndvi

combinePointValue1=cbind(pts,ndvi)

write.table(combinePointValue1,file='combinedPointValue1.csv',
            append=FALSE, sep= ",", row.names = FALSE,
            col.names=TRUE)

########################################################

This didn't work so I used the kml files from the FINCA folders instead of the csv files:
#########################################

F1 <- st_read("FINCA01/FINCA01.kml")
plot(F1)

#open your raster
r<- raster("NDVI.tif")
r

#For example:
ndvi= extract(r, F1)
ndvi

combinePointValue1=cbind(F1,ndvi)

write.table(combinePointValue1,file='combinedPointValue1.csv',
            append=FALSE, sep= ",", row.names = FALSE,
            col.names=TRUE)

########################################################

The extract function gave this error but it produced a csv file with the point values.
Warning messages:
1: In proj4string(x) : CRS object has comment, which is lost in output
2: In .local(x, y, ...) :
  Transforming SpatialPoints to the CRS of the Raster

Can anyone tell me how I can put all these combinedpointvalues .csv files together as one file or if there is any other non repetitive way to get the values from all farms in one code because right now I'm writing the same code over and over for different farms.
Regarding the datosPatia. csv file, is it giving NA values Is it something to do with the structure of the .csv file and the column headers?


texts: shortest_path is not working on large source/target

texts:
I tried to find a route using shortest_path function of pgrouting, and it took to the server (postgres) a lot time to process  my query (more then 5 minutes!) even I have only 10 rows in my table.

After playing with my table I found out that the IDs of source & target must be short (up to 6-7 digits) in order to get the query results fast. But these IDs I get them from the map so I don't want to change them.

For Example, this is the table:

id   source     target   cost
--   -------    -------  -----
1    123456     2        1
2    3          2        2
3    3          4        3
4    4          5        4
5    5          6        3.3

This query takes ~30ms:

SELECT * FROM shortest_path('SELECT * FROM temp_table_2',123456, 3, false, false);

But the same query when I change the source 123456 to 1234567 I get the result after ~400ms! (and of course, If I try to increase the source to 12345678 it will take more than few seconds...)

How can I use pgrouting to get the shortest path (and not just wait...) from my table?

Thanks!


texts:
The DataPlotly plugin is not loading correctly on flatpak QGIS 3.12 for linux mint 19.3. I am running python 3.6.9 but the plugin seems to be seeking another version. Can anyone point me in the direction of a fix that does not involve this?
A pip install dataplotly terminal command is also yielding an error

Couldn't load plugin 'DataPlotly'

ModuleNotFoundError: No module named 'plotly'  Traceback (most recent call last):   File "/app/share/qgis/python/qgis/utils.py", line 312, in loadPlugin
    __import__(packageName)   File "/app/share/qgis/python/qgis/utils.py", line 743, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)   File "/home/eoin/.var/app/org.qgis.qgis/data/QGIS/QGIS3/profiles/default/python/plugins/DataPlotly/__init__.py", line 25, in
    from .data_plotly import DataPlotly   File "/app/share/qgis/python/qgis/utils.py", line 743, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)   File "/home/eoin/.var/app/org.qgis.qgis/data/QGIS/QGIS3/profiles/default/python/plugins/DataPlotly/data_plotly.py", line 31, in
    from DataPlotly.gui.dock import DataPlotlyDock   File "/app/share/qgis/python/qgis/utils.py", line 743, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)   File "/home/eoin/.var/app/org.qgis.qgis/data/QGIS/QGIS3/profiles/default/python/plugins/DataPlotly/gui/dock.py", line 14, in
    from DataPlotly.gui.plot_settings_widget import DataPlotlyPanelWidget   File "/app/share/qgis/python/qgis/utils.py", line 743, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)   File "/home/eoin/.var/app/org.qgis.qgis/data/QGIS/QGIS3/profiles/default/python/plugins/DataPlotly/gui/plot_settings_widget.py", line 75, in
    from DataPlotly.core.plot_factory import PlotFactory   File "/app/share/qgis/python/qgis/utils.py", line 743, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level)   File "/home/eoin/.var/app/org.qgis.qgis/data/QGIS/QGIS3/profiles/default/python/plugins/DataPlotly/core/plot_factory.py", line 14, in
    import plotly   File "/app/share/qgis/python/qgis/utils.py", line 743, in _import
    mod = _builtin_import(name, globals, locals, fromlist, level) ModuleNotFoundError: No module named 'plotly'

Python version: 3.7.4 (default, Nov 10 2011, 15:00:00) [GCC 9.2.0]  QGIS version: 3.12.3-București București, 82342615

Python Path: /app/share/qgis/python /home/eoin/.var/app/org.qgis.qgis/data/QGIS/QGIS3/profiles/default/python /home/eoin/.var/app/org.qgis.qgis/data/QGIS/QGIS3/profiles/default/python/plugins /app/share/qgis/python/plugins /usr/lib/python37.zip /usr/lib/python3.7 /usr/lib/python3.7/lib-dynload /usr/lib/python3.7/site-packages /app/lib/python3.7/site-packages /app/lib/python3.7/site-packages/numpy-1.14.3-py3.7-linux-x86_64.egg /app/lib/python3.7/site-packages/MarkupSafe-1.0-py3.7-linux-x86_64.egg


texts: Obtaining latitude and longitude on Earth from Satellite location and pointing angles

texts:
Considering the Earth as a perfect sphere and knowing the position (latitude,longitude and height) of a satellite.

I would like to transform the pointing angle from the satellite (defined by two angles) to latitude and longitude over the Earth surface.

I could try to compute intersections of lines and spheres to get it, but I think there are some packages in Python such as Pyproj that compute coordinate transformations (faster and more accurate), but I haven't managed to get something as I request.


texts: Set print composer to map canvas extent using python

texts:
I have a QGIS python script which successfully loads a specified print composer template.

At the end of this script I have attempted to set the composer to the map canvas extent, but I receive the following error message:

AttributeError: 'QgsComposerView' object has no attribute 'getComposerItemById'

I'm running Windows 10, QGIS 2.14.7

The code is as follows:

import os,sys
from qgis.core import*
from PyQt4.QtCore import*
from qgis.utils import*
from qgis.gui import *
import qgis
from PyQt4.QtXml import *

canvas = qgis.utils.iface.mapCanvas()
layers = canvas.layers()

myFile = os.path.join('P:\Office Admin\GIS\Templates\Plan templates\QGIS template plans\Title deed templates', 'A0_landscape_title deed.qpt')
myTemplateFile = file(myFile, 'rt')
myTemplateContent = myTemplateFile.read()
myTemplateFile.close()
myDocument = QDomDocument()
myDocument.setContent(myTemplateContent, False)
comp = iface.createNewComposer()
comp.composition().loadFromTemplate(myDocument)

map_item = comp.getComposerItemById('Map 1')
map_item.setMapCanvas(canvas)
map_item.zoomToExtent(canvas.extent())


texts: How can I mosaic images in ArcMap without it affecting the colour?

texts:
This is 4 images (NE, SE, NW, SW) sections of imagery. They have some overlap between them.

Working with 1m 3-band orthoimagery.
Different parts are sometime flown at different times/sources hence the difference in times of the year.
The images are in .tif format.
Imagery properties says it's 8 bit pixel Depth and the pixel Type is unsigned integer

When I mosaic them the colours change.

Looks the same as below
Number of bands = 3
Mosaic Operator = Last
Mosaic Colormap Mode = Last

Looks the same as below
Mosaic Operator = Last
Mosaic Colormap Mode = Match

Mosaic Operator = Last
Mosaic Colormap Mode = First

This image looks little better, features are better defined. Less Haze to.

How can I mosaic these images without it affecting the colour?
The brightness of the image has increased..

Is there alternative software to Arcmap that merges raster's better for free?

EDIT
Images are 8 Colour bit. Updated better Mosaic image. Still using Raster Mosaic still alters the image colour.

I used the Mosaic to New Raster tool.

Is there a way to remove the black parts?


texts: How can I convert SRTM data to DEM or DTED?

texts:
i am working on a project which is not commercial and i need some digital elevation data for visualization. in net i have found a link and i have downloaded .asc data from following link..

SRTM Data Link

with arcmap in toolbox, Conversion Tools > To Raster > ASCII to Raster , with floating options i can convert it to what i want. but when i export it as tiff and taking info with gdalinfo i am getting no coordinate system information and three band (RGB)...

in normal dted2 file with gdalinfo i am getting this result with other coordinate info:

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"]]

Pixel Size = (0.000277777777778,-0.000277777777778)

Image Structure Metadata:
  INTERLEAVE=BAND

Band 1 Block=3601x1 Type=Int16, ColorInterp=Gray
  NoData Value=-32767
  Unit Type: m

Can you help me how can i convert srtm asc data to dem or dted as below with arcmap or gdal? and is there any batch conversion for this?


texts: Assigning object in shapefile to location based on another shapefile?

texts:
I have two different shapefiles to work with: one of them has polygons representing tree canopy data and the other is a shapefile for census tracts in the area I am working on. I am trying to figure out what percent of each census tract's area is covered by tree canopy.

Is there a way that I can assign the polygons to the census tract they are in based on their location?

The attribute table for the polygons only has the object id, shape length, and shape area.


texts: Transforming a shapefile from EPSG:3338 to EPSG:4326?

texts:
I'm trying to transform this shapefile from EPSG:3338 (Alaska Albers) to EPSG:4326.

When I display the file in EPSG:3338 in QGIS it looks as you'd expect:

But when I display it in EPSG:4326 I see all kinds of problems:

And when I export it as EPSG:4326 I see the same problem. This happens both whether I export it in QGIS, or use ogr2ogr:

ogr2ogr -f GeoJSON -s_srs EPSG:3338 -t_srs EPSG:4326 .
  mv_oil_and_gas_basin_py.geojson mv_oil_and_gas_basin_py.shp

How can I transform this file to EPSG:4326 without getting problems around the international dateline?


texts:
I'm working on mobile app and I want to use OSM as its free solution. Whats the purpose of application?

User is geotargeted and the map like this: http://mc.bbbike.org/mc/ appears. No POI, buildings number etc. All I need is map with roads and buildings. Then if user is moving, we can track that and next tiles are loading when needed. Nothing fancy i think.

I came across few problems and items I don't understand. I know there are tiles (graphic) and .osm or .map file. I'm trying to understand how its working together.

What is the purpose of .osm file? From what I saw its just XML file with roads, POI, etc. If I need only streets and their name can do I even need this? Maybe it could be served by tiles?

Can someone explain me how tiles and xml file are working with each other?


texts:
I've been running a code that searches a desired field in a feature class and performs an specified SQL query based on whether the field is string or double.  The variables are just temporary fillers to test the code, I'll be shifting them to GetParameterasText() as soon as I know that all of the code works.  I need the code to work so that it can be turned into a ArcGIS toolbox script.

import arcpy
from arcpy import env
env.overwriteOutput = True

#Create working variables
fc = "C:/gisclass/assignment5/Results/CITIES.shp" #vector feature class
affectedField = "CAPITAL"
oldValue = -1
newValue = 2
integer = 0

fieldlist = arcpy.ListFields(fc)

for field in fieldlist:
    if field.type == "String":
        queryString = """"%s" LIKE '%s'"""%(affectedField, oldValue) or """%s" = '%s'"""%(affectedField, oldValue)
    elif field.type == "Double":
        queryString = """"%d" >= '%d'"""%(oldValue, integer) or """"%d" <= '%d'"""%(oldValue, integer)

with arcpy.da.UpdateCursor(fc, (affectedField,), queryString) as cursor:
    for row in cursor:
        row[0] = newValue
        cursor.updateRow(row)

mess = arcpy.GetMessageCount()
print arcpy.GetMessage(0)
print arcpy.GetMessage(mess-1)

I have gotten the string where clause to work, but the integer clause is still producing the same result - Invalid SQL query.  Any idea on what the solution to the problem is?


texts: GDAL merge of DEM data resulting in large black tiles, elevation min/max problem?

texts:
I am new to QGIS (former Arc user) and having trouble with the raster merge tool.

I am merging 7 tiles of SRTM arc 1 data to cover a province in eastern Ecuador. Data is resulting in large areas of total black areas. I've done quite a bit of research on GIS stack exchange to troubleshoot:
- I've made sure my output data type is consisent with input (int16)
- I've tried merging with both projected and unprojected data
- I've tried some general tinkering around with the style settings on each DEM
- I've tried the tool in both QGIS 2.18 and 3.0

At this point, I'm fairly positive it has to do with the large variations in my minimum and maximum values in each tile -- I am working with an area that goes from high elevation of the Sierra mountain range down into the Amazon basin, so my real range is between 127 and 6218. My output comes out 0-6218, which is close enough to the real range (127-6218) that adjustment of the min/max of the output DEM to reflect real range doesn't do much. Extremely changing the min/max values in the final output DEM allows me to see data where there is black, but simultaneously "erases" data at the other end of the spectrum (If the high elevation data appears, the low elevation data disappears, and vice versa).

I do not see any parameter to adjust in the merge tool to fix this. Anyone have any ideas?? Below are some screenshots of original data, merge parameters, and output in order:


texts: ArcGIS API for Javascript QueryTask not working

texts:
I'm trying to use query tasks to extract attributes of a feature based on what was selected from a filtering select list(setVal = "Chain Lakes Provincial Recreation Area"). But I get the error

TypeError: Cannot read property 'length' of undefined(…)

from the line of code

for (var i = 0, il = resultFeatures.length; i < il; i++) {

in the following:

        queryTask = new       QueryTask("http://www.arcgis.com/sharing/rest/content/items/1d890f37f9f049d1b499141511235b5b/data");
        //initialize query
        query = new Query();
        query.returnGeometry = true;
        query.outFields = ["OC_NAME", "TYPE", "ACRES"];

        //create symbol for selected features
        symbol = new SimpleMarkerSymbol();
        symbol.setStyle(SimpleMarkerSymbol.STYLE_SQUARE);
        symbol.setSize(10);
        symbol.setColor(new Color([255, 255, 0, 0.5]));

        function showResults(featureSet) {
            //remove all graphics on the maps graphics layer
            map.graphics.clear();

            //Performance enhancer - assign featureSet array to a single variable.
            var resultFeatures = featureSet.features;
            //Loop through each feature returned
            for (var i = 0, il = resultFeatures.length; i < il; i++) {
                //Get the current feature from the featureSet.
                //Feature is a graphic
                var graphic = resultFeatures[i];
                graphic.setSymbol(symbol);

                //Set the infoTemplate.
                //graphic.setInfoTemplate(infoTemplate);

                //Add graphic to the map graphics layer.
                map.graphics.add(graphic);
            }
        };

        function executeQueryTask() {
            //set query based on what user typed in for population;
            query.where = "OC_NAME = '" + selVal + "'";
            console.log(query);
            //execute query
            queryTask.execute(query, showResults);
        };


texts: Coordinate System changes after reprojection

texts:
I have changed the projections of some Natural Earth Shapefile to "WGS 84/Pseduo Mercator", that works fine (after cutting out the artic), but now the co ordinates of the shapefile are all out of wack.

I can display the shapefile(s) in my application (using MapWinGIS), but the lat longs that i'm passing for points etc are in decimal degrees and the shape file seems to be in something else - meters maybe?

What can i do to keep the shapefile in degrees, or is there a conversation between to two.


texts: Remove/Delete a layer or layergroup with mouse click

texts:
I add a ol.layer.Vector OR ol.layer.Group with a mouse click in the map. By a new mouse click should the layer removed/deleted AND a new layer should be add. etc.

I have search and search and don't found a simple methode. There is all complicated or don't work, or with Ajax or JQuery (i don't use it).

Then i use this and that is the best solution:

if (map.getLayers().array_.length > 1) {
map.getLayers().array_.pop();
}

That works absolutelly perfect and is absolutelly simple. OpenLayers does not seem to have such a simple solution.


texts:
I'm creating a flexviewer application (arcgis API for flex 3.4) with application builder. I'm new with flexviewer. I want to add thematic widget. I was able to add the widget. But when I want to run the widget, an error occurred saying I need to have dynamic layer. I went googling all day..but I could not find a way to do it in arcgis 10.0. Most of them are guides on 10.1.

So, is it possible to enable dynamic layer using arcgis server manager 10 or arccatalog 10?
If it is possible, can someone direct me how? maybe some links?

Thanks


texts:
I've hit a wall trying to loop through multiple fields with arcpy.UpdateCursor.  I have a list of files (2) with multiple fields (10-15'ish) that have its own error check in SQL, and an error value to be written.  I am wanting to loop through and update each field with its programmed value.

    import arcpy
    FILE = 'X:\\SOME_FILE'
    FIELDS = [ "FIELD_X", "...", "FIELD_N"]
    SQL = "SOME_FIELD = SOME_VALUE" # or less than, greater, etc
    for FIELD in FIELDS:
        rows = arcpy.UpdateCursor(FILE, SQL, "", FIELD)
        row.FIELD = 'FIELD_ERROR_VALUE'
        rows.updateRow(row)
    del rows

There's a bit more to it, but the short is that row.FIELD is read as row.'FIELD_NAME', with FIELD_NAME as a string value.  I'm rather new to the backside of Python programming and got thrown in just to learn.  I've worked through many issues, but this one has got me stumped.  If row.FIELD would work, it would be easier than going through each single field with it's NAME, SQL, UPDATE_VALUE, etc manually.  It works if I type in row.WHATEVER_FIELD = ..., but calling the variable (from string) adds a ' to both ends as a string value.  The only thing not working in the loop is row.FIELD.

I'm working with python2.6, arc10.0, mostly win7.

Any help would be greatly appreciated.

GM


texts: How to load JQueryUI into Web AppBuilder?

texts:
There is an example widget using JQuery in {WAB root}/client/stemapp/widgets/samplewidgets/UseJQuery and the code can be easily adapted for other versions of JQuery. However, I want to use JQueryUI as well. I tried to duplicate jimu/loaderplugins/jquery-loader within the require array in the define statement:

'jimu/loaderplugins/jquery-loader!http://code.jquery.com/jquery-1.9.1.min.js',
'jimu/loaderplugins/jquery-loader!http://code.jquery.com/ui/1.10.3/jquery-ui.min.js'

But since doing this, the widget doesn't open in webappviewer and its code is no longer visible to FireBug, so this clearly is not the way. Seems I need some jqueryUI-loader, not existing in the standard jimu/loaderplugins directory. Is there any way to enable JQueryUI on the widget level, or do I need to configure init.js, hack the WAB code or even write my own jqueryUI-loader?

EDIT: I found an example of including the library in the head of Widget.html file. AFAIK widget templates should be short snippets not containing the head section so this seems too hacky and I would expect some trouble after doing this. Is this the only way to do it on the widget level?

I plan to publish the widgets to more apps, including those owned by those outside our organization, so tweaking the widget is preferred to doing it anywhere else. Changing WAB is not completely out of question though, if it's cleaner than this one.


texts: Snapping "Point to Line" Results to Correct Streets in ArcGIS for Desktop?

texts:
Here is the problem I am facing. After geocoding some addresses, I used the "point to line" tool. Some long streets with minimal matching address points created long lines that intersect the road layer. What I need is for the dark blue line to "snap" as close as possible to the light blue road layer.

What I've tried: Densifying the dark blue line and then edge-snapping to the roads layer. The line then snaps to the wrong road sometimes. The vertices of the dark blue line are the points that matched to the light blue line.

Is there anyway I can give the "points to line" result the same street name and then snap to the matching street?

Using ArcMap Desktop Advanced 10.2.


texts: Finding xy-coordinates inside a polygon, which are the closest to the centroid of that specific polygon

texts:
I have the following problem:

I have zip code regions, which basically have one single value attached to them most of the time. Then I have a grid consisting of grid boxes, which you can put over these zip code regions of interest.

What I want is the following: For every zip code region of interest, I want to attach the value from the zip code region to the grid coordinate, which is the closest to the centroid (see first figure below) of that zip code region. That grid box must be inside of that zip code region. It's enough actually, when I know which grid box in the zip code area is actually the closest one to the centroid of the zip code region. I can use these coordinates to index the grid boxes of relevance to me from the pixels plotted in figure 2 (see figure 2, where I didnt select the grid boxes closest to the centroids of a zip code region [figure 1]).

trueCentroids.zip_code_regions <- gCentroid(SpatialPolyDf.zip_code, byid=TRUE)  % Gives me the centroids of my zip code regions

Initially, I didnt care and just gave every grid box, which was inside of a specific zip code region the value, which is attached to this region (In the case a value existed for a specific region). The following figure and code should make that clearer:

> class(Spatialpixels_grid)
[1] "SpatialPixels"
attr(,"package")
[1] "sp"

Spatialpixels_grid basically just contains coordinates, which I attach values to from the SpatialPolygonsDataframe object SpatialPolyDf.zip_code (see figure 2). For some regions, there are no values, so you see white pixels. The points represent the grid boxes, which are saved in Spatialpixels_grid.

Spatialpixels_grid$value_zip_code_region <- over(Spatialpixels_grid, SpatialPolyDf.zip_code)$value

Example for using the function over

Is there the possibility to modify over, so I can achieve my goal?


texts: Entity Framework & Sql Server 2012 spatial

texts:
I'm struggling to find any concrete performance of Entity Framework paired with SQL server 2012.

We're looking at storing nearly a million points in a database an right now the most common scenario would be to find all points within radius of N.

I've mocked up a small prototype using EF5 and it seems fine. However, a senior developer has mentioned that in his opinion EF is bloated and I should consider using a micro-orm or even pure ADO.NET as EF is a bloated monster (lol).

Performance is crucial on this project and I would like to know if anyone has experience to prove this guy otherwise (or right for that matter)?


texts: Different timevalues in tool make service area ModelBuilder

texts:
In ModelBuilder i want to make a lot of service area layers.
The input data would be the same but the day of time would change every time.
So from 6:00:00 untill 22:00:00
i created already a table with the different times (in format date) (see figure 2)
and add this table into the model (see figure 3)

now when i add %Value% in the field day of time (see figure 4) and click next on the  button ok.. nothing happend, the field of day of the time is again empty..

So how could i run the tool 'make service layer' for different time steps?

And how does it works to change the output name automatic?


texts:
I have a GeoServer running, I have made a process that simply returns the requested tiff.
Now here's the deal : when I send my request from the integrated WPS demo in GeoServer, the GUI gives me a list of available maps to choose from.
However, when I try the same thing with my (working) HSLayers client (a simple HTML page with an embedded js WPS client), the GUI asks me for an URL and I have no list of maps.
I would prefer having them shown.

Additionnaly, if I take any random process, the GUI will work fine. If I switch between processes, everything works fine (updating input fields, etc.). But once I switch to my AddGeoTIFF process, the GUI generator bugs and doesn't switch anymore (I have to reload the page).

Is this due to a small error in the HSLayers client or is it a problem in the prototyping of my process ?

Here's my process' code :

@DescribeProcess(title = "AddGeoTIFF", description = "Returns a GeoTIFF\n")
public class AddGeoTIFF implements GSProcess {

private static final CoverageProcessor PROCESSOR = CoverageProcessor.getInstance();

@DescribeResult(name = "result", description = "The resulting GeoTIFF")
public GridCoverage2D execute(
        @DescribeParameter(name = "coverage", description = "The GeoTIFF you wish to obtain") GridCoverage2D coverageA,
        ProgressListener progressListener) throws ProcessException {
            return (GridCoverage2D) coverageA;
    }
}

I am not good with Java and based it on an existing process, so there might be some horrors.

I have been looking for the cause of the issue, and it is (what I believe to be) a parsing error. The (on the fly) GUI generator's part which troubles me can be found in [HSLayers install file]\source\HSLayers\WPSClient\complexField.js, on lines 77 to 82. If someone experienced could have a look at it, It'd help.

Thanks


texts: Selecting features from another shapefile using ArcPy cursors and SelectLayerByLocation?

texts:
Having the code below im trying to select all the table features from an school_enroll.shp that are within each bigger zones from the hs_bndy.shpwhen i use the arcpy.SelectLayerByLocation_management for some reason the selection isn't happening. Am i doing something wrong here? I like to be able to sum some of the fields from the selected school_enroll.shp  features

#hs
fields = {i:str(field.baseName) for i,field in enumerate(input)}
input = arcpy.ListFields("HS_Bndy1516.shp")
school_ampa = "C:/school_enroll.shp"
rows = arcpy.UpdateCursor("HS_Bndy1516.shp")

for row in rows:
    total = 0
    row_name =  row.getValue("NAME")
    sql = "{} = '{}'".format(fields[3], row_name)
    arcpy.MakeFeatureLayer_management("HS_Bndy1516.shp", 'school_lyr', sql)
    arcpy.SelectLayerByLocation_management('school_lyr',"HAVE_THEIR_CENTER_IN",school_ampa)
    with arcpy.da.SearchCursor('school_lyr',"HSadd")as cursor:
        for row in cursor:
            total = row[0] + total


texts:
I would like to be able to change or remove the tooltip when the mouse is over a layer in TOC.
I've seen in source code that i must to use the QgsLayerTreeViewIndicator class to perform it, but it doesn't work (crash).
layerNodes = QgsProject.instance().layerTreeRoot().findLayers()
ltw = self.iface.layerTreeView()
for node in layerNodes :
    if (node.layer() is not None) :
        indicator = QgsLayerTreeViewIndicator(node.layer())
        ltw.addIndicator(node, indicator)

If is it possible to change the default publicSource method to display a better information, i take as solution too.


texts: Loading huge amount of data in OpenLayers 3

texts:
I have been working on a PoC to load a huge amount of data ~ 350k-500k but I haven't been successful so far. I tried loading JSON then converting them to ol3 features (failed at 130k) then I used GeoJSON where I could load 150k data but zoom-in , zoom-out makes the map crash. Also, I have been using Clusters. How do I get around the problem? I can't use WMS since my only source of data is a JSON whose format can be modified.

INFO: The GeoJSON I prepared is roughly 40MB. I was thinking of using gzip compression. Would it help me in anyway?


texts: Copy shapefile to SDE feature class using arcpy

texts: Is it possible to copy shapefiles to SDE geodatabases using ArcPy? I have tried using Copy_management, but it says that it cannot copy between different workspaces. Is there a different function that should be used, or perhaps a workaround?


texts: QGiS: Problems converting vector to raster to asc wrt Maxent (java)

texts:
I'm trying to create a "bias file" to feed into Maxent (a java program). To do this, I use QGiS to turn a vector shapefile into a raster with the "heatmap kernel density estimation" package, and then turn that into an asc. However, Maxent gives me the following error when I try to give it the asc:

java.lang.NumberFormatException: For input string: "NODATA_value"

I've tried checking for this problem, but my asc's NODATA_value seems to pretty much match the one that's initialized:
I've already tried adding/removing the decimals from both the initial NODATA_value and the one in the top left.
My prof thinks it might be something that popped up when QGiS converted the shapefile to a raster. However, I don't know how or why that might be.
Thoughts?


texts: Snap mean coordinate to layer in QGIS and creating MCs for multiple shapefiles

texts:
1: I wish to create a mean coordinate for a group of points. However, I would like the new point to be snapped to the nearest point from the set of points used to generate it.

How might one go about this, if at all possible?

2: I have come to the stage where I have many shapefiles I wish to create mean coordinates for. Is it possible to make QGIS generate these all in one go? As opposed to manually generating the coordinate through inputting one shapefile at a time.


texts:
I'm trying to add a NOAA hosted layer to an OpenLayers map.  I've seen some examples of doing it with an older (pre 3) version of OpenLayers and of services no longer offered by NOAA.  Has anyone tried to present any of the current NOAA services (such as the NEXRAD [ https://nowcoast.noaa.gov/arcgis/rest/services/nowcoast/radar_meteo_imagery_nexrad_time/MapServer/0 ] image) to a map with the current version of OpenLayers (5.2.0 in my case)?

Below is the code I was attempting to use, but it is definitely wrong as it ends up trying to retrieve the HTML instruction page instead of an image.

var source = new ol.source.ImageWMS({
    attributions: ['NOAA'],
    url: 'https://nowcoast.noaa.gov/arcgis/rest/services/nowcoast/radar_meteo_imagery_nexrad_time/MapServer/0',
    projection: 'EPSG:3857'
});

var layer = new ol.layer.Image({
    title: 'NOAA Radar',
    zIndex: 1,
    visible: true,
    source: source,
    opacity: 0.7
});

map.addLayer(layer);


texts: Divide polygons based on internal differences in shape using PyQGIS

texts:
For my research, I'm trying to assess where we could potentially plant small, medium and large trees (with respectively 5m², 20m² and 45m² soil area needed). Besides patch size, also patch shape is an important factor. It could be that we have a patch of 20m² (=large enough for 1 medium tree), but it has a complete strip-like shape (= dimensions 1m x 20m) and is thus unsuitable for a medium tree because we assume that the roots won't be able to reach the end of the patch. Therefore, a patch may have a maximal value of 1.2 for the patch fractal dimension index (this is a value between 1 and 2, where 1 is a very simple shape like a square and 2 a very complex shape like a line; ).
The goal is to do this in an urban territory and as you can imagine, plant-able patches have very complex shapes. An example of such a shape is shown in the figure below.

The main issue is that this polygon provides space for 1-2 medium size trees in the green circle, and for some small trees in the black circle. However, I'm struggling in developing an algorithm that divides the polygon automatically into (in this case) two patches: 1 for the medium trees and 1 for the small trees...
I work in QGIS 3, so I'm rather looking for an answer that can be executed by PyQGIS


texts:
I am following the ArcGis10.0 desktop help to bring cad data into ArcMap but I'm having trouble accessing the Autocad layers set up in the Autocad file. I am:

using the add data to reference in the cad file - in the ArcMap layers 4 layers appear (points, polyline, polygon and multipatch).
I select the polyline layer (as the drawing only contains polylines), right click and select 'convert cad feature layer'.
This brings in the linework but all three layers in the Acad drawing are collapsed into 1 layer in ArcMap including layers which were off in Acad.

What I really want is to have my Acad layers available in Arcmap - I'm more familiar with Acad than ArcMap can you help?


texts: Overpass not returning data within my given bounding box

texts:
I have been using this endpoint to query data from the highway tags. But I am receiving data way outside the bbox? Is there something wrong with my query?
http://overpass-api.de/api/interpreter?data=(node["highway"](34.21092368471588,-118.6453487197585,34.212720315284116,-118.64752128024149););out center;