texts: Find orientation of rectangle in PostGIS for QGIS composer

texts:
I am storing rectangle in my DB in order to automatize the creation of technical maps of a networks. I am creating these maps in QGIS using the print composer. Many of these rectangles are oriented (they do not follow a strict W/E orientation).

I have a field in my table to store the orientation value, with the idea that I would be able to give this value in the print composer to generate an atlas with all the maps at once, properly oriented.

I however do not manage to find the proper formula in PostGIS to fill the value of an orientation that works in QGIS. The left side of the map is always the smaller side of the rectangle, and I do not want to have the North arrow pointing further South than East or West (i.e. the orientation in the QGIS composer has to be between -90° and 90°).

I found two queries that seems to do the job, but nonetheless the angle is systematically reversed (I get 45° when I need -45°), they also don't take into account the "North arrow never pointing down" parameter.

Here they are :

SELECT ((degrees(ST_Azimuth(ST_PointN(ST_ExteriorRing(geom),1),ST_PointN(ST_ExteriorRing(geom),2))) + 180)::numeric % 360 ) - 180
FROM table;

SELECT CASE WHEN degrees(ST_Azimuth(ST_PointN(ST_ExteriorRing(geom),1),ST_PointN(ST_ExteriorRing(geom),2))) < -180 THEN degrees(ST_Azimuth(ST_PointN(ST_ExteriorRing(geom),1),ST_PointN(ST_ExteriorRing(geom),2))) + 360
ELSE degrees(ST_Azimuth(ST_PointN(ST_ExteriorRing(geom),1),ST_PointN(ST_ExteriorRing(geom),2))) END
FROM table;


texts: Graduated symbology makes some polygons disappear

texts:
I am working with QGIS 3.8 and I have layer of the world countries which I did Join in order to create data about each country (e.g population, corona cases and more).

Whenever I try to create graduated symbology based on one of the columns, some of the countries polygons dissapear, seems like only the ones that has no data.
So as you can see, when I have sungle symbol the map looks like this:

But when I try to create graduted symbology based on column where not all the coountries has data I get this:

I believe the countries that dissapere are the ones that have no values, meaning, if it has 0 it shows the country but if in the table it's empty it has nothing.

Is there any easy way to fix it, which doesen't require me to go to Excel and write zero in any country that has no value?

Is there any way to make those values null and also add them with grey color to the map?


texts: Finding out which U.S. state a tweet was generated in

texts:
I have a database of records with latitude and longitude information (imagine this is a CSV file containing tweets and their geo-coordinates). I also have shapefiles of U.S. states. I want to add another column which indicates the U.S. state where a tweet was generated.

In R, this would be very simple: just use the over function in the rgdal package, which is really simple and doesn't require the use of any for loops.

However, I'm less sure about how to do this in Python. It looks easy to check whether any point lies within any U.S. state, (though you won't be able to tell which state) and one can adapt the code to loop over each U.S. states.

However, I'm wondering if there's a simple function which wouldn't require the use of a for loop? If not, would someone be willing to demonstrate a relatively easy way to do this in Python?


texts: Counting event points by categories and matching to polygons, normalized by population

texts:
I'm new to CartoDB, but I'm a seasoned user of ArcGIS. My university was sorely lacking in GIS training; I only learned how to create static maps, not web content. After graduating 3 years ago, I'm playing major catch-up.

Anyhow, I have two primary datasets: Events (points) and categories (table).

The events are points in New York City of events such as "Yoga Meetup, Location: Central Park, 10 am on 6/1/2014" or "Python Pals, Location: Tech Center, 6 pm on 6/15/2014."

Each unique ID from the event table is listed in the concept table. The concept table tags every event by category - for example "Sports and Recreation" (category 1) or "Technology" (category 2) and so on. Some events have null categories because they haven't been tagged yet.

I want to link the categories to the events, and then count how many events fall under each category.

Next, I want to link these points to 5 different polygon datasets: NYC boroughs, City Council districts, Census tracts, and so on. There is a column in each table for the population of each polygon.

The goal is to count how many event points of each category are within each type of polygon -- normalized by population.

The end result is to see which neighborhoods are hosting the most events for each issue. Which borough has the most "Sports" events by population? Breaking it down further, which Census tracts have the most events by population?

Thanks for any help you can provide!

-S.G.


texts: How can I easily detect gaps in a polyline using arcmap 10.0 (arcview license)?

texts:
I'm making routes by selecting links from a streets layer which I export as a new temporarely layer in arcmap 10. This layer which I call route_temp unfortunately can be very long and contain gaps, for instance when i forget to select very small links which are difficult to see because they are overlaid by other layers. But sometimes the route contains gaps which have to be there because a cluster is in between it. What I would like to be able of is to visiualize all the gaps in the route, so that I can judge mysel which gaps have to be filled and which gaps have to stay in the route. Can anybody help me with this question, knows a tool for this, maybe a fieldcalculator-script or something else?

Thank you in advance


texts: Reduce memory used when reading LAS/LAZ using lidR

texts:
I want to save multiple height normalized .laz files as variables to my global environment in R, but the readLAS() function starts decompressing the .laz files, and there comes a point where my laptop can no longer allocate more RAM to R for storing variables. I tried using reducing the size of a laz file variable by  using the select function select = "xyz0" according to this question, which has reduced the size to 428 MB which was initially increased to 1.6 GB (without the select function) from 78.4 MB (original .laz file size).
Can I use the filter function to exclude anything else from the .laz file in order to further reduce the variable size, when the purpose is to create convex hull polygons? Basically, I want to reduce memory allocation size due to the way R stores LAS/LAZ files as variables.


texts: Printing from Feature Layer via ArcGIS API for JavaScript?

texts:
I am trying to print a map using a print service. I have followed the example on the ArcGIS site namely:

https://developers.arcgis.com/javascript/jssamples/widget_print_esri_request.html

However I can't get it to work using a local feature layer.

Can someone please give some advice?

This is what I have done so far:

    parser.parse();

    app.printUrl = "http://local_server/arcgis/rest/services/Layout/WSC_LAYOUT1/GPServer/Export%20Web%20Map";

    esri.config.defaults.io.corsEnabledServers.push("local_server");

    app.map = new Map("map");

    //Base map
    var basemap = new ArcGISDynamicMapServiceLayer("http://local_server/arcgis/rest/services/NAME/Basemap/MapServer/");
    basemap.setVisibleLayers([0, 1, 2, 3, 4, 5, 6, 7]);

     app.map.addLayer(basemap);

    // add graphics for feature
    var url = "http://local_server/arcgis/rest/services/NAME/NAME/FeatureServer/n";
    var pipeFeatureLayer = new FeatureLayer(url, {
      "mode": FeatureLayer.MODE_SNAPSHOT
    });
    app.map.addLayer(pipeFeatureLayer);


texts: Auto-generate Grid Reference for Point and Polygon features within existing Shapefile, QGIS

texts:
I have various point and polygon layers (shapefiles) which I use regularly in QGIS (3.6.3), and often need to work out their location as an 8 figure grid reference (British co-ordinates, EPSG-code 27700, same code used by each project). For polygons this means working out the central point of the shape (centroid?) first and then converting this to a grid ref.

At the moment I am using the FSC plugin tools, 'OSGR Tool' and 'Add GRs to Layers Tool' to do this, but this means I have to go through the process every time I have a new/edited feature which I need the info for.

I would prefer to be able to automatically generate the grid ref every time a point/polygon is added/edited, within the grid ref field of the existing layer (see below). Is this possible?

I presume it could be done using the 'default value' option for the existing field, but I cannot work out what expression I would need to do this.


texts: Drawing connecting lines between points to nearest within layer using ArcMap

texts:
I need to connect points within one layer to the nearest with polyline. Tool to the nearest would be ideal, but there is nothing in the attribute table.

I tried to find tool in X Y coordinate, but nothing works. It requires another attribute to run it to. Layer ID is not in any order. Tried XY tool. Tried Creating best fit line from point data using ArcGIS for Desktop?

Only connects some, and not always to the nearest.
None seem to work.
Any Suggestions?


texts: Is it possible to change interactivity in Carto.js map when called from Carto

texts:
I am making a carto map using carto.js and a map/data from my Carto account. I know I can set the interactivity there for clicking on polygons, etc. in the wizard, but am not able to alter that in my carto.js code.

Example:

cartodb.createLayer(map, vizjson)
                .addTo(map)
                .on('done', function(layer) {
                    layer.setInteraction(true);

                    var subLayerOptions = {
                        sql: "SELECT * FROM final_1",
                        cartocss: "#final_1{polygon-fill: #F84F40; line-color: #000;}",
                        interactivity: "name"
                        }

                    var sublayer = layer.getSubLayer(0);

                    sublayer.set(subLayerOptions);

                    sublayers.push(sublayer);

                    alert("this is a map");

                });

My subLayerOptions interactivity line does not seem to have any effect on the data when I include that line. Instead, it just displays a popup that is forever "loading". How can I set the interactivity in the carto.js code instead of on the Carto builder side?


texts: Missing Processing Tools in Vector Menu of QGIS?

texts: I'm missing all the usual processing tools in the Vector menu in QGIS 2.16 (eg vector overlay, analysis, etc). They do appear in the processing toolbox but there's no sign of them in the Vector menu. Has anyone come across this issue?


texts: Error connecting to a PostgreSQL database from QGIS using SSL certificates

texts:
I'm trying to connect QGIS (v3.14) to a PostgreSQL database hosted on Google Cloud.  I've successfully connected to the database through pgAdmin (v4.23) so my assumption is that there isn't a problem with the database itself, more likely a QGIS related issue.
I have also successfully connected to another PostgreSQL database that doesn't use SSL authentication, just a username and password.  So I suspect its related to the certificates / key and username / password combination.

.
The only other things I can think of are 1. there isn't a root.crt file (do I need one if I already have the certificate and key specified / referenced elsewhere?) or 2. maybe I need to make a change to the the pg_hba.conf (it currently just contains the standard info).
I couldn't find anything else on here that specifically answered the problem, but if I've missed it please point me to the link.


texts:
I'm a recent convert to QGIS but have been confused about the lack of support for drawing curved shapefiles, especially when Arc has almost too many. A lot of the work that our consultancy do is based around natural features such as lakes/rivers so there aren't many straight lines and being able to accurately represent or work with these features requires smooth curves as opposed to the tools I've been able to find on QGIS so far.

Other threads I've found relating to this have stated that being able to draw curves would require a fundamental redesign of QGIS, how is this though?

Other threads also recommend to use a rough shapefile first then simplify/smooth over to produce your curve but i would be worried about using this simply because you can't see the finished lay of the line when you draw it, which could be frustrating if your boundary moves after smoothing.


texts: Is there a method in SQL Server (2008 or 2012) that will return the equivalent of the Esri LabelPoint for geometry?

texts:
I'm looking for a method that is equivalent to the Esri LabelPoint for SQL Server Geometry (or fairly close).  Does it exist by a different name?   Or is there an after-market implementation?

I found this site, but it deals with actual centroids:
http://www.spatialdbadvisor.com/sql_server_blog/263/alternate-centroid-functions-for-sql-server-spatial

Related:
Algorithm for finding irrregular polygon centroid (label point)

Is there a way to get centroid of geographic object in MSSQL?


texts: Iterating an ArcGIS tool over multiple Fields with selection boxes

texts:
I have one shapefile that contains many fields and I'd like to write a simple script tool that allows me to run IDW for the fields I select. Is this possible?

For example.....

Input: 1 shapefile

Sample  field1  field2  field3  field4 .....
1       26.5    27.5    22.2    32.3
2       19.1    12.2    45.3    33.3
3       11.1    15.6    44.4    34.4

Options:

field 1 ☑

field 2 ☑

field 3 ☑
.....

Output: field1_IDW, field2_IDW, field3_IDW .....

EDIT: The following code works properly

# Import arcpy module
import arcpy
from arcpy import env
from arcpy.sa import *

# Check out any necessary licenses
arcpy.CheckOutExtension("Spatial")

inPointFeatures = arcpy.GetParameterAsText(0)

fieldList = arcpy.GetParameterAsText(1)

# Define local variables
cellSize = 250.0
power = 2
searchRadius = RadiusVariable(12)

for fieldName in fieldList.split(';'):
    outIDW = Idw(inPointFeatures, fieldName, cellSize, power, searchRadius)
    # Save the output
    outIDW.save("C:/temp/idwout_{0}.tif".format(fieldName))


texts: Retaining Feature class Information Using Construct Polygons Modify Feature in ArcGIS Pro

texts:
I have a polyline feature I need to turn into polygons but I don't have the license which would allow me to use the feature to polygon tool. The workaround as suggested on the Construct Polygons From Features help page is to use the construct polygons modify feature. The problem is that the fields information doesn't carry over. I've tried using a spatial join to get this information from the polyline to the polygon in various ways but nothing turned out quite right and I'm at a loss if this would even be the correct thing to do.
How can I get this information from my polyline feature to the newly created polygon feature?
Edit with Whitebox Tools Lines to Polygons
The tool only partly worked. See the screenshot below which outlines sections that broke during the polygon creation process.


texts: How to quickly change EPSG in QGIS (Openlayers; 4326 vs. 3857)

texts:
I'm using Qgis 2.12.2-Lyon on Ubuntu 14.04 LTS.

With the OpenLayers plugin I added some OpenLayers to my Canvas.
I usually use the WGS 84 (EPSG:4326) projection for my project, but I have to switch on the OpenStreetMap or BingAerial layers sometimes.

When I forget to change the projection to WGS 84 / Pseudo Mercator (3857) a priori then a python error occurs:

(Traceback (most recent call last):
          File "/home/.../.qgis2/python/plugins/openlayers_plugin/openlayers_layer.py", line 120, in request
            self.setup_map()
          File "/home/.../.qgis2/python/plugins/openlayers_plugin/openlayers_layer.py", line 185, in setup_map
            olSize = QSize(int(olWidth), int(olHeight))
        OverflowError: cannot convert float infinity to integer).

And the project is loading in perpetuity (sometimes crashes running on Windows).

Is there a way to avoid this behaviour?

Is it possible to create a short-key to change from any EPSG to WGS 84 / Pseudo Mercator (3857)?

Is it possible to have a drop down menu at the bottom right of the Canvas where you can change the EPSG a bit faster?


texts: Data for depth contours/DEM for Mediterranean sea

texts: I am searching for data for depth contours for the Mediterranean sea, especially the Aegean sea in any format (raster, shp). Does anyone know if there are any free available data?


texts:
I have added all data in PostgreSQL with the help of PostGIS.  I made a store in GeoServer, but after adding data when I click on the OpenLayers option it does not open, and a file downloaded is shown at the bottom (showing WMS service).

What is the problem?


texts: How to create a single contour defining the boundary between sea and land

texts:
Hello!
I am still learning so forgive me if my question is very simple... I am trying to create a map with ArcGIS. The area I am working on is a seasite and I found that in some parts the cells appear as white. Do you have any idea how I can fix it? Moreover, I would like to create a single contour along the seacoast just to highlight the boundary between the land and the sea. Does anybody know how I can do it? Thank you in advance! :)


texts: Prevent Data Driven Pages Strip Index tool from rotating polygons

texts:
When using Strip Map Index Features Tool the pages rotate with the linear feature, is there a way to prevent this?

So at the moment I get

or

and I want..

No rotation of polygons.

Using 10.6
Is this possible? I have tried changing between Vertical and Horizontal and the orientations.


texts:
Summary: How to use the MDB file (link below) as a locator to reverse geocode (point->address) with ArcMap?

I've implemented this sample:
http://help.arcgis.com/en/sdk/10.0/arcobjects_net/conceptualhelp/index.html#/How_to_find_the_address_closest_to_a_point_using_reverse_geocoding/0001000001n9000000/

The sample works with a file ESRI provided, Atlanta.mdb, and with DataMaps10 US Streets locator.

I am trying to get my own MDB files to work.  I can do the following with them:

//Address -> Point
//Address -> Closest Address Match
//Point -> Nearest intersection
//Intersection -> Point

I can't get Point -> Address to work.
I don't know how to identify the problem.  I project my IPoint coordinates to the projection of the MDB, but that doesn't work.  I've tried everything I can think of.  Can you give me some ideas of things to try to isolate the problem?   I have a temporary license for ArcMap for another week, so I can use that to help.  I'm using ArcEngine 10, C#, & VS 2010.

Edit:

I've uploaded a zip file here: http://www.filedropper.com/roads

It has an MDB, a shapefile and a projection.  I used the shapefile to create the MDB file.

I need to know how to reverse geocode this MDB file.  I need to just one more operation:

//Point -> Nearest Address

I would like to do it programmatically, but if you tell me the steps to do it in ArcMap, then I will have a decent idea how to implement it programatically, and that will be good enough.

To summarize:  How to use the MDB file as a locator to reverse geocode with ArcMap?

Edit #2:
When I implement the sample here: http://help.arcgis.com/en/sdk/10.0/arcobjects_net/conceptualhelp/index.html#/d/0001000001n9000000.htm

There is this line of code:

IPropertySet addressProperties = reverseGeocoding.ReverseGeocode(point, false /*return intersection */);

If the second parameter, return intersection == true, then it works.  But if the second parameter is equal to false, it throws a COM error,  0x80040217, which corresponds to
GEOCODING_E_REVERSE_GEOCODE_NO_ADDRESS_FOUND .  I've tried a lot of different coordinates and I have my search radius equal to 1 km, so something should be found.

Edit #3:

I found this
http://resources.arcgis.com/content/kbase?fa=articleShow&d=38314

Which points to this:  9.3.1 Locator Styles
http://resources.arcgis.com/gallery/file/geocoding/details?entryID=12D8D400-1422-2418-34B0-4FE1CC06C0ED

I am hoping that I can get these 9.3.1 styles in the same way that the empirical programmer got it to work, and that might be a solution.  If I can I'll post it, and give you the points empirical.


texts: Does gdalinfo report tifftag_datetime capture date?

texts:
Is gdalinfo meant to report the date a geotiff was captured? When I run it in qgis the tifftag_datetime seems to report the date it was saved to disk, rather that capture/creation date.

Is there a way I can see the capture date?


texts: Need to Extract by Mask on Features in a Feature Class

texts:
I have a polygon feature class of a country that has water basin polygons as individual features. I am trying iterate each feature and Extract by Mask using a set of rasters - for instance on Basin_ID I need to extract by mask from crop_a, crop_b, crop_c. I was able to create this sequence in Model Builder using the "Iterate Feature Selection" tool but I want to script this so that I have more control and because I need to collect a sum from zonal statistics down the line.

I have "tried" opening a cursor on the country, creating a for loop (row in cursor:) and then running a select query on the feature to make a new feature layer to use as a mask for the extract by mask tool. I assume that is what the model builder is doing but I can't seem to recreate this myself.

The beginning of my model essentially looks like this:


texts:
I have an entire Spain agronomic parcels layer in a PostgreDB and I´m triying to serve it through GeoServer.

I split by villages to improve the speed of the petitions, but how can I serve all the territory with one client petition (and quickly)?

My machine is an ubuntu 18 google machine. In postgres, "spain table" is divided by village (+5000 tables). All of these are in geoserver gruped by province (state) (50 groups). My first idea was group the 50 groups in one and call it, but is very low (minutes to load in the client web.

How can I do that?

I'm going to try with geowebcache but I'm not sure it can pregenerate the enormous amount of tiles.


texts: PostGIS - Calculate of the average distance between points for a delimited area

texts:
I am working with PostGreSQL

I have 2 data sets:
1) A points layer where each point represents a building.
2) A Polygons layer created from the points layer using convexhull and clusterwithin fonctions to aggregate points in areas.

Now what I am trying to do is calculate the average distance between all points (by area) within a radius of 2 km (my projection is EPSG2950 so i work with meters).

The result should be a new column in my polygon layer with a number that indicates what i stated above.

Any idea how such a request would look like?


texts: Whether ArcGIS can use the IIS published image as the base image and not use ArcGIS Server?

texts:
I posted the tile base map on IIS. I wanted to load it through ArcGIS for JS, but I could not use ArcGIS Server.
Just like OpenLayers
var center = ol.proj.transform([112.233778,22.688963], 'EPSG:4326', 'EPSG:3857');
var map = new ol.Map({
    view: new ol.View({
        center: center,
        zoom: 1,
    }),
    target: 'map'
});
var offlineMapLayer = new ol.layer.Tile({
    source: new ol.source.XYZ({
        url: 'http://localhost:2674/roadmap/{z}/{x}/{y}.png'
    })
});
map.addLayer(offlineMapLayer);

Or MapBox
var map = new mapboxgl.Map({
    container: 'map', // container id
    style: {
        'version': 8,
        'sources': {
            'raster-tiles': {
                'type': 'raster',
                'tiles': [
                    'http://localhost:2674/roadmap/{z}/{x}/{y}.png'
                ],
                'tileSize': 256,
                'attribution':
                    'Map tiles by Stamen Design, under CC BY 3.0. Data by OpenStreetMap, under CC BY SA'
            }
        },
        'layers': [
            {
                'id': 'simple-tiles',
                'type': 'raster',
                'source': 'raster-tiles',
                'minzoom': 0,
                'maxzoom': 22
            }
        ]
    },
    center: [-74.5, 40], // starting position
    zoom: 2 // starting zoom
});


texts: Batch processing near feature distance in ArcGIS Desktop?

texts:
I have 2000+ point shape files. These are all GPS tracks. I want to calculate the distance from each point to the nearest feature (road, train line, restaurant etc) within a 30 meter radius. I am using Near feature from proximity analysis in ArcGIS. I am facing two problems.

  Using batch processing one can add as many input rows in Near feature tool box as s/he wants. But this requires selecting each file
  one by one which is time consuming especially when there are thousands
  of files. So I tried using a python scripting as follows. But after
  processing the result is written in the same file. I want to write the
  table separately. How can I do that?
  In case I want to add different columns for each features say road, train network for each point record in a file using near feature
  how can I do that?

Currently, I am using the following python code but all the output are generated on the current file. I want to write it separately (I mean I want to generate a new table/file after near_analysis)

import arcpy # Import arcpy module

# path where all my point shp files are kept
arcpy.env.workspace = r'C:\ArcGIS sample data\'

#feature to which distance will be calculated
nearFeature_shp = 'C:\\GIS route network\\road.shp'

# looping through all the files
for file in arcpy.ListFeatureClasses ():

    # Process: Near
    arcpy.Near_analysis(file, nearFeature_shp , "30 Meters", "LOCATION", "ANGLE")


texts: Create Materialized View to Flatten ArcGIS JavaScript API Database

texts:
I'm doing a project using ArcGIS JavaScript API v3.5 hosted on ArcServer v10.1. While reading the question below, remember that I do not have control over the datasource and cannot just create a materialized view in the database to solve the issue.

Problem:

The issue I'm having is that parts of the database I'm querying are very very normalized -- one table has 11 Foreign Keys that point to parent tables that are even more normalized.

You can imagine that to join all the information for a single record for all these relationships is cumbersome. It requires many requests and lots of looping and checking which is slow and error prone.

Questions:
I have 2 questions:

What other options are there in the JS API to make querying easier from the client? I see layer types like QueryDataSource in the API reference that seem to be what i want. But there are no examples of how that works and no one online seems to be using that solution. Has anyone here used QueryDataSource?
I'm open to options on how to request that the datasources be modified. The most simple solution is to make them create a flattened table joined with all the data. What other options are there? Since I don't know the ArcGIS Server aspect as well I'm a little blind in what to request. For example, is there a way in ArcMap to create materialized-like views for a normalized datasource? Another possible solution is that the API should gather all related records up for me and spew out one big JSON object when i request a child table with many Foreign Keys. This is how I create API's in Django. All related records become nested hashes of a record. Is there someway to construct your API like this with ArcServer?

I can't wait to hear what thoughts/solutions others have for this problem


texts:
I have used 'ogr2ogr' for ESRI shape files conversion to Mapinfo tab files, but issue is that as soon as I open those converted .tab files over Mapinfo, the software stop working with a message

  MapInfo Professional has stopped working, Close the program.

Please help me out with the solution for the same.

I need a software which could convert n-number of ESRI shape files at one go into MapInfo tab files.


texts: Unable to run v.clean in QGIS 3.0.1 in Python console

texts:
I am trying to use v.clean algorithm in QGIS 3 python console. Code is as follows.

    lyr = iface.activeLayer()
    new_layer = lyr.materialize(QgsFeatureRequest().setFilterFids(lyr.selectedFeatureIds()))
    QgsProject.instance().addMapLayer(new_layer)

    alg = 'grass7:v.clean'

    extent = new_layer.extent()
    xmin = extent.xMinimum()
    xmax = extent.xMaximum()
    ymin = extent.yMinimum()
    ymax = extent.yMaximum()

    params = {"input": new_layer, "type": 1, "tool": 0, "-b": True, "-c": True,
              "output": "E:/out.shp", "error": "E:/errLayer.shp",
              "GRASS_REGION_PARAMETER": "%f,%f,%f,%f" % (xmin, xmax, ymin, ymax),
              "GRASS_SNAP_TOLERANCE_PARAMETER": 5, "GRASS_MIN_AREA_PARAMETER": 5,
              "GRASS_OUTPUT_TYPE_PARAMETER": 2}

    test = processing.run(alg, params)

When I run this code I get following error:

    Traceback (most recent call last):
    File "C:\PROGRA~1\QGIS3~1.0\apps\qgis\python\plugins\splitroads\Split_Roads.py", line 211, in run
    test = processing.run(alg, params)
    File "C:\PROGRA~1\QGIS3~1.0\apps\qgis\python\plugins\processing\tools\general.py", line 84, in run
    return Processing.runAlgorithm(algOrName, parameters, onFinish, feedback, context)
    File "C:\PROGRA~1\QGIS3~1.0\apps\qgis\python\plugins\processing\core\Processing.py", line 189, in runAlgorithm
    feedback.reportError(msg)
    TypeError: reportError() missing 1 required positional argument: 'fatal_error'

I have also checked input parameters of v.clean algorithm by using

   processing.algorithmHelp("grass7:v.clean")

and got following output

   Input parameters:

   input:  <QgsProcessingParameterVectorLayer>
           Layer to clean
   type:  <QgsProcessingParameterEnum>
          Input feature type
          0 - point
          1 - line
          2 - boundary
          3 - centroid
          4 - area
          5 - face
          6 - kernel
    tool: <QgsProcessingParameterEnum>
          Cleaning tool
          0 - break
          1 - snap
          2 - rmdangle
          3 - chdangle
          4 - rmbridge
          5 - chbridge
          6 - rmdupl
          7 - rmdac
          8 - bpol
          9 - prune
          10 - rmarea
          11 - rmline
          12 - rmsa
   threshold: <QgsProcessingParameterString>
              Threshold (comma separated for each tool)
   -b:  <QgsProcessingParameterBoolean>
        Do not build topology for the output vector
   -c:  <QgsProcessingParameterBoolean>
        Combine tools with recommended follow-up tools
   output:  <QgsProcessingParameterVectorDestination>
            Cleaned
   error:  <QgsProcessingParameterVectorDestination>
    Errors
   GRASS_REGION_PARAMETER:  <QgsProcessingParameterExtent>
    GRASS GIS 7 region extent
   GRASS_SNAP_TOLERANCE_PARAMETER:  <QgsProcessingParameterNumber>
    v.in.ogr snap tolerance (-1 = no snap)
   GRASS_MIN_AREA_PARAMETER:  <QgsProcessingParameterNumber>
    v.in.ogr min area
   GRASS_OUTPUT_TYPE_PARAMETER:  <QgsProcessingParameterEnum>
   v.out.ogr output type
    0 - auto
    1 - point
    2 - line
    3 - area

   Outputs:

   output:  <QgsProcessingOutputVectorLayer>
     Cleaned
   error:  <QgsProcessingOutputVectorLayer>
     Errors

Everything seems fine but I am unable to run it.Can you please let me know how to run v.clean successfully in python console using QGIS 3?


texts: Math formula to calculate GPS location for weighted GPS data points

texts:
I'm working to improve the location and accuracy of a map / data-point driven application.
I have a collection of locations of interest. Each location of interest has a collection of data points that include lat/long and are weighted from 1 to 100. 100 means the exact location should be used and 1 means this location should be considered the least in the output of whatever formula to find the center of the location of interest.
How can I find the center of the location of interest while including the weights to help improve accuracy?
I found a code sample to calculate this but it doesn't consider weights in the data points.
/**
 * Calculate the center/average of multiple GeoLocation coordinates
 * Expects an array of objects with .latitude and .longitude properties
 *
 * @url http://stackoverflow.com/a/14231286/538646
 */
function averageGeolocation(coords) {
  if (coords.length === 1) {
    return coords[0];
  }

  let x = 0.0;
  let y = 0.0;
  let z = 0.0;

  for (let coord of coords) {
    let latitude = coord.latitude * Math.PI / 180;
    let longitude = coord.longitude * Math.PI / 180;

    x += Math.cos(latitude) * Math.cos(longitude);
    y += Math.cos(latitude) * Math.sin(longitude);
    z += Math.sin(latitude);
  }

  let total = coords.length;

  x = x / total;
  y = y / total;
  z = z / total;

  let centralLongitude = Math.atan2(y, x);
  let centralSquareRoot = Math.sqrt(x * x + y * y);
  let centralLatitude = Math.atan2(z, centralSquareRoot);

  return {
    latitude: centralLatitude * 180 / Math.PI,
    longitude: centralLongitude * 180 / Math.PI
  };
}

// expect ~ 37.790831, -122.407169
const sf = [{
  latitude: 37.797749,
  longitude: -122.412147
}, {
  latitude: 37.789068,
  longitude: -122.390604
}, {
  latitude: 37.785269,
  longitude: -122.421975
}];

console.log(averageGeolocation(sf));

// expect ~ 8.670552, -173.207864
const globe = [{ // Japan
  latitude: 37.928969,
  longitude: 138.979637
}, { // Nevada
  latitude: 39.029788,
  longitude: -119.594585
}, { // New Zealand
  latitude: -39.298237,
  longitude: 175.717917
}];

console.log(averageGeolocation(globe));


texts:
In Noosa 3.6.3 when I run a tool (for example: Regular Points) I save the output to file and choose SHP files. But the name I choose do not appear in the table of content. Instead, when the tool finish running, it adds a generic name of "Regular Points" and not the name I picked for the SHP file.
Is there a way to force QGIS tools to add the output file actual name and not the name of the tool like its doing now? See snapshots.


texts: Understanding proj4, epsg and esri *.prj in depth?

texts:
Where can I find a good resource that teaches me

1) about the difference between proj4, EPSG and ESRI *.prj? (Why do several coordinate system definitions exist? When should I use which? How are they defined?)

2) how the conversion between these systems works in detail.

I came across these questions while encountering a problem with EPSG 2056 which I defined using the proj4 string (I wrongly assumed that the EPSG is an ID for the proj4 string). When exporting this coordinate system to a shapefile using the function st_write from the R package sf I realized that the ESRI *.prj file did not translate correctly. I assume it has something to do with the towgs84 paramter which ESRI does not recognize.


texts: Enabled GeoServer Layer has Yellow Triangle with Exclamation Mark Warning

texts:
I am having an issue with one of my layers in GeoServer (version 2.1.3). The layer does not appear on the Layer Previews page or in GetCapabilities for either WMS or WFS; however, it does appear using OpenLayers on my website.

When I look at my Layers and Cached Layers pages, the layer has a yellow exclamation mark in the "Enabled?" field. What does this warning mean? Why would a layer not appear in layer preview or GetCapabilities, but appear on the website?

Edit: It is a large vector layer with points pulled from a directory of shapefiles. I do not see any error messages in the log related to this layer. There is the following error for a different layer, not sure if it's related:
2019-11-08 12:57:49,078 ERROR [geoserver.ows] - org.geoserver.wfs.WFSException: Feature type MHC:pa_areas unknown

The "Enabled" and "Advertised" boxes are checked on the Edit Layer page. Updating to a current version of GeoServer is not an option at the moment.

I am new to using GeoServer.


texts: Create a line string from street geometry

texts:
I'm currently facing some hard times trying to create a simplified line from a street geometry.
The data that I'm currently using comes from OpenStreetMap, using osm2pgsql.(Which means PostGIS, but I'm open to solutions that don't use it, or use external code, like C, JS, or whatever).

The problem:
I need to create a simplified line string from given streets, however, most streets have multiple segments, such as:

In this example, the street is forked and segmented, the GeoJSON is available here.

Dumping the points of the street is useful which could be converted to a line string, however, their order is messed up, causing bizzare line formations.
I could partially solve this issue ordering by lon/lat or lat/lon, but in several cases of points in the same axis or due geometry angle, it does not order them properly, and thus, does not work for all streets.

Few solutions I came up with:

Using ST_Buffer on every segment, joining, and ST_ApproximateMedialAxis on the resulting geometry. It didn't work because some streets, even with a generous buffer have all it's segments touching each other and also ST_ApproximateMedialAxis doesn't always result in a line string (Not a function fault, it's the expected behavior).
Ordering by LAT/LON ~ LON/LAT and creating a line from the result: Also didn't work as explained above.
Using ST_LineMerge: Nope, in regular streets it does work, but on forked or segmented streets, it fails, as explained here.
Iterating over points to create a line from every point nearest neighbor: Fails due the inability to find the real first and last points.

The hardest task so far is getting the farthest points from given set, to define the first and last point, getting this informations would solve the problem for me, but I'm clueless so far?

I thought about using ST_VoronoiPolygons, but I have no idea how to extract some result from that?


texts: QGIS 3.10.7 - Subscript/superscript in the legend text

texts:
I'm trying to add subscript/superscript in the legend text of a layout in QGIS. I can do it in a simple text box just click on View as HTML and then using for example <sub> in order to write a subscript.
How could I do this in the legend text? I don't find anything to enable HTML on the legend property panel.


texts: selecting a percentage (75%) of a group of points based on the distance from one separate point in ArcGIS?

texts:
This is ArcGIS specific.

I have 2 point shapefiles A and B, the first one (A) is a single point that contains a lat long, the second (B) is a multitude of points (over 12k) that each contain their lat and long.  What I am trying to do is to automate the selection of 75% of shapefile B's points based on the distance from shapefile A.  In other words I want to select the nearest 75% of shapefile B points to shapefile A's one point.


texts: How do I get statistics out of a table using Python in Arc 10?

texts:
I am very new to programming (I started learning this semester) and would very much appreciate helop with one piece of a project. I have the following table, and I need to write a script compatible with Arc10 that will extract some statistics from it. Specifically, I need to obtain a total of all the areas, and then a percentage describing how much of the total area each class of land takes up. If someone knows of a tool or command I should use, that would be a lot of help, even without a script. Thank you in advance to anyone who wants to help an extremely raw beginner.

Id_1    LandClass   F_AREA

    0   Mixed Forest    15612315.78
    0   Softwood Forest 2161121.455
    0   Wetland         2255656.822
    0   Grassy          27009516.24
    0   Hardwood Forest 10345393.4


texts: Styling polyline to view which part is above another in self-intersection in QGIS?

texts:
I need to view a polyline constructed from waypoints of a gps track that has self intersections in curves; like an exit loop from a highway.

I need to visualize it in a way that shows which part is above the other in the intersection.

Option 1: Style with solid color

This does not work, as the intersection is filled without delimiting the side of the path.

Option 2: Style with some predefined road from QGIS:

Almost doesn't work, as each segment is surrounded by the line, not only at both sides but also in the head and tail of the segment.

Option 3: Enabling symbol levels:

This does not fit my need as it merges the intersections, thus preventing from seeing the precedence.

I've managed to do the some styling after finding this article in http://planet.qgis.org: http://nathanw.net/2012/07/22/qgis-style-tricks-using-styles-to-help-fix-kerb-line-directions/

The key for me was as easy as using a style of a marker line, with a single triangle marker, and a solid line. Then, I applied transparency so that street names are visible.

This is the result:

What other option could I try?


texts: Why does tool take longer to run when in Python script?

texts:
I have written a python script that does some euclidean allocation. For that I use arcpy.gp.EucAllocation_sa() method. Now, my python script takes 30-40% extra time than spatial analyst euclidean allocation tool.

    arcpy.gp.EucAllocation_sa(in_raster_name, out_allocation_raster_name, "", "", cellSize, "VALUE", out_distance_raster_name, "")

Question:

Is there any process that help me to improve the running time of this method? Is there any environment setting issue?


texts: Discovering gravimetric station data, including location and measurement?

texts:
I would like to refine an precision estimate of local gravity based upon current GPS position. Ideally, similar to the National Weather Service API, I would like to interrogate a server with a REST interface, providing location and accepting gravimetric data in return. If no such service is available, I would consider hosting mgal data and referencing its value by location, possibly with interpolation.

Together with other natural measurements, this gravity data would be an input to a developing Android application. (No, I don't deem the device's native sensors sufficiently accurate for this; thus, I am seeking a benchmark reference.)


texts: How to define map extents based on a feature's geometry?

texts:
There is a polygon geometry in my QGIS project which represents the extents of my map. Does anybody know how to use the data defined override to get the X/Y min and X/Y max coordinates of the polygon in the print composer?

Unfortunately this isn't working:

xmin(attribute(getFeature('layername', 'gid', 1), 'geom'))


texts:
In the previous version of QGIS when you loaded a multiband color image into raster calculator the bands were separated.

In version 3.6 I am only seeing 1 band how do I get the other bands displayed ?


texts: Is it possible with TileMills and CartoCSS?


texts:
I am following a tutorial at https://courses.spatialthoughts.com/pyqgis-in-a-day.html
that gives some instructions on how to set up a new toolbar item with a custom icon in a QGIS map using PyQGIS in the integrated python console
import os
from datetime import datetime

icon = 'question.svg'
data_dir = os.path.join(os.path.expanduser('~'), 'Downloads/pyqgis_in_a_day/')
icon_path = os.path.join(data_dir, icon)
print(icon_path)

def show_time():
    now = datetime.now()
    current_time = now.strftime("%H:%M:%S")
    iface.messageBar().pushMessage('Time is {}'.format(current_time))

action = QAction('Show Time')
action.triggered.connect(show_time)
action.setIcon(QIcon(icon_path))
iface.addToolBarIcon(action)

The code runs without error but no icon gets placed in the toolbar area of QGIS
I have checked the path to the icon and tested that also by supplying the full path, copy and paste, to avoid typos but no success. I am using QGIS version 3.14 and checked over the documentation and the methods seem to be correctly
Can anyone see what I might be doing wrong?


texts: Unresolved externals on simple mapnik c++ application

texts:
I'm attempting to write a simple win32 x64 console application in c++ using VS 2017 that uses the Mapnik v3.0.9 SDK to render a png from a pgsql database.

All I've done so far is include a reference to the mapnik/map.hpp header file and attempting to build produces the following unresolved externals:

  unresolved external symbol "private: static enum mapnik::logger::severity_type mapnik::logger::severity_level_" (?severity_level_@logger@mapnik@@0W4severity_type@12@A)

  unresolved external symbol "private: static char const * * mapnik::enumeration::our_strings_" (?our_strings_@?$enumeration@W4well_known_srs_enum@mapnik@@$01@mapnik@@0PEAPEBDEA)

  unresolved external symbol "private: static class std::basic_string,class std::allocator > mapnik::enumeration::our_name_" (?our_name_@?$enumeration@W4well_known_srs_enum@mapnik@@$01@mapnik@@0V?$basic_string@DU?$char_traits@D@std@@V?$allocator@D@2@@std@@A)

  unresolved external symbol "private: static char const * * mapnik::enumeration::our_strings_" (?our_strings_@?$enumeration@W4aspect_fix_mode@Map@mapnik@@$08@mapnik@@0PEAPEBDEA)

I've included the mapnik.lib and all boost lib files for the linker, I'm not sure how to proceed here. Any ideas on what I could be missing?

#include <mapnik/map.hpp>

int main()
{
    return 0;
}


texts:
I had a geodatabase file that contained two geodataset, which are line and area. I clicked on the topology to area and then instead of deleting the topology, all the information and edits were deleted from the  ArcCatalog.
Is there any chance to restore my data and files? I lost too many fields and too much information.


texts: Finding equidistance between two overlapping polygon shapefiles in ArcGIS Desktop?

texts:
I have 2 polygons that represent districts from neighboring cites. The Polygons extend 1.5 miles out from the city border. There is an area between the 2 cities where these polygons overlap each other.

Is there a tool that i can input the 2 polygons and have it output them with the borders meeting in the middle between the two cities instead of overlapping?


texts: Deleting points which are close to each other

texts:
I have a vector file containing thousands of points. Where points are within x distance from each other, I wish to delete one of the points and keep the other (either by snapping to one point, or simply erasing one of the two points)

However, in Whitebox the vector cleaning tool only works with polygons and lines.

Is there some other way I can do this with my shapefiles?


texts: Reclassified slope raster attribute table in QGIS

texts:
I have reclassified a slope raster in Qgis 3.4, but how to open or generate attribute table of that reclassified slope raster?

Slope ranges are on below picture.

and my desired attribute table is as below:


texts: Edit tilemill legend to add circle symbol

texts:
I am trying to edit my tilemill legend so that instead of displaying a red rectangle as below, it displays a red circle (I'm making a proportional symbol map).

What line would I need to add in my legend code (see below) to make this happen?

On the mapbox website it says : "The advantage with html/css is that you can quickly make edits to the legend directly in TileMill, and maintain the ability to manipulate the legend styling with css even after the map has been exported to MBTile format. However, you are limited to right angles and solid colors, and may have to write many lines of code to create a relatively simple design." https://www.mapbox.com/tilemill/docs/guides/advanced-legends/

Hopefully "limited to right angles" does not mean that I cannot create a circle?


texts: cannot add second interaction, openlayers 3

texts:
This is my code in Openlayers 3.9.0

I just want to pass some extra Interactions to the map after I create it. I get no errors in my console. The weird part is that I can dragRotate, this works fine, but the keyboardZoom does not work. I fail to understand why keyboardZoom fails, or how to fix this.


texts: Click to Openlayers shows blank pop up - All Formats for Images Works

texts:
I am totally new to GeoServer and all this WebServer stuff, so please understand that all the gray hairs on my head over the last 2 weeks are from trying to get this all set up and working.

When I select my layer to preview with OpenLayers, I get a pop up that opens and is blank, except for the bounding box and the text to the bottom left and right. Not sure how to fix this so I click and have it show the map item, which in this case is a shapefile.

I should note that I can select the 'All Formats' form the Layer Preview and if it is an image, I can see and/or download the shape as a jpg, pdf, tiff, gif, so I know the data is being read, just not displayed.

Anyone else have this issue and resolved it?


texts: ArcPy Insert Cursor giving RuntimeError: Error in InsertRow

texts:
When I run the code snippet below, to create a test table, add a field to it, and populate that field with a value in its first row, I get an error:

Traceback (most recent call last):
  File "C:\Temp\test.py", line 9, in <module>
    cursor.insertRow(row)
  File "C:\Program Files\ArcGIS\Pro\Resources\ArcPy\arcpy\arcobjects\arcobjects.py", line 230, in insertRow
    return convertArcObjectToPythonObject(self._arc_object.InsertRow(*gp_fixargs(args)))
RuntimeError: Error in InsertRow

import arcpy
arcpy.env.overwriteOutput = True
arcpy.env.workspace = r"C:\Temp\Projects\TestProject\TestProject.gdb"
arcpy.CreateTable_management(arcpy.env.workspace,"TestTable")
arcpy.AddField_management('TestTable',"ID","SHORT")
cursor = arcpy.InsertCursor('TestTable',["ID"])
row=[50]
cursor.insertRow(row)
del cursor


texts: Georeferencing Raster automatically fails

texts: While trying to auto register a raster image in georeferencing - using a wms raster layer as reference - I get the simple error message "Failed to auto register". Thereby the Link Table remains empty. I am currently using ArcGIS Desktop 10.3.1 and followed all the steps in the related ESRI Article (up to #7). The WMS Layer is projected in WGS84, data frame rotation is set to zero.


texts: GDAL VRT Pixel Function no output with signed or float data type

texts:
I am trying to average multiple partially overlapping rasters using a VRT and a Python based Pixel Function (https://gdal.org/drivers/raster/vrt.html#using-derived-bands-with-pixel-functions-in-python). The rasters have cells with no data.
The function on this page was my starting point: https://gis.stackexchange.com/a/324680/169329. The problem I am running into is that I am unable to change the dataType to one that allows me to have fractions and negative values without having useless output.
I have tried multiplying by 1000 and rounding the values before building the VRT and averaging the values using a Pixel Function, but this works only for positive values with dataType UInt32.
My expectation is that it should be possible to use other dataTypes in Pixel Functions besides those that only contain positive and rounded numbers. How can I do this?
I have used the following code after first multiplying the rasters by 1000 to try and get an average before dividing by 1000 again:
<VRTRasterBand dataType="Int32" band="1" subClass="VRTDerivedRasterBand">
<PixelFunctionType>average</PixelFunctionType>
<PixelFunctionLanguage>Python</PixelFunctionLanguage>
<PixelFunctionCode><![CDATA[
    import numpy as np

    def average(in_ar, out_ar, xoff, yoff, xsize, ysize, raster_xsize,raster_ysize, buf_radius, gt, **kwargs):
        div = np.zeros(in_ar[0].shape, dtype = 'Int32')
        for i in range(len(in_ar)):
            div += (in_ar[i] != 0)

        div[div == 0] = 1
        y = np.sum(in_ar, axis = 0, dtype = 'Int32')
        y = y / div
        out_ar[:] = np.round(np.clip(y,-100000,100000))
]]>
</PixelFunctionCode>


texts: Mapserver tiles - trouble with projection

texts:
I'm trying to render tiles with shp2img tool from Mapserver suite and my tiles aren't matching. Tiles are numbered like in openstreetmaps.org. This is an attempt to draw tile (level 1, x=1, y=0) which has extent [xmin=0, ymin=0, xmax=180, ymax=90]:
http://img69.imageshack.us/img69/4112/z0x1y0.png
This picture is rendered with extent [xmin=0, ymin=0, xmax=180, ymax=90].
See the problem? Ymin is clearly not zero here.
And this one with extent [xmin=0, ymin=0, xmax=90, ymax=90]:
http://img834.imageshack.us/img834/5489/z0x1y0v2.png
Now ymin looks zero as it should. It seems like Mapserver tries to keep aspect ratio and is extending rendered area to do that. But this ruins my tile service.
How do I setup projection in my MAP file to get matching tiles? Or what else to do to get tiles right.

Update: this is how tile z=1, x=1, y=0 looks in openstreetmaps.org, I guess I need reproject my map in such projection?

Some details:
ogrinfo output:

Layer SRS WKT:
GEOGCS["Unknown",
    DATUM["MIF 0",
        SPHEROID["WGS_84_MAPINFO_Datum_0",6378137.01,298.257223563]],
    PRIMEM["Greenwich",0],
    UNIT["Degree",0.017453292519943295]]

Projection section from MAP file:

PROJECTION
    "init=epsg:4326"
END


texts: Using --config flag with the Python GDAL gdal.Translate()

texts:
So I've recently learned that GDAL and OGR command line utilities can be used as functions now.  See this site for more info.  For example:

from osgeo import gdal
tmp_ds = gdal.Warp('temp', 'in.tif', format = 'MEM', dstSRS = 'EPSG:3857')
gdal.Translate('out.png', tmp_ds, format = 'PNG')

This works amazing!  But one thing I can't figure out is how to pass in the --config flag.  In my case, I am using gdal_translate to convert a US Topo GeoPDF to GeoTiff, exporting out certain PDF layers:

gdal_translate -of “GTiff” --config GDAL_PDF_LAYER “Map_Collar,Map_Frame.Terrain” ustopo.pdf ustopo_terrain.tif

I've tried adding to the gdal.Translate in 2 ways:

gdal.Translate('out.tif', tmp_ds, format = 'GTiff',options='--config GDAL_PDF_LAYER “Map_Collar,Map_Frame.Terrain”')

and

gdal.Translate('out.tif', tmp_ds, format = 'GTiff',config='GDAL_PDF_LAYER “Map_Collar,Map_Frame.Terrain”')

Neither work.  I am totally stumped.


texts: Is there some kind of limitations of number of features when we work with Network Analyst in SDE ? My net has 170.000 features, is it a problem ? Because it take at least 6 hours to finish the "Build Network", and some times with errors. Is it normal?


texts: Mapbox use point result from Geocoding/Search to queryRenderedFeatures on polygon layer that intersects

texts:
I've been quite stumped on this and while I have seen it used in other online examples I can't find anything on stack or in the documentation.

Right now, I am trying to use the search bar result (a point) to identify which polygon from my existing layers added through paint it lies within.

Right now, the code below functions, but it always returns the wrong polygon within the fed_searched variable.  In tests it is often within a short distance, but not reliably correct.  I also found that without the moveend function, if the point is close to the edge of the polygon layer it doesn't find anything for the fed_searched variable.  I am wondering if this is because my source layer for the polygons is an .mbtiles source.

[This example from mapbox][1] is what I pulled a lot of my logic from, though it is based off of a clicked point.

I am struggling because I have been able to run this same time of process with click points without any issues in selecting the correct polygon. See below:

    map.on('click', function (e) {
        // query which fed has been clicked
        var fed_clicked = map.queryRenderedFeatures(e.point, {
            layers: ['fed_click_fill']
        });

        console.log(fed_clicked);
    });

See my current code below:

    // geocoder/searchbar
    var geocoder = new MapboxGeocoder({ // Initialize the geocoder
        accessToken: mapbox_token, // Set the access token
        mapboxgl: mapboxgl, // Set the mapbox-gl instance
        marker: false,
        flyTo: {
            zoom: 10,
            speed: 1, // make the flying slow
            curve: 1, // change the speed at which it zooms out
            easing: function (t) { return t; }
        },
        countries: 'ca'
    });

    // Add the geocoder to the map
    map.addControl(geocoder);

    // Add the source(empty) and point layer for geocoding results to the map
    map.addSource('search', {
        type: 'geojson',
        data: {
            type: 'FeatureCollection',
            features: []
        }
    });

    map.addLayer({
        id: 'search_point',
        source: 'search',
        type: 'circle',
        paint: {
            'circle-radius': 3,
            'circle-color': '#448ee4'
        }
    });

    geocoder.on('result', function (e) {
        map.once('moveend', function () {
  // add the result as a point in the 'search_point' layer to show up as marker
            var geocoder_result = e.result.geometry;

            console.log(geocoder_result);

            map.getSource('search').setData(geocoder_result);

            var coordinates = e.result.geometry.coordinates;

            var fed_searched = map.queryRenderedFeatures(coordinates, {
                layers: ['fed_click_fill']
            });

            console.log(fed_searched);
// should be polygon feature containing the searched point

        });
    });

  [1]: https://docs.mapbox.com/mapbox-gl-js/example/queryrenderedfeatures-around-point/


texts:
I have a large set of points and need to determine the individual point with the maximum value for a given attribute within 12 km quadrants and extract these to a separate shp. Currently, I have nothing but the points so the quadrants would need to be generated.

I am relatively new to GIS so my efforts have not shown any progress. How can I accomplish this?


texts: How can I draw points on my mobile and then pass it to my computer?

texts: we are doing a survey to know how many people still use outside bathrooms in the Caribbean, I want to go house by house using a shapefile image on mobile and point the location of points of interest in different houses and them merge all this info on the PC, is there a way to do this?


texts: Sandbars feature extraction from satellite image

texts:
I want to know if we can analyse
the permanancy of sandbars from a satellite image? Whether the formation of the sandbar is permanent or temporary.


texts: How to install additional driver to GDAL installed with apt-get

texts:
I have installed GDAL 2.2.3 on Ubuntu 18.04 using apt-get. I'd like to add the GeoJSONSeq driver.

How can this be done?


texts: What projection is being used to display in QGIS for Geographic Coordinates

texts:
Google/Bing Maps etc. or web mapping applications use WGS84 coordinates. When you moves the mouse on map, you see WGS84 latitude and longitude. Although coordinate system is WGS84, map is projected in Web Mercator on the screen. Spatial data with projected coordinates is in 2D and it's easy to show on screen, because it's in 2D and already a plane coordinates.

I wonder what is projection system QGIS uses for spatial data that has just geographic coordinates.


texts: How to access a deleted features attributes (before committing)

texts:
I have a signal setup for a feature deletion event on a layer:

self.layer.featureDeleted.connect(self.feature_deleted)

When I try and grab an attribute from the deleted feature on the slot, a StopIteration error is thrown (i.e it cannot retrieve the feature).

I would have thought I  should still be able to retrieve attributes from a deleted feature before the changes to the layer are committed? Otherwise I don't see what use the fid parameter is.

def feature_deleted(self, fid):
    feature = self.layer.getFeatures(QgsFeatureRequest().setFilterFid(int(fid))).next()
    idx = self.layer.fieldNameIndex('my_attribute')
    print feature.attributes()[idx]

How can I access the attributes of a deleted feature (before the changes to the layer are committed)?


texts: Automatically adding increasing transparency values to gradient stops in QGIS

texts:
I have a point layer in QGIS 3.14 on which I applied a heatmap symbology based on an intensity attribute (here a price value):

But the heatmap colors extends all over the basemap layer, covering it all with the lowest values, here represented as a yellowish hue.
Is there a simple and automated way to apply a linearly increasing transparency from lowest to highest values of the heatmap colors so that the basemap layer stays visible in regions where the intensity is low?
E.g. something like this:

In order to achieve the above result, I actually had to tweak the transparency values for each of the 9 first gradient stops using a +10% increment, as represented by the italic number under the first stops on this figure (0=full transparency, 100=full opacity):

This can be tedious if you want to expand it linearly on all the gradient stops as it can be represented by the green line over the color ramp, or let's say, half of them, which I think would be enough. You would have to calculate the amount of transparency to add for each step, and doing it manually, with the (high) risk of accidentally moving a gradient stop when clicking on it.
That's why I'm searching for a way to do that simply and automatically, and of course being able to save the results so that it's reusable.
Please notice that I'm not searching to apply a global transparency.


texts: Specify KML pop-up window appearance when exporting shapefiles from R to KML

texts:
I am trying to export shapefile data from GIS to KML with R and want to specify the appearance of the pop-up window of the KML file when using it in Google Earth. I know of course that it's easily possible within the GIS programs but want to do it with R.

The example below works fine but the pop-up windows don't really look good.

library(raster)
library(plotKML)
library(colorRamps)

col <- getData(country = "COL", level = 1)

kml_open("colombia.kml")
kml_layer(obj = col,
      folder.name = "testData",
      altitude = 100,
      plot.labpt = F,
      html.table = paste("Info 1: ",as.character(col@data$NAME_1),"  ||  ",
                          "Info 2: ",as.character(col@data$HASC_1),"  ||  ",
                          "Info 3: ",as.character(col@data$NAME_0),"  ||  ",
                          "Info 4: ",as.character(col@data$TYPE_1),"  ||  ",
                          "Info 5: ",as.character(col@data$ENGTYPE_1)))

kml_close("colombia.kml")

I would like specify line breaks e.g. "\n" instead of "  ||  " but that doesn't work with paste. The function cat, where "\n" works fine, doesn't seem to work within the html.table command.

I didn't find any further information or examples regarding the html.table command in the kml_layer function and also tried to solve my problem with the htmlTable-package but didn't have success.

Does anyone have a tip or link with information to solve it?


texts: Linear Referencing Methods: Station Equation origins

texts:
I am working in the pipelines field, and I get in touch with the concept of the linear referencing methods, specifically with station equation method.

I had to do research about this method, I'm quite confident about my understanding of the topic and the concept associated to it, but a set questions came to my mind:

Where did this topic originate?
Who proposed it?
When was it proposed?

After an intensive Google search I have not found a clue about the history of station equation method. I got a lot of technical documents that describe the concepts and how this concept works, but not its origins.

I am seeking any book or paper suggestions


texts: Dice a polygon with over 21 million vertices in ArcMap?

texts:
Based on posts like this I suspect I'm bumping up against issues when trying to work with a feature class that contains features with a single very complicated polygon. These errors include 000426, 000072, 99999, and 99998. I'm attempting to use the Project tool, Cut Polygons Tool, Add Geometry Attributes tool, and Feature Class to Feature Class tool. My operations complete successfully on other feature classes with features that have far fewer vertices.

I isolated the feature in the problematic feature class with the greatest number of vertices (over 21 million) and put it in its own shapefile. I've tried using the Dice tool to break this one polygon feature into several smaller features. I've tried setting the vertex limit both low (100) and high (1,000,000), as I'm not sure which would be better for minimizing memory requirements. The tool most often completes with this same error:

ERROR 000426: Out Of Memory

But I've also gotten it to intermittently complete successfully. I've only had success with Dice under these conditions:

Recently restarted computer.
Freshly restart ArcMap, but not necessarily the first start since the computer has been turned on.
Allow the Dice tool's default recommended output name and location, on the first run of the tool in this ArcMap session (saves to [user]\Documents\ArcGIS\input_file_name_Dice.shp).
Vertex limit set to 1,000,000 (no commas).

After completing successfully under these conditions, I can run Dice again with a different output name and/or location (shapefile or feature class, at a different directory or not) but cannot change the number of vertices. It will re-run successfully, some unpredictable number of times (0-4, perhaps)

I've already run this feature class through the Repair tool -- which crashes ArcMap, but completes without error when run as an arcpy standalone script. Conversely, when I run the Dice tool through an arcpy standalone script, it errors out the same as above, but with an additional error prepended:

argisscripting.ExecuteError: ERROR 000072: Cannot process feature with OID 0
ERROR 000426: Out Of Memory

I have been able to perform some operations on the original feature class, now with the most complicated polygons removed (such as Project[ion]) which I couldn't do before. But some of the remaining features in the original feature class also have millions of vertices and are causing similar problems.

My system has 16gb of RAM, and literally has nothing installed on it other than what come with Windows 10 and the ArcMap 10.6 suite of software.

So here's my question: Am I missing something that would empower me to consistently run the ArcMap Dice tool on these large-number-of-vertices feature classes? Does anyone have workarounds other than restarting the computer, opening ArcMap, immediately running Dice and hoping it works on this attempt? Perhaps I've misidentified the root cause of these errors (but I don't think so).

I should mention that the organization I'm working with requires the exclusive use of ESRI/ArcGIS software tools to achieve this. If any tool was available to me, I'd definitely use PostGIS... but this is not an acceptable solution to my question.


texts: Modelling zero-inflated spatially autocorrelated data with extreme values

texts:
I have spatially explicit data (above; Longitude, Latitude) with measurements on a random variable (z). This data was not collected in a random or systematic fashion, but instead was recorded incidentally during other survey activities.

The data (see above; left hand figure is the raw data, right hand figure are log-transformed non-zero observations) is zero-inflated (approximately 90% of observations are zero), with a small number of extreme values thrown in for good measure.

I want to model this data in space so I can interpolate/predict values for the entire range of the data, but I am unsure what is my best approach. There appears to be so many options (e.g. geostatistics, GAMs, thin plate splines), but I am unsure where to start. R is my preferred analysis platform.


texts: Temporal segmentation algorithms using Google Earth Engine

texts:
I'm using the temporal segmentation algorithm Vegetation Regeneration and Disturbance Estimates - Verdet - for land cover change detection. The output layer is created but I can't display or export it (GEE returns processing errors). Does anyone know how can I work with this tool?

// NDVI image collection from Landsat TM
var collection5 = ee.ImageCollection('LANDSAT/LT05/C01/T1_SR')
    .filter(ee.Filter.eq('WRS_PATH', 219))
    .filter(ee.Filter.eq('WRS_ROW', 72))
    .filterDate('2008-01-01', '2011-12-31')
    .filterMetadata('CLOUD_COVER', 'equals', 0);

  var segIndex = function(img){
      var index = img.normalizedDifference(['B4', 'B5']);
      return index;
  };
  var collection5 = collection5.map(segIndex);

// NDVI image collection from Landsat OLI
var collection8 = ee.ImageCollection('LANDSAT/LC08/C01/T1_SR')
    .filter(ee.Filter.eq('WRS_PATH', 219))
    .filter(ee.Filter.eq('WRS_ROW', 72))
    .filterDate('2013-01-01', '2017-12-31')
    .filterMetadata('CLOUD_COVER', 'equals', 0);

  var collection8 = collection8.map(segIndex);

// Merge collections
var merge = ee.ImageCollection(collection5.merge(collection8));

// Run VerDET
var verdet = ee.Algorithms.TemporalSegmentation.Verdet(merge);


texts: How to connect CartoDB.js toggle buttons to a dropdown menu?

texts: I have gone through the CartoDB.js toggle map view tutorial, and I have been successful in linking up SQL queries to a handful of buttons. Now I would like to instead control these functions with a Bootstrap (or similar) dropdown menu, to reduce the clutter on the screen. However, I am having trouble modifying the code properly. All I have managed to do is get the dropdown button show up on the CartoDB map. How do I make the items in the dropdown menu show up when I click on the dropdown button? How do I link up the items in the dropdown menu to the SQL query functions that I have written?


texts: Arcpy FeatureClassToFeatureClass_conversion is causing ArcDesktop to Crash

texts:
I am new to coding in python and specifically Arcpy. Through trial and error I was able to piece together this code that was working fine yesterday until I installed service pack 5 for ArcGIs Desktop 10.0, on a 64 bit Windows machine. The reason I decided to install this service pack is because the code was incapable of doing more than ~60 iterations of the loop without stopping unless everything on my machine was closed except ArcCatalog. I saw on a discussion board that updating to a new service pack might solve my problem. Now my problems is the machine crashing halfway through one loop.

Through trial and error I have determined that using the line of my code which calls on the
FeatureClassToFeatureClass_conversion crashes ArcCatalog. The only thing that I can trace the problem back to is upgrading to service pack 5.

import arcpy
arcpy.env.scratchWorkspace=r'C:\Users\fk128474\Documents\Freds_Bitterroot_Model\wrf_data\Scratch'
arcpy.env.workspace=r'C:\Users\fk128474\Documents\Freds_Bitterroot_Model\wrf_data\2001_LW'

import os

xcelFiles = arcpy.ListFiles('*.xls')
count = 1

arcpy.env.workspace=r'C:\Users\fk128474\Documents\Freds_Bitterroot_Model\wrf_data\2001_LW'

for count in range(0, 140):
    fileName = xcelFiles[count]
    tableName = fileName[:-4] + '_Table_View'
    layerName = fileName[:-4] + 'Point_Layer'
    xyName = fileName[:-4] + '_pts'
    tempName = 'LW' + fileName[11:-4]

    shpPath = r'C:\Users\fk128474\Documents\Freds_Bitterroot_Model\wrf_data\2001_LW\Shape'
    spRef = r"Coordinate Systems\Geographic Coordinate Systems\World\WGS 1984.prj"'

    points = arcpy.MakeTableView_management(fileName + '\\sheet1$', tableName)

    b = arcpy.MakeXYEventLayer_management(points,"F2","F3",tempName, spRef, "F4")
    arcpy.FeatureClassToFeatureClass_conversion (b, shpPath, tempName)

I have not included the rest of my code as this last line causes the program to crash.

Thanks

Fred


texts: Displaying single symbol in legend for proportional symbol layer

texts:
I'd like to be able to show just one line in the legend for my layer of multiple proportional symbols. Previously, my workaround was to create a secondary layer not displayed on the map, and to bring that into the legend to serve as a "dummy" symbol. But now, the client wants a .lyr file with symbology preserved, so I must be able to replicate a single-line legend feature.

Is there any way to do this in ArcMap?


texts: How to access a table in a personal geodatabase using C#

texts:
I've a personal geodatabase which has some shape files included. To get a "complex" relation I've created a single table which will be have 1:1 or 1:N references which the different shape files (atm. 3).

Now I want to know how I can:

access the table in the personal geodatabase,
add new entries to this table,
remove/modify existing entries

using C#. Does anybody done this before and can help me to solve this question?

Thanks a lot!


texts: Is it possible to show markers from MySQL database with Javascript?

texts: I'm interested if it's possible to show markers (latit, longit) from MySQL database with Javascript libraries like Leaflet?


texts: SQLite Driver available with R's rgdal:writeOGR?

texts:
Does anybody know how I can get OGR's SQLite Driver to work with R's rgdal package?
When listing the drivers in R with the ogrDrivers() command the SQLITE driver is missing..


texts:
I try to use Intersects filter but I get an empty geometry in return. What did I miss in the syntax?
https://data.gov.au/geoserver/brimbank-stormwater-pits/wfs?request=GetFeature&version=1.0.0&typename=brimbank-stormwater-pits:ckan_3f07bec5_99c2_491f_8d39_b6061891ebb1&outputFormat=GML2&FILTER=<Filter xmlns="http://www.opengis.net/ogc" xmlns:gml="http://www.opengis.net/gml"><Intersects><PropertyName>geom</PropertyName><gml:Polygon xmlns:gml="http://www.opengis.net/gml" srsName="EPSG:4326"><gml:outerBoundaryIs><gml:LinearRing><gml:coordinates>144.7769072860001,-37.733680591086539 144.77697968000007,-37.733955687086535 144.78617369000006,-37.734882328086506 144.78611577600009,-37.734592753086503 144.7769072860001,-37.733680591086539</gml:coordinates></gml:LinearRing></gml:outerBoundaryIs></gml:Polygon></Intersects></Filter>


texts: Simplifying code to get each feature class record from shapefile?

texts:
How to simplify this code and do some improvement to get each feature class record from the shapefile?

import arcpy

featureClass = r"C:\Users\TL\Desktop\812test\SplitAllRoutes_FacID\FacID1655.shp"

count = int(arcpy.GetCount_management(featureClass).getOutput(0))
# In this case, count = 6

input = featureClass
output0 = r"C:\Users\TL\Desktop\812test\test\Route0"
arcpy.Select_analysis(input, output0, "{0} = 0".format(arcpy.AddFieldDelimiters(input,"FID")))

input = featureClass
output1 = r"C:\Users\TL\Desktop\812test\test\Route1"
arcpy.Select_analysis(input, output1, "{0} = 1".format(arcpy.AddFieldDelimiters(input,"FID")))

input = featureClass
output2 = r"C:\Users\TL\Desktop\812test\test\Route2"
arcpy.Select_analysis(input, output2, "{0} = 2".format(arcpy.AddFieldDelimiters(input,"FID")))

input = featureClass
output3 = r"C:\Users\TL\Desktop\812test\test\Route3"
arcpy.Select_analysis(input, output3, "{0} = 3".format(arcpy.AddFieldDelimiters(input,"FID")))

input = featureClass
output4 = r"C:\Users\TL\Desktop\812test\test\Route4"
arcpy.Select_analysis(input, output4, "{0} = 4".format(arcpy.AddFieldDelimiters(input,"FID")))

input = featureClass
output5 = r"C:\Users\TL\Desktop\812test\test\Route5"
arcpy.Select_analysis(input, output5, "{0} = 5".format(arcpy.AddFieldDelimiters(input,"FID")))

My original plan is to change the code like this, but it doesn't work. Any suggestions?

import arcpy

featureClass = r"C:\Users\TL\Desktop\812test\SplitAllRoutes_FacID\FacID1655.shp"

count = int(arcpy.GetCount_management(featureClass).getOutput(0)) # In this case, count = 6

for i in range(count):
    input = featureClass
    output%s = r"C:\Users\TL\Desktop\812test\test\Route%s" % (i, i)
    arcpy.Select_analysis(input, output%s % i, "{0} = %s".format(arcpy.AddFieldDelimiters(input,"FID") % i))


texts: How to create a Value Map widget type to select multiple values (checklist) from the dropdown?

texts:
I have a polygon layer that represents tree/shrub planting sites that were done. As part of the attribute table, I have fields representing things like tree density, species, soil work, tree size, etc. Taking tree sapling size (height) as an example field, my value map looks like:

< 1m

1m

,>1m
(ignore comma)

Currently, in the Attribute Form, I simply have a Value Map widget type, making it possible to select 1 of the 3 value options.

However sometimes, we plant a mix of tree sizes. Is there a way to make the value map a checklist, so that we'd be able to choose 2, or all 3, of the above value options?

I am using QGIS 3.4.


texts: What is the difference between mosaicking and merging?

texts: What is the difference betweem image Mosaicking and Merging?


texts: Area calculations using a layer in a projected data frame versus a projected feature class in ArcGIS

texts:
So I've been tasked with automating a formerly manual process in model builder in which I'm trying to calculate land area.

In the former manual process the feature class was brought into ArcMap and the dataframe's projection was switched from a Lambert conformal conic projection to a Albers equal area conic projection. The area was then calculated through the Calculate Geometry dialog in the attribute table.

To replicate this in a model I use the Project tool which creates a new feature class in the Albers projection and then the Calculate Field tool using !shape.area@SQUAREKILOMETERS! as the expression.

Generally speaking this all works fine. The issue that comes up is that in comparing previous results from the manual process to the automated process there are slight differences in the calculated area. Most of the time this is on the order of a fraction of a square meter in some features with the difference for the whole feature class hovering around 1 square meter total.

Is there a reason for this difference? Is one method more correct than the other? My thought is that it has something to do with tolerances and vertices snapping to different locations during projection of the feature class but I'm not sure.

Ultimately the difference is miniscule but there is a difference and I'm wondering why?


texts: Intersection of 2 Linestrings (not interpolated)

texts:
I am struggling with a query in POSTGIS

I want to find linestring geometries where the geometries are containing the same point.

I don't want to calculate the intersection points. I only want to get a point, if the point is a real member of both geometries.

Example

  SELECT ST_AsText(ST_Function(ST_GeomFromText('LINESTRING(1 2, 2 2, 0 0)'),
          ST_GeomFromText('LINESTRING(0 1, 2 2, 3 0)') ) )

should return

  POINT(2 2)

and

  SELECT ST_AsText(ST_Function(ST_GeomFromText('LINESTRING(1 2, 10 10)'),
          ST_GeomFromText('LINESTRING(1 1, 2 2, 3 3)') ) )

should return

  EMPTY

or something similiar for empty.

How can I achieve that in a single query?

I tried to use ST_Intersection(), but this one calculates the point (which is wrong to my requirements).


texts: Requesting Styled Google Maps XYZ tiles with QGIS seems to ignore colors

texts:
At the time of posting this question unauthenticated requests for Google Maps tiles work.

E.g. the following URL scheme would correctly load the tile shown below

https://mt1.google.com/vt/lyrs=r&x={x}&y={y}&z={z}&
apistyle=s.e%3Al.i%7Cp.v%3Aoff&s=Galileo&style=api%7Csmartmaps

Note that this is a styled tile as s.e:l.i|p.v:off style directive turns icon visibility off. Without this directive a highway number sign would be visible.

This works as an XYZ tile layer in QGIS 3.

However setting custom colors would not work. For example to make the roads grey I would add an additional style directive s.t:3|s.e:g|p.c:#ffe0e0e0.

The modified URL scheme:

https://mt1.google.com/vt/lyrs=r&x={x}&y={y}&z={z}&apistyle=s.e%3Al.i%7Cp.v%3Aoff
%2Cs.t%3A3%7Cs.e%3Ag%7Cp.c%3A%23ffe0e0e0&s=Galileo&style=api%7Csmartmaps

Produces the expected modified result.

However the corresponding XYZ Tile layer in QGIS seems to ignore this directive and still shows the roads white.

Is there a way to customize Google Maps tile colors in QGIS?


texts:
I'm using QGIS 3.0 and I have a problem with checking
validity of the geometries of a shapefile.

I imported a DXF, converted the polylines into polygons, manually corrected the wrong or intersecting polygons. Now I would like to cut out a raster using my building shapefile as a mask, but continue to
generate the message in figure 1:

I then made the topological check, then the geometric one, but nothing
seems to solve the error as seen in figure 2:

I'm almost thinking of doing all the polygons manually, but I would like to understand how to solve the problem in case I find myself developing projects with thousands of polygons in the future.


texts: Adding new attribute field to processing output algorithm in PyQGIS?

texts:
I have a processing script that will carry out a selection by expression of the input layer, then do a join attributes by location with a join layer.

The next stage will be to add a new field and calculate the field, I got a python script that runs this but I am struggling to adapt it into the processing tool.

The code is below:

def initAlgorithm(self, config=None):
    """
    Here we define the inputs and output of the algorithm, along
    with some other properties.
    """

    # We add the input vector features source. It can have any kind of
    # geometry.
    self.addParameter(
        QgsProcessingParameterFeatureSource(
            self.INPUT,
            self.tr('NRD Layer'),
            [QgsProcessing.TypeVectorPoint]
        )
    )
    self.addParameter(
        QgsProcessingParameterFeatureSource(
            self.JOIN,
            self.tr('Parliamentary Constituency Layer'),
            [QgsProcessing.TypeVectorPolygon]
        )
    )
    self.addParameter(
        QgsProcessingParameterField(
            self.JOIN_FIELDS,
            self.tr('Choose Join Field'),
            parentLayerParameterName = self.JOIN,
            allowMultiple=False,
            optional=False
        )
    )
    self.addParameter(
        QgsProcessingParameterVectorDestination(
            self.OUTPUT,
            self.tr('Output layer')
        )
    )
def processAlgorithm(self, parameters, context, feedback):
    selectData = processing.run(
        'qgis:selectbyexpression',
        {
            "INPUT":parameters['INPUT'],
            "EXPRESSION":'\"mcmcode\" !=\'1\' AND \"mcmcode\" !=\'9\' AND \"mcmcode\" !=\'960\' AND \"mcmcode\" !=\'999\'',
            "METHOD":0,
            'OUTPUT':parameters['OUTPUT']
        },
        context=context,
        feedback=feedback
    )
    joinoutput = processing.run(
    'qgis:joinattributesbylocation',
    {
         'INPUT':QgsProcessingFeatureSourceDefinition(selectData['OUTPUT'], True),
         'JOIN':parameters['JOIN'] ,
         'PREDICATE':[0],
         'JOIN_FIELDS':parameters['join_field'],
         'METHOD':1,
         'DISCARD_NONMATCHING':True,
         'PREFIX':'',
         'OUTPUT':parameters['OUTPUT']
    },
    context = context,
    feedback = feedback
    )
    NRDLayer = joinoutput['OUTPUT']
    F_Onset = QgsField('F_Onset',QVariant.Double)
    NRDLayer.dataProvider().addAttributes([F_Onset])
    NRDLayer.updateFields()

    return {'OUTPUT':joinoutput['OUTPUT']}

The error suggests you can't use dataprovider on joinoutput['OUTPUT'].

I thought I might need to use ParameterasSink to get my output in order to be able to use things like .featureCount(), .dataProvider() with it etc.

However not understanding this very well I would rather avoid that.
The tool works up to that point.

Any ideas on what to do in QGIS 3?


texts: Problem opening Natural earth quick start for QGIS

texts: Trying to open the file "Natural_Earth_quick_start_for_QGIS.qgs" a window opens, listing 20 "bad layers". The window is entitled "Handle bad layers". Consequently I cannot even get started with Natural Earth in QGIS, although all advice I find assumes that opening this project is a piece of cake. what should I do?


texts: Grouping latitude and longitude coordinates by city neighborhood

texts:
I'm working on a model to predict whether craigslist apartment listings are over or under market value. I'm a complete novice with GIS.

I have the lat and long as well as what neighborhood the author says the apartment is in. But seeing how some neighborhoods are nicer or better known than others, the apartment's actual address is often in other neighborhoods.

The picture is an example of apartments listed under the neighborhood SOMA but many are outside the neighborhood borders.

Is there any way in Python to group lat and longs by neighborhood?


texts: Determine the range of values in one raster associated with certain values in another raster

texts:
I have two rasters, one with thematic data, and the other with continuous data (floats).  I am trying to calculate/determine the the range of float values associated with the different classes in my thematic raster.  For instance, in my thematic raster, I have a class A.  If there are 5 cells in my continuous raster data that overlay thematic cells of class A, and these 5 cell have values of 1.2, 3.7, -4.5, 6.23, -0.8, I would like the end result of my calculation to show:

Class A:  Range [-4.5 to 6.23]
Class B:   .....and so on.

Does anyone know of a way to do this?  I've been trying to determine if the raster overlay tools in ArcMap can provide this, but I haven't come across a solution yet.  I'm hoping to avoid a very complicated equation in the raster calculator.

Another way to produce the analysis I am looking for would be to generate a histogram for each thematic class in one raster based on the associated (overlaid) values in another.  Is this possible?


texts: Select by Attribute within a specific area (using Select by Location?)

texts:
I have a map with four layers. Two of these layers cover an entire state, and the other two focus on the area I'm actually interested in within the state.

I'm trying to Select by Attribute with one of the larger focused layers.

Is there a way I can Select by Attribute but only within the bounds of the smaller layers? I think there is a way, possibly using Select by Location, but I can't seem to figure out what it is.

Hoping there's a way to do this within the GUI of ArcMap rather than any Python scripting, though I'm pretty sure using a Python script would make this so much easier.


texts: QGIS points disappear after changing CSV file location, changing data source doesn't help

texts:
I moved my database and shape files, and wanted to relink them to my QGIS project. The shape file worked fine, but the CSV file wasn't reloading properly no matter what I tried. I could see the fields and the attribute table all fine, but the points just wouldn't appear, and there wasn't even a checkbox to the left of the layer in the layers panel like there usually is.
I don't want to remove the layer and re-add it, although that did work in making it visible. But then I would lose all the properties that I've already set to that layer.
I'm using QGIS v3.14.
I'm new to QGIS and mapping in general.


texts:
I'm using QGIS 2.18 (Las Palmas)

I have hourly Data (of a whole year) of 11 Wind Measuring Points in an area.

How would i best structure my Data (8760 hours for a whole year) for 11 measuring points in order to have the possibilty to interpolate within the points?

Or can anyone perhaps name some scientfic papers i could consult?


texts: Calculating statistics for defined grid cells in ArcMap

texts:
I have a vector that defines the grid cells (black lines in attached image) and a raster dataset that contains landscape metrics for over 50,000 patches (in different columns in the attribute table). One of the landscape metrics is patch area. I want to calculate summary statistics (e.g. mean, majority/mode, median) for that specific landscape metric within each grid cell.

I've tried using zonal statistics but it doesn't seem to work, as there is no where I can select the 'Area' column for which the stats will be calculated.


texts: where can I get land use data for Beijing

texts: Iam looking for land use data for Beijing of 1985 and 2010


texts: How can I check that a PostGIS database (schema) is correctly backed up?

texts:
I'm something of a novice user of PostGIS. I have some valuable data in a single schema on a PostGIS database. I believe that I've saved a backup of the data etc to a local file using pgAdmin - but many of the options available when doing this are unfamiliar/confusing. I'd like to have some way to check that the data is correctly backed up (without risking the existing live data). Can this be done simply?

I have access to pgAdmin and QGIS.

....Additional....

I believe there to be longer term backups of the whole database (etc) as part of my company's ordinary IT processes - what I'm wanting to have to hand is a process where (as one of the supervising GIS staff) I can take steps to protect the data - this might be before making bigger changes to large amounts of data, or before bigger changes to table structure.


texts: How to upgrade PostGIS 2.0 to 2.1 on Windows 7?

texts:
PostGIS 2.1 was recently released and I would really like to upgrade from 2.0, but can't quite figure it out.  I had previously installed PostgreSQL 9.2 and PostGIS 2.0 using StackBuilder.  Running SELECT PostGIS_full_version();, shows I currently have the following install base:

  "POSTGIS="2.0.3 r11132" GEOS="3.3.8-CAPI-1.7.8" PROJ="Rel. 4.8.0, 6
  March 2012" GDAL="GDAL 1.9.2, released 2012/10/08" LIBXML="2.7.8"
  LIBJSON="UNKNOWN" TOPOLOGY RASTER"

So I downloaded the tar.gz file, and after uncompressing it, where do I copy the new 2.1 folders/files so that I can do a "soft upgrade" using the ALTER EXTENSION postgis UPDATE TO "2.1.0"; method?

Or am I missing additional intermediate steps, too?


texts: Clip global data by polygon using rioxarray fails (off by 180 longitude)

texts:
I want to extract CMIP6 climate data within polygons defined by a selection of the Large Marine Ecosystems of the world. The data(CMIP6) are found in ranges from 0 to 360 degrees longitude, while my Shapefile containing the polygons of the LMEs are -180 to 180 degrees longitude. I use Cartopy to do my plotting and I define the projection as ccrs.PlateCarree(central_longitude=-180). When I plot the polygons of interest the projection of my selected polygons are shown correctly (see image)

Next I read in my CMIP6 data (open archive on Google Cloud) and reassign their coordinates to be rolled 180 degrees longitude. I then try to extract data residing within the selected polygons. My extraction works but the data found are 180 degrees off (data for California are from the Mediterranean), suggesting the data are still not rolled correctly 180, or that I am doing something else wrong.

I wonder if anyone can shed some light on to what I am doing wrong here in terms of rolling the data 180. Is there some option I have to pass to rioarray? This should be trivial....

projection=ccrs.PlateCarree(central_longitude=-180)
ax = plt.subplot(1, 1, 1, projection=projection)
ax.coastlines()

extent = [-248, -100, 20, 85]
ax.set_extent(extent)
ax.add_feature(cfeature.LAND)
ax.add_feature(cfeature.COASTLINE)

records = pices.get_LME_records()

for key in config_pices_obj.dset_dict:
    ds=config_pices_obj.dset_dict[key]

    ds = ds.assign_coords(lon=(((ds.lon + 180) % 360) - 180)).sortby('lon')
    ds = pices.preprocessing_wrapper(ds)

    # Set the CRS of the CMIP6 data
    ds = ds.rio.write_crs(4326)
    print("Geometric projection of dataset : {}".format(ds.rio.crs))

    counter=0
    all_lmes=[]
    for LME_GEO, record in zip(records.geometries(), records.records()):
       LME_NAME = record.attributes['LME_NAME']
       LME_NUMBER = int(record.attributes['LME_NUMBER'])

       if (LME_NAME in config_pices_obj.LMES):
           print("=> LME_NAME {} ({})".format(LME_NAME,LME_NUMBER))
           ax.add_geometries([LME_GEO], projection,
               facecolor="None", edgecolor='k')

           geometries = [get_geometry_for_LME(LME_NAME)]

           # Clip the data inside the LME:
           clipped = ds.rio.clip(geometries=geometries, crs=ds.rio.crs)
           all_lmes.append(clipped)

      #  ds.tos.isel(time=0).plot(ax=ax, cmap='coolwarm',transform=ccrs.PlateCarree())

      #  plt.show()

   ds_lme=xr.merge(all_lmes)
   ds_lme.tos.isel(time=0).plot(ax=ax, transform=ccrs.PlateCarree(), cmap='coolwarm')

   plt.show()


texts:
I have a spatial table (polygons,polylines) in a database. These are digitized by multiple people. For example, a state is divided into multiple Grids and given to multiple people to digitize using QGIS; which is connected to a postgres database. Now, I have to merge all those polygons and polylines.
I don't know how to do this via ST_Union.
I Used the following snippet to merge the polygons and polylines
CREATE TABLE schema_name.polygon_un2 AS
SELECT upg.un_geom,npg.Column1,npg.Column2,npg.Column3,npg.Column4,npg.Column5
FROM (
 SELECT st_union(the_geom) As un_geom,(
  SELECT oid FROM schema_name.polygon
  WHERE oid = 100748191
 ) FROM schema_name.polygon WHERE oid IN (100748191,100763519)
) upg
INNER JOIN schema_name.polygon npg ON npg.oid = upg.oid;

DELETE FROM schema_name.polygon WHERE oid IN (100748191,100763519);

I have the interesecting polygon oid in a Seperate Table, I want to do this process in an update query and that should be a dynamic plpgsql function.

Green : Grid , Other : Polygon
But This snippet did not work !!
The snippet which i expect is

It should create a new row for the merged polygons and polylines (eg oid - 1 and 2 merged to oid 3).

The polygons which are not merged (eg 1 and 2 should be deleted)


texts: ArcGIS Server and Microsoft SQL Server database connection

texts: I have an ArcGIS Server and I need to register a Microsoft SQL Server database as a data storage. When I connect to database directly from ArcMap everything is fine, but when i'm trying to register a database I have a message, that 'connection properties lost parameter or it cannot be defined. Trying to connect with older version of sql server client, that is not supported' (I have a localized version). Where can be a problem? I connect to SQL Server DB using domain authorization.


texts: Extract representative slope profiles from hillslope polygons

texts:
With the help of the Whitebox Geospatial Analysis Tools I was able to easily extract hill slopes from a river network in combination with a flow direction map.

Now I would like to extract representative slope profiles from each of this hill slopes. Representative slope profiles are 2D plots giving me the distance and elevation to the river for a representative flow line in a single hill slope.

I was wondering if there is already a open source tool for this since erosion models like 'WEPP' or 'Catflow' need this kind of input for modeling?

There are some toolboxes for ArcGIS but I would like to use open-source tools for several different reason.


texts:
It's been a while since I was here. I have this problem now.
I have electricity line feature class and point feature class intersecting line features. What I need to do is simplify this line feature class but effects must not effect those vertices which intersect with points. Using ArcGIS 10.2.2.
Any ideas?


texts: Custom marker colours in leaflet based on attribute

texts:
I'm aiming to customize my markers (leaflet maps) into 3 different colours based on different values in an attribute column.
Each marker is assigned a value A-C and I would like to represent each value as a different colour marker.
I have been trying to use a switch statement but cannot seem to get it to work, any tips?

here is my current code:

    var customMarker = "A";

    switch(Data[i].value){
        case "A":
            markercolor (orange);
        break;
        case "B":
            markercolor (blue);
        break;
        case "C":
            markercolor (green);
        break;
        default :
            markercolor (pink);
        }

     var markerOptions = {
          color: customMarker,
          opacity: 1,
      };

     var marker = new L.circleMarker(markerLocation,markerOptions).addTo(map);


texts: How does QGIS calculate the area of a polygon?

texts: I am working on a project that requires the calculation of the area of a polygon by two different methods.  I need to know what algorithm QGIS uses to calculate the area of a polygon.  I know ESRI ArcGIS uses the shoelace (Gaussian) method and am looking for GIS software that uses a different method.


texts: How do you link QGIS legends in different print composers?

texts:
I have a long linear QGIS project with about 20 separate print composer layouts.  I'd like to link each legend to the original one, so I don't have to copy and paste 20 times every time I need to update the legend.

What's the best way of doing this?


texts: Remove item from legend on a several mxd's

texts:
I used this code from this question here: Remove item from legend using ArcPy?

However I keep getting an error at the following line:

legend = arcpy.mapping.ListLayoutElements(mxd,"LEGEND_ELEMENT")[0]

that says:

  AttributeError: 'list' object has no attribute 'listLegendItemLayers'

However, I clearly have a legend in the mxd's.

import arcpy
import arcpy.mapping
import os
import sys
from arcpy import env
import string

env.workspace = r"C:\Project"
for mxdname in arcpy.ListFiles("*.mxd"):
    print mxdname
    mxd = arcpy.mapping.MapDocument(r"C:\Project\\" + mxdname)
    legend = arcpy.mapping.ListLayoutElements(mxd,"LEGEND_ELEMENT")[0]
    for lyr in legend.listLegendItemLayers():
        if lyr.name == "atikot":
            legend.removeItem(lyr)
            print 'remove'
    mxd.save()
del mxd


texts: Projection change in QGIS after import & export in Grass

texts:
I'm working on an application written with Grass command in python and I have a projection problem.
First I create a location with belgian lambert 72 (EPSG:31370) in Grass 7.0.

After importing all the environment variables I run this code:

gscript.run_command("r.in.gdal", input=path+"initialDEM.adf",flags='o', output="DEM")
gscript.run_command("g.region",rats="DEM",res=10)
gscript.run_command("r.out.gdal",input="DEM",type='Byte',output=pathOutput+'finalDEM',overwrite=True,flags='f')

Now I want to import the output raster (finalDEM) in QGIS. And I don't know why the projection change as you can see on the picture. (in red initDEM and in blue finalDEM)

The SCR of the raster finalDEM is USER:100000, *SCR generated ....
Do you know why I have this change of projection? Do you know how I can avoid it in my python code?