texts: Can't install GitBook 3.2.3 from Geoext3 workshop

texts:
I am trying to setup the Geoext3 workspace based on this tutorial for web mapping aplication:
https://github.com/geoext/geoext3-ws
But i get this npm error after starting 'npm start', it can't install GitBook 3.2.3 and generate the log below:
 0 info it worked if it ends with ok
1 verbose cli [
1 verbose cli   'C:\\Program Files\\nodejs\\node.exe',
1 verbose cli   'C:\\geoext\\geoext3-ws-master\\node_modules\\npm\\bin\\npm-cli.js',
1 verbose cli   'run',
1 verbose cli   'doc:serve'
1 verbose cli ]
2 info using npm@5.1.0
3 info using node@v14.16.0
4 verbose run-script [ 'predoc:serve', 'doc:serve', 'postdoc:serve' ]
5 info lifecycle geoext3-ws@1.0.0~predoc:serve: geoext3-ws@1.0.0
6 info lifecycle geoext3-ws@1.0.0~doc:serve: geoext3-ws@1.0.0
7 verbose lifecycle geoext3-ws@1.0.0~doc:serve: unsafe-perm in lifecycle true
8 verbose lifecycle geoext3-ws@1.0.0~doc:serve: PATH: C:\geoext\geoext3-ws-master\node_modules\npm\bin\node-gyp-bin;C:\geoext\geoext3-ws-master\node_modules\.bin;C:\Program Files\nodejs\node_modules\npm\node_modules\npm-lifecycle\node-gyp-bin;C:\geoext\geoext3-ws-master\node_modules\.bin;C:\Users\Administrator\AppData\Roaming\npm;C:\Program Files\nodejs\;C:\Python39\Scripts\;C:\Python39\;C:\Program Files (x86)\Common Files\Oracle\Java\javapath;C:\app\Administrator\product\11.2.0\dbhome_2\bin;C:\Windows\system32;C:\Windows;C:\Windows\System32\Wbem;C:\Windows\System32\WindowsPowerShell\v1.0\;C:\Python39\Scripts\;C:\Python39\;C:\Program Files (x86)\Common Files\Oracle\Java\javapath;C:\app\Administrator\product\11.2.0\dbhome_2\bin;C:\Windows\system32;C:\Windows;C:\Windows\System32\Wbem;C:\Windows\System32\WindowsPowerShell\v1.0\;C:\Python39\Scripts\;C:\Python39\;C:\Program Files (x86)\Common Files\Oracle\Java\javapath;C:\app\Administrator\product\11.2.0\dbhome_2\bin;C:\Windows\system32;C:\Windows;C:\Windows\System32\Wbem;C:\Windows\System32\WindowsPowerShell\v1.0\;C:\Windows\System32;C:\gdal1\bin;C:\gdal1\bin\gdal\java;C:\Program Files\nodejs\;C:\ProgramData\chocolatey\bin;;C:\gdal1\bin;C:\gdal1\bin\gdal\java;C:\Program Files\nodejs\;C:\ProgramData\chocolatey\bin;;C:\gdal1\bin;C:\gdal1\bin\gdal\java;C:\Program Files\nodejs\;C:\ProgramData\chocolatey\bin;C:\Users\Administrator\AppData\Local\Microsoft\WindowsApps;C:\Users\Administrator\AppData\Roaming\npm
9 verbose lifecycle geoext3-ws@1.0.0~doc:serve: CWD: C:\geoext\geoext3-ws-master
10 silly lifecycle geoext3-ws@1.0.0~doc:serve: Args: [ '/d /s /c', 'gitbook serve src' ]
11 silly lifecycle geoext3-ws@1.0.0~doc:serve: Returned: code: 1  signal: null
12 info lifecycle geoext3-ws@1.0.0~doc:serve: Failed to exec doc:serve script
13 verbose stack Error: geoext3-ws@1.0.0 doc:serve: `gitbook serve src`
13 verbose stack Exit status 1
13 verbose stack     at EventEmitter.<anonymous> (C:\geoext\geoext3-ws-master\node_modules\npm\lib\utils\lifecycle.js:289:16)
13 verbose stack     at EventEmitter.emit (events.js:315:20)
13 verbose stack     at ChildProcess.<anonymous> (C:\geoext\geoext3-ws-master\node_modules\npm\lib\utils\spawn.js:40:14)
13 verbose stack     at ChildProcess.emit (events.js:315:20)
13 verbose stack     at maybeClose (internal/child_process.js:1048:16)
13 verbose stack     at Process.ChildProcess._handle.onexit (internal/child_process.js:288:5)
14 verbose pkgid geoext3-ws@1.0.0
15 verbose cwd C:\geoext\geoext3-ws-master
16 verbose Windows_NT 10.0.14393
17 verbose argv "C:\\Program Files\\nodejs\\node.exe" "C:\\geoext\\geoext3-ws-master\\node_modules\\npm\\bin\\npm-cli.js" "run" "doc:serve"
18 verbose node v14.16.0
19 verbose npm  v5.1.0
20 error code ELIFECYCLE
21 error errno 1
22 error geoext3-ws@1.0.0 doc:serve: `gitbook serve src`
22 error Exit status 1
23 error Failed at the geoext3-ws@1.0.0 doc:serve script.
23 error This is probably not a problem with npm. There is likely additional logging output above.
24 verbose exit [ 1, true ]


texts: Is it possilbe to change SLD from client side?

texts: How to change predefile SLD in client side? The reason I want this is I want to turn label on and off on client side.


texts:
I'm trying to import geometry lines from Oracle Spatial database to ArcGIS enterprise geodatabase.

and some data will cause the error.

Its says not enough points for feature class, but when I look at offending object.

it clearly has 2 or more coordinates stored, but I notice the start and end point is the same.

I have check the following

a) Using geometry filter to take out lines only

b) Use FME inspector to view the lines.

b) Used geometry validator

4) Used duplicate coordinate -> Valid dataset will cause the same error

How can I solve this problem?

May i know is it possible to set FME to skip those data


texts: GDAL crashes with segmentation fault when trying to open a NetCDF subdataset with time dimension

texts:
I am trying to create a NetCDF file which is compliant with the NetCDF Conventions CF-1.4.

The dimensions of the file are (TIME, Y, X) where Y and X are the polar stereographic grid.

The file that I generate contains three variables: runoff (TIME, Y, X), lat (Y, X) and lon (Y, X). With gdalinfo I can access lat and lon, but not runoff:

$ gdalinfo NETCDF:"myfile.nc":runoff
Segmentation fault (core dumped)

In other words it looks like the TIME dimension is causing problems. I've tried writing all the data in the runoff attribute as zeroes but that didn't help.

The file was created using xarray, with output set to NETCDF4. The ncdump of the file is as follows. Is there anything strange about it that might be causing GDAL to have problems?

netcdf myfile {
dimensions:
    X = 756 ;
    TIME = 12 ;
    Y = 786 ;
variables:
    int64 X(X) ;
        X:units = "meters" ;
        X:standard_name = "projection_x_coordinate" ;
        X:point_spacing = "even" ;
        X:axis = "X" ;
    int64 TIME(TIME) ;
        TIME:standard_name = "time" ;
        TIME:axis = "TIME" ;
        TIME:units = "days since 1958-01-01 00:00:00" ;
        TIME:calendar = "proleptic_gregorian" ;
    int64 Y(Y) ;
        Y:units = "meters" ;
        Y:standard_name = "projection_y_coordinate" ;
        Y:point_spacing = "even" ;
        Y:axis = "Y" ;
    double runoff(TIME, Y, X) ;
        runoff:_FillValue = -9999. ;
        runoff:long_name = "runoff" ;
        runoff:units = "mmWE" ;
        runoff:standard_name = "runoff" ;
        runoff:grid_mapping = "polar_stereographic" ;
    double lon(Y, X) ;
        lon:_FillValue = -9999. ;
        lon:grid_mapping = "polar_stereographic" ;
        lon:units = "degrees_west" ;
        lon:standard_name = "longitude" ;
    byte polar_stereographic ;
        polar_stereographic:grid_mapping_name = "polar_stereographic" ;
        polar_stereographic:scale_factor_at_central_origin = 0. ;
        polar_stereographic:standard_parallel = 70. ;
        polar_stereographic:straight_vertical_longitude_from_pole = -45. ;
        polar_stereographic:latitude_of_projection_origin = 90. ;
        polar_stereographic:false_easting = 0. ;
        polar_stereographic:false_northing = 0. ;
    double lat(Y, X) ;
        lat:_FillValue = -9999. ;
        lat:grid_mapping = "polar_stereographic" ;
        lat:units = "degrees_north" ;
        lat:standard_name = "latitude" ;

// global attributes:
        :Conventions = "CF-1.4" ;
        :history = "redacted" ;
        :institution = "redacted" ;
        :title = "redacted" ;
        :source = "redacted" ;
        :proj4 = "+init=epsg:3413" ;
        :nx = 756. ;
        :ny = 786. ;
        :xmin = -1790000. ;
        :ymax = -65000. ;
        :spacing = 5000.


texts:
I want to assign values to a polygon shape, according
to another polygon shape taking the closest polygon of the last shape.
I´d prefer to use Arcgis but any software is ok.


texts: Adding new nodes or moving segments with Vertex Tool in QGIS 3.10

texts:
I'm used to work with simple polygons and the vertex tool. But now the Vertex Tool doesn't show me the capability of adding new nodes in centerline neither move polygon's segments. The function just disappeared.

Now the vertex tool only gives me the option to move the existent vertices. There's no more the "+" symbol giving me the option to add it. But, otherwise, it still gives me the option to add a new vertex while at the end of a line:


texts: How to create lines between pairs of two coordinates?

texts:
I have a table with pairs of points. Out of that I want to create a line layer (that keeps all attributes of the table). One line for each raw from point to point. I want to calculate the length of these lines later.

A simplified example table would look like that:

Can I make this in QGIS?

I am using QGIS 1.8.0 on a Mac 10.6.8


texts: How to keep color changes / corrections of all GeoTIFF files during merge with gdal or in QGIS?

texts:
I am trying to merge two GeoTIFF files in QGIS or w/ gdal, but when I merge them I lose the color changes / corrections I made to the second GeoTIFF.

Is there any way to keep the color changes I made to the GeoTIFFs after merge?

I clipped a region of the first GeoTIFF that needed color correction and then made the color changes, then when I merge them back together I lose the changes.

I have tried using QGIS merge tool as well as gdal_merge.py & gdal warp to achieve the desired result, but always the same issue is present.

The clipped color corrected raster is last in the merge order, so I am sure it is being merged over the top of the original.

Any ideas?

Original:

Clipped & color corrected:

Both together, clipped on top:

After merge

Here's a link to download the files I am working on.

https://drive.google.com/open?id=17wyRKbkDE1twWBSMDsY2upY7cuoHmOq1

Here is a similar question, but I don't understand how I could apply the answer to this question to my issue....

QGIS: Is it possible to merge two stacked raster layers keeping color rendering and transparency?


texts: Layer does not show in both versions of QGIS

texts: I was originally using QGIS 3 and i recently download QGIS Desktop 3.8.2 with Grass 7.6.1. All of the layers can be seen in both except a layer called Masspoint. This layer is not showing up in QGIS Desktop 3.8.2 with Grass. when I check the properties, it shows the feature count but does not show the points. Why does this show up in QGIS3 but not in QGIS Desktop with Grass? All other layers are visible.


texts: Adding large numbers of labels using QGIS

texts:
In 2007, when I was young and foolish and before I knew about OpenStreetMap, I started an urban historical map project. I was working in Illustrator, it was going to be an interactive Flash piece, and my process was to draw the maps first, with the thought that I'd label some, but not all, of the street later on.

As we know Flash was began to die about 2010 and I put the project away for a number of years. I picked it up again a couple years ago and continued my earlier practice of just drawing streets and water features, this time with the intention of making it a conventional web map. Now I'm pretty close to finishing the drawing of a five-layer (1871, 1903, 1932, 1952 and 2016) historical map of a medium-sized city, though it still lacks labels.

My problem now is how to add large numbers of labels, many of them duplicates. There could be as many as 10,000 for all five layers, though as a practical matter I may have to settle for a smallish fraction of that number. Based on web searches I gather my workflow is unusual and that mine is therefore an unusual problem.

I've exported my maps and brought them into QGIS and played with the software a little. The process of adding labels to objects doesn't seem terribly efficient or user-friendly, but that's probably due to my unfamiliarity with the program.

Are there any tricks to speed up the painful process of adding large numbers of duplicate labels in QGIS?

Since so many of the streets exist in all five layers, functionality like the ability to select multiple objects in different layers and edit their attributes simultaneously in the Attribute Table would be a godsend. (Doesn't seem possible.) So would the ability to copy the attributes from one object and paste them onto other objects. Or the ability to do either of these things in Illustrator via a plugin and then export the data along with the shapes to a GIS program.


texts: ogr2ogr shapefile to PostgreSQL table | cannot change geom field

texts:
ogr2ogr -lco GEOMETRY_NAME=geom -f "PostgreSQL" "PG:host=localhost user=postgres dbname=lidar password=postgres" "C:/Users/Rose/Documents/faulty_rasters_test/raster_summary/gdaltindex/1.shp" -lco FID=gid -lco PRECISION=yes -nlt PROMOTE_TO_MULTI -nln lidar.test_dem

ogr2ogr -f "PostgreSQL" "PG:host=localhost user=postgres dbname=lidar password=postgres" "C:/Users/Rose/Documents/faulty_rasters_test/raster_summary/gdaltindex/1.shp" -lco GEOMETRY_NAME=geom -lco FID=gid -lco PRECISION=yes -nlt PROMOTE_TO_MULTI -nln lidar.test_dem

Used above to import shapefile into PostgreSQL database, using python.subprocess() and batch shell.
No sterr messages were outputted e.g.

The PostgreSQL table was created but the geom column was incorrect.
The geom column was still the default (wkb_geometry), even though it was specified in command.
Using PostgreSQL 11.8, pgAdmin 4.20, Python 3.8 e.g.

How do I change this in this command? Is this possible?


texts: Apply color ramp to (QgsSingleSymbolRenderer) Simple marker

texts:
I'm using a QgsSingleSymbolRenderer to display various points and lines etc for a Line geometry.

I have a Single symbol symbology to display the last vertex as a: Line > Marker line > Marker > Simple marker which I'd like to colour according to one of the pre-defined colour ramps (e.g. Blues, with range between 0 - 100, based on an attribute value "count"), using PyQGIS.

I can do this manually in QGIS by opening the layer's properties and going to > Symbology > Line > Marker line > Marker > Simple marker where I can use the colour assistant to set the symbols fill colour based on the value of "count".  Is it possible to do this using PyQGIS?  A snippet of my current code looks something like:

# Define symbol layers
symbol_layers = QgsSymbol.defaultSymbol(vector_layer.geometryType())
symbol_layers.changeSymbolLayer(0, simple_line_symbol_layer)
symbol_layers.insertSymbolLayer(1, last_vertex_simple_marker_symbol_layer)

# Create renderer based on symbol layers, and apply to vector layer
renderer = QgsSingleSymbolRenderer(symbol_layers)
vector_layer.setRenderer(renderer)


texts:
I am tracing building plans in adobe acrobat. There may be a time in the future I will be asked to take the polygons and overlay them in a GIS for a borough-wide map. Are there any steps I can take from the start in preparation for such a task.

The image below is a very simple plan I have traced. Many others are far more complex, and the colour scheme (and label within polygon) are important information I would like to retain.


texts: Creating and importing marker symbols (or 'Style References') in ArcMap?

texts:
I work with ArcMap 10 and many different federal lands all of which can be associated with iconic logos (think National Park Service, National Wildlife Refuge System, etc.).  I'd like to create a 'Style Reference' that has all of the different logos to denote point data.

This will be a working style reference as I encounter new features from time to time.

How can I maintain or edit a style reference?

I'm assuming the preferred format would be .emf (vector) to maintain transparency, but I'm not at all familiar with the style format structure or file type.


texts: How do I create an animation in google earth using png overlays with kml timestamps?

texts:
I have the following kml files, each with a different timestamp and png:

<?xml version="1.0" encoding="UTF-8"?>
 <kml xmlns="http://www.opengis.net/kml/2.2"
 xmlns:gx="http://www.google.com/kml/ext/2.2">
 <GroundOverlay>
 <name>gx:LatLonQuad Example</name>
 <Icon>
<href>
/Users/asd/Dropbox/0mylayerExport.png
</href>
 </Icon>
<gx:LatLonQuad>
<coordinates>
-71.007186998722,42.2869344810739 -70.9807964189347,42.360133347777 -71.1410513483055,42.3916360561188 -71.1680083772048,42.3178859041619
</coordinates>
<TimeStamp>
<when>1944-06-06T06:2010-10-26T06:30:08
></when></TimeStamp><rotation>100</rotation>
</gx:LatLonQuad>
</GroundOverlay>
</kml>

What I'm trying to do is to animate the sequence of png overlays, but I can't figure out how. I can view all of the kml files in google earth, and click through them manually, but I was wondering if there was a neater way of animating this, instead of ticking them on and off.

I have tried manually adding timestamps, and time intervals in google-earth, and am still failing to automatically animate overlays, so I would appreciate any advice on how to do this too.

Thanks!


texts: How can I import .tab file from mapinfo to QGIS with same color on symbol?

texts:
I am apprentice Drivetest Analyzer. I want to see color in QGIS like Mapinfo. This is my step.

1.I'm open log file from Drivetest Engineer by Actix

2.I'm export this file to .tab and open in Mapinfo

3.And add vector layer from .tab file to see in QGIS but the color dosen't same like Mapinfo

This is my question : How can I do to see .tab in QGIS with same color in Mapinfo


texts: Problem with showing WFS layer with OpenLayers

texts:
I'm using OpenLayers and GeoServer in combination with Apache for creating GIS-based WEB application.

I have created layers and styled them in GeoServer, but when I try to load those layers with OpenLayers as a WFS I don't get any result over the map. I read somewhere that HTML file should be added to the same host as GeoServer (for me it is localhost:82/geoserver).

I tried to put HTML file to "Tomcat 8.5\webapps\geoserver\data\www" but when I try to open http://localhost:82/geoserver/www/index.html I get HTTP Status 404 - description The requested resource is not available. I can open HTML files that were there before ol-demo.html and wfs-t.html.

Is there any suggestion what I should do?


texts:
I have received files in .qgd and .qgs extensions from a friend who works with QGIS.
Is it possible to use that in ArcMap 10.8?


texts: Has TIN been superseded by Terrain Dataset in ArcGIS?

texts:
I've had version 10 for a few months but it wasn't until today that I really tested the 3D Analyst extension. (Other then visualizing data)

I already noticed a lot of changes; a number of tools have been re-named or re-worked all together which also impacts my usual work-flow.

I will not list all the details but the question is whether the 'TIN format' has been superseded by the improved 'Terrain Dataset'? Should I be using terrain datasets now even for simple surfaces?

Example: Hardclip and Softclip functions do now work in version 10 (as far as I can tell) and I noticed that defining the edges of a TIN dataset should be accomplished by using the 'Deliniate TIN Data Area' tool instead. However, to do this Average Point Spacing is required. Strangely this needs to be obtained from the input XYZ file not from the TIN itself. Alternatively, knowledge of the data is needed to take some manual measurements and to estimate the average spacing. This information is also needed to create a Terrain Dataset so it might be easier to just migrate to Terrain Datasets all together. But, di functions such as 'Extrude Between' support Terrain Datasets?
TIN not processed with 'Delineate TIN Data Area' yields very strange contour results so it looks like things have to be done this way.


texts: How to correlate Image DN with Spectrometer reflectance?

texts:
I am trying to develop a system to measure NDVI from the commercial grade digital camera.

I am trying to get two bands NIR and RED from the image and apply it for NDVI measurement using diffent filters. I also want to apply the Spectrometer reading on field to validate the NDVI measured from the Images.

Can anyone help me, how to correlate spectrometer reading(Normally Reflectance of the surface) and the image DN(Three Layer Image) ??


texts: Print composer messes up layout of different maps

texts:
I made a map including 3 levels of zoom as shown in the screenshot. The problem is that whenever I save the map to a PDF, the layout of these maps is changed seemingly random.

I have a hunch that it has to do something with 'Lock layers' or 'Lock styles for layers' functions, but they don't really seem to work. Right now, I have the lock layers box checked on all maps and Lock styles for layers unchecked.


texts: Is it possible to generate Vector Tiles using the 2.2.0 Mapnik release? If so, please give a sample configuration or point to an example.


texts:
I have a dataset on rivers and a dataset on the location of dams. I want to calculate the distance from the dams to the river’s mouth at the coast along the river (not the shortest distance) for those dams located next to one of the rivers in the dataset.
These are the datasets:

var rivers = ee.FeatureCollection('users/basins/ne_10m_rivers_lake_centerlines_scale_rank');
var coast = ee.FeatureCollection('users/basins/GSHHS_i_L1');
var grand_dams = ee.FeatureCollection('users/basins/GRanD_dams_v1_3');

I guess the spatial join of the rivers and dams looks something like this:

// Define a spatial filter, with distance of 2 km.
var distFilter = ee.Filter.withinDistance({
  distance: 2000,
  leftField: '.geo',
  rightField: '.geo',
  maxError: 10
});

var saveAllJoin = ee.Join.saveAll({
  matchesKey: 'scenes',
});

var dam_river= saveAllJoin.apply(grand_dams, rivers, distFilter);

I don’t know where to start to solve the problem.


texts: GDAL/OGR SQL Can't select correct records (Python)

texts:
I have two attributes that contain a value of ' ', among others, which are not selected when I apply an SQL statement in my python script. Trying the same statement in ArcMap does result in them being selected.

I have an attribute "BYGNSTYPE" which has as distinct values:

' ' <-
'Bygning' <-
'Drivhus' <-
'Ikke tildelt'
'Pelsdyr Opdræt' <-
'Tank/Silo'

The values with an arrow next to them are the ones I want. When I use the WHERE clause: BYGNSTYPE <> 'Ikke tildelt' AND BYGNSTYPE <> 'Tank/Silo', ArcMap returns: 4,517,327 records.

When I apply the clause in my script using GDAL it returns 4,514,343 records.

sql =  "SELECT * FROM BYGNING WHERE BYGNSTYPE <> 'Ikke tildelt' AND BYGNSTYPE <> 'Tank/Silo'"
selectData = inData.ExecuteSQL(sql)
print sql
if selectData is None:
    print 'No data is selected...'
    sys.exit(1)
else:
    print 'Finished selecting data'
print "Selected features: %s " % selectData.GetFeatureCount()

When providing the oposite statement it also does not select records where BYGNSTYPE = ' '. How can I select ' ' values in Python GDAL using SQL?

FIXED:

Adding 'OR BYGNSTYPE IS NULL' to the SQL statement in Python/GDAL results in the same records being selected.


texts: Finding all neighbors for all polygons in Spatialite

texts:
Excuse me if this is a basic question but I'm new to Spatialite...

I have a dataset with zip code polygons and I would like to create a list where for each zip code we have listed all the neighboring (i.e. adjacent) zip codes.

If I do the following it works but only for 1 zip code. What I would like to do is to repeat this automatically for all the zip codes in the data (since I have thousands of codes).

SELECT p1.ROWID, p1.ZCTA5CE10 FROM "US_zcta5_2010" p1, "US_zcta5_2010" p2
WHERE p2.ROWID=1 AND
Touches (p1.Geometry, p2.Geometry)=1;

I'm using the spatialite-gui on a Mac.


texts:
Using QGIS 3 I'm having difficulty with what seems to be a simple problem: for each US state and Canadian province with a coastline (on Shapefile maps I downloaded), replace the coastline segments with a minimum number of segments that will encompass any belonging islands.
Essentially, I want to delete all the coastline vertices except one on each end of the coastline, then extend the coastline by drawing a few line segments out into the ocean (or lake) that will surround the belonging islands, and delete the island vertices from the feature.
I have found very laborious and error-prone ways to do it, but nothing efficient.  Ideally, I'd like to select two points on a coastline and delete all the points in between with a command.  Selecting the in-between points is clumsy the way I do it, using the vertex tool, trying to "box" parts of the coastline and delete in pieces, zooming in and out, etc.  Is there a simple way to do this?
Once I have replaced the coastline with a single line segment, I'd like to stretch it out around the islands as a couple of additional segments.  I visualize this as dropping a couple of points into the coastline segment, then dragging those points outward.  Alaska has a bundle of islands.  It would be great to extend the coastline around the islands, then with some command deletes all the separate polygons (islands ) contained inside the big polygon. In any case, there are probably a number of ways to get the end result.


texts:
this my first question in this forum and I don't have much experience on QGIS

I have a layer with lines representing streets in the city and I have to show the distances of intersection points of that line/street from the beginning of it

How can I do this?


texts: Converting GPX track_points layer to CSV with lat/lng with ogr2ogr

texts: What's the correct command to convert a GPX track_points layer to CSV with lat/lng coordinates using ogr2ogr?


texts: GeoServer: Cross Layer Filtering vs Layer Groups

texts:
To pull features using WFS from several layers in one request, I have two options: One, create a layer group, and then get features from the layer group, as opposed to individual layers. Two, use cross layer filtering to obtain features for each desired layer.

Now, I'm not adept at cross layer filtering. But it appears to me as the more flexible and robust solution. What I wish to achieve is obtain points within a selected administrative boundary. The administrative boundaries and the points/population layers are separate.

I'll appreciate any advice on the way to go.

My environment: GeoServer 2.1.3, OpenLayers 2.11, PostGIS 1.5.


texts: Is it possible to save an OpenLayers map as a shape or tab file in QGIS. I would like to have a part of an area from Google Maps as a shapefile for later use?


texts:
I selected my climate data from wordclim ( tile 29 ) and then used QGis to convert the data
to .asc format . My occurrence data is in degree decimal.But when I run Maxent on any layer
it says that the occurrence points are out of environment bound.
What should I do ??


texts: How to georeference a tile downloaded from Google Maps

texts:
I'm trying to georeference a set of image tiles downloaded from a Google Maps site. Based only on each image's Google Maps grid reference and zoom level, how can I calculate the real-world coordinates?

The page at http://www.maptiler.org/google-maps-coordinates-tile-bounds-projection shows the coordinates and Google Maps grid reference for web mercator tiles.

For example, the centre tile is x=49, y=30, z=6 and has the corner coordinates 10644926.307106785, 626172.1357121654, 11271098.44281895, 1252344.27142432:

Is it possible to calculate the corner coordinates based solely on the grid reference and zoom level (49, 30, 6)?

This is similar to How to georeference a web mercator tile correctly using gdal? but I need a completely programmatic solution, without needing to look anything up online. (Or, if I do need to look up values online it needs to be programmatic).

Edit: this screenshot shows the location of a sample tile (zoom level 11, x = 1593, y = 933). How would I georeference that tile?


texts:
In asked question I have attached image which is some separator line on globe.

I  googled/stackoveflowed but did not find answer for following questions (for kid's home work).

What is the name of that line?
What is the significance of that line?
Why line shaped like hook? History behind it if any, any special consideration for shape of line.


texts: WFS layer has no features in OpenLayers?

texts:
I have a layer with many features and when I add it to a map as WMS layer it works.
But if I add the layer as WFS it seems there are no features.

The code I use to load WFS layer is:

          var layerName = '...';
          var workspace = '...';
          var vectorLayer = new OpenLayers.Layer.Vector('questolayer', {
            projection: new OpenLayers.Projection('EPSG:4326'),
            strategies: [new OpenLayers.Strategy.BBOX()],
            protocol: new OpenLayers.Protocol.WFS({
                version: "1.1.0",
                url: "/geoserver/wfs",
                featureType: layerName,
                featureNS: workspace,
                srsName: "EPSG:4326"
            })
          });
          vectorLayer.events.on({
            'loadend': function(event) {
                console.warn('00 features: ', vectorLayer.features);
            }
          });
          openlayersmap.addLayer(vectorLayer);

If I check the response with Firebug I can see an XML with gml:featureMembers tags and many features into of it.
But the warn message '00 features' in the log contains an empty vector.

Do you know why?


texts: convert ordnance survey easting and northing data to WGS84 latitude longitude values

texts: The C# code published here is not very accurate as also mentioned by tanoshimi. I think this has to do with the lack of 'approximate Helmurt transformation'. Just curious is anyone aware of some accurate C# code or is it even possible to achieve using SQL Server's spatial functionalities?


texts: Using PostGIS geometry for ArcSDE database?

texts:
I just enabled my PostgreSQL/PostGIS database with ArcSDE. Since I did this, when I copy and paste a dataset into my database via ArcCatalog, the dataset is not visible in QGIS.

This is because there is no geometry according to QGIS. I want to use ArcSDE database with all its benefits in ArcMap, but still be able to add this dataset in QGIS via a PostGIS connection.

How do I make sure that ArcGIS uploads the table with PostGIS geometry?

I almost have a solution by changing the ArcSDE table sde_dbtune to contain Default Geometry_Storage of PG_GEOMETRY, but if I upload a new dataset via ArcCatalog, the table is still not visible in QGIS.

It does work if I change the config.keyword to PG_GEOMETRY in Catalog (see image), but shouldnt it have to work on DEFAULTS too since I changed the settings?

I use PostgreSQL version 9.4, PostGIS 2.1.8 and ArcGIS 10.2.


texts: Where to break the loop to make several polyline instead of one?

texts:
I am trying to solve a problem, I am pasting part of the code:

from arcpy import env
arcpy.env.overwriteOutput = True

spatialRef = arcpy.SpatialReference("WGS 1984")
polylineFC = arcpy.CreateFeatureclass_management( r"H:\GIS5044\GIS5044_Nov 06\GIS5044\homework\hw#3\ ", "RhinoParsing.shp", "POLYLINE", "", "", "", spatialRef)
fieldList=["Name"]
arcpy.AddField_management(polylineFC,fieldList[0],"Text","#","#",8)
point = arcpy.Point()
array = arcpy.Array()
featureList = []
for rhino in nameCord_Dict:
    CordList= [nameCord_Dict[rhino]]
    print CordList
    print len(CordList)
    for feature in CordList:
            for cordpair in feature:
                point.X = cordpair[0]
                point.Y = cordpair[1]
                array.add(point)
    pointPolyline = arcpy.Polyline(array, spatialRef)
    with arcpy.da.InsertCursor(polylineFC,("SHAPE@")) as iCursor:
        array.removeAll()
    iCursor.insertRow((pointPolyline,))

Here I am trying to make one polyline for each [rino] of the dictionary "nameCord_Dict". but when I run the code it creates one line feature instead of 5 as there are five keys in the dictionary. I am not understanding where should i do to end polyline once all the coordinate related to one element of "nameCord_Dict" is finish and start a new line for a new element.

This is the polyline i got:


texts: Querying points inside and outside of polygons?

texts:
In PostGIS I have two tables:

Table that contains polygons:

table a.provinces

gid | geom       |   name    |
------------------------------
1   | 010604 ... | Champagne |

Table that contains points:

table a.stations

gid | geom     | province_short |
---------------------------------
1   | 0132 ... | ch             |

They are not connected with foreign key!

I would like to get:

Only points that are inside of the polygon
Only points that are outside of the polygon

SRID: 3857

Something like this:

Can it be done just like this or I have to modify database?


texts:
How do i create a linestring from points sort by date in PostgreSQL?

I tried:

SELECT observations.id,
   ST_MakeLine(observations.geom ORDER BY observations.date) AS newgeom
   FROM   observations
   GROUP BY
   observations.id;

This is going wrong. It returns me a linestring with same coordinates.
So, instead of a line from a to b, it gives me a to a.

What am I doing wrong?


texts: ArcPy: Drop fields from list of fields for arcpy.da.SearchCursor

texts:
I have a feature class with several fields in the attribute table, and I want to export a few but not all of those fields out to a tab delimited text file.  The code below sort of works, half of the fields in the uneeded_fields variable are dropped, but three (SHAPE, long, and utm_y) still end up in my text file.  I'm not sure why some of the fields are deleted while the other are passed on to the search cursor.  I have tried copying the names from the field set up in Catalog and pasting them into my script so it isn't a typo issue.  I'm not sure what I'm missing here.

As a side note, will the SHAPE field always be the first field returned unless specified otherwise?  If so del fields[0] can take care of that, but the other will still be a problem.

fields = [f.name for f in arcpy.ListFields(fc)]

uneeded_fields = ('SHAPE', 'OBJECTID', 'lat', 'long','utm_x', 'utm_y')

for i,f in enumerate(fields):
    if f in uneeded_fields:
        del fields[i]

with open(outfile, 'w') as f:
    f.write('\t'.join(fields)+'\n') #column headers
    with arcpy.da.SearchCursor(fc, fields) as cursor:
        for row in cursor:
            f.write('\t'.join([str(r) for r in row])+'\n')


texts: pyqgis: could not load WFS with bbox parameter

texts:
I try to load a WFS including a bbox parameter via pyqgis. The WFS is loaded without problems, but unfortunately the bbox parameter does not work. If the bbox parameter is added to the query i get all the results, including the results that are outside of the bbox range. When I start the request in the browser or using the QGIS GUI, i get the correctly filtered features. It looks like pyqgis does not take into account the bbox parameter. My code looks like this:

vl = QgsVectorLayer("http://geoserver_url/ows?style=&request=GetFeature&service=WFS&version=1.0.0&typeName=dataset:water&bbox=3562017.31020000018179417,5679733.8194363797083497,3563253.31000480009242892,5681088.4972945200279355&srsname='EPSG:31467'", 'water', 'WFS')
QgsMapLayerRegistry.instance().addMapLayer(vl)

Is this a known bug or something missing in my request? As already said, in the browser or via the QGIS GUI the request works, but not with pyqgis.


texts: GeoServer WFS layer loading but not being displayed (OpenLayers 3)

texts:
I am trying to load a WFS layer comprised of points. Fiddler shows that data is being returned but it is not being displayed.

Does anyone know what I am missing? Is it something to do with it being in OSGB projection?

I have created a fiddle to show my code (with the base layer removed):
https://jsfiddle.net/7mm1bges/

It seems that the method of getting WFS data through OpenLayers has evolved over time, so I have opted for the methodology shown here:

http://openlayers.org/en/latest/examples/vector-wfs.html

Here is a sample of the JSON data that Fiddler shows is being returned:

{"type":"FeatureCollection",
"totalFeatures":1803,
"features":[{"type":"Feature","id":"EAFloodAlertsWarningCentroids.1",
"geometry":{"type":"Point","coordinates":[339095.07019032,345438.84249476]}

...


texts:
I have like 6 .xyz files that are to be included in a single map (they sort of overlap). However, I'd like to save time and automatically merge them together  through python code into one .XYZ before processing them through my script.

Is there a seamless way to do this in ArcGIS Desktop?

BOTTOM LINE: MERGE XYZ's into one XYZ file. Any suggestions or alternatives to this?


texts:
Does anyone know how to calculate the "rank" of a field but within 1-100?
For example i have a field called "ValSumBasin" with all these values. I am wanting to rank these values from 1-100.

I can rank it from 1-3084 but am wanting 1-100.
This is what I have so far:

import arcpy

from itertools import count, groupby

fc = "my shapefile"
sort_fields = ["ValSumBasin"]
rank_field = "Rank"
sql_orderby = "ORDER BY {}".format(", ".join(sort_fields))
with arcpy.da.UpdateCursor(fc,
                           sort_fields + [rank_field],
                           sql_clause=(None, sql_orderby)) as cur:
    cntr = count(1)
    for key, group in groupby(cur, lambda x: x[:-1]):
        dense_rank = next(cntr)
        for row in group:
            cur.updateRow(row[:-1] + [dense_rank])

del cur


texts: Raster to polygon conversion's raster field takes only "Value"

texts:
Recently I am working with raster data and arcobjects.net. Here I find something new for me. And here some question arises in my mind. I tried to search some document that can satisfy my mind, but I couldn't find. Any pointer or document or link would be helpful.

When I perform Raster to polygon conversion's, fourth parameter is raster_field. Sometimes I can give either "Value" field or "Count" field in it, sometimes, its only "Value" field. when and why? When does it take only "Value" field?


texts: Pixel values change when opening in other version

texts: I have a raster layer with pixel values ranging from 0-0.346, which is a result after calculations in QGIS raster calculator (QGIS version 2.18.26). I want to open the same layer in QGIS 3 on another computer to continue my work there. However, when I open in in QGIS 3, suddenly the values ranges from 0-0.404. This is after having compressed the raster together with others to which the same thing happened), uploaded them to my drive, downloaded them and unzipped them on the computer with QGIS 3. Does anyone know what might be the problem?


texts: set default value of the first entry made in attribute table in qgis

texts: I will make an example so it's easier to understand what I am trying to do. So let's say I have an Attribute table with 3 columns(district-name, building, district-number). What I want to do is the following: When the first entry is made, I want QGIS to set the district-number of the first entry to be the default number for every other entry that is going to be made with the same district-name. With entry I mean when you create a feature. I can't really provide any code, because I don't even know if I have to do this with python or if I can directly do this in QGIS.


texts: Change the reclassify method in reclassify tool in ArcGIS Pro

texts:
I'm trying to get use to ArcGIS Pro and I have tried the reclassify tool.
The thing is that I know that in ArcMap there is an option to see histogram of the raster values and also choose the classify method (e.g equal count etc) but I can't find this option in the Pro.

Is there any method to see the histogram and change the reclassification method?


texts: Styling GeoJSON by Attributes with Cesium

texts:
Situation:

I have 30 3D heat maps (for 30 years) stored in GeoServer, my goal is to create an animation showing the change year by year. I found Cesium and was really impressed. Right now I am trying to style the GeoJSON layers according to the heat map symbol ID (a field in my layer) and I am having trouble getting it to work. Cesium has a sandbox showing an example styling GeoJSON by random colors but that's not exactly what I'm trying to do and I can't figure out what I am doing wrong. I think it has something to do with the "if" that is within "if (!color)", but again I don't know.

JS Code:

        var promise = Cesium.GeoJsonDataSource.load('My_GeoServer_GeoJSON');
        promise.then(function (dataSource) {
            viewer.dataSources.add(dataSource);
            var entities = dataSource.entities.values;
            var colorhash = {};

            for (var i = 0; i < entities.length; i++) {
                var entity = entities[i];
                var Symbol = entity.SymbolID;
                var color = colorhash[Symbol];

                if (!color) {
                    if (entity.properties.SymbolID === 0) {
                        color = Cesium.Color.FORESTGREEN.withAlpha(0.4);
                        colorhash[name] = color;
                    } else if (entity.properties.SymbolID === 1) {
                        color = Cesium.Color.GOLD.withAlpha(0.4);
                        colorhash[name] = color;
                    } else if (entity.properties.SymbolID === 2) {
                        color = Cesium.Color.ORANGERED.withAlpha(0.4);
                        colorhash[name] = color;
                    } else {
                        color = Cesium.Color.MAROON.withAlpha(0.4);
                        colorhash[name] = color;
                    }
                }
                entity.polygon.material = color;
                entity.polygon.extrudedHeight = entity.properties.Count * 20;
            }
        }).otherwise(function (error) {
            window.alert(error);
        });

Question:

Can anyone help me figure out how to style a 3D GeoJSON by attributes in Cesium? I'm 100% open to suggestions and am eager to learn!

If you need more information, let me know and I will try my best to explain.

Edit: I changed around my code a bit and I think I'm getting close but I'm not there yet!


texts: Virtual layer from query with special characters in field name using PyQGIS

texts:
I have a layer loaded from a WFS server that contains fields with æøå in some of the names. I am making SQL queries to make virtual layers from this layer, but whenever I query a field with æøå in its name, it only returns an empty layer. I am mainly trying to do this with pyqgis, but have tried from the QGIS GUI as well with the same results. Note that there is no problem if the value contains æøå, just the field name. Does anyone know how to fix this, or a workaround?

example of my SQL query:

 "SELECT * FROM TettstedInngangBygg WHERE dørtype = 'skyvedør'"

I got the layer TettstedInngangBygg from the WFS: https://wfs.geonorge.no/skwms1/wfs.tilgjengelighettettsted?request=GetCapabilities&service=WFS

EDIT:
By writing the following code in python console:

field_names = [field.name() for field in vlayer.pendingFields() ]
print field_names

I get the name u'd\xf8rtype'. I have tryed to write the name like that, but withe the same result as before.

UPDATE:
It seems that if i add a where statement before a statement with special characters, it works fine. I don't know why this is. By adding something that is always true, i make it work, but it does not look very pretty, and i would like to know why it only works like this.

eksample of working query:

"SELECT * FROM TettstedInngangBygg WHERE lokalId > 0 AND dørtype = 'skyvedør'"


texts:
I am using GeoServer 2.5.

I am having some trouble setting my data directory on Windows server 2012. currently the default is "
C:\ProgramData\Boundless\OpenGeo\geoserver " I have gone through the steps here which tell you to set an environment variable which I have followed and used both relative and absolute paths. I have restarted both geoserver and the server after each attempt with no luck.

what I want to do is change the data dir from what it currently is to D:\geoserver where i have copied the schema from the c drive.

Has anyone else experienced a similar problem and if so can you send me in the right direction?

Here is an image of what the env. variable looks like


texts: Adding field to table using QGIS graphic modeler?

texts:
I was wondering if anyone could help as I'm struggling to get anything out of the QGIS graphic modeler.  I need to be able to add a new field to an attribute table and update that field with a value that the user will enter.

Is there any way to do that with the graphic modeler?

I have no programming skills whatsoever.


texts: Adding combo box to QGIS attributes form

texts:
I want to customize attributes form in QGIS 3.6 with a combo box which will fill with a predefined value.

I found on the web examples but with the version 2.x or with Python code. This last solution doesn't seem to work very fine.

Is there a simple way to associate a field with a combo box?


texts: Calling canvas.saveAsImage from script in QGis -Qstring problem

texts:
I want to create a processing script that saves the QGIS MapCanvas exactly as it is seen on the screen and then Zip up the tifw and a created PRJ. The existing QgsMapCanvas.saveAsImage()-method does exactly what I want. Saving to PNG works, but I want TIF. The third argument takes a QString which I'm supposed to provide like QString("TIF").

However I can't seem to work out how to create a QString. It seems the present wrappers around Qt has removed this class. Just sending in a string or a unicode string ( eg u"TIF") does not work.

Nor does the solution in this thread work:
QGIS Map tools and Map canvas--QString errors

I'm clueless as to how to get around this problem without rewriting the logic in saveAsImage().

UPDATE: I got it to work by passing None as per gcarillo's solution, however I have since moved on by using the code from Martin Dobias in this thread: http://osgeo-org.1560.x6.nabble.com/How-to-use-QPixmap-argument-in-saveAsImage-td4112041.html. This is a better solution (although more complicated) as you can set the output resolution.


texts: Using JavaScript created select options as Leaflet map setContent

texts:
Working on a creating a custom CRS Leaflet map.  On map rightclick, show a popup with select boxes for users interaction.

var submit_Loc_popup = L.popup();
map.on('contextmenu', onPopupBoxSubmit);
var featureGroup = L.featureGroup().addTo(map);

I have the select box's hidden, until something is chosen, shown here Fiddle#1: https://jsfiddle.net/ay9pso8L/

I'd like to have those select options inserted into my maps popupContent. Unsuccessfully tried different ways but nothing produces the JSfiddles clean selects/hidden content...
I manually created the Div + select(s) with a gigantic HTML string, but I'm unable to hide content as they are embedded in quotes since it's a string... JSFiddle #2 Here: https://jsfiddle.net/0x16zqkc/

Excerpt from JSfiddle#2;

function onPopupBoxSubmit(eb) {

       var xy = map.project([eb.latlng.lat,eb.latlng.lng],6);
       console.log(xy);
       var coords3 = convertToInGameCords(xy);
       console.log(coords3);
       var popupContent = <seeJSFiddle#2 for HTML>
submit_Loc_popup
                .setLatLng(eb.latlng)
                .setContent(popupContent + '"x": ' +coords3.y+', "y": '+coords3.x)
                .openOn(map);
          };

Reading the documentation: https://leafletjs.com/reference-1.4.0.html#popup-setcontent
Maybe I've tried a function to create the DIV and pass that in but its still not modifiable.

I've been researching for a week and this Answer is somewhat similar to what I want... Kind of, but its for a marker and when I tried to do that with the popUpcontent it didn't work.
Dynamically create leaflet popup via JavaScript Object

Long story short;

How can I merge JS fiddle#1 hidden_select_(#node_resource)
into JS Fiddle#2's popup.setContent()?

EDIT#1 Based on Toms reply, JS Fiddle updated
https://jsfiddle.net/q5cdhw3e/1/
Drop downs are working when including the handlers INSIDE a mymap.on('popupopen',
However, my onSubmit is not working, I will dig around to see what I can find.

EDIT#2 JS Fiddle updated
https://jsfiddle.net/Ls1zax5c/

Placing the eventlistener first solved the on_submit problem

  document.querySelector("#node_submission_form2").addEventListener("submit", function(e) {


texts:
I have old topographic maps that I have scanned and converted into PDF in Acrobat XI 11. I need this PDF to be text-searchable (not mere image). To do this I tried text recognition in Acrobat but it fails, maybe this map is so cluttered so it fails.

I think both manual and automatic OCR may be the solution i.e. first OCR and correct or edit this OCR. Sample of map is attached

I need to generate text-searchable pdf, and asking for help on the correct workflow.


texts: Can OpenLayers 3 run without browser or is it designed to be used only in browser?

texts: I need to run OpenLayers app which will render some map into canvas and then I need url of that canvas ( Canvas.toDataUrl() ), I would use NodeJS and Node Canvas to get it done. I wonder can OpenLayers run without browser?


texts: REtrieve field info from web map services in ArcMap

texts: I am trying to work out what fields are associated with a map service in ArcMap.  The service is very slow in drawing and I wondered if it was possible to retrieve attributes associated with a layer another way?


texts:
I have created an app using the CMV and I've incorporated the Search and Attribute widgets as well.

I would like to pass parameters via the URL to achieve similar results as in my old ArcIMS application here:  http://ims.wvgs.wvnet.edu/wvog/viewer.htm?ActiveLayer=12&Query=API%3D4704971001&QueryZoom=Yes&Layers=0000111100001000

My new CMV map (in testing) is here: http://www.wvgs.wvnet.edu/gis/og/wvog/index.html

So for instance, I am looking at a well with an API number (identifier) of 4704971001.  It works when I do a search for that number in the Search box for my CMV map.   I'm not sure how to pass that number into the URL to achieve similar results as in my IMS map.  I was assuming appending something like ?Search=API%3D4704971001 would work but it doesn't.


texts: Optimising map cache creation performance in ArcGIS for Server 10.2.1

texts:
I'm relatively new to ArcGIS for Server, so I hope that someone is able to point me in the right direction in case what I've been doing isn't good practice.

I have 2 boxes with ArcGIS for server 10.2.1, both on the same site.
Both boxes have 4 processors and 16gb RAM.  Both boxes run on Windows Server 2008.

The site is used to both provide a couple of base map services to a small amount of users (<5) and to generate cache tiles for future services.

I am currently generating cache tiles for a mapping service (~50GB).
I would have expected to see the CPU usage on the 2 boxes running quite high.  But it tends to be sitting between 15% and 30% on each box.

The maximum instances for the caching tools is set to 6.

The maximum number of instances per machine is set to 3.

Am I wrong in assuming that I should be seeing a higher usage of the CPU?

Have I not put in the correct figures?

Or is my set-up not best practice? i.e. should I be using one site just for serving maps and another site just for caching?

I think I've followed the guidelines mentioned here and here.  But I'm pretty sure that the caching is running slower than it should.  After 19 hours it has only cached 1.17% of all my tiles.

Any best practice suggestions much welcome.

UPDATE:  After 21 hours the CPU usage on both machines is down to nothing:

machine 1:

machine 2:

The cache status "in progress" bar on the server is still moving, but the cache % hasn't increased in the last 2 hours.


texts: Alias's in the Attribute Table's Column Preview/Column Filters

texts:
We've spent some time setting up alias's in Layers/Properties.../Fields (Fig. 1) for database fields displayed as the column headers in the Attribute Table (Fig. 2).

Switching the layout in the Attribute Table from table view to form view and selecting between the Column Preview, the column names are now no longer aliased (Fig. 3). (This is also the case for the Column Filters, but is not shown here)

Is there any way to get the filters to use the aliased names? If there is a PyQGIS solution I'd be happy to hear about it.


texts: Rasterize vector data with the extent and resolution of reference raster

texts:
I have long tried to rasterize a point shapefile with the same resolution (32.461 km) as another raster layer I have loaded into the QGIS canvas too.

I have seen multiple cases, here and here where they suggest the GDAL rasterize tool in QGIS. However, in the tool they describe, I can set an input raster. This is nowhere to be found in the QGIS 3.4 Madeira version.

So how can I then make sure to rasterize my shapefile with the same resolution as the underlying raster layer?
EDIT:
I have tried to set the pixels and the extent of the new raster layer equal to the loaded raster file. So I get this:
Pixel Size
32.46300147045617734,-32.46300031136775033
Extent
-5645.5752507352281100,-4626.0776915619335341 : 5684.0122624539781100,4366.1733946869335341

But there is no raster layer appearing on the canvas, even after reprojecting it.
EDIT EDIT:
Okay okay(!) I get something now, but when I put the pixels the size of the raster layer, I end up with a lot larger pixels than first hoped, see below:

Why?


texts: Looking for common PyQGIS functions for QGIS 3

texts:
I would like to rewrite some of my QGIS plugins to use more common functions used by all of my plugins.

There exist some older functions (GIS2) like

https://github.com/NathanW2/parfait

https://github.com/qgis/pyqgis_wrappers

and something newer like

https://github.com/boundlessgeo/lib-qgis-commons
https://pypi.org/project/qgiscommons/
or
https://github.com/inasafe/inasafe/tree/master/safe/utilities

Are there any other known collections of PyQGIS3 common functions that can be used to do basic tasks like load vectorlayers, show comboboxes, create toolbar-actions or something like that without having to write these code-snippets everytime on your own?

Example:

If I want to lad a vectorlayer I can make use of the function load_layer from https://github.com/qgis/pyqgis_wrappers/blob/master/parfait/layer_wrappers.py

 def load_layer(filename, name = None):
    '''
    Tries to load a layer from the given file
    :param filename: the path to the file to load.
    :param name: the name to use for adding the layer to the current project.
    If not passed or None, it will use the filename basename
    '''
    name = name or os.path.splitext(os.path.basename(filename))[0]
    qgslayer = QgsVectorLayer(filename, name, 'ogr')
    if not qgslayer.isValid():
        qgslayer = QgsRasterLayer(filename, name)
        if not qgslayer.isValid():
            raise RuntimeError('Could not load layer: ' + unicode(filename))

return qgslayer


texts:
I am completely new to GitHub, I've never had a need to use it, but now I would like to use rio-color plugin from rasterio, mapbox. I am on Ubuntu.

I've already managed to start using rasterio, load my tiffs and so on, but now I am stacked on installing rio-color.

How should I install it?

Should I download all directories and run setup.py or is there something I am missing?


texts: Cannot connect to SQL Server Express 2014 from ArcMap 10.1?

texts:
I have ArcMap 10.1 and SQL Server Express 2014. I am trying to add a database server from ArcMap by following Arcmap>Catalog>Add database server>.

But the connection is not being established. as far as i know the server\sqlinstance name are correct. Is there something else I am missing?


texts: Set 500 subtypes in ArcGIS 10.0 Desktop from a csv list

texts: Does anyone know of a way I could import a long list of Subtypes from a csv list and if there is a maximum number of subtypes that one can import? My list consists of biotops and contains 500 entries. For data entry purposes in ArcGIS, it would be advantageous to have these a subtypes. But in reality maybe 500 is too many for ArcGIS. Has anyone had a similar case?


texts: Seeking software for least cost path analysis?

texts:
Background: running a least cost path analysis across mountains in middle east for archaeology project. Ideally, pick several points across an arbitrary interval distance (e.g. 1 km apart) and same for destination on other side of mountains.

Not sure which software to use. ArcGIS doesn't seem like it would do well with this kind of thing.  I ruled out ArcGIS mainly due to favoring the r.walk module in GRASS and the ability to employ K-moves. The accumulated cost surface isn't all that relevant in this project because the scope is the Middle/Late Bronze Age and things like land use and vegetative cover (NDVI) just don't apply.

Considering the scope involves several origin/destination points in the foothills of a mountain range on either side, the cost surface involves slope and aspect purely (would've involved hydrology features but the data isn't there).

Moving forward, the project is broadly established as an LCA involving a set of paths that minimizes time, and another set that minimizes energy. The time/cost function seems ideally suited in GRASS r.walk; however I'm lost as to the specific workflow. The energy/cost seems like it could involve one of two processes. 1) manage a workflow with r.walk; or (since I'm much more familiar with ArcGIS) reclassify a slope raster with energy consumption data (i.e. Minetti 2002) and run the cost path tool while maintaining the aspect raster as backlink.

Thoughts?


texts: pop() method returns error with feature class list in ArcGIS python toolbox

texts:
I am trying to do something very simple, but I'm pulling my hair out with this one.

In an ArcGIS python toolbox, for one of my tools I ask the user for a list of feature layers, using the following code in the .pyt file:

    def getParameterInfo(self):
    """Define parameter definitions"""
    param0 = arcpy.Parameter(
        displayName="Input Line Networks",
        name="InputFCList",
        datatype="GPFeatureLayer",
        parameterType="Required",
        direction="Input",
        multiValue=True)
    param0.filter.list = ["Polyline"]

The key parameter here is multiValue=True.  Then, in my function (where InputFCList is sent), I use this code to pop off the first item in the list of feature classes:

def main(fcList, fcBndPoly, boolSeg, fcOutCombine):
    fcToLine = fcList.pop
    # rest of the code...

This code runs fine during testing, when I hard-coded the list of feature classes.  However, when running the tool in ArcMap, I get the following error:

AttributeError: 'unicode' object has not attribute 'pop'

According to the ESRI documentation, setting the multiValue parameter to True should allow for the creation of a list of input feature classes, but there's obviously something else going on if it thinks the list is actually a unicode object.  Any ideas?


texts: Pan sharpening Landsat 8 composites in Google Earth Engine (Python API)

texts:
I am currently trying to use pan sharpening on Landsat 8 composites with Google Earth Engine. What I did is simply create a composite from an image collection as described in the documentation:
collection = ee.ImageCollection('LANDSAT/LC08/C01/T1').filterDate('2015-01-01', '2015-12-31')
custom_composite = ee.Algorithms.Landsat.simpleComposite(collection=collection,
                                                             percentile=75,
                                                             cloudScoreRange=5)

map_sharpened = folium.Map(location=[37.76664, -122.44829], zoom_start=13)

map_sharpened.add_ee_layer(custom_composite, {
    'bands': ['B4', 'B3', 'B2'], 'max': 128, 'gamma': [1.1, 1.1, 1]}, 'rgb')
display(map_sharpened)

This works as expected. Next, I want to use pan sharpening on my composite. I am basically following the example in the Earth Engine documentation (https://developers.google.com/earth-engine/guides/image_transforms) with the exception that I use a composite instead of one specific image. However, when I try to display the new layer, I only get a white map.
hsv = custom_composite.select(['B4', 'B3', 'B2']).rgbToHsv()
sharpened = ee.Image.cat(
    [hsv.select('hue'),
     hsv.select('saturation'),
     custom_composite.select('B8')]).hsvToRgb()

map_sharpened.add_ee_layer(sharpened, {
    'min': 0, 'max': 0.25, 'gamma': [1.3, 1.3, 1.3]}, 'pan-sharpened')
map_sharpened.add_child(folium.LayerControl())
map_sharpened


texts:
I am new to WMS.

I requested some data (I assumed it'd be in vector format).  I received a WMS link, which I opened and I am not totally sure how.

Could someone shed some light on how I opened this?

I am referring to arcgis on vtransmap01.aot.state.vt.us (in the picture below).

Why is the link not shown under Add WMS Server?


texts: I'm not able to add SQL file to postgresql database

texts:
I started with basics and did these steps:

createdb -U postgres -E UTF8 -O Peete gis -OK
createlang -U postgres plpgsql gis (was already installed) -OK
psql -U postgres -d gis -f PATH_TO_POSTGRES/share/contrib/postgis-2.1/postgis.sql -OK
psql -U postgres -d gis -f PATH_TO_FILE/4326.sql (I downloaded file from Spatialreference and chose "PostGIS spatial_ref_sys INSERT statement") - OK
Then I ran this command shp2pgsql -c -s 4326 -I roads.shp > roads.sql-OK
Now I wanted to add this SQL file to my database created above so I ran this command: psql -d gis -f roads.sql-FAILED

I got this error message:

ERROR:  AddGeometryColumn() - invalid SRID

CONTEXT:  SQL statement "SELECT AddGeometryColumn('',$1,$2,$3,$4,$5,$6,$7)"

PL/pgSQL function addgeometrycolumn(character varying,character varying,character varying,integer,character varying,integer,boolean) line 5 at SQL statement
********** Error **********

The shapefile I'm using is in WGS84 coordinate system EPSG:4326

Any recommendations?


texts: Python and ArcObjects- Stand Alone open MXD

texts:
I wrote a Python function using comtypes to access ArcObjects to programmatically change the size of layer symbols in an MXD.  The function works great in the Python window of ArcMap.  However, I need to do this stand alone to do this for many MXD's.

Using arcpy I can open an MXD in memory like this:

mxd = arcpy.mapping.MapDocument(r'C:\some_folder\map.mxd)

I need to figure out how to do the above using ArcObjects (except return an instance of Application/IMxApplication).  In my ArcObjects function, I am starting with a reference to the current open document.  However, I cannot figure out how to do this without having a map document open.

I thought I could do something like this (using my own version of snippets):

import arcobjects
import comtypes.gen.esriFramework as esriFramework
import comtypes.gen.esriArcMap as esriArcMap
arcobjects.InitStandalone()
app = arcobjects.NewObj(esriArcMap.Application, esriFramework.IApplication) #this isn't doing the trick

And then calling app.OpenDocument() to point to the file path to actually open the MXD. This is not working and I do not know enough about initializing things in ArcObjects to know what I'm doing wrong.  I know I cannot access anything from the AppROT to get a reference to the current running application because ArcMap/Catalog are not open.

In order to be able to call my function to set the symbol sizes, I need to pass in a reference to an IApplication instance (which I normally get this interface by using GetApp() or GetCurrentApp()).  I cannot seem to find code samples for this anywhere.  Does anyone know how to do this stand alone in Python?

Here is what I'm using to set the symbol sizes:

def setSymbolSize(pApp, layer_names=[], pointSize=12, lineWidth=1, autoSave=True):
    """sets the size/width of simple esri Marker or Line symbols

    Required:
        pApp -- open mxd instance, must be of type esriArcMapUI.IMxApplication

    Optional:
        layer_names -- list of layer names to set size of as they appear in TOC.
        pointSize -- size for all point layers specified.  Default is 12
        lineWidth -- width for all line layers specified.  Default is 1
        autoSave -- option to save document automatically after making changes.
            Default is True.
    """
    import comtypes.gen.esriDisplay as esriDisplay
    import comtypes.gen.esriArcMapUI as esriArcMapUI
    import comtypes.gen.esriCarto as esriCarto

    if isinstance(layer_names, (basestring)):
        layer_names = [layer_names]

    # defaults for point (shapeType 1) is 12, default for line (shapeType 3) is 2
    defaultSymbol = {1: esriDisplay.IMarkerSymbol,
                     3: esriDisplay.ILineSymbol}

    # iterate through layer list
    for lyr in iterLayers(pApp):

        if lyr.Name in layer_names:

            # cast to esriCarto.IGeoFeatureLayer to get renderer
            flyr = CType(lyr, esriCarto.IFeatureLayer2)
            iGeo = CType(lyr, esriCarto.IGeoFeatureLayer)
            renderer = iGeo.Renderer

            # cast to simple renderer and then to simple symbol
            simple = CType(renderer, esriCarto.ISimpleRenderer).Symbol

            # cast to marker symbol interface and set size or width
            marker = CType(simple, defaultSymbol[flyr.ShapeType])
            if flyr.ShapeType == 1 and pointSize is not None:
                marker.Size = pointSize
            elif flyr.ShapeType == 3 and lineWidth is not None:
                marker.Width = lineWidth

    # refresh Map
    pDoc = pApp.Document
    pMxDoc = CType(pDoc, esriArcMapUI.IMxDocument)
    pMxDoc.UpdateContents()
    pMxDoc.ActiveView.Refresh()

    if autoSave in (True, 1):
        pApp.SaveDocument()
    return


texts: NDVI index based on QA bands for MCD43A4 , in GEE

texts:
I want to calculate the NDVI index based on the MODIS/006/MCD43A4 collection and based on quality bands, but as a new GEE user I have no idea If I am on the rght way to do this.
So far I have the following :

  var imageCollection = ee.ImageCollection("MODIS/006/MCD43A4")
        //.select(['Nadir_Reflectance_Band2', 'Nadir_Reflectance_Band1'], ['NIR', 'RED']).multiply(0.0001)
        .filterBounds(geometry)
        .map(function(image){return image.clip(geometry)})
        .filter(ee.Filter.calendarRange(10,11,'month')),
    geometry = /* color: #d63000 */ee.Geometry.Polygon(
          [[[74.31655930070895, 56.84897198026647],
            [71.72378586320895, 53.54030739149663],
            [81.04019211320895, 51.563412328672136],
            [87.80777023820895, 55.429013452403986]]]);;

var modNDVI = imageCollection.select("Nadir_Reflectance_Band2","Nadir_Reflectance_Band1","BRDF_Albedo_Band_Mandatory_Quality_Band1","BRDF_Albedo_Band_Mandatory_Quality_Band2");

       var filterqa = function(feature){
     var image = ee.Image(feature);
     var QAImage = image.select(['BRDF_Albedo_Band_Mandatory_Quality_Band1']);
    var QAImage2 = image.select(['BRDF_Albedo_Band_Mandatory_Quality_Band2']);

     var maskedImage = ee.Image(0)
                         .where((QAImage.eq(0),1)&&(QAImage2.eq(0),1));

     return image.mask(maskedImage);
};
var md_masked= modNDVI.map(filterqa);

        var ndvi = (function(image) {
          return image.expression('(NIR - RED) / (NIR + RED)', {
            'NIR': image.select('Nadir_Reflectance_Band2').multiply(0.0001).float(),
            'RED': image.select('Nadir_Reflectance_Band1').multiply(0.0001).float(),
          }).clamp(-1, 1);
        });

    print(imageCollection)

        //var ndviParams = {min: -1, max: 1, palette: ['blue', 'white', 'green']};
        //Map.addLayer(, ndviParams, 'NDVI');

From the imageCollection, I see that each image has 14 bands from which, 7 (7-13) have quality information. How can I calculate the NDVI index based on the good quality observations?
I am really confused as there is no a SummaryQA band simillar to the one from the MOD13Q1.


texts: Exporting query to identify the exact points that fall into a polygon using ArcGIS

texts:
I am using ArcGIS 10.6.1. I have a shapefile consisting of points (facilities names and adresses). I also a shapefile consisting of polygons made up of enumeration areas as well as the Local Government Areas and States they fall into.

I want to find out how many points (that is the facilities) falls into each polygon as well as the exact points. I have tried this using the spatial join query. The points are 16,000 while the polygon are 687. The spatial join query is only giving me the number of points that fall into each polygon. But i also need to know the details of the point that fall into each polygon. Please how do I go about this?


texts: Adding WMS or WFS services to Google Earth Engine

texts:
I'm new to using Google Earth Engine. I would like to know if data from WMS or WFS services can be added to this platform. I need this information to do a data validation

The information is on the servers:

http://idecor-ws.mapascordoba.gob.ar/geoserver/idecor/wms

http://idecor-ws.mapascordoba.gob.ar/geoserver/idecor/wfs


texts: How to color lat/lon points in given RGB data?

texts:
I am using QGIS on a Mac.

I am working with hydroacoustic data for mapping and interpreting the seafloor.
Data are collected, processed and interpreted. The result is a txt file with Lat and Lon in the first columns followed by R, G and B data, e.g. [54.20,7.78,0.94,0.62,0.37]. The files can have 100,000 lines or more.
I am not really good in GIS but I need to use it for some reasons. I seem to be unable to plot these data in GIS. There is obviously no symbology/style way to tell the system to take the color values from a file (i.e. the next 3 values in the same line of the file). It is so easy in MatLab: scatter(lat, lon, symbolsize, [R G B]) but I can't work easily with multiple layers in Matlab.


texts: Cause of ERROR 000539: Error running expression?

texts:
The following errors in ArcMap are preventing me from moving forward:

1) In Raster Calculator using "Aspect":  "%DEM_Aspect%" <= 90 or "%DEM_Aspect%" >= 270

  ERROR 000539: Error running expression: rcexec()
  Traceback (most recent call last):
    File "", line 1, in
    File "", line 5, in rcexec
  ValueError: The truth value of a raster is ambiguous. Invalid use of raster with Boolean operator or function. Check the use of parentheses where applicable.

I'm using SRTM1N34W087V2, 11-FEB-00 (Void Filled) raster.

I'm attempting to use ArcMap to accomplish the same results done in the QGIS tutorial Lesson 3:  Terrain Analysis.

The lesson is found here:  http://manual.linfiniti.com/en/rasters/terrain_analysis.html


texts: Why no more error messages after importing QGIS modules into python?

texts:
I have a problem when writing python code whilst including the qgis modules. As can be seen below, when I usually call a variable that does not exist (for example), I receive the usual error message. Just what you'd expect to happen. But after having imported qgis.core and qgis.utils into python, all error messages disappear (see below). This obviously makes code debugging extremely tedious, as I never get any indication, why and where a code failed.

>>> print this_variable_does_not_exist
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
NameError: name 'this_variable_does_not_exist' is not defined
>>> from qgis.core import *
>>> from qgis.utils import *
>>> print this_variable_does_not_exist
>>>                          (<< why is there no error here?)

How can I stop this "error-suppressing" behavior?


texts: Assigning Values from Raster Layer to Vector Layer in ArcGIS Desktop?

texts:
I'm using ArcGIS Desktop 10.3.

I have two layers. The first is a raster layer representing density of birds over a given area. The second layer is a polygon layer representing property boundaries, of which there are roughly 14,000.

I'm attempting to prioritize the polygon features based on the corresponding density of shorebirds. The nature of the polygon layer is such that some polygons cover multiple cells in the raster, while some are small and in some cases cover less than a single cell in the raster.

Ideally what I would like to do is where a polygon overlaps with the raster layer, I would like to assign the highest corresponding value from the raster layer to an attribute field of the given polygon feature.  In a perfect world, where the polygon is large enough that it covers multiple cells in the raster, I would like to either take the mean value of those cells, or alternatively assign the highest value to the polygon.

I can't figure out how or if this can be done.


texts: Selecting only polygons with certain variable in them using ArcMap?

texts:
I have green space polygons for London and also illegal pollution areas in polygons in ArcMap. How would I go about selecting only the green space polygons that have pollution areas in them (so the green spaces with no pollution in are not selected). I know this could be done by individually selecting by attributes but I have around 1,400 green spaces to consider so this would take a long time.

Does anyone know what tool/s I’d need to use to select just these green space polygons?


texts:
I want to create a new MultiPolygon shapefile, I tried this code

Mp_Layer = QgsVectorLayer('MultiPolygon?crs=EPSG:26191', 'layer', 'memory')
QgsProject.instance.addMapLayer(Mp_Layer)

But when I check the wkt_geom if the created shapefile its gives a polygon type.

It’s possible to create a Multipolygon directly by PyQGIS?


texts: Adding a security layer to my tile server

texts:
I have a tile server (apache2 + mod_tile serving OpenStreetMap data) that I'm going to use via Leaflet. Everything's working fine but before putting this public, I have a big concern about the security of this set up.

If I just make my server public, anyone will be able to point their own Leaflet / OpenLayers / whatever to my tile server and use it inadvertently and without authorization.

What's the best approach to make it safer?

I thought of a few possible solutions:

I could use Apache to restrict the origin of the request

Pro: easy to implement

Con: easy to forge and bypass
I could issue some sort of key and use it to validate the requests

Pro: relatively easy to implement

Con: also easy to steal
I could build my own authorization service to be between the clients and the tile server with more sophisticated techniques to avoid unauthorized uses

Pro: relatively safe

Con: harder to implement and may also add some significant latency when validating every single tile request

What do you guys suggest? I have the impression that if I just make the tile server public it's a matter of time until somebody starts abusing it. Is there any best practice when hosting your own OpenStreetMap solution?

The solution is hosted on AWS so any solution involving AWS stack is valid


texts: Finding the point closest to a point in another feature class and the paths between them

texts: There are two point feature classes (let's call them Sinks and Sources). Working with ArcGIS Desktop, I need to find, for each Sink, its nearest Source. Normally I'd do a spatial join, but this time I need to find the path between each Sink and its nearest Source along roads, with the path itself to become a line feature, to be used in subsequent analysis. Am I right that the Network Analyst extension is absolutely necessary for this, or are there any workarounds?


texts:
I'm new to QGIS and working with GIS in general (student of architecture), so please bear with me.

I have obtained shapefiles with topographic curves of an ocean floor. My intention is to export these as CAD files and render as 3D models of the ocean floor in Revit. I have done this procedure earlier with terrain shapefiles. My problem now is that these files mapping depth does not contain z values, so when exporting them into CAD they are completely flat rather than 3D models.

The depth information in meter exists as attributes in the file (see attached image). What I want is to transform this depth value into real z coordinates for each line so that the curves have different heights when importing as CAD file in Revit.

Thankful for any pointers in the right direction.


texts:
I've seen this question was asked before but I couldn't find any solution to the problem.

Since today, my GeoServer does not reply on any requests anymore, I can't even access the GeoServer UI. I didn't change any settings since the last time I used GeoServer and it worked fine then. This is the HTTP error message:

HTTP ERROR: 503
SERVICE_UNAVAILABLE
RequestURI=/geoserver/web
Powered by jetty://

Going further, access to my PostgreSQL/PostGIS Database is denied as well, telling me

  could not connect to server: Connection refused (0x0000274D/10061)  Is
  the server running on host "127.0.0.1" and accepting TCP/IP
  connections on port 54321?

Has anybody experienced something similar? If lines from the log are needed, let me know which ones.

This is probably why it fails:

  ERROR [context.ContextLoader] - Context initialization failed
  org.springframework.beans.factory.BeanCreationException: Error
  creating bean with name 'gwcCatalogConfiguration' defined in URL
  [jar:file:/C:/Program%20Files%20(x86)/OpenGeo/OpenGeo%20Suite/webapps/geoserver/WEB-INF/lib/gwc-2.1-SNAPSHOT.jar!/geowebcache-geoserver-context.xml]:
  Cannot resolve reference to bean 'rawCatalog' while setting
  constructor argument; nested exception is
  org.springframework.beans.factory.BeanCreationException: Error
  creating bean with name 'rawCatalog' defined in URL
  [jar:file:/C:/Program%20Files%20(x86)/OpenGeo/OpenGeo%20Suite/webapps/geoserver/WEB-INF/lib/main-2.1-SNAPSHOT.jar!/applicationContext.xml]:
  Initialization of bean failed; nested exception is
  java.lang.RuntimeException: java.lang.IllegalArgumentException:
  Resource must be part of a store


texts: Where is the QGIS API action for "deselect"?

texts:
Where is the QGIS action for the “deselect” toolbutton in the 'Attributes" toolbar? I do not see this action in the QgisInterface class or do not recognize it by name.  Is it missing from iface?

I'd like to use this action in a plugin directly instead of the button on the toolbar.


texts: Web Application with data editing capabilities by users stored in Geoserver

texts:
I have been using GeoExt/Extjs/Openlayers/Geoserver for some time.I have created a web application for mapping features in a county.I wanted to know whether there is a way I can add (edting/Addnew feature) capabilities on my application with authentication which reads the users on Geoserver.

This is how I expect it to work; All viewers can see the features but when they click on addNew button , an authentication window pops up for them to login.The people allowed to edit are the users of the office whose usernames and passwords are on Geoserver.

Is this possible and if yes, how?


texts: Just curious if anyone has used the ASP MVC 3 framework for creating web mapping applications.  I have one coming up and am planning on using the ArcGIS javascript API.  Would be interested in hearing if anyone has had success or their pitfalls using this framework.  Anyone try connecting to SDE with MVC3?


texts: Check for licensing level in ArcGIS 10.2.2

texts:
I am in process of migrating some legacy code from 10.0 to 10.2.2. As part of application, I need to make sure user is using ArcEditor (now called Standard) level license. Current logic looks at the registry setting to establish the license level. I wanted to get some feedback from developer community if there is a better way of doing this without using registry setting.

        pRegSubKey = Registry.LocalMachine.OpenSubKey("SOFTWARE\ESRI\License10.0")
        pRegLicenseValue = pRegSubKey.GetValue("SOFTWARE_CLASS")

        If Not pRegLicenseValue Is Nothing Then
            Select Case pRegLicenseValue
                Case "Viewer"
                    pArcMapLicense = ArcMapLicenseLevel.ArcView
                Case "Editor"
                    pArcMapLicense = ArcMapLicenseLevel.ArcEditor
                Case "Professional"
                    pArcMapLicense = ArcMapLicenseLevel.ArcInfo
            End Select
        End If


texts: Dealing with large sets of data on web map without using Table of Contents?

texts:
I'm working on a map design that integrates many (around 300 separate map layers) datasets that need to be made available to users within the same map (just for viewing). I've looked at creating a table of contents with high level and sub categories but it just seems very clunky.

We are harvesting the datasets as web map services. The WMS URLs are stored in an oracle database where they are configured into 9 main categories to provide a top level summary table of contents. We also have an advanced table of contents (also configured in Oracle) which arranges the dataset layers within the main category list. The whole thing is becoming unmanageable and not at all user friendly but we still need to make all the datasets available.

Does anyone have any good guidelines/best practice/inspirational ways of dealing with large sets of data on a map?


texts: Shortest path (point to point) in iteration

texts:
I am using Network Analysis-> Shortest path (point to point) to compute the distance between two points, following a network from OSM. It works like magic for single requests, but I do not seem to be able to make it work in an iteration. The script is as follows:

import processing
import csv

with open('.../Locations.csv', "r") as Tele_points:
    TeleLoc=list(csv.reader(Tele_points, delimiter=','))
    for line in TeleLoc:
        x=line[2]
        y=line[3]
        location=y+','+x
        x1=line[10]
        y1=line[11]
        location1=y1+','+x1

        processing.run("native:shortestpathpointtopoint", {'INPUT':'...\\OSM_Net.shp|layername=OSM_Net','STRATEGY':0,'DIRECTION_FIELD':None,'VALUE_FORWARD':'','VALUE_BACKWARD':'','VALUE_BOTH':'','DEFAULT_DIRECTION':2,'SPEED_FIELD':None,'DEFAULT_SPEED':50,'TOLERANCE':0,'START_POINT':''+location+'','END_POINT':''+location1+'','OUTPUT':'.../Path.shp'})

This script overwrites the output shapefile and I would like to have it appended instead so that the results of each pair of locations can be on a different line in the attribute table.


texts: Allow WFS searches but no total data download

texts: I am using GeoServer and I published some data. I want that users can extract features information by predefined queries but I do not want that they can also download the entire dataset with georeferencing. Is it possible?


texts:
As I can read in many posts, it is necessary to convert Sentinel 2 L1C pixel values from Digital Count to TOA reflectance by dividing with the scaling factor.

However when I open my image in SNAP, the values seem to already be in TOA reflectance, or at least between 0-1. I have added my screenshot, showing the histogram and a random pixel value.

Is there an explanation why I see these values, and are they truly in TOA reflectance?


texts: Volume calculation using Mosaic & Cut Tools

texts:
I'm new to ArcGIS and I'm taking an introductory course for it. I have a problem trying to quantify the volume of a lake. The lake is enclosed by two dams which are to be extended from 1909 to 1932 m.a.s.l. I calculated the volume of the original elevation (1909m).

Now we are told to use the "Mosaic to raster" tool to merge our DTM with the raster of the dam with the new height and raster for the lake. We then have to fill the basin to the new elevation of the dam and then use the "Cut Fill" tool to quantify the volume.

I have tried this but I never get the right result. I have been at it for hours, but I do not seem to find where my mistakes lie. I tried building the mosaic with different input sequences and then filled the lake up, but to no avail. The lake either doesn't get filled fully or the dam doesn't have 1932m as height value and remains at 1909m.

Does anyone have any advice or has an idea for a different approach?


texts: Debugging ERROR 000358: Invalid expression from SelectLayerByAttribute?

texts:
I'm creating a script tool and using the SelectLayerBy Attribute function (below),

arcpy.SelectLayerByAttribute_management("cities_lyr", "NEW_SELECTION", '"CRIME_INDX" <= 0.02')

However, I keep getting this error - ERROR 000358: Invalid expression.

Can anyone please offer some clarity here.


texts: Using Arcade, AGOL, and Attribute Expressions to Create a Text field with Quotes

texts:
I am working on creating a Attribute Expression that would display HTML with quotes in the field. I tried putting double quotes for the inner quotes, but I keep getting "Parse Error: Unexpected String". Here is how I am writing out the expression which I know works with Python.

"PHONE NUMBER"

What it should produce:

PHONE NUMBER

Does anyone know how to print out quotes within a text field using Attribute Expression?


texts:
Can someone please tell me what's happening here. I'm still new to Python for ArcGIS and I really need to understand this piece of code for my project.

How does dct_max work and how is tpl created?

At first dct_max seems to me like an empty dictionary. Later its changed to a list?. What does this line mean - if from_g in dct_max:. I can't find out what in does. Also the line tpl = dct_max[from_g] seems to create a list called tpl with a single item. Then how can there be any meaning to tpl[1]???

def main():
    import arcpy
    arcpy.env.overwriteOutput = True

# input fc
fc = r"C:\Forum\Gemeinde\NRW_Gemeinde_SpatiaJoin.shp"
fld_from = "GEN"
fld_to = "GEN_1"
fld_val = "Value"

# output field (will be added to input fc)
fld_out = "Gemeindes"

# output fc
fc_out = r"C:\Forum\Gemeinde\NRW_Gemeinde_diss01.shp"

# detemine the max commutors per from gemeinde
dct_max = {}
with arcpy.da.SearchCursor(fc, (fld_from, fld_to, fld_val)) as curs:
    for row in curs:
        from_g = row[0]
        to_g = row[1]
        val = row[2]
        # check that from_Gemeinde is not equal to to_Gemeinde
        if from_g != to_g:
            if from_g in dct_max:
                # check if value is higher
                tpl = dct_max[from_g]
                max_to_g = tpl[0]
                max_val = tpl[1]
                if val > max_val:
                    # update entry in dct
                    dct_max[from_g] = (to_g, val)
            else:
                # insert value
                dct_max[from_g] = (to_g, val)


texts: How to calculate the average velocity over the flowpath length?

texts:
I’m working on building a GIS based groundwater travel time distribution (TTD) model similar to the one developed by Schilling and Wolter (2007)
http://lilt.ilstu.edu/ewpeter/GEO435/Readings/2012/Schilling_Kelly.pdf.

The conceptual model of cell TTD is illustrated below (Schilling and Wolter, 2007).

As you know
Velocity = -k (dh/dl)/n

I have calculated the velocity using raster calculator. I used the permeability values  from the soil map as a substitute for hydraulic conductivity, the slope from the DEM as substitute for hydraulic gradient (dh/dl) and porosity (n) one value only from literature. I calculated the distance till the nearest river using Euclidean Distance.

After getting the velocity and the distance, I can easily calculate the time using raster calculator. However, it is not straightforward like this. Schilling & Wolter (2007) described how to calculate the time as follows “The time needed for groundwater in each 5-m cell to travel from the cell to the stream network was then determined by averaging all the individual cell velocities over the ﬂowpath length from the cell to the stream”

So, in order to calculate the time, I have to calculate the average velocities over the flowpath length which I don’t know how to do. It would be much appreciated if you can help me in this regard.


texts: Associating critical zones using ArcGIS Pro

texts:
I am new to ArcGIS. I have 5 layers in ArcGIS Pro. Four of them are about the Zones (named as Zone 1 to Zone 4) and are polygon in shape. Zone 1 is considered as the most critical zone. Zone 2 is considered a critical zone so on and so forth.  The fifth layer is about the conductor which is polyline in shape. The conductor layer is passing through a different zone and has shape length as one of the attributes. Some conductors are located in more than one zones. I am trying to find a zone of those conductors who are located in different zones in such a way that if a small length of those conductors comes in critical zones and large length comes in the less critical zone then still we say that conductor is associated with the critical zone. For example, if a conductor length is 50 meters and 5 meters comes in zone 1 and 45 meters comes in zone 3, we would assign zone 1 to the conductor because zone 1 is more critical than zone 3.
I tried the intersection tool(as shown in the figure)

but it did not help.
Could anyone guide me?


texts: Exporting spatial population density data to Excel

texts:
I have spatial population density data that I would like to export to Excel (or similar) so that I can use it in analysis with another dataset. Specifically, I want the XY coordinates and the values in a table. The pop. dens data is a raster dataset, and I have tried converting it into integer form, then converting that to points, joining the points to the integer raster data, then getting the XY coordinates. I run into problems with the join, because the object ids in the raster file don't match those in the vector/point file...which was created from the raster file. I understand why not, but I'm not sure how to fix this.

Any suggestions on how best to go about this? I suspect my dataset is too large for the sample command; I tried running it for 48 hours with no output.

I am using ArcMap 10.3.1 for desktop.

Edited to add that I have this problem regardless of whether I convert the raster data to points or to polygons (and then get the XY coordinates from the centroid of each polygon)


texts:
I used the following code from this github repo I know it take a point coodinates and return area of 10 sequare Km around it by transforming Km to degree, but I need to make this area to 5 sequare Km for example. how can I achieve that? The only number I found is (3.75) but I'm not sure what it represents and how to change it to get the area I need.
I know this might sound obvious for some, but I have zero background in GIS and GEO data and I played around with the code but the results I had didn't make any sense at least for me.
import math

# Distances are measured in kilometers.
# Longitudes and latitudes are measured in degrees.
# Earth is assumed to be perfectly spherical.

earth_radius = 6271.0
degrees_to_radians = math.pi/180.0
radians_to_degrees = 180.0/math.pi

def change_in_latitude(kms):
    "Given a distance north, return the change in latitude."
    return (kms/earth_radius)*radians_to_degrees

def change_in_longitude(latitude, kms):
    "Given a latitude and a distance west, return the change in longitude."
    # Find the radius of a circle around the earth at given latitude.
    r = earth_radius*math.cos(latitude*degrees_to_radians)
    return (kms/r)*radians_to_degrees

def ten_km_square(latitude, longitude):
    slat, nlat = latitude+change_in_latitude(-3.75), latitude+change_in_latitude(3.75)
    wlon = longitude+change_in_longitude(latitude,-3.75)
    elon = longitude+change_in_longitude(latitude, 3.75)
    return(nlat, wlon, slat, elon)

def main(lon, lat):
    '''First argument degrees longitude (E is positive, W negative)
        of the landslide location,
        second argument latitude (N positive, S negative),
        in decimal format(not minutes etc.)'''
    nlat, wlon, slat, elon = ten_km_square(lat,lon)
    #print("(NLat:{:.4f},WLon:{:.4f}),(SLat:{:.4f},ELon:{:.4f});".format(nlat, wlon, slat, elon))
    print("var region = '[[{:.4f},{:.4f}], [{:.4f},{:.4f}], [{:.4f},{:.4f}], [{:.4f},{:.4f}]]';".format(wlon,nlat,elon,nlat,wlon,slat,elon,slat))
    print("var rectangle = [{:.4f},{:.4f},{:.4f},{:.4f}];".format(wlon,slat,elon,nlat))

#change these to longitude, latitude
main(18.07769,44.14354)


texts:
I am writing a VB.net console application that calls ArcGIS Desktop 10.1 ArcToolbox geoprocessing tools from my code. Some of these tools use Python syntax and I am encountering a problem that I suspect is caused by conflicting architectures. Is there a way to control which version of Python that ArcGIS uses?

System: 64 bit Windows 8.1
ArcGIS: 32 bit ArcGIS 10.1 Service Pack 1
Python: Both 32 and 64 bit versions installed (latter from ESRI's background processing download)
My App: 32 bit VB.net (not Any CPU) console application.


texts: QGIS development version(2.99) and python api questions

texts:
Can I use QGIS 2.99 PyQGIS in standalone apps like the others QGIS version?

Do I need python 3.6 to use this python api?


texts: QGIS Convert MTM coordinate reference system

texts:
I have X,Y coordinates that imported from csv file into QGIS and assigned WGS84 to the layer, when I open the attribute table and select Geometery $X it just shows the state X value and doesn't convert to lat/lng. Below is value of X,Y trying to convert.

294556.96,5039387.933

GEOGCS["GCS_WGS_1984",
 DATUM["D_WGS_1984",
 SPHEROID["WGS_1984",6378137.0,298.257223563]],
 PRIMEM["Greenwich",0.0],
 UNIT["Degree",0.0174532925199433]]
PROJCS["NAD_1983_MTM_8",
 GEOGCS["GCS_North_American_1983",
 DATUM["D_North_American_1983",
 SPHEROID["GRS_1980",6378137.0,298.257222101]],
 PRIMEM["Greenwich",0.0],
 UNIT["Degree",0.0174532925199433]],
 PROJECTION["Transverse_Mercator"],
 PARAMETER["False_Easting",304800.0],
 PARAMETER["False_Northing",0.0],
 PARAMETER["Central_Meridian",-73.5],
 PARAMETER["Scale_Factor",0.9999],
 PARAMETER["Latitude_Of_Origin",0.0],
 UNIT["Meter",1.0]]


texts: Odd colors rendered when using Bicubic interpolation in GeoServer

texts:
Just upgraded GeoServer to latest version (2.15). I used to apply bicubic interpolation method to render some WMS basemaps in the previous install (v 2.10) but now when displaying the same maps, odd colors appear because of the interpolation.

Here some examples:

The colors are Ok when I use other interpolation methods (e.g. Nearest Neighbor). However, I still want to use Bicubic because results were better when I used it in the previous install.

There have been other changes because the host server has been upgraded as well, so for example the new server runs Java 11 while the old one ran Java 8.

Are there any modifications that I could apply whether in Geoserver or in Java (or anywhere else) so bicubic interpolation provides good results?


texts: No output when rendering GeoJSON using Leaflet() package in RShiny

texts:
I am creating a Rshiny interface on that can take in the user input when a user click on a point on the map, query the user's selection from postgresql and return another map layer. The package I am using is leaflet. (https://rstudio.github.io/leaflet/). I have followed the examples, but am unable to generate any map layers.

The following are my ui.r and server.r codes

ui.r

shinyUI(fluidPage(
  titlePanel("title here"),

  sidebarLayout(position = "right",
     sidebarPanel(
       selectInput("time",
                   label = "Choose a time period",
                   choices = list("1", "2", "3")
     )
   ),

  mainPanel(leafletOutput('map'))
))

server.r

shinyServer(function(input,output) {

  drv <- dbDriver("PostgreSQL")
  con <- dbConnect(drv, dbname="postgres",host="localhost",user=***, password=***, port=5432)
  geojson_fromsql <- dbGetQuery(con, "SELECT ST_ASGeoJson(geom) FROM table WHERE time='input';)

  geojson <- fromJSON("file path/file name.geojson")

  map <- leaflet() %>%
    addTiles() %>%
    addGeoJson(geojson)

  output$map <- renderLeaflet(map)

}

geojson gives me the following:

$type
[1] "Feature Collection"

$crs
$crs$type
[1] "name"

$crs$properties
                         name
"urn:ogc:def:crs:OGC:1.3:CR84"

 $features
 $features[[1]]
 $features[[1]]$type
 [1] "Feature"

 $features[[1]]$properties

 year
 "2012"

 $features[[1]]$geometry
 $features[[1]]$geometry$type
 [1] "Point"

 $features[[1]]$geometry$coordinates
 [1] 103.77749    1.28983

 ...

 $features
 $features[[17]]
 $features[[17]]$type
 [1] "Feature"

 $features[[17]]$properties

 year
 "2012"

 $features[[17]]$geometry
 $features[[17]]$geometry$type
 [1] "Point"

 $features[[17]]$geometry$coordinates
 [1] 103.67749    1.32983

while geojson_sql has the following structure:

(data.frame)
                                        st_asgeojson
1{"type":"Point","coordinates":[103.77748, 1.28983]}
...
17{"type":"Point","coordinates":[103.67749, 1.32983]}

addGeoJson(geojson) gives me nothing in the main panel (fig 1). *addGeoJson(geojson_sql) gives a leaflet widget, with no tilelayers (fig 2). *addGeoJson(geojson_sql[1,]) gives a leaflet widget, with no tilelayers as well (fig 2). I believe I have some syntax errors in setting the tilelayer as well.

But more importantly, I would like to know what format/structure does the leaflet package accept for geojson, and how I can derive it in R?

PS: Suggestions for alternative packages or methods are welcome too. I have tried rCharts and have hit a wall here as well: https://gis.stackexchange.com/questions/151536/how-to-convert-a-geojson-from-postgresql-to-geojson-format-readable-in-rcharts-i
And I am not certain that rCharts takes in the user's input when users click on points on the map.

fig 1 - reading geojson directly into R

fig 2 - reading geojson from postgresql


texts: OSM Tools plugin is broken for QGIS 2.8.6-Wien

texts:
When installing osm tools for Wien 2.8.6 I get an error: the plugin is broken python said no module named yaml.

Yaml is on my linux mint box 18 and I can't find documentation with google or on the developer website.

Any insights or answers?


texts: Changing coordinate of point shapefile in ArcGIS using ArcPy?

texts:
I'm working on coding a function in ArcGIS trying to measure the number of pixels between a point (feature class) and its nearest building horizontally within a UDEM raster, then return the number.

Here is my idea:

Execute Zonal Statistics as Table and create a cursor to see whether the value in output table equals 0(the elevation of ground in my UDEM are set to 0 while elevation of building remain the real value)
If so(value=0), make x coordinate plus one, which means move the point 1 meter to the right(the resolution is 1×1) and refresh the coordinate of the point.
Iterate over step2 until the value doesn't equal 0,thus showing the location of point has moved to the "boundary" and I'll know the distance.

My problem now is having no idea how to perform Step2.

I tried some code to change its location  but doesn't work.


texts: Read (GEO)JSON file from web service using FME Desktop

texts:
I have a webmap of bike stations and there is a json service for these stations behind the map.

Is it possible to provide this URL as the source dataset, using FME Desktop?

Does FME Desktop support real-time data applications?

I need to be able to access the URL and write a workbench that would capture the bike station status in 5 minutes increments and put the results in a .csv ot .txt file.

Any ideas?