texts:
I am doing hydrological analysis for a river network in a DTM.
I have produced a graph like (b).
Now, I want to take it one step further and plot a log slope-area graph like graph (a).

How do I calculate slope and area for the catchment?
I assume the (+) on the graph show slope of diferent rivers?  If so, how do I calculate the slope from the DTM?

I am working with ArcGIS Pro latest version and MATLAB.


texts: Overview window in basemap viewer in ArcMap hides toolbar and can't be moved

texts: When adding a basemap to ArcMap project an overview window appears that hides the toolbar and can't be moved. How can this window be adjusted or removed?


texts: Is there an opensource software to route openstreetMap data based on database?

texts:
I have a project that consists in aggregate traffic information on OpenStreetMap for routing proposes. My first idea is to use the OSM data from the region on a database, so I can aggregate the traffic information on the table, probably adjusting the cost column.

The problem is that all the routing software that I have already downloaded and tested creates a bin or obj file from the OSM instead of working with databases, so does anyone know an opensource software that creates and draw routes on OpenStreemaps based on database?


texts:
SOLUTION: For everyone interested: turned out the solution was moving all files that use DoJo beneath the jQuery UI files, so that jQuery UI is fully loaded before any DoJo stuff loads. Moved them all in the footer in that order (jQuery, then jQuery UI and then alle DoJo related files).

Edit: to keep things organized:

This one IS working (WordPress 4.0.1): https://topografieindeklas.nl/topotrainer/provincies/
This one IS NOT working (WordPress 4.1): https://tpgrf.nl/testserver/alpha/topotrainer/provincies/

I'm running a Wordpress website which uses both jQuery and DoJo (part of the ArcGIS JS API, loaded from http://js.arcgis.com/3/12/).
After recently upgrading from Wordpress 4.0.1 to 4.1, things broke.
Dojo is returning multiple Error: multipleDefine errors.

Basically it looks like DoJo tries to load a new instance of jQuery, which it shouldn't, since WordPress has loaded jQuery already and there isn't a jQuery module available in the DoJo map. Furthermore, this issues seems to only exists in WordPress, not in a jsfiddle I made including exactly the same scripts.

For debugging purposes I first disabled all WordPress plugins, without any result; the problem kept existing. I went on and found a couple of interesting things. But I'm really stuck in the further debugging proccess. I know the scope of this question is rather broad.

What I know so far/that i'm wondering about:

-Things stopped working after upgrading to WordPress 4.1

-In the info of the multipleDefine error it states https://js.arcgis.com/3.11/jquery.js, which is weird: it shouldn't be looking for jQuery there, since it's already included in Wordpress

-Things break as soon as I included the ArcGIS Javascript API, which includes DoJO, but only in WordPress; on a jsfiddle things work fine http://jsfiddle.net/dhunink/mry5vn3s/7/

<script>
  var dojoConfig = {
    isDebug: true,
    async: false
  };
</script>

<script src="http://js.arcgis.com/3.12/"></script>

The things mentioned above leads me to the conclusion that something changed in the last WordPress release that effects my script, although I can't find anything in the release notes that could be linked to my issue.

Feel free to see all in action at https://tpgrf.nl/testserver/alpha/topotrainer/provincies/.


texts:
I wondered if you could share some of your experience with plot files in ArcGIS. I've been trying to make these, as I have a good expeirence with them from autoCAD (every drawing is exported to a PLT file, which stands by if another copy is wanted later on).

In arcGIS, I saw there is only options to export to PRN or RTL, but no PLT (which is a HP designed file, if I'm not mistaking).

Both PRN and RTL save the maps (containing vector and raster features) horridly, making the maps produce from them very discolored and messy.

From what I've deduced from the Internet, PRN is usally for documents, so it's probably not meant for rasters, hence the poor capabilities.

(for all those who want to ask: Saving only as a MXD isn't much good, as the layers can be changed over time, moved to different locations, and anyway, handling MXD can be done only by GIS trained personnal while printing plot files can be done by any given secretary).

Any thoughts? Has something in this field advanced in ArcGIS 10 (I'm in 9.3.1)??


texts: GeoServer integration with GeoWebCache. Empty miss reason

texts:
I enabled the direct integration with GeoWebCache within GeoServer for a set of static maps. But currently all the requests are returning with the Header geowebcache-cache-result: MISS and the geowebcache-miss-reason: with and empty string.

I understand that for the request to be served with the cache it needs to fulfill all the requirements:

WMS Direct integration is enabled (you can set this on the Caching defaults page)
tiled=true is included in the request
The request only references a single layer
Caching is enabled for that layer
The image requested is of the same height and width as the size saved in the layer configuration
The requested CRS matches one of the available tile layer gridsets
The image requested lines up with the existing grid bounds
A parameter is included for which there is a corresponding Parameter Filter

The requests when sent trough the demo site in the GeoWebCache configuration return a geowebcache-cache-result: HIT so the cache is working. But its sending the requests directly to the GeoWebCache

The requests that I make trough the WMS integration have the following parameters

https://GEOSERVERHOST.com/geoserver/wms?SERVICE=WMS
&VERSION=1.3.0
&REQUEST=GetMap
&FORMAT=image/png
&TRANSPARENT=true
&LAYERS=test:flood_riverbed
&TILED=true
&WIDTH=256
&HEIGHT=256
&CRS=EPSG:3857
&STYLES=
&BBOX=-215246.6716510579%2C4520180.104672182%2C-195678.7924100528%2C4539747.983913188

The Gridset available for the layers is the default EPSG:900913

I cannot find where the issue is, more with the "miss-reason" beeing empty. Any ideas?


texts: QGIS algorithms do not work on copies of shapefiles?

texts:
I'm using QGIS 2.6.0. I'm trying to run various algorithms such as "SAGA/Shapes-points/remove duplicates" & "SAGA/Shapes-Grid/Add grid values to points" through the sextant toolbox.

When I run these algorithms on a shapefile (initially created via the QGIS option "Create a Layer from a Delimited Text File" and then saving the output as a shapefile), it works fine. However, if I create a copy of the shapefile through QGIS (i.e. using Save As) even with no changes (e.g. same coordinate system) and try to run these algorithms using the same input parameters, it fails. An example of the error that appears is as follows:

"Oooops! The following output layers could not be open
Result: C:\Users\stephen\AppData\Local\Temp\processing\bf494593f4a44905bb4ec0bbe0a7f847\RESULT.shp
The above files could not be opened, which probably indicates that they were not correctly produced by the executed algorithm
Checking the log information might help you see why those layers were not created as expected
This algorithm requires SAGA to be run. A test to check if SAGA is correctly installed and configured in your system has been performed, with the following result:
SAGA seems to be correctly installed and configured"

The log output is as follows:

"SAGA execution console output

C:\PROGRA~1\QGISBR~1\bin>set SAGA=C:/PROGRA~1/QGISBR~1/apps\saga

C:\PROGRA~1\QGISBR~1\bin>set SAGA_MLB=C:/PROGRA~1/QGISBR~1/apps\saga\modules

C:\PROGRA~1\QGISBR~1\bin>PATH=PATH;C:/PROGRA~1/QGISBR~1/apps\saga;C:/PROGRA~1/QGISBR~1/apps\saga\modules

C:\PROGRA~1\QGISBR~1\bin>saga_cmd shapes_points "Remove Duplicate Points" -POINTS "E:/Studies/Data/Consolidated/Records_06092014_Copy.shp" -FIELD "Precision" -METHOD 0 -NUMERIC 1 -RESULT "C:\Users\stephen\AppData\Local\Temp\processing\bf494593f4a44905bb4ec0bbe0a7f847\RESULT.shp"
_____________________________________________
#####   ##   #####    ##
###     ###  ##       ###
###   # ## ##  #### # ##
### ##### ##    # #####
##### #   ##  ##### #   ##
_____________________________________________

library path:   C:/PROGRA~1/QGISBR~1/apps\saga\modules
library name:   shapes_points
module name :   Remove Duplicate Points
author      :   O.Conrad (c) 2008
_____________________________________________
go...

Load shapes: E:/Studies/Data/Consolidated/Records_06092014_Copy.shp...

ready

okay

Remove Duplicate Points: Invalid parameters!

Shapes: Points

error: executing module [Remove Duplicate Points]

C:\PROGRA~1\QGISBR~1\bin>exit"

Just wondering whether anyone has come across a similar problem and if they have found a solution?


texts: Editing the configuration of a join in QGIS

texts:
Is there any way to edit the conditions of an existing join in QGIS? I have a .csv joined to a (polygon) .shp layer; because I only needed a few fields from the .csv, I only joined those I needed.

Now, I'm realizing I need a few more fields--is there any way to edit the current join instead of starting from scratch?

I could also see this being relevant if I decided to change the field by which the tables are joined.

This should be simple but I don't see any options for it.


texts: Rasterize vector file using gdal.RasterizeLayer() without losing attribute information in python

texts:
Is it possible to rasterize a vector file with gdal.RasterizeLayer(), saving attribute as Raster value ? Like gdal_rasterize can do with :

gdal_rasterize -a attribute -ts rasterXsize rasterYsize -l layername shapePath rasterPath

As suggested here, we can rasterize the attribute in color levels (0-255) on the 3 bands, but my objective is to keep the attribute value intact on a single raster band. I found many other examples here but None fit my goal.


texts: Create tracks from GPS points honoring polygon barriers

texts: I have fish tracking data from a river. The data include positions, times and tag ids. Time interval between positions may be several days.  I would like to create a polyline track for each fish tag that stays within a polygon that defines the shoreline of the river.  I am looking for a procedure or extension or add-in to use in ArcGIS or one of several open source GIS.


texts:
I can't open this TIFF data in GIS (can´t see the Attribute Table or convert them to polygon)

https://gitlab.com/hotmaps/vol_res_curr_density/blob/master/data/vol_res_curr_density.tif


texts: ST_Buffer in PostGIS produces different results for the same set of lines

texts:
We have a single LineString (a route on a map) made of multiple points (blue line) which we are then buffering to produce an area around it. We have found that everything looks as expected unless the route intersects itself, at which point the buffered area then gains considerable longitudinal curvature.

An example shows the differing results: http://geojson.io/#id=gist:jgwconsulting/1e2a6e8bad9f018f2c6321016a527bef&map=6/55.621/-3.618

Grey: entire route is buffered as a single LineString, resulting in curvature of the buffered area.
SELECT ST_AsGeoJSON(
ST_Buffer(ST_GeomFromGeoJSON('{"type":"LineString","coordinates":[[-14.466666666666667,55.25],[7.233333333333333,54],[0,58.86666666666667],[-11.7,52.06666666666667],[-12.6,56.55]]}')::geography, 9260, 'quad_segs=30 endcap=round join=round')
);

Red / Orange: the LineString is split into individual straight lines made of two points, each of which is then buffered and the resulting Polygons are merged using ST_Union.
SELECT ST_AsGeoJSON(
ST_Union(
ARRAY[
ST_Buffer(ST_GeomFromGeoJSON('{"type":"LineString","coordinates":[[-14.466666666666667,55.25],[7.233333333333333,54]]}')::geography, 9260, 'quad_segs=30 endcap=round join=round')::geometry,
ST_Buffer(ST_GeomFromGeoJSON('{"type":"LineString","coordinates":[[7.233333333333333,54],[0,58.86666666666667]]}')::geography, 9260, 'quad_segs=30 endcap=round join=round')::geometry,
ST_Buffer(ST_GeomFromGeoJSON('{"type":"LineString","coordinates":[[0,58.86666666666667],[-11.7,52.06666666666667]]}')::geography, 9260, 'quad_segs=30 endcap=round join=round')::geometry,
ST_Buffer(ST_GeomFromGeoJSON('{"type":"LineString","coordinates":[[-11.7,52.06666666666667],[-12.6,56.55]]}')::geography, 9260, 'quad_segs=30 endcap=round join=round')::geometry
]
));
Yellow: the LineString is split into individual straight lines made of two points, each of which is then wrapped into a MultiPolygon using ST_Collect.
SELECT ST_AsGeoJSON(
ST_Collect(
ARRAY[
ST_Buffer(ST_GeomFromGeoJSON('{"type":"LineString","coordinates":[[-14.466666666666667,55.25],[7.233333333333333,54]]}')::geography, 9260, 'quad_segs=30 endcap=round join=round')::geometry,
ST_Buffer(ST_GeomFromGeoJSON('{"type":"LineString","coordinates":[[7.233333333333333,54],[0,58.86666666666667]]}')::geography, 9260, 'quad_segs=30 endcap=round join=round')::geometry,
ST_Buffer(ST_GeomFromGeoJSON('{"type":"LineString","coordinates":[[0,58.86666666666667],[-11.7,52.06666666666667]]}')::geography, 9260, 'quad_segs=30 endcap=round join=round')::geometry,
ST_Buffer(ST_GeomFromGeoJSON('{"type":"LineString","coordinates":[[-11.7,52.06666666666667],[-12.6,56.55]]}')::geography, 9260, 'quad_segs=30 endcap=round join=round')::geometry
]
));

Is anyone able to explain the differing results we're seeing, ostensibly with the same set of coordinates and particularly why an intersecting line changes the buffer so significantly?

This is an example of creating a (different) route, showing how the buffer geometry changes considerably as soon as a line is added which intersects another part of the route.

Click this link: https://vimeo.com/320203046


texts:
In the list of raster pixel data types that are available with GDAL, the "smallest" types are Bytes, UInt16 and Int16. I would like to store integers from 0 to 100 in the smallest image type, but with the constraint of using -1 and -2 for NoData value and other specific value. GDAL Bytes are unsigned, so the smallest type would be Int16.
http://www.gdal.org/gdal_8h.html#a22e22ce0a55036a96f652765793fb7a4

I use this type to create rasters, with the Rasterize function in Python bindings.

Do you know if there is a trick to use signed Bytes ?


texts:
I have a python script which works with Grass 6.4. At some point during the project I decided to split this big python script into several scripts. I have made all necessary variables global, created methods, imported scripts in other scripts, etc.

The problem now is that during execution, I get an error on importing grass. More specifically, on "import grass.script as grass". This happens in any of the new scripts that have the import line. The strange thing is that I had the same import line in the old script (the big one) and it never failed.

Anyone has any ideas? Is there anything I am missing here?


texts:
I'm following these instructions to setup osm2psql . I need it to import osm data file inside a PostGIS DB.

The following steps worked fine:

createdb -E UTF8 -O user mydb
createlang plpgsql mydb

I successfully added PostGIS extensions too:

psql mydb < postgis.sql
psql mydb < spatial_ref_sys.sql

I start getting problems when trying to alter table's owner:

psql -d mydb -c "ALTER TABLE geometry_columns OWNER TO user"
psql -d mydb -c "ALTER TABLE spatial_ref_sys OWNER TO user"

both above commands produce the following syntax error:

  ERROR:  syntax error at or near "user" LINE 1: ALTER TABLE
  geometry_columns OWNER TO  user

I've also tried to add a semicolon at the end of the query, but I still get the same error:

psql -d mydb -c "ALTER TABLE geometry_columns OWNER TO user;"
psql -d mydb -c "ALTER TABLE spatial_ref_sys OWNER TO user;"

I'm not so confident with SQL nor with Postgres.
Could anyone tell me where's my mistake?


texts:
I've a scenario that I need to display many lines on the map (for instance all road links of Manchester in UK). I am using Leaflet as map plugin. I use geoJSON as source of data, when I try to load more than 1000 features the map and browser become unresponsive.
What do you propose to improve performance when number of features to be drawn is larger than 1000 (even more..)?

I was thinking of a way to just display only roads which are in the shown map extent. But still if user zooms out number of features will be exceeded.


texts:
So I have a feature class : City Limits

I have set it to have a fill of 50% gray. This places a gray shade over the entire city.

I have a 'Map Page' polygon feature class called 'Map Page' which is a square mile boundary outlining the basic contents necessary for each map book page. I am building a data driven mapbook based on the 'Map Page' feature. The Data Driven pages are showing 110% of each Map Page, meaning that for each Data Driven page there is roughly a 1/5th of a mile border around the Map Page feature made visible additionally.  I have clipped the data frame to the data driven index feature (Map Page) and now it is showing the 50% gray fill (from the City Limits feature) inside of each Map Page boundary.

I want to have the 50% gray fill of the City Limits ONLY SHOW UP OUTSIDE of the Map Page boundaries (which would put it ONLY visible between the map page boundary and the printable area of each data driven page. This would mean that each page that gets printed would show a gray border around each Map Page boundary that runs to the edge of the print area. Instead, when I setup clipping, it shows the 50% gray INSIDE the Map Page boundary and it's blank outside of it.

Can anyone help me out? I basically want to INVERT the clip so that the City Limits feature does NOT clipped inside the MAP PAGE boundary, but shows up outside of the boundary.


texts:
Someboby please help me with understanding of the following question.
I have many postgreSQL databases with osm data on a server. Each of them contains the map of the certain country. I want to look through the maps on my notebook via offline WMS client. Can I convert databases to .mbtiles for tileserver or is there another way for a best perfomance?


texts: Dzetsala plugin on QGIS, missing scikit-learn

texts:
I want to run random forest classification from Dzetsaka plugin on QGIS 3.4.7 but I can't as a message appears as following:

  Scikit-learn library is missing on your computer. You must use Gaussian Mixture Model, or consult dzetsaka homepage to learn on to install the missing library.

I installed scikit-learn in the following way:
Open OsGeo shell, then :

py3_env.bat
python3 -m pip install scikit-learn --user

The following messages appears:

  Requirement already satisfied: scikit-learn in c:\users\gugir\appdata\roaming\python\python37\site-packages (0.22)
  Requirement already satisfied: joblib>=0.11 in c:\users\gugir\appdata\roaming\python\python37\site-packages (from scikit-learn) (0.14.0)
  Requirement already satisfied: numpy>=1.11.0 in c:\progra~1\qgis3~1.4\apps\python37\lib\site-packages (from scikit-learn) (1.15.4)
  Requirement already satisfied: scipy>=0.17.0 in c:\progra~1\qgis3~1.4\apps\python37\lib\site-packages (from scikit-learn) (1.1.0)

So it seems that scikit-learn is present on my computer, but still I can't run random forest classification form dzetsaka plugin. What can I do?


texts: Getting accurate length for custom dash pattern in QGIS 3.4

texts:
I try to set a custom dash pattern for simple line using QGIS 3.4. Both my project and layer are in a projected CRS using meter as unit (EPSG:2154, layer are stored as Multilinestring in a geopackage).
This dash pattern should display a repeating 1.5 meter dash followed by a 5 meter blank space. I set my symbology unit as "Meters at Scale" (I also did try with "Meters" as unit with the exact same result) and use a custom dash pattern with [1.5;5] value.
As you can see in the image below the dash are rendered as 0.9m and the space between the dash mesure around 3.3m (as I'm measuring 3 dash and 2 space the total length should be 14,5m the measuring tool give a total of 9.5m). I try changing Cap style between flat and square, checking/unchecking the "Clip Feature to Canvas Extent" option and using different ellipsoid for measurement (the one associated with my CRS and none/Planimetric) but none of that give me a dash pattern close to 1.5m/5m.

Is this a bug still present in QGIS 3.4 (I find some related bug for the old QGIS version 2) or is there a way to get accurate dash pattern ?


texts: Parsing WKB with python geos: ParseException: Unknown WKB type 16

texts:
I am pulling a WKB strings from the MYSQL table in pandas / geopandas, Python 2.7. As it is MySQL, the only way to read geometry is to parse it on the client.

However, when parsing, geos raises an Exception, and I am not sure what is wrong with it.

ParseException: Unknown WKB type 16. At the same time, I know that production team uses geos binding for Ruby and everything works fine.

In particular, this is what is happening:

from shapely import geos
import binascii
geos.WKBWriter.defaults['include_srid'] = True
from shapely import wkb

connection = pymysql.connect(host=***,
                         user=***,
                         password=***,
                         db=***
                         )

df = pd.read_sql_query('SELECT geometry, building_count, created_at  FROM boundaries LIMIT 10', con=connection)
df.head(3)

    geometry    building_count  created_at
0   �\t+~��}R�)�Ǻ\D@�ʡE�}R�'�UH�]...    2892    2011-03-24 01:59:15
1   �.�o�}R�R�8�cD@�'�X}R��ǵ�bd...  1317    2011-03-30 15:02:45
2   �\r-��燎R�����]D@{�/L��R��?�Ŋ^D@...  5810    2011-03-30 23:05:02

z = df.loc[0, 'geometry']
zh = z.encode('hex')
zh

'e6100000010300000001000000090000002b137ea99f7d52c029cb10c7ba5c4440c1caa145b67d52c027da5548f95d444081b22957787d52c0d1e80e62675e44403a92cb7f487d52c0ad86c43d965e4440c959d8d30e7d52c0ca54c1a8a45e4440bc3fdeab567c52c075594c6c3e5e4440aa436e861b7c52c07767edb60b5d444094a46b26df7c52c09bc937dbdc5c44402b137ea99f7d52c029cb10c7ba5c4440'

wkb.loads(zh, hex=True)

-> ParseException: Unknown WKB type 16

I am sure that there is something with encoding or WKB specification, but can't figure out how exactly to resolve this.


texts: A way to find out when credits were purchased in ArcGIS Online

texts:
Is there a way to tell when ArcGIS Online credits were purchased?  Our organization keeps everything behind the firewall in the Portal but one of our groups bought ArcGIS Online credits sometime in Fall/Summer 2017 (my predecessors documentation is not good).  We don't do much in there and these credits have not been touched and will not get used up before they expire.  Is there a way to find out when they were purchased so I can notify the group before they vanish in a few months?
I email our Esri account rep but no reply yet.


texts:
I have a world map divided into zones and need to compute the area for each polygon. I tried using Vector > Geometry Tools > Export/Add geometry columns and successfully get the area and perimeter but I don't know the units of this measurement.

I assumed this was in square meters because the Layer CRS is WGS84 and I have 'Units for area measurement' in Project > Project Properties set to square meters but it seemed off by couple orders of magnitude. If I assumed these measurements were in square kilometers the sum of all the areas is still off by about 20% from the earths surface area.


texts: Using ArcPy to Generate Near Table in memory

texts:
I have a point SHP in WGS84 (EPSG 4326) and I am trying to tabulate the distances to the nearest point in meters. I need to do the analysis using ArcPy, but I want to get the results in a Pandas dataframe to do other things.
I'm stuck on getting the GenerateNearTable_analysis() in an in-memory format (throws the error below). And then from there, I actually don't know how to pass it into a dataframe.
import arcpy
import pandas as pd

filein = [
    r'D:\mypath\test1.shp',
    r'D:\mypath\test2.shp',
]

for f in filein:
    fout = f.replace(".shp", "")
    TempTab = 'in_memory\\{}_Table'.format(fout)
    arcpy.GenerateNearTable_analysis(f, f, TempTab, method='GEODESIC')

Above throws this error:
arcgisscripting.ExecuteError: ERROR 000210: Cannot create output in_memory\D:\mypath\test_Table
Failed to execute (GenerateNearTable).

I also tried arcpy.GenerateNearTable_analysis(f, TempTab, "memory", method='GEODESIC') as suggested in this answer, and got the following error: arcgisscripting.ExecuteError: Failed to execute. Parameters are not valid. ERROR 000732: Near Features: Dataset 'in_memory\D:\mypath\test_Table' does not exist or is not supported Failed to execute (GenerateNearTable).
Using ArcMap 10.7


texts: Getting layers to html with mouseover for all layers?

texts:
I'm a casual user of Qgis for my job.  This latest project has taken me to unfamiliar waters and I need advice on last remaining steps to get me to my goal.  My goal is to have a choropleth map of the lower 48 states of America with another layer of cities. Using data from the attribute tables imported by csv, I want the mouse over (tool tip) feature to work for both layers simultaneously rather than just the one selected.  Then I want to image map this to html so I can use it in web form.

I've gotten all the way to creating html tool tips in the layer properties screen but this is where I'm stuck.

How do I get the tool tips to work on two layers simultaneously in an html format for web use?

With new projection suggested by user3338197.


texts: How can I use gdal_retile on RGBI (4 band) GeoTIFFs while preserving band information?

texts:
I want to retile large GeoTIFFs with 4 bands: R = Red, G = Green, B = Blue, NIR = Near Infra Read; in short "RGBI" if I understand correctly.

gdalinfo on the input file prints as band information:

Band 1 Block=512x512 Type=Byte, ColorInterp=Red
Band 2 Block=512x512 Type=Byte, ColorInterp=Green
Band 3 Block=512x512 Type=Byte, ColorInterp=Blue
Band 4 Block=512x512 Type=Byte, ColorInterp=Undefined

After running gdal_retile.py, the fourth (NIR) channel is converted to an Alpha-Channel.
gdalinfo on the output file prints as band information:

Band 1 Block=1024x2 Type=Byte, ColorInterp=Red
  Overviews: 512x512, 256x256, 128x128
  Mask Flags: PER_DATASET ALPHA
  Overviews of mask band: 512x512, 256x256, 128x128
Band 2 Block=1024x2 Type=Byte, ColorInterp=Green
  Overviews: 512x512, 256x256, 128x128
  Mask Flags: PER_DATASET ALPHA
  Overviews of mask band: 512x512, 256x256, 128x128
Band 3 Block=1024x2 Type=Byte, ColorInterp=Blue
  Overviews: 512x512, 256x256, 128x128
  Mask Flags: PER_DATASET ALPHA
  Overviews of mask band: 512x512, 256x256, 128x128
Band 4 Block=1024x2 Type=Byte, ColorInterp=Alpha
  Overviews: 512x512, 256x256, 128x128

The conversion of the Infrared-band into an Alpha mask is causing me troubles. Can I prevent this behavior (of gdal_retile.py)? Or alternatively: Can I reverse the effect by changing the Alpha band back to an "Undefined" band?

The gdal_retile.py parameters I ran with (on GDAL version 1.11) are:

gdal_retile.py
-s_srs EPSG:31464 [georeference]
-ps 1024 1024 [pixel width/height)
-tileIndex "tileIndex.shp"
-csv "tileGeoRef.csv"
-targetDir "...\outputfolder"
"...\inputfile.tif"

Purpose of the operation is to cut the large GeoTIFF into small 1024x1024 images.


texts:
I have point(picket?) between 2 contour lines with even number altitudes.
Now i got these new 5 meters steps contour in between. As you see, this point between red lines was OK, but now that this 175 apeared it's wrong.
I want to find faulty points between contours.
I would love to avoid creating new polygons from contour lines, because there's alot of features. Any ideas?


texts: QGIS composer - How to create a table with attributes and style of legend?

texts:
Is it possible to create a table in QGIS composer in which I can display organized additional information about items from one layer as well as the symbol of the categorized symbology?

For example:
Layer1

species | regional redlist-state | nationwide rl-state | etc..
bird1 | 1 | 2 | ..
bird2 | 3 | * | ..

The symbols for Layer1 are categorized by the species-column.

I know the attribute table, which is almost perfect. The problem is that I didn't found a solution to show the symbol of the categorized symbols for the layer as well.

It would be nice if the symbol would be shown in the first column. I'm pretty sure that there might be a way, but don't know how to achieve it.


texts: How do I set the type of geometry when creating a SQL table?

texts:
The table I am trying to recreate has these properties
Data Type      File Geodatabase Feature Class
Feature Type   Simple
Geometry Type  Polygon    <<<< WHAT I WANT TO RECREATE

My target SQL table is this:
CREATE TABLE [dbo].[MyAggregateTable](
    [OBJECTID] [int] IDENTITY(1,1) NOT NULL,
    [SHAPE] [geometry] NULL,
    [Status] [nvarchar](25) NULL,
    ...

And the resulting properties are:
Data Type      Enterprise Database Feature Class
Feature Type   Simple
Geometry type  Any         <<<< NOT THE SAME!

FURTHER PROGRESS / EDIT
Now I have discovered that I can manipulate the geometry type by using the ArcCatalog App. I can change the Type to "Polygon Features" manually on the FeatureClass properties - but it's not permanent. As soon as I exit, it 'forgets'.
Then I dicovered that if I manually "Load Data" in the ArcCatalog App, the FeatureClass properties are set, the data is loaded, and we are good to go.
However, if I flush the table, (not DROP), the geometry type reverts back.
So basically a SearchCursor / InsertCursor approach does not work >>> it has to be done manually.


texts: nested iterator, 2 iterations in one model, collect values, cell statistics, batch processing dialogue box

texts:
(ArcGIS10, SP2, ArcInfo)

This is my model, and it works perfectly, but...

I do have a large amount of data to process. The interim output files (from raster calculator, collect value) would reach the capacity of my hard drive. Therefore, I cut my input shape file into tiles and want to run it on all tiles and just overwrite the interim output files which are create for each feature.

However, as we all know, the model builder does not allow using two iterators in one model.

I used to do this with the batch processing function (http://www.lib.uoguelph.ca/resources..._arcgis_10.pdf).

However, the cell statistics tool does not allow me edit the output location within the batch dialogue box (after right clicking on the model, batch).

The fields are greyed in the table and I cannot add further output locations. Thus, I cannot run the model on several input files.

I guess it is a buck of the combination of cell statistics and collect values. Only the combination of those two tools resulted in this problem.

Potential solution for my problem are:

• Fixing the problem with the output location within the batch dialogue box.

• Substituting the combination of cell statistics and collect values by anything else resulting in the same outcome.

• Using any other possible way to run a feature iterator on several feature classes (something like a nested iterator)

Any suggestions?


texts: Seeking public database to input date and city to determine country at time and present-day?

texts:
I'm looking to automate a common task for us.  For example, Jaromir Jagr was born in Kladno in 1972.  I want to be able to take that and output a string.

Born in Kladno, Czechoslovakia (present-day Czech Republic).

I'm a little bit familiar with geonames.org, but I'm not sure if it has the historical component.  Any places that have that?


texts:
I am trying to perform a batch clip using arcpy. I cannot seem to get the script to save where I want them to, with a specific name. I remove the "Ohio_" and rename them with "Summit_" what is in the directory ex. "Summit_CensusTracts".

Any input?

import arcpy
import os
from arcpy import env

arcpy.env.workspace = "H:/Python_GIS/Summit_County_Project/Shapefiles"
# Variables
fcInput = "H:/Python_GIS/Summit_County_Project/Shapefiles/Ohio_Counties.shp"
fcOutput = "H:/Python_GIS/Summit_County_Project/Shapefiles/Summit_County_Updates/Summit.shp"

# Expression
whereClause = "COUNTY_NAM = 'SUMMIT'"

# Select
arcpy.Select_analysis(fcInput, fcOutput, whereClause)

# Processing

outFolder = "H:/Python_GIS/Summit_County_Project/Shapefiles/Summit_County_Updates"
clipFeature = "H:/Python_GIS/Summit_County_Project/Shapefiles/Summit_County_Updates/Summit.shp"
shpList = ["Ohio_CensusTracts.shp", "Ohio_Cities.shp", "Ohio_Hydro.shp", "Ohio_Roads.shp"]

for root, dirs, files in os.walk(outFolder):
   for name in files:
      if name.startswith("Ohio_"):
         # shapefile name without extension
         shpName = os.path.splitext("Summit")[0]
         # absolute file path
         absFile = os.path.abspath(os.path.join(root,name))
         # output file path
         output_layer = outFolder
         # buffer
         arcpy.Clip_analysis(shpList, clipFeature, output_layer)


texts: Update feature class' data via Python Script Tool without removing layer from .mxd?

texts:
I have an mxd with 20+ Feature Classes all with varying degrees of 'formatting' (i.e. most are organized in groups, some are filtered with definition queries, and all have their own unique symbologies).

Occasionally these layers must be updated with new data.  I have written a script to automate the process of overwriting an outdated FC with new CSV data.

Here's the twist: The code works perfectly fine.  A new FC is indeed produced from the information in the new CSV.

Again there are no errors with my code, but just to give you a look at what I'm doing:

arcpy.MakeXYEventLayer_management(Source_CSV, "Lng", "Lat", Temp_Layer, SR, "")
# Following line removes FC from TOC...
arcpy.FeatureClassToFeatureClass_conversion(Temp_Layer, Target_GDB, Target_FC_name, "", Field_Mapping)
# ...so I re-add it here...
addLayer = arcpy.mapping.Layer(Target_file)
arcpy.mapping.AddLayer(df, addLayer)

The problem in terms of my desired workflow is that FeatureClassToFeatureClass removes the target FC from the map document.  I want it to remain there, hence the last two lines of code.  Of course when my script re-adds the FC it lacks all of its previous formatting (i.e. symbology, grouping, filtering, etc.) that it had before the tool was run.

How can I update a FC within an open mxd via a Python Script tool without having the FC be removed and thus lose all of its formatting?

Perhaps there's a tool that won't remove the target FC?  XYTableToPoint sounds like a good candidate, unfortunately it does not exist in my version of Arc (10.2.2).

Perhaps a SearchCursor / UpdateCursor combination will provide the solution I need?  Just seems like a lot more code and processing time than with the above script, so wanted to check with online community before starting over.


texts: "Driver not loaded Driver not loaded" when trying to connect to MSSQL database in QGIS 3.8 on MacOS

texts:
When trying to create a new MSSQL Connection, I get an error saying: Error opening connection: Driver not loaded Driver not loaded.

I have tried both QGIS 3.8.0 and the newest pre-compiled QGIS 3.8.1 on MacOS Mojave. I am not sure how to run QGIS in a mode that might provide a better error message.

Anyone have ideas on what to change to make this work?


texts:
I'm having trouble editing text in a pop up window for a layer I have in CARTO.

The data in the pop up window comes from an Excel spreadsheet. I would like to know if there is a way for me to create line breaks in the window so that the separate terms in there remain separated on their own lines. For example, the first term is 'genesis electric skateboard' and the second term is 'abigail ratchford skateboard.' I need those terms on their own lines, but as is they are running together.

I tried creating line breaks in the Excel cell using Alt+Enter but that doesn't make a difference when I import that sheet to this layer.

I would also like to know how to adjust the width of the pop up window so I can fit the entire term like 'genesis electric skateboard' on one line without it wrapping to a line below.

I would also like to know how to adjust the title of the pop up window to get rid of the underscores between the words TOP 5 GOOGLE SEARCHES.


texts:
I'm currently running on a VPN for work from home and I have access to several network drives.  I'm trying to create a new ArcGIS Pro project on one of the network drives and it seems to hang the software (spins and spins and never opens resulting in a task manager visit).  It works fine if I save my project to the local drives but not the network drives.  I  have full access to the network drives (read and write) and know issues working with files there.

Any thoughts on why this might be a problem.  Again it's not slow it just doesn't work...


texts: Changing legend window in Carto new editor?

texts:
I have created a map with an automatic legend window.

I'd like to move the window to the lower right corner.

I tried to edit the CSS but I cannot figure out how to do with the new editor.


texts: Calculating Distance and Bearing North Azimuth (0 degrees pointing North) in ArcGIS Near tool?

texts:
I am trying to calculate distance and bearing in ArcMap using the tool Near, and although it works I need the results to be calculated using North Azimuth (0 degrees pointing North) and not Polar which is the default. I can change these settings for when I am editing but not for toolbox tools.

How do I calculate distance and bearing in ArcMap using the tool Near and  North Azimuth (0 degrees pointing North)?

To be honest I think I need to understand the results the tool is giving me before running any code. According to the Near tool help, when the GEODESIC method is used, the angle is within the range of -180° to 180°, with 0° to the north, 90° to the east, 180° (or -180°) to the south, and -90° to the west. However this is what I am getting which doesn't match that. I calculate the bearing from the blue points to the green square. Values are in decimal degrees.

So even if I run the python code I don't think I will get the correct results.

Am I missing something?


texts: Configure Geoprocessing Service after publishing?

texts: I want to use some variables in a geoprocessing service like a URL, that should be configurable from the arcgis ardmin or manager. Is it possible have a tool input, which value will be specified by the admin, after the service have been published?


texts: What is the inverse of QgsGeometry.asGeometryCollection() in qgis python?

texts:
What is the inverse operation of QgsGeometry.asGeometryCollection() as documented here: http://qgis.org/api/classQgsGeometry.html#a7cb1cbf8cfa441456067a66ad3a25a94

I'd like to do something like:

from qgis.core import (QgsFeature, QgsGeometry,
                   QgsVectorLayer, QgsMapLayerRegistry,
                   QgsField)
from PyQt4.QtCore import QVariant
from qgis.utils import iface
import math

def myFunction(arg):
    return arg

layer=iface.mapCanvas().currentLayer()
feature = layer.selectedFeatures()[0]
geom = feature.geometry()
if geom.isMultipart():
    newgeom=QgsGeometry.INVERSE([myFunction(part) for part in geom.asGeometryCollection ()])
else:
    newgeom=myFunction(geom)
changeGeometry(feature.id(),newgeom)

Things that do not work:

# doesn't work:
newgeom=QgsGeometry.combine([part for part in geom.asGeometryCollection ()])
# doesn't work:
newgeom=([part for part in geom.asGeometryCollection ()])
# doesn't work:
newgeom = [part for part in geom.asGeometryCollection()]
# doesn't work: argument 1 has unexpected type 'list'
newgeom = QgsGeometry.fromPolyline([part for part in geom.asGeometryCollection()])
# QgsGeometry has no attribute 'MultiPolyline':
newgeom = QgsGeometry.MultiPolyline([part for part in geom.asGeometryCollection()])
# MuliPolyline is not defined:
newgeom = MultiPolyline([part for part in geom.asGeometryCollection()])
newgeom = QgsGeometry.fromMultiPolyline([part for part in geom.asGeometryCollection()])
newgeom = QgsGeometry.fromMultiPolyline([part.asPolyline for part in geom.asGeometryCollection()])

newgeom = QgsGeometry.fromMultiPolyline([part.asPolyline for part in geom])
newgeom = QgsGeometry.fromMultiPolyline([part.asPolyline for part in geom.asMultiPolyline()])

Works with Polylines / LineStrings:

newgeom = QgsGeometry.fromMultiPolyline([part for part in geom.asMultiPolyline()])

newgeom = QgsGeometry.fromMultiPolyline([userFunction(part,*args).asPolyline() for part in geom.asGeometryCollection()])

http://docs.qgis.org/testing/en/docs/pyqgis_developer_cookbook/geometry.html has some documentation on accessing and creating the various types.  http://qgis.org/api/classQGis.html lists the enumerations of WKB* types and plain geometry types.

It appears that each multipart type must be treated separately, and the multi-* types can't be disassembled and reassembled generally.  In other words, you need to spell out and process the individual geometries separately like:

parts = ([userFunction(part,*args) for part in geom.asGeometryCollection()])

if geom.wkbType() == Qgis.WKBMultiLineString:
    newgeom= QgsGeometry.fromMultiPolyline([x.asPolyline() for x in parts])
elif geom.wkbType() == QGis.WKBMultiPolygon:
    newgeom= QgsGeometry.fromMultiPolygon([x.asPolygon() for x in parts])
elif geom.wkbType() == QGis.WKBLineString25D:
    newgeom= QgsGeometry.fromMultiPolyline25D([x.asPolyline25D() for x in parts])
elif geom.wkbType() == QGis.WKBMultiPolygon25D:
    newgeom= QgsGeometry.fromMultiPolygon25D([x.asPolygon25D() for x in parts])
...

layer.changeGeometry(feature.id(),newgeom)

Is there a function or combination of functions that would simplify the above?  Something like:

newgeom = CAST_AS_GEOMETRY_TYPE([part for part in geom.asGeometryCollection()],geom.WKBType())


texts: Is it possible to conduct a flood model analysis with QGIS using a DEM?  Need to simulate rising water levels on a 250k map.


texts:
I need to export raster (results of r.walk and r.drain) to ArcMap for comparison.

I did it on data from guatemala - WGS84, UTM16N as r.out.gdal and it works normally.

But when I did this on Slovak data (S-JTSK / Krovak East North, EPSG 5514) and add exported layer into ArcMap, it says:

I don´t understand, why it´s not working,I set projection at the start.
Any ideas,what to do?
I´m rookie in Grass.


texts: Registration of a panchromatic satellite image to a multispectral image

texts:
I have a multispectral satellite image (RGB and NIR band) which will serve as the reference image. I also have several panchromatic (PAN) satellite images which I need to register to the multispectral (MX) images. The geographical extent covered by MX image is very large and all the PAN images fit somewhere into the bigger MX image.

The MX and PAN images are in different co-ordinate systems, however I am converting the PAN images to the system of the MX images by gdal library.

I have worked previously on registering one PAN image to another PAN image by feature extraction and feature matching using SURF features, and I tried using this with the new problem statement as well (MX with PAN) but the results were poor.

MX images have poorer resolution hence I guess SURF algorithm is not able to detect some good features in those images.

Am I approaching the problem the right way or are the ground principles different when working with satellite images with different modalities?


texts: Displaying direction of movement with points2paths in QGIS?

texts:
I have a data set containing ship movements. I used the Points2paths plugin in QGIS to connect the points with the ship's names as factor and the time the signal was picked up as sorting value. Now I would like to display the direction of movement in the form of little arrows on the lines.

Not sure if the plugin has this function or if there is any other plugin I might use?


texts: Adding output to existing shapefile using Modeler

texts:
now I try to automatize the creation of variable distance buffers for a line vector-layer with the help of the Graphic Modeler in QGIS. I don't want to save the buffer results as a new shapefile but instead adding the new features automatically to an existing shapefile (just adding not overwriting).

The problem thereby is that I can't specify the output in this way with the Graphical Modeler. Is there any possibilty to achieve this?


texts: Spatial Join features of File GDB to shape - problems

texts:
I tried to build a Model in ArcGis 10.1. It should be very simple. At one hand  I have the features of a File Geodatabase as Feature to join on the other I have the target shape. The features from the GDB are Clips from former layers with the target shape.

I am using a Iterator (feature classes), in the hope, that the three columns of each feature class would be spatial joined to the target shape. But Model aborts after it joined the first feature class from the GDB. And even this Join is not correct. Instead of the field values that should be joined, it is always the value "0" represented in the output.

I really tried every option of the Spatial Join Tool, with same or similar unuseful outputs.

Does somebody have an idea?

Here is the model:


texts: Export the values in featurecollection using Google Earth Engine

texts:
I want to export the values of NDVI and NDBI for the feature collection.

var data= L8TOA1.filterDate('2019-02-01', '2019-03-28')
              .filter(ee.Filter.equals('WRS_PATH', 146))
              .filter(ee.Filter.eq('WRS_ROW', 40))
              .sort('CLOUD_COVER').first();
print(data);
//Get cloud cover information
var cloud= data.get('CLOUD_COVER');
print('cloud_cover is',cloud);
//Cip with boundary
var Study_region = data.clipToCollection(aoi);
// Define the visualization parameters.
var colour = {bands: ['B5', 'B4', 'B3'],min: 0,max: 0.5,gamma: 1.4,}
// Add layer to the map
Map.addLayer(Study_region,colour,'FCC');
Map.centerObject(aoi, 9);

//Indice
var NDBI = ee.Image (Study_region.normalizedDifference(['B6', 'B5']).rename('NDBI'));
var NDVI = ee.Image (Study_region.normalizedDifference(['B5', 'B4']).rename('NDVI'));
var SI= NDBI.addBands(NDVI)
print("Stacked_Image",SI);
var viz = [
'#FFFFFF','#808080','#000000','#FF0000', '#800000', '#FFFF00', '#008000', '#00FFFF', '#0000FF', '#FF00FF', '#800080'];
Map.addLayer(NDBI, {min: -1, max: 1, palette: viz},'NDBI');
// Get the Indices values to the point
var VTP = ee.FeatureCollection(SI.select(['NDBI','NDVI'])
  .reduceRegion({
  reducer: ee.Reducer.mean(),
  geometry: Urban.geometry(),
  scale: 30
}));
var Urban = Urban.set(VTP) ;
print('a',Urban);
print("Indices values to the point",VTP);
Export.table.toDrive({
collection: Urban,
description :"Indcies",
folder:"GEE",
fileNamePrefix: "Indcies",
fileFormat:"CSV",
selectors:['NDBI','NDVI']
});

When I exported the data, the value of NDBI and NDVI are not exported. Only the names were there exported.


texts: QGIS 3.6 Noosa on Windows 10 cannot start

texts:
Recently downloaded QGIS 3.6 and I cannot open it from the desktop shortcut. It doesn't work when I run as an administrator either. The only way I am able to open the program is if I go to open an existing QGIS file saved on my old version of QGIS (2.18 Las Palmas).

Is there a setting that is making this happen?


texts:
How do I crop a raster layer in QGIS?

I want to make contours from a section of a layer, and including the whole layer will create far more data than I need and slow down the rendering.

I have found that I can do it by using the Raster calculator and specifying a new layer extent, but how do I know what values to enter for the rows and columns?

I also tried cropping the geotiff in Photoshop, but it seemed to lose the coordinate data.

Is there any way to create contours from part of a raster layer without having to crop it first?


texts:
I am trying to calculate NDVI from rasters that I've downloaded with download rasters, clipping them with some features (selection_feature) that I'd split from an original archive. I am a begginer in Python scripting, and I can't understand what is going wrong with this script!

import arcpy
import os
from arcpy import env
from arcpy.sa import *

arcpy.CheckOutExtension("Spatial")

NDVI_tempo = 'Z:\\04 - Temp\\Charles\\NDVI_tempo\\NAL_PL'

arcpy.env.workspace = 'Z:\\04 - Temp\\Charles\\NDVI_tempo\\NAL\\NAL.gdb'

arcpy.SplitByAttributes_analysis("Area_produtiva", NDVI_tempo, ['ID_TALHAO'])

down = 'Z:\\04 - Temp\\Charles\\NDVI_tempo\\NAL_PL'
image_server = 'Y:\\E - IMAGEM SERVER\\IMAGEM_SERVER.gdb\\F01_PLANETLAB'

arcpy.env.workspace = r'Z:\\04 - Temp\\Charles\\NDVI_tempo\\NAL_PL'
fcs = arcpy.ListFeatureClasses()

for fc in fcs:

    arcpy.DownloadRasters_management (image_server, down, "", fc, "CLIPPING", "CONVERT_AS_REQUIRED", "TIFF", "NONE", "", "NO_MAINTAIN_FOLDER")

rasters = arcpy.ListRasters("*", "TIF")
env.workspace = 'Z:\\04 - Temp\\Charles\\NDVI_tempo\\NAL'
for tiff in rasters:

    arcpy.MakeRasterLayer_management(tiff, "Z:\\04 - Temp\\Charles\\NDVI_tempo\\NAL\\temp_rast1.tif", "", "", "3")
    arcpy.MakeRasterLayer_management(tiff, "Z:\\04 - Temp\\Charles\\NDVI_tempo\\NAL\\temp_rast2.tif", "", "", "4")

    arcpy.CopyRaster_management("Z:\\04 - Temp\\Charles\\NDVI_tempo\\NAL\\temp_rast1.tif", "copy_rast1", "","", "", "", "", "", "", "", "", "")
    arcpy.CopyRaster_management("Z:\\04 - Temp\\Charles\\NDVI_tempo\\NAL\\temp_rast2.tif", "copy_rast2", "","", "", "", "", "", "", "", "", "")

    Property1Result = arcpy.GetRasterProperties_management ("copy_rast1", "MINIMUM")
    Property2Result = arcpy.GetRasterProperties_management ("copy_rast2", "MINIMUM")

    Property1 = Property1Result.getOutput
    Property2 = Property2Result.getOutput

    Correcao1 = Minus ("copy_rast1", "Property1")
    Correcao2 = Minus ("copy_rast2", "Property2")

    Rastdiv1 = Minus ("Correcao2", "Correcao1")
    Rastdiv2 = Plus  ("Correcao2", "Correcao1")

    NDVI_tiff = Divide (Rastdiv1, Rastdiv2)

The script is running pretty well until make raster layer, I don't understand what is happening there, it seems like it doesn't recognize the downloaded image! The error message is that:

Traceback (most recent call last):
  File "Z:\04 - Temp\Charles\NDVI_tempo\NAL\NAL.gdb\NDVI_RF#Script_NDVIRF.py", line 38, in <module>
  File "c:\program files\arcgis\pro\Resources\arcpy\arcpy\management.py", line 6219, in MakeRasterLayer
    raise e
  File "c:\program files\arcgis\pro\Resources\arcpy\arcpy\management.py", line 6216, in MakeRasterLayer
    retval = convertArcObjectToPythonObject(gp.MakeRasterLayer_management(*gp_fixargs((in_raster, out_rasterlayer, where_clause, envelope, band_index), True)))
  File "c:\program files\arcgis\pro\Resources\arcpy\arcpy\geoprocessing\_base.py", line 506, in <lambda>
    return lambda *args: val(*gp_fixargs(args, True))
arcgisscripting.ExecuteError: Failed to execute. Parameters are not valid.
 ERROR 000865: Input raster: 20170312_124653_0e1f_381_s17.tif does not exist.
Failed to execute (MakeRasterLayer).
Failed to execute (Script).


texts: Calculate time duration from military time with ArcMap Field Calculator

texts:
I have two text fields (STRTTIME, ENDTIME) which are are formatted as military time with leading zeros. For example, 0540, 1655.

I want to calculate a new field with the time duration as integers (in minutes).

How can I go about this using the Field Calculator in ArcMap?


texts: What to do with old GIS books/manuals/etc?

texts:
We're having an office clear-out and have come across a number of old GIS paper materials. Random examples (to give a taste - there are about 20-30 in all) include:

Postcodes: The New Geography (Book - 1992)
Intermap (Product Handbook and Quick Start Guide (2004))
ARC Macro Language (Self-study Workbook 1993)
Understanding GIS; The ARC/INFO Method (ESRI Book, 1990)
Institute of British Geographers; Transactions New Series; Volume 16 Number 1 1991
arc/info; Computer mapping and geographica information management system (Folder, probably training materials; Nov 1985!)

So I'm wondering - what should we do with them? Is there anywhere that would take "historic" GIS materials (museum, library etc)? They're all a bit too specialised for a regular charity shop but it seems a shame to just dump them.

(Note: I'm in the UK if that makes any difference)


texts: Searching an ArcGIS Server 10.1 search service with ArcGIS for Sharepoint

texts:
Could someone kindly tell me if it is possible to search an ArcGIS Server  10.1 Search Service index using ArcGIS for SharePoint?

I am basically looking for a comparable search functionality to the Enterprise Search in ArcMap 10.1.


texts:
I am trying to add attributes to a layer. I have a select function that finds all points within 10 miles of a feature. If the point is within the radius distance I add it to an array and make another "Selected Set" out of it for a map and report. However they now want to sort the selected point based on distance from.  How do I push this vale to the selected point?

var selPts = [];
    function SelectPoints(lat,lon){
        var dist = 10; //Get Radius distance
        xy = [lat,lon];  //center point of circle
        var theRadius = dist * 1609.34  //1609.34 meters in a mile

        selPts.length =0;  //Reset the array if selecting new points

        sites.eachLayer(function (layer) {
            // Lat, long of current point as it loops through.
            layer_lat_long = layer.getLatLng();

            // Distance from our marker To current point in meters
            distance_from_centerPoint = layer_lat_long.distanceTo(xy);

            // See if pt is within radius, add the to array
            if (distance_from_centerPoint <= theRadius) {
                 selPts.push(layer.feature);
                 // how do I push distance_from_centerPoint into the feature properties?
            }
        });


texts: Exporting styles from ArcGIS Desktop to MapServer *.map file?

texts: How to export styles from ArcGis 9.3 (or 10) to a .map file. Is it possible?


texts: Convert DWG to Shapefile Using FME and R

texts:
I converted a DWG file to an ESRI Shapefile using FME Workbench. Now I wish to create an elevation raster. I'm not sure how to extract the elevation info.

library(sp)
library(raster)

# https://www.dropbox.com/s/n3nvssl3swa8ryn/Massings.shp?dl=0
shapename <- read_sf('Massings.shp')
plot(shapename, axes=TRUE)

# Data has no assigned CRS. Set it to what website says it is.
st_crs(shapename) <- 2019

# Now can transform CRS to lat/long OK
shp_4326 <- st_transform(shapename, crs=4326)
plot(shp_4326, axes=TRUE)

# Create empty raster
ext <- extent(c(xmin = -79.395, xmax = -79.385, ymin = 43.638, ymax = 43.644))
gridsize <- 5e-5
r <- raster(ext, res = gridsize)

# Fill with elevation data from shapefile ???
# dwntwn_raster <- rasterize(shp_crs, r, field = ???)

Shapefile 24.8MB


texts: Loading table from Postgres database using PyQGIS

texts:
I'm trying to load a table in QGIS layers tree using PyQGIS (see here for example). Said table has got a primary key with no null/duplicate entries; a geometry column (which is empty for now); and works perfectly fine when I load it using QGIS database manager.

But, when I try to load it using PyQGIS, the layer does not load correctly: there's a /!\ next to its name in QGIS layers tree, and print(myLayer.isValid()) returns False. And I can't figure out what is wrong with my code...

Here is my code so far (host, bdd, password etc.. have been defined as global variables earlier in the code):

uri = QgsDataSourceUri()
uri.setConnection(host, port, bdd, user, password, QgsDataSourceUri.SslDisable)
uri.setDataSource("ref", "grds_types_milieux","empty_geom","")
myLayer = QgsVectorLayer(uri.uri(), "grds_types_milieux", "postgres")
QgsProject.instance().addMapLayer(myLayer)


texts: Python Calculate Cumulative Polygon Areas (Basins)

texts:
I've been asked to automate the calculation of cumulative polygon areas within a basins shapefile. The calculation will take the area of the upstream basin and add it to the area of the downstream basin, and continue adding to that cumulative value until reaching the farthest downstream basin. It also needs to take into account the areas from two or more stream basins due to stream confluences. This leads to two problems:

How to add the value from one row to the value in another row (assuming with an update cursor).
How to calculate the cumulative area moving from upstream to downstream.

The screenshot below shows the test shapefile I was given with the cumulative square miles manually calculated. As far as I know, the Basin_ID should always be unique. The table has some sample rows selected for Stream04 and Squaw Creek. For this example, BASIN17 needs to be the cumulative square miles of Basins 3 and 8 from Stream04, and Basins 4, 11, 14, 16, and 17 on Squaw Creek. In that order.

The "US" in Fields US_1, US_2, US_3, and US_4 stands for upstream (basin). "DS" is downstream basin. If the US_1 field is blank, then that basin is the most upstream basin on the stream and the starting point for the cumulative area calculation. I think a field may need to be added to flag this basin as the upstream starting point, and is used in the example below (but isn't in the table).

Initially it seems like it should be:

shapefile = path
fields = ["Basin_ID","Cum_SqMi","AreaSqMi","US_1","US_2","US_Start"]
with arcpy.da.UpdateCursor(shapefile,fields) as uc:
  for row in uc:
    if row[4] == " ":
      row[6] = 1 # calc as 1 if true
    else:
      row[6] = 0

Now I run out of logic. For the Stream04 example, we know that BASIN3 is the upstream most point. To calculate the cumulative area for BASIN8, it's simple math: BASIN3 AreaSqMi + BASIN8 AreaSqMi. These are separate rows, and as far as I know, this can't be done with a cursor since it's only looking at one row at a time.

So, what needs to be done to add the value from a row that meets a condition to the value in another row? What about a dictionary or an array? I'm not very familiar with those, hence the reason why I keep going back to the cursor. I know what's below won't work, but it shows what should be done. It also doesn't handle the multiple starting basins that will equal 1, so I may need to make that equal to the Basin_ID instead.

if row[6] == 1:
  startBasin = row[0]
  startBasinValue = row[3]
if row[4] == startBasin:
  row[2] = row[3] + startBasinValue

And then to how to continue with the calculation based on if the upstream cumulative area has been calculated yet? It seems like it would be possible to cascade the calculations from the starting point. I could manually add a field and calculate a sequential order for calculating the rows, but there are around 50 shapefiles that would have to be done for, and the user wants it as automated as possible.

Am I in over my head? Going about it entirely wrong? Should they stick to calculating the areas manually in a spreadsheet?


texts: UpdateCursor - how to pass over certain rows based on value using conditional statement

texts:
I have a script that is supposed to update several fields based on the values of several other associated fields.  My script is almost working as it should.  However, it is overwriting my "High" and "Medium" values with "Low" despite my conditional statement.  Any suggestions of where I am going wrong?

# fill polygons with "Medium", but do not overwrite those already filled with "High"
for species, primary in zip(SGCN_list, primary_habitat):
    if species not in (f.name for f in arcpy.ListFields(SGCN_FC)):
        pass
    elif primary is None:
        pass
    else:
        with arcpy.da.UpdateCursor(SGCN_FC, [species, primary]) as cursor:
            for row in cursor:
                if row[1] is not None: # if record is already populated with "High"
                    pass
                elif row[1] >= 0.10:
                    row[0] = "Medium"
                else:
                    pass
                cursor.updateRow(row)

# fill polygons with Low
for species, secondary in zip(SGCN_list, secondary_habitat):
    if species not in (f.name for f in arcpy.ListFields(SGCN_FC)):
        pass
    elif secondary is None:
        pass
    else:
        with arcpy.da.UpdateCursor(SGCN_FC, [species, secondary]) as cursor:
            for row in cursor:
                if row[1] is not None: # if record is already populated with "High" or "Medium"
                    pass
                elif row[1] >= 0.10:
                    row[0] = "Low"
                else:
                    pass
                cursor.updateRow(row)


texts: How to host ArcGIS 9.3, 10.0 & 10.1 licenses on same server?

texts:
We have 3 license servers for ESRI products.

Server 1: 10.1 License Server
Server 2: 10.0 License Server
Server 3: 9.3 License Server

We are in the middle of a massive server migration/merge and would like to combine our license servers into one and store it on a single centrally-located machine. But we have some specific applications that require 9.3, a few that require 10.0, and our more common use-cases are migrating to 10.1 so we can't just upgrade the licenses. But keeping them on different servers won't work because we are losing some of our servers due to a reshuffle.

So, with 9.3 Concurrent use licenses being utilized, can we store those alongside 10.0 and 10.1 licenses on the same server without losing them?


texts: Getting max/min value from list of rasters in ModelBuilder?

texts:
I have 3 Rasters and would like to have the max/min value from all of these raster using ModelBuilder.

Unfortunately the "Get Raster properties" tool gives me the value of each Raster separately.


texts: Create from Vector Tiles a custom svg section

texts:
I'm planing to enclude some mapping features to my wordpress website . Therefore I'm asking you whether my plan makes sense and is practical. And also how much a freelancer would cost (one payment).

I created this style.json above with Mapbox to style my Vector Tiles.
The style.json should be applied on the OSM2VectorTiles download. Then somone need to help me by installing them on my A2 Hosting.
Most complicated would be probably a custom wordpress plugin which create from the Vector Tiles a svg (like on the openstreetmap hompage - just with custom styles). I want to use this svg maps in the Mapplic Plugin (probably).


texts: Verifying that it is ok to use multiple band combinations as I create polygons for my supervised classification?

texts:
I am sorry if this is a naive question, but I am very new to remote sensing analysis.

I am conducting a supervised classification in Google Earth Engine. I want to verify that as I develop my training samples that it is ok for me to toggle between spectral band combinations (as the different combinations help me see more clearly what to identify as which landcover type). I think that doing so is standard practice, but I want to verify!

i.e. I want to confirm that when I am viewing the different layers, represented by different spectral combinations, and drawing polygons that I identify as different landcover types, that even though the visualization of the pixels are changing when I toggle between spectrum combinations, when I actually run the algorithm that trains the classificatoin based on my samples that the algorithm isn't running it based off of the premise that I used just one specific band combination.

I hope my question makes sense. thank you!


texts: Reproject multiple vector layers in QGIS

texts:
How can I reproject multiple vector  layers in QGIS? I have a thousand layers (around 100 in each folder) without projection info. All need to be set in WGS84.
I know about the batch Processing tool but they are so many and in different folders. Is there a command or script in python?


texts: Can I use .crd and .dwg files as QGIS data sources?

texts: I  want to use my autocad files(.dwg) and my field files (.crd), and want to know if they are compatible with qgis.


texts: QGIS server: use of FILTER clause with boolean values

texts:
I have set up a QIS server (3.4) to serve WMS layers with a PostGIS database.
The WMS queries include a dynamically created FILTER clause.
However, the official documentation is quite lacking explanations for the FILTER clause and I have troubles with filtering based on boolean attributes.
I want to filter based on  zico IS true where, in the PostGIS database, zico is a boolean field of the table locations.
I have tried:
http://qgis.localhost/cgi-bin/qgis_mapserv.fcgi?MAP=/home/qgis/projects/my_qgis_project.qgs&service=WMS&request=GetMap&layers=locations&styles=&format=image/png&transparent=true&version=1.3.0&srs=EPSG:3857&width=256&height=256&crs=EPSG:3857&bbox=2504688.5428486555,3757032.8142729844,3757032.814272984,5009377.085697314&FILTER=locations:%22zico%22%20IS%20True

i.e. including  FILTER:locations:'zico' IS True in the query to the WMS server,
for which the server response was the standard error message

The filter string "zico" IS False has been rejected because of
security reasons. Note: Text strings have to be enclosed in single or
double quotes. A space between each word / special character is
mandatory. Allowed Keywords and special characters are
IS,NOT,NULL,AND,OR,IN,=,<,>=,>,>=,!=,',',(,),DMETAPHONE,SOUNDEX. Not
allowed are semicolons in the filter expression.

What is the proper syntax for FILTER with boolean attributes?


texts: Automatically populate r.reclass rules based on table values

texts:
In QGIS I have a raster map showing land, snow and cloud zones, having resp. pixel values 1, 2, and 3. I extracted the mean altitude values (in m) of each zone on the basis of an underlying DEM and Zonal Statistics of the LecoS (Landscape Ecology Statistics) plugin. The results show in a table like:

On the basis of the mean values of zone 1 and 2, using r.rclass from the Processing Toolbox, I reclassified my zones using this reclass rules text file:
2588 thru 4497  = 1
4498 thru 5122  = 12
5123 thru 7533  = 2
(2588 is the minimum altitude of the DEM, 7533 the maximum altitude).
The result is a new raster map with zones land (1), snow (2) and mixed land-snow pixels (12). This works fine in separate steps and manually preparing the reclass rules text file.

However I have thousands of files to process!
How would it be possible to use the output of the table mean values to automatically populate a reclass rules text file, that can then be directly used in r.reclass?
r.reclass doesn't have an option to import  values from an outside table file. I guess I need to write a script for this. I have no experience with scripts so far.


texts:
I'm trying to apply symbology to my dataset after it is populated to the map display, but I am getting the following error. Have I placed the ApplySymbologyFromLayer code incorrectly?

Traceback (most recent call last):
  File "C:\script\script.py", line 61, in <module>
    arcpy.ApplySymbologyFromLayer_management (Output_Feature_Class, in_symbology_layer)
  File "c:\program files (x86)\arcgis\desktop10.2\arcpy\arcpy\management.py", line 6004, in ApplySymbologyFromLayer
    raise e
ExecuteError: Failed to execute. Parameters are not valid.
ERROR 000732: Input Layer: Dataset C:/data/cluster_outlier_fieldname.shp does not exist or is not supported
Failed to execute (ApplySymbologyFromLayer).

My script:

for fieldName in fieldList.split(';'):

    # Iterate
    Output_Feature_Class = "C:/data/cluster_outlier_" + fieldName + ".shp"
    outCluster = arcpy.ClustersOutliers_stats(inPointFeatures, fieldName, Output_Feature_Class,
                                              Conceptualization_of_Spatial_Relationships, Distance_Method,
                                              Standardization, Distance_Band_or_Threshold_Distance,
                                              Weights_Matrix_File, Apply_False_Discovery_Rate_FDR_Correction)

    # Add layers to map
    mxd = arcpy.mapping.MapDocument("CURRENT")
    df = arcpy.mapping.ListDataFrames(mxd, "*")[0]
    addLayer = arcpy.mapping.Layer(Output_Feature_Class)
    arcpy.mapping.AddLayer(df, addLayer, "BOTTOM")

    # Symbolize
    in_symbology_layer = "C:/data/reference/design.lyr"
    arcpy.ApplySymbologyFromLayer_management (Output_Feature_Class, in_symbology_layer)
    arcpy.RefreshActiveView()
    arcpy.RefreshTOC()
    del mxd, df, addLayer


texts: Importing WFS Feature in ArcGIS Desktop?

texts:
I am trying to Import an open source  WFS Feature in ArcGIS. I have filled all the Options correctly but the Toolbox does not work and I receive error that feature is empty.

How can I correct it?

Here is the link:

https://security.geobasis-bb.de/wss/service/VG-WFS/guest?SERVICE=WFS&REQUEST=GetCapabilities&

I have tried different other open source WFS Features too but none of them work. Even I tried the sample from ESRI but it behaves the same.

http://sampleserver6.arcgisonline.com/arcgis/services/SampleWorldCities/MapServer/WFSServer?request=GetCapabilities&service=WFS&


texts: GDAL ImportError: DLL load failed: The specified module could not be found

texts:
I'm running Python 3.5.4 on Anaconda and i'm trying to install geopandas.  I've followed the instructions to install geopandas here -
 http://geoffboeing.com/2014/09/using-geopandas-windows/ - but whichever way i try i get the same error.
I've looked at many pages on gis.stackexchange and other sites that seem to have similar error messages, but none of them had the exact same error message, and no solutions worked for me.  I've set up the Path environment variable, C:\Users\joseph.mabe\AppData\Local\Continuum\anaconda3\Lib\site-packages\osgeo; but to no avail.

I wondered if my use of Jupyter was part of the issue?


texts: Making Shapely box larger using GeoPandas?

texts:
I'm creating a box that covers the the area of a GeoPandas GeoDataFrame,
so that I can use it as an inverted map to later cover up unwanted spill-over data.

import geopandas as gpd
import geoplot as gplt
from shapely.geometry import box

denmark = world[world.name == 'Denmark']
denmark_box = gpd.GeoDataFrame(
    [box(*denmark.total_bounds)],
    columns = ['geometry'],
    geometry='geometry',
    crs = denmark.crs)

denmark_cover = gpd.tools.overlay(denmark, w3, how="symmetric_difference")

fig, ax = plt.subplots(figsize=(20,20))
ax.set_aspect('equal')

# shp is some other data
# that plots a bit over the country borders
shp.plot(ax=ax, column='sum', edgecolor='none', cmap='hot', alpha=1, legend=True)

denmark.plot(ax=ax, color="none", edgecolor='white', facecolor="None")
denmark_cover.plot(ax=ax, edgecolor='white', alpha=1)

I need to make the box about 5% bigger to cover all unwanted data.

What is a good aproach to make a shapely box bigger?


texts: How can I find the average bearing using a set of zig-zaggy coordinates, knowing the user travelled in a straight line during that period?

texts:
Some background:

I'm developing an app which uses augmented reality and GPS data in order to improve the True North accuracy of the device.

Currently, iPhones have a True North accuracy of about 15º at best, which is fine for displaying a heading on maps, but not good enough when you want to overlay things onto the world using Augmented Reality. As the AR "scene" can be aligned to gravity and True North, I thought I could use the GPS updates as the user moves through the world in order to improve the True North accuracy, down to perhaps within 1º of accuracy.

Each time a new coordinate is received, about once per second, I store it along with the user's current position within the AR scene (measured on a x/y/z axis). After 10 seconds, if the user has moved in a straight line, I'm able to compare the user's position 10 seconds ago with their position now, in order to determine their bearing, and I can compare data points across that span of time to ensure they walked in a straight line.

Knowing that they did walk in a straight line for 10 seconds, and their bearing, I can theoretically take their GPS coordinates from the readings across the same period and take that bearing. Thus, I can then measure the angle of difference between the user's bearing and the GPS bearing, and know how inaccurate the user's bearing is.

The problem is that the GPS data is a lot more zig-zaggy. Here's a screenshot from a period when I was walking in a straight line, along a pavement.

The green circles on the map represent each recent location estimate and their range of accuracy (in this case about 25m). The yellow line is the route between each of those GPS points. The other data shown can be ignored - green icon is the most accurate GPS estimate, and the orange icon is where the user is presumed to be based on that estimate and their movement since then, with the blue line just being a line between those points.

I've used a few techniques to try and get the true bearing of the GPS path, but because of inaccurate data points, like the point where it suddenly jumps north incorrectly, that throws off the average bearing.

The question:

So, the data I have:

I know the user has moved in a straight line for the past 10 seconds (and indeed in this example for more than 15 seconds before that)
I have a coordinate for the user's location, taken every second

The 25 most recent coordinates:

lat: 51.481508456031484, long: -0.014432226002130824
lat: 51.481507477162651, long: -0.0144049132597354
lat: 51.481506884611207, long: -0.014389807458603288
lat: 51.481509960771142, long: -0.014373092399400777
lat: 51.481507855339153, long: -0.014348176495456588
lat: 51.481509945317967, long: -0.014317286536546627
lat: 51.481505189244366, long: -0.014294481275542248
lat: 51.481506703609959, long: -0.014272413434476617
lat: 51.481505235807262, long: -0.01425324957115618
lat: 51.481498779846952, long: -0.014239295878236219
lat: 51.481501712942894, long: -0.014214523322743366
lat: 51.481508305052301, long: -0.014194804404044298
lat: 51.481512575463469, long: -0.014166769148728857
lat: 51.481519965341164, long: -0.014157078200848361
lat: 51.481521638215781, long: -0.01413473757481863
lat: 51.481520698352988, long: -0.014114195453429512
lat: 51.481518763912348, long: -0.014092665529632791
lat: 51.481516373462448, long: -0.014073392691598534
lat: 51.481517130888442, long: -0.014052945595763552
lat: 51.481522471260533, long: -0.014031560984244714
lat: 51.481519806277156, long: -0.014018366235482322
lat: 51.481521076237108, long: -0.0140037919026139
lat: 51.481522270228304, long: -0.013984366007679964
lat: 51.48152335719741, long: -0.013964033303506903
lat: 51.481522384989098, long: -0.013942519179059738

Goal: Use 10 of these points to find an average bearing of the route, which isn't made inaccurate by the outlying results.


texts: To check if the area within circle radius has any feature of layers from GeoServer

texts:
I have used google map to show the WMS layers from the GeoServer. Now what I am trying to do is, I have drawn a circle on a Google map. Whenever I change the radius of the circle I want to check if within that radius have any feature of layers from the GeoServer.
I am doing this project in angular 9 using agm.
In order to achieve this I have done the following steps.

On changing the radius I got circle bound

    function a(){
            circle.getBounds().then(bounds => {
            var bbox = {
                x1 : bounds.getNorthEast().lat(),
                y1 : bounds.getNorthEast().lng(),
                x2 : bounds.getSouthWest().lat(),
                y2 : bounds.getSouthWest().lng()
            }
            })
        }

I send them to server, there I run the below URL to get the list of layers on the server
http://localhost:8080/geoserver/wms?service=wms&version=1.1.1&request=GetCapabilities
By this URL I got bbox for each layers.
Then I use below code to check if the two bbox overlaps each other
function overlaps(a, b) {
    // no horizontal overlap
    if (a.x1 >= b.x2 || b.x1 >= a.x2) return false;

    // no vertical overlap
    if (a.y1 >= b.y2 || b.y1 >= a.y2) return false;

    return true;
}

This is what I got in capabilities URL
<EX_GeographicBoundingBox>
<westBoundLongitude>72.246193</westBoundLongitude>
<eastBoundLongitude>167.972839</eastBoundLongitude>
<southBoundLatitude>-67.013055</southBoundLatitude>
<northBoundLatitude>-9.799550</northBoundLatitude>
</EX_GeographicBoundingBox>
<BoundingBox CRS="CRS:84" minx="72.246193" miny="-67.013055" maxx="167.972839" maxy="-9.799550"/>
<BoundingBox CRS="EPSG:4326" minx="-67.013055" miny="72.246193" maxx="-9.799550" maxy="167.972839"/>
<BoundingBox CRS="EPSG:3857" minx="8042409.436216" miny="-10159779.218526" maxx="18698650.892001" maxy="-1096238.744138"/>
<BoundingBox CRS="EPSG:4283" minx="-67.013055" miny="72.246193" maxx="-9.799550" maxy="167.972839"/>

But still this code is not working. What I am doing wrong or there is any other way to do this?


texts: Table Join only returns first record from CSV table in QGIS 3.0.3

texts:
I have a .csv file that contains instances of hospitalization with diagnosis codes and zip codes and I want to join it to the corresponding polygons in the TIGER 2017 zcta510 shapefile. Every time I attempt to do this by selecting the Joins button in Layer Properties...then +...Join layer is my .csv import...join field is Zip...Target field is GEOID10.

Anyway, it seems to be only importing the initial instance when each zip code should have multiple instances. I am trying to use the symbology tools to create a shaded view of where patients are coming from.

I have a feeling that I'm missing something simple...but I've over-thought it now.


texts: Getting the departure time based on a desired arrival time for driving

texts:
Context:
I'm working on a logical problem as part of my Msc. thesis in GIS. I want to calculate the departure time (output) based on a desired arrival time (input). I have chosen to use an API to outsource much of the coding and maintenance.

Current Problem:
Tried to use the Google Distance Matrix and Google Diretions API's. Both do not allow a travel mode of driving combined with an arrival time as input.

Question:
What API can be used with the (batch size 1 - 50) input of a desired arrival time (/datetime) and a transport mode of driving to give an estimated departure time/travel duration taking into account traffic conditions such as rush-hour?

Related Posts:

API - Travel time based on arrival time | No clear answer is given here

Google Maps Directions API - Arrival time? | 2500+ views on stackoverflow since 2.5 years, no clear answer yet. One of the comments mentions the HERE API, but their documentation suggests it is not possible. See the note in the below picture.

After Taras Comment: Looking for free web service that calculates driving distance between 2 addresses | Question is about driving time, not driving time based on desired arrival time per se. Bing maps, google maps, HERE maps, and arcgis are recommended. As far as I have been able to find, these do not allow for specifying the 'desired arrival time' (except for public transport).

REST calculateroute want to supply arrival time to determine when to depart | Question from 2012 about the same problem, but specific for the HERE maps API. Answered on 4 oct 2018 (6 year later?) that only public transport modality is supported for now.


texts: Creating a layer from another layer's extent

texts:
Using QGIS, I have to map the general area of my project. However, this area depends on another layer representing sub-areas.

I know that I can use the dissolve tool on the second layer to get the extent, which gives me the exact result I am looking for.

However, as the second layer will undergo further modifications which might impact its extent, I was wondering if there was a way to create some kind of a "virtual layer" that represents the extent of the sub areas from the second layer? I wish to avoid dissolving the layers over and over again.

I can't use a simple style trick, as my main extent has to be surrounded by a dash line (management requirement).


texts: Where to find documentation on algorithm parameters when writing scripts in QGIS 3?

texts: For older versions of QGIS, the algorithm parameters were documented in the 'Console usage' section in the QGIS documentation (e.g. the reproject layer documentation). In the QGIS 3 documentation (e.g. reproject layer documentation) there is no 'Console usage' section. Where can I find out how to use the algorithms in my scripts, especially since many of them have changed?


texts: ArcMap 10.6: Is it possible to extract point values from points?

texts:
I've got three point shapes with different fields that I cannot seem to merge as I get ERROR 000229: cannot open 'name of shapefile'. Instead of solving this error, I wanted to try something different: in edit mode, I selected and copied all points into a new shapefile. However, all attributes got lost.

Now I would like to somehow migrate the values from the field 'TOP', which is present in all source shapefiles, to the new shapefile.

I know there is a possibility to extract values from a raster dataset ('Extract Values to Points'). Is it also possible to extract values FROM points, to points in a different shp if they are in the exact same location?

EDIT: In a comment below, I was advised to try Spatial Join which by the look of it would suit my needs perfectly, except that I get the same ERROR 000229.
I guess I just have to solve my initial problem.


texts: I can get a image to overlay and display just fine, but can't seem to get a popup on hover or click. Here is my code for version 1.6.0:


texts: Vector and Processing toolbars not available

texts:
Using Mac OS Yosemite 10.10.5.  Installing QGIS 3.2.2 Bonn

The vector and processing toolbars are not available after installation. Also, because processing is a core plugin, I am not able to add it via the plugins toolbar.

I imagine this has something to do with installation settings or possibly my system.


texts: Excluding some columns in the attribute table for the public use

texts:
I am trying to post a feature service for an Open Data Hub site but want to exclude some columns in the attribute table for the public use.

Is there a way to use a feature service but alternate the attribute table to exclude information such as 'moduser' 'moddate' or 'idnumber'?

These are things that we don't want the public to see, but want to be publishing our data to our open data hub site. Without creating a new service, is there a way to edit the attribute table without deleting it from the original source?

Note: Hiding columns does not work for when the dataset is being downloaded by a public user.


texts: Making selection with ModelBuilder using parameter

texts:
I am using ArcMap 10.5.
I am trying to calculate the volume and area of a pit using ModelBuilder (I need a non-python solution). The pit will by preselection have the largest volume of any pits in the selected area, but not necessary the largest area (the no-change area will usually be larger). I can isolate the volume of the pit in the CutFill raster using Summary Statistics (yielding a max_volume value).
I am having trouble configuring a selection to use the max_volume result of Summary Statistics to select the entry in the CutFill attribute table with that value. I cannot figure out how to get Select to recognize the max_volume variable for an SQL query. I tried %max_volume%, but no go.
In the attached image, I would be happy using either the extract by attribute or select tool to select from the cut_fill_raster raster using the max_volume value from Max_volume_value. Question marks indicate where I am adrift.


texts: With the same code the Python script won't work, while Python window will

texts:
I'm just trying to create a tool that zooms in on a layer and then exports the image as PDF, I'll add more functionality later. The code below will run and produce a PDF, but won't zoom to the layer. While running the same basic code in the Python window inexplicably works just fine. I've tried using zoomToSelectedFeatures and definitionQuery, both yield the same result. So now I'm rather lost.

import arcpy, csv, os
pole=arcpy.GetParameterAsText(0)
worksheet=arcpy.GetParameterAsText(1)
out=arcpy.GetParameterAsText(2)+"\_"
mxd=arcpy.mapping.MapDocument("CURRENT")
df=arcpy.mapping.ListDataFrames(mxd)[0]
i=0
worksheet=r'C:\Users\Samuel\worksheet.csv'
with open(worksheet) as f:
    for row in csv.reader(f):
        i=row[0]
        text=row[1]
        arcpy.MakeFeatureLayer_management("pole","lyr","FID="+str(i))
        arcpy.SelectLayerByAttribute_management("lyr","NEW_SELECTION")
        df.zoomToSelectedFeatures()
        arcpy.SelectLayerByAttribute_management("lyr","CLEAR_SELECTION")
        arcpy.RefreshActiveView()
        arcpy.mapping.ExportToPDF(mxd, out+"_pole_"+i,"PAGE_LAYOUT")
        break
arcpy.Delete_management("lyr")


texts: I need to align several images using FOSS. There is Correlator tool in GDAL (source) that seems to be suitable for the job. But I was unable to find documentation for it or a tutorial. Can someone explain what does it actually do and how to use this tool? The only help I could find is in this document


texts: Problem getting value from arcpy script tool with validation code driven drop-down: code window popup on tool execution

texts:
I am attempting to provide the user a drop down list of field values based on a predefined shapefile or table. The string value users choose then needs to be taken by the script tool for use later in the script.

I have the 'drop down' value list dynamically populated using the following validation code:

def updateParameters(self):
   self.params[2].filter.list = [str(val) for val in sorted(set(row.getValue('FC_field') for row in arcpy.SearchCursor("FC_path",None, None, 'FC_field')))]
   return

I then have the following script tool code which attempts to retrieve and store that value:

#---Import system modules.---
import arcpy

ProjectName = arcpy.GetParameterAsText(0)
ProjectBoundaries_FC = arcpy.GetParameterAsText(1)
Airport = arcpy.GetParameterAsText(2)

arcpy.AddMessage("INPUT: Project Name (0): <" + ProjectName + ">.")
arcpy.AddMessage("INPUT: Project boundaries (1): <" + ProjectBoundaries_FC + ">.")
arcpy.AddMessage("INPUT: Project Airport (2): <" + Airport + ">.")

The tool Parameters are set according to the image:

When the script tool is executed the text editor window pops up (illustrated below). Once I close it the tool processing bar stops and the rest of the code which would utilize that value never runs. In other words, when I close the editor window the script tool reports processing completed without errors, but never processes the "AddMessage" code (illustrated in the following image).

I suspect this is largely due to my unfamiliarity with using validation code.

As requested, here is a screenshot of the tool window showing the tool parameter options in use.

I am using ArcGIS v10.1 for Desktop - Advanced (ArcInfo) License

I think the most interesting thing here is why the tools script pops up in the text pad editor window at all?


texts: QGIS: GRASS GIS plugin not accepting UTM point coordinates

texts:
I am using QGIS 3.16 with the GRASS plugin.
I imported a UTM point cloud (x,y,z) and successfully generated a TIN surface based on it (using QGIS own functionality).
But when I want to use GRASS v.surf.rst (using the same input points layer) I get the error:
C:\Users\mk\Documents>g.region n=2132823.133 s=2132590.681 e=439329.924 w=439144.2887 res=100.0
ERROR: Illegal latitude for North: 2.13282e+006
...
WARNING: Some points outside of region (ignored)
WARNING: There are points outside specified 2D/3D region - 663 points ignored
WARNING: Zero points in the given region
ERROR: Input failed

Apparently GRASS has a problem with the coordinate system (UTM).
What can I do?


texts: Integrate Leaflet / OpenStreetMap autocomplete functionality to an existing search form on a WordPress template

texts:
I am trying to add the autocomplete functionality similar to what you have here https://community.algolia.com/places/examples.html#advanced to the search form in my website. Unfortunately there are no pre-existing plugins for Leaflet https://leafletjs.com/ & OpenStreetMap maps integrations like you have for Google Places autocomplete for WordPress templates.

Is there an easy way of integrating in existing template searchforms like you have in the examples here?

Leaflet address search/geocoding
http://esri.github.io/esri-leaflet/examples/geocoding-control.html

Looking for some suggestions.

The page I need to change: http://demo.ohmyproperty.com/


texts:
So I was wondering if it's possible to convert the KML files to an CSV attribute table (of those geometric features).

So I've found this code that converts the KML to CSV, but really just prints the tags there in a CSV :

import csv
from lxml import etree

out = raw_input("Name for output file:  ")
if out.strip() is "":
  out = "trader-joes-all-locations.csv"

out_data = []

parser = etree.XMLParser(recover=True, remove_blank_text=True)

file_name = "trader-joes-all-locations.xml"

root = etree.parse(file_name, parser)

tag_list = [ "name", "address1", "address2", "beer", "city", "comingsoon", "hours", "latitude", "longitude", "phone", "postalcode", "spirits", "state", "wine" ]

out_data.append(tag_list[:])

def missing_location(p):
  lat = p.find("latitude")
  lon = p.find("longitude")
  if lat is None or lon is None:
    return True
  else:
    return False

def get_poi_info(p):
  # if latitude or longitude doesn't exist, skip
  if missing_location(p):
    print "tMissing location for %s" % p.find("name").text
    return None
  info = []
  for tag in tag_list:
    # if tag == "name":
    #   print "%s" % p.find(tag).text
    node = p.find(tag)
    if node is not None and node.text:
      if tag == "latitude" or tag == "longitude":
        info.append(round(float(node.text), 5))
      else:
        info.append(node.text.encode("utf-8"))
        # info.append(node.text.encode("ascii", "ignore"))
    else:
      info.append("")
  return info

print "nreading xml..."

# get all <poi> elements
pois = root.findall(".//poi")
for p in pois:
  poi_info = get_poi_info(p)
  # print "%s" % (poiInfo)
  if poi_info:
    out_data.append(poi_info)

print "finished xml, writing file..."

out_file  = open(out, "wb")
csv_writer = csv.writer(out_file, quoting=csv.QUOTE_MINIMAL)
for row in out_data:
  csv_writer.writerow(row)

out_file.close()

print "wrote %sn" % out

There is a copy of the KML here ( KML ), and what I'm wondering is how to extract this metadata info and convert it to an attribute table, in csv, similar to this [Image][2]. I've searched a little bit and it seems the site http://www.gpsvisualizer.com/ does create a text table, containing only geometric data (long/lat).

Could you guys help me out here?

Thanks.


texts: How to change the field value of a specific feature in a Shapefile using GDAL/OGR?

texts:
I am trying to change the field value of a shapefile, but only in some specific features. However the code I have changes the value for all features in the shapefile. Is there any way of accessing the desired feature(s) in a simple way? This is what I have so far:

import ogr

driver = ogr.GetDriverByName('ESRI Shapefile')
file = '/home/src.shp'
dataSource = driver.Open(file, 1)
layer = dataSource.GetLayer()
feature = layer.GetNextFeature()
los = 1

while feature:
    feature.SetField("LOS", los)
    layer.SetFeature(feature)
    feature = layer.GetNextFeature()

dataSource.Destroy()


texts: Storing points from multiple GDBs as object

texts:
I often download data that comes in 30 minutes to 1-hour windows and then have to compile that data manually into one feature class, from multiple GDBs. I am trying to iterate through the separate GDBs and store each point feature class in an object. Here's what I have so far (Python 2.7.14):

import arcpy
import os
from arcpy import env
outputOverwrite = True

print('imported')

arcpy.env.workspace = r'C:\Vector\20190902'
print('ws set')

workspaces = arcpy.ListWorkspaces(workspace_type = 'FileGDB')

for workspace in workspaces:
   datasets = arcpy.ListDatasets(feature_type = 'All')

for ds in datasets:
   arcpy.env.workspace = ds
   fc = arcpy.ListFeatureClasses(feature_type = 'Point')
   print(fc)

I am unable to post pictures or copy and paste my code, as the system is on a totally different network.

Edited to reflect changes made by @BERA.

This code outputs

imported
ws set
>>>

It doesn't print the feature classes like i expected it to.

I traced this back to the line of code

for workspace in workspaces:
    datasets = arcpy.ListDatasets(feature_type='All')

When I tried to do

print(datasets)

It gave me blank lists

[]
[]
[]
[]
...

Hope this clears up the ambiguity on my end.


texts: Using QGIS, how can I click on a route and obtain the measure or distance for that clicked point?

texts:
I'm looking for the QGIS-equivalent of ArcMap's Identify Route Locations tool

I do use the LRS plugin, but I don't see this functionality in that plugin. Thank you.


texts: Updating a field with something like "if feature_1 is in field_1 then write category_1 to field_2" in QGIS Field Calculator?

texts:
I'm grouping different geological time scales into more general time scales.

An example of what I would want is

"if Oligocene-miocene is in descript_title, write Oligocene-Miocene in descript_category".

How would I write this as an expression?

As suggested by @ahmadhanb, I could extract a specific part of the feature to the target field.

However, how should I write this code properly? The goal is to do what @ahmadhanb but there is a specific line for each feature.

IF 'Basement Complex (pre-jurassic)' IN descript_title
substr("descript_title",0,17)


texts: MapServer: How to display and style a Text layer?

texts:
I'm trying to display a text layer using MapServer, but I'm kind of confused how to get it to display properly.

When opening the .TAB file in a GIS application (in this case, MapInfo Profressional) the labels show up correctly. They all show up in the correct location with varying sizes and angles.

When looking at the output of MapServer for this layer, there are the following problems:

The layer is blank unless I set STYLEITEM "AUTO" for the layer. If I try to set the style manually, the settings seem to be ignored.

When displaying the layer using STYLEITEM "AUTO", all of the labels are rendered in the same font size and with no rotation.

Excerpt from the .map file is below (I've excluded other layers):

        MAP
        NAME    mymap
        IMAGETYPE   PNG
        EXTENT          259997 100002 399997 160002
        SIZE        550 450
        SHAPEPATH   "./"
        CONFIG "MS_ERRORFILE" "/var/log/mapserver/mymap_error.log"
    WEB
            LOG '/var/log/mapserver/mymap.log'
            METADATA
            WMS_TITLE   "My WMS Title"
            wms_srs     "EPSG:27700"
            wfs_srs     "EPSG:27700"
            wfs_title   "WFS Test"
            wfs_enable_request  "*"
            ows_enable_request  "*"
            wfs_onlineresource  "http://[AN IP]/maps/vectors/mymap&"
            ows_onlineresource  "http://[AN IP]/maps/vectors/mymap&"
            wfs_maxfeatures "50"
        END #METADATA
    END #WEB

    PROJECTION
        "init=epsg:27700"
    END #PROJECTION

LAYER
        NAME "Mastermap_Text"
        DEBUG 2
        METADATA
                "wfs_typename"          "Mastermap_Text"
                wfs_title               "Mastermap_Text"
                wfs_srs                 "EPSG:27700"
                gml_include_items       "all"
                gml_featureid           "PRN"
                wfs_enable_request      "*"
                ows_enable_request      "*"
        END #METADATA

        TYPE ANNOTATION
        STATUS OFF
        DUMP    TRUE
        CONNECTIONTYPE OGR
        CONNECTION "/path/to/maps/vectors/Mastermap_Text.tab"
        STYLEITEM "AUTO"

        PROJECTION
            "init=epsg:27700"
        END #PROJECTION

        CLASS
            NAME "Mastermap_Text"
        END #CLASS
END #LAYER

When trying to view this layer, I'm using an URL similar to:

http://[AN IP ADDRESS]/maps/vectors/mymap?LAYERS=Mastermap_Text&FORMAT=image%2Fpng&SERVICE=WMS&VERSION=1.1.1&REQUEST=GetMap&SRS=EPSG%3A27700&BBOX=328800,146642,329680,147522&WIDTH=256&HEIGHT=256

As states, the above .map file results in an image of the text labels being generated, but they are all the same font size and they have no rotation.

If I remove the STYLEITEM "AUTO" line and add a STYLE block to the CLASS (with COLOR and SIZE attributes), I get a blank image.

Similarly, adding a LABEL block with the respective styles also doesn't help.

What am I doing wrong?

Edit to add images:

Mastermap_Text layer viewed using MapInfo Professional

Mastermap_Text layer rendered by MapServer when using the above URL

Mastermap_Text layer generated by shp2img: shp2img -m mymap.map -o test.png -l Mastermap_Text


texts: What is the best way to display point data (with standard deviation) using QGIS/SAGA/Manifold?

texts:
I would like to display point data values (indicating biomass at a point) but also to show the standard deviation value for each point. What would be the best method to do this in QGIS (I'm a basic user!). I do not have ±100 points. I am imagining a point with size = biomass and a buffer showing SDev, but I don't know how to do this.
The next step would be to show the breakdown of values in each point by species (about 5 species biomass values making up the total biomass at each point). Would the best way be with a pie graph per point? And how to show each StDev?

I am using QGIS, Manifold and (less often) SAGA. Any suggestions welcome, many thanks


texts: Getting local incidence angle map from SAR-Imagery?

texts: I have a dual pol TanDEM-X imagery with a bi-static acquisition mode. To calculate the vertical wavenumber, I need a local incidence angle map of the Master and Slave image. As a first step I co-registered both images and then I proceed a Multilooking. The next step was a terrain correction (with SNAP) or the geocoding and radiometric calibration procedures with SarSCAPE to get the local incidence angle maps. After I check both results (retrieved from SNAP and SarSCAPE) I found that both local incidence angle maps from the Master and Slave imagery are identical what surprised me very much. Do you have an Idea what could be the reason for this result? What should I do?


texts: Using Color Ramp Manager Plugin in QGIS?

texts:
I followed instructions mentioned here and installed the Color Ramp Manager by choosing install experimental plugins as well.

Now, when I go to Properties > Style > Single Band PsuedoColor and I choose New Color Ramp > Cpt-city > This empty box opens up.

How do I hillshade my region?