texts: Points to raster (image) in Google Earth Engine

texts:
I have a bunch of points representing bird observations, and I need to convert this to an image in GEE where each pixel value is the total number of observations in each 1 km pixel. This is a very common task with ecological data, and in ArcGIS I would use Point to Raster. But I can't figure out what the equivalent would be in GEE.

As an example, here are some points and a bounding polygon:

var points = ee.Feature(ee.Geometry.MultiPoint(
    [[57.43610672115551, -20.298790475544028],
     [57.566569367639886, -20.221490947834383],
     [57.695658723108636, -20.319397176387103],
     [57.676432648889886, -20.171225637115814],
     [57.58991531490551, -20.216336278108052],
     [57.57892898678051, -20.225356838032646],
     [57.572062531702386, -20.227934044777072]]),{});

var poly =
  ee.Feature(
    ee.Geometry.Polygon(
    [[[57.36469558834301, -20.087414125113316],
      [57.36469558834301, -20.499588506612334],
      [57.85908035396801, -20.499588506612334],
      [57.85908035396801, -20.087414125113316]]], null, false), {});

So the resulting image should have the extent of poly, a scale of 1 km, and pixel values corresponding to the number of points in each pixel. I thought reduceToImage would do this, but that seems to be more for summarising properties of polygon geometries.


texts:
I would like to insert the current position of the device at the first waypoint. With a function onLocationFound I get the position but I can't add it to the waypoint. Does not work.

is there a best solution?

function onLocationFound(e) {
        var radius = e.accuracy / 2;
        var curr_latitude = e.latitude;
        var curr_longitude = e.longitude
        L.marker(e.latlng).addTo(map).bindPopup("You are within " + lat + " meters from this point").openPopup();
        L.circle(e.latlng, radius).addTo(map);
    }

function onLocationError(e) {
        alert(e.message);
    }

    map.on('locationfound', onLocationFound);
    map.on('locationerror', onLocationError);

    map.locate({setView: false, maxZoom: 16});

map.on('locationfound', onLocationFound);
    map.on('locationerror', onLocationError);

    map.locate({setView: false, maxZoom: 16});

L.Routing.control({

routeWhileDragging: false,
    draggableWaypoints: false,
    addWaypoints: false,
    waypoints: [
        L.latLng(curr_latitude, curr_longitude),
        L.latLng(45.5062563,9.1864641)
    ],
    language: 'it',
    showAlternatives: false,
    lineOptions: {
    styles:
            [
                {color: 'black', opacity: 0.15, weight: 9}, //sombra
                {color: 'white', opacity: 0.8, weight: 6}, // Contorno
                {color: 'blue', opacity: 1, weight: 4}] // Centro
   },
          summaryTemplate:  '
Percorso: {name}

Distanza: {distance}
',
          show: false,

}).addTo(map);


texts: Removing newline from AddMessage() with ArcPy

texts:
I am looking for a way to add text in the same line in the progress dialog box.
The AddMessage() adds a newline and this is not convenient.

For example:
I want to have one line with the text:

arcpy.AddMessage("Sorting data...Done")

where "Sorting data..." should appear before the execution of the tool and "Done" after it.

Is there a way to do this?


texts: Landsat Level-2 products vs Landsat ARD products?

texts:
I understand that they both pre-process to give surface reflectance, but other than that I do not understand the difference.

Can anyone please explain it to me?


texts: Spatial Analysis: use of buffer function inside GeoDjango's LineString model

texts:
I'm trying to make a buffer around a line. I've a model with a LineString field and I would to use the buffer function like this to make my buffer. The aim of this is render the buffer on a webmap using MapBox.

I have this simple model:

class LinesSource(models.Model):
    geom = models.LineStringField()

def __int__(self):
    return self.pk

@property
def buffer(self):
    """pseudo-code"""
    distance = 200 ##I need the distance in meters and my model SRID is in degree but this is not a problem
    mybuffer = geom.buffer(width=distance, quadsegs=8) ##something like GEOSGeometry.buffer(width, quadsegs=8)
    return mybuffer

My problem is that is not clear to me how I can change the pseudo-code to real code because for me is not clear how I can use the spatial analysis function inside GeoDjango.

Someone can made a simple example?


texts: How to access a table row in CARTO.js histogram dataview

texts:
How do I access a certain row in a CARTO source?

For instance, I have Chart.js reading the ticketcount column using this:

var source = new carto.source.SQL('select * from ticketspercounty_copy');
var histogram = new carto.dataview.Histogram(source, 'ticketcount',{bins: 46});
histogram.on('dataChanged', (data) => {
       //console.log(data)

         let chart = new Chart(histo, {
           type: "bar",
           data: {
             labels: ['Abbeville',
                      'Charleston',
                      'York',
                    ],
             datasets: [
               {
                 data: data.bins.map(element => element.freq)
               }
             ]
           }
         });
       });
      // Add the histogram to the client
      client.addDataview(histogram);

but this generates an object that does not resemble the data like so:

I've pointed it to the freq field for testing but ultimately, I would like it to just read the singular value in each of the 46 rows.

Here is a snapshot of the table:

All I want to do is create a histogram of the ticketcount per ticketcounty and that information is not in the bin array.


texts: How to acess and store ImageWMS layers information?

texts:
I want to click on the Map and store some information from a set of layers at the click position.

I thought that I could do it with this code:

map.on('singleclick', function(evt) {

    document.getElementById('info').innerHTML = '';

    var viewResolution = view.getResolution();

    var url = [];

    for(var i=0; i<layers.length; i++){
        (function(i){
            url[i] = layers[i].getSource().getGetFeatureInfoUrl(
                evt.coordinate, viewResolution, "EPSG:3857",{'INFO_FORMAT': 'text/html'}
                );
        })(i);
    }

    var v = [];

    for(var j=0; j<url.length;j++){
        (function(){
            var m = j;
            document.getElementById('info').innerHTML = '<iframe id="frame" seamless src="' + url[m] + '"></iframe>';

            document.getElementById('frame').onload = function(){

                v[m] = document.getElementById('frame').contentWindow.document.getElementById('tpl').innerHTML;

                if(m === url.length-1){
                    alert(v);
                }
            };
        })();
    }

});

The output from this function shows that v was modified just in the last iteration eg. v=[,,,,12].

It probably has something to do with the asynchronous nature of the onloadevent and I thougth that the looping clousure would deal with that.

I am new to OpenLayers.


texts: Unable to Merge Multiple .DXF Files and Output .DXF Using Merge Vector Layers

texts:
I'm using QGIS 3.8.3 with GRASS 7.6.1 to merge multiple .dxf files into a single .dxf file.

When I attempt to use the Merge Vector Layers via Vector>Data Management tools>Merge Vector Layers, I get the following error "Creation of field path failed (OGR error: DXF layer does not support arbitrary field creation, field 'path' not created.)" when I attempt
1. Add number of .dxf files as the input layer.
2. Output it as a .dxf file.

When I output it as a temporary file, it doesn't show anything even though it looks like the attribute table is populated with features.

I don't have any issues viewing each individual .dxf file when I import them into qGIS.

I also don't have any issues when I convert them via GRASS v.patch but I was hoping that Merge Vector Layers would be faster as it looks like v.patch imports each file via v.in.ogr and processes it as a GRASS vector layer, merges each of the GRASS vector layers via v.patch and before using v.out.ogr to output it.

I've been told that doing the merging these .dxf files via Global Mapper into a single .dxf file takes 15 to 20 mins but using QGIS's GRASS v.patch takes a bit over an hour so I was wondering whether there's a better way to do it in QGIS/GRASS.


texts: Zooming in to area using Copernicus Open Access Hub?

texts:
I can pan but I cannot zoom in when using the Copernicus Open Access Hub.
This makes it difficult to order data.

I need to know how to zoom in on an area.

I know how to zoom to a product.

Copernicus Open Access Hub.

https://scihub.copernicus.eu/userguide/GraphicalUserInterface

This interface. I see Pan but not zoom.


texts: Detecting invalid polygons (self-intersecting or touching)

texts:
An invalid polygon can be defined as one that intersects itself, or one where the lines make contact, ie, the two shapes at the top of this picture.

Given the XYPoints of a polygon, what logic (using C# preferably, but as long as I see logic, I can translate) could be implemented to detect these situations?

I've looked at the Bentley-Ottmann_algorithm, however this detects if lines intersect and not if they touch. Not to mention, it says in the wikipedia article that it's slow.


texts: Feature Geometry/shape changed event in arcObjects

texts:
I have a feature class(polyline). there are 5 features in the feature class. I want to get notification when user changed property any of feature of feature class. How it can be achieved ? I am using arcmap 10.0.
say, user got to arcmap edit mood. then make a straight line feature to curved feature. How can I get this notification that user changed its shape??

python and .net is my development technology.


texts: I have a date feild in my data and when I add it to an info window, it shows in zulu time (2017-11-06T00:00:00.000). Is there a way to format the cartoCSS and HTML so that only the date, and not the time, shows?


texts: Exporting raster data from QGIS in chosen projection

texts:
EDIT: it is seeming like some of the time my last method will lead to a correct projection, some of the time it will not.  Unsure what's happening there but at least trial and error will get me through what I'm doing now.  Any suggestions for a more robust method would be appreciated.

apologies if this has been answered somewhere, but I've not found a solution I could get working.

I have a DEM file that I need to change the projection of, and then export for use in another program (accepts gdal).  I am using QGIS.  I can import it to QGIS and change the projection no problem, but cannot figure out how to export the result for rasters (fine for vectors).

I have tried using the clipper to save a copy, however this changes the projection to WGS84, defeating the purpose.  I have set the project to my projection of choice AGD66, and in Settings/Options have set all options to this projection.  The clipper always converts to EPSG:4326 - WGS 84. I have tried both with and without on the fly reprojections.

I have tried to use the GRASS plugin and use r.out.gdal.gtiff and r.out.gdal, however output files are all blank.

Through the Sextante toolbox I have tried GDAL Conversion/Translate however this too puts everything in EPSG:4326 - WGS 84.  The option "Override the projection for the output file" with AGD66 selected led to the error

An error has occured while executing Python code:

Traceback (most recent call last):
  File "C:/Users/.../.qgis//python/plugins\sextante\gui\CrsSelectionDialog.py", line 47, in okPressed
    self.authid = self.selector.selectedAuthId()
AttributeError: 'QgsProjectionSelector' object has no attribute 'selectedAuthId'

Python version:
2.7.2 (default, Jun 12 2011, 15:08:59) [MSC v.1500 32 bit (Intel)]

QGIS version:
1.8.0-Lisboa Lisboa, f440b60

Python path: ['C:/Users/.../.qgis//python/plugins\sextante', 'C:/PROGRA~2/Quantum GIS Lisboa/apps/qgis/./python', 'C:/Users/.../.qgis//python', 'C:/Users/.../.qgis//python/plugins', 'C:/PROGRA~2/Quantum GIS Lisboa/apps/qgis/./python/plugins', 'C:\PROGRA~2\Quantum GIS Lisboa\bin\python27.zip', 'C:\PROGRA~2\Quantum GIS Lisboa\apps\Python27\DLLs', 'C:\PROGRA~2\Quantum GIS Lisboa\apps\Python27\lib', 'C:\PROGRA~2\Quantum GIS Lisboa\apps\Python27\lib\plat-win', 'C:\PROGRA~2\Quantum GIS Lisboa\apps\Python27\lib\lib-tk', 'C:\PROGRA~2\Quantum GIS Lisboa\apps\qgis\bin', 'C:\PROGRA~2\Quantum GIS Lisboa\apps\Python27', 'C:\PROGRA~2\Quantum GIS Lisboa\apps\Python27\lib\site-packages', 'C:\PROGRA~2\Quantum GIS Lisboa\apps\Python27\lib\site-packages\PIL', 'C:\PROGRA~2\Quantum GIS Lisboa\apps\Python27\lib\site-packages\win32', 'C:\PROGRA~2\Quantum GIS Lisboa\apps\Python27\lib\site-packages\win32\lib', 'C:\PROGRA~2\Quantum GIS Lisboa\apps\Python27\lib\site-packages\Pythonwin', 'C:\PROGRA~2\Quantum GIS Lisboa\apps\Python27\lib\site-packages\wx-2.8-msw-unicode', 'C:\PROGRA~2\Quantum GIS Lisboa\apps\qgis\python\plugins\fTools\tools']

I thought I had a solution using Raster/Converstion/Translate (convert format) as this runs and lets me set the Target SRS, however the output is still EPSG:4326 - WGS 84.  I do get the following error

Warning 6: A dataset opened by GDALOpenShared should have the same filename (C:\Users...\demh\w001001x.adf) and description (C:...\demh)

Any suggestions for how I can get raster data OUT of QGIS would be very much appreciated.  Or - if you can't get data out of QGIS, can you suggest an alternative free GIS system that works nicely with projections?  I have tried SAGA, however have not had much luck getting reprojections to work properly.

Thanks


texts: Is there a way to locally add links to an OSM data set

texts:
Here's what I'm trying to do:

I'm modeling some transportation logistics data
I want to test the effects of a "new highway" on existing travel demand
The new road will be created from an arbitrary point (usually not on the existing roadway) to any point on the road network (road line layer)
The new road should acquire the same fields / attributes as the road it connected too

Is this possible with OSM data?


texts: openstreetmap->download data... does not download data

texts:
I cannot get the OpenStreetMap -> download data... to actually download.

It acts like there is no problem, says the file is 0MB (it is actually 1KB) and it can be opened as a vector layer, however, the vector layer says there is no data.

I have ensured that it is in EPSG 3857 and the "on-the-fly CRS transformation" is enabled. The scale of the area is 1:500,000 as said by the QGIS website to test this.

The lecturer demonstrated it using 2.14.3 and the file he saved was about 17.5MB so around 17MB is the expected size for this.

Is there any way to work around this? I am using version 2.14.16 64bit


texts: Qgis exported SLD file is not valid on Geoserver

texts:
i'm trying to export my symbols from Qgis as an SLD file to use it as a style for the layer in geoserver here is an example of my SLD:

<?xml version="1.0" encoding="UTF-8"?>
<StyledLayerDescriptor xmlns="http://www.opengis.net/sld" xmlns:ogc="http://www.opengis.net/ogc" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" version="1.1.0" xmlns:xlink="http://www.w3.org/1999/xlink" units="mm" xsi:schemaLocation="http://www.opengis.net/sld http://schemas.opengis.net/sld/1.1.0/StyledLayerDescriptor.xsd" xmlns:se="http://www.opengis.net/se">
  <NamedLayer>
    <se:Name>stclines_streets</se:Name>
    <UserStyle>
      <se:Name>stclines_streets</se:Name>
      <se:FeatureTypeStyle>
        <se:Rule>
          <se:Name>Single symbol</se:Name>
          <se:LineSymbolizer>
            <se:Stroke>
              <se:SvgParameter name="stroke">#d575ef</se:SvgParameter>
              <se:SvgParameter name="stroke-width">0.26000000000000001</se:SvgParameter>
              <se:SvgParameter name="stroke-linejoin">bevel</se:SvgParameter>
              <se:SvgParameter name="stroke-linecap">square</se:SvgParameter>
            </se:Stroke>
          </se:LineSymbolizer>
          <se:LineSymbolizer>
            <se:Stroke>
              <se:GraphicStroke>
                <se:Graphic>
                  <se:Mark>
                    <se:OnlineResource xlink:type="simple" xlink:href="ttf://Dingbats"/>
                    <se:Format>ttf</se:Format>
                    <se:MarkIndex>65</se:MarkIndex>
                    <se:Fill>
                      <se:SvgParameter name="fill">#ff0000</se:SvgParameter>
                    </se:Fill>
                  </se:Mark>
                  <se:Size>2</se:Size>
                </se:Graphic>
                <se:Gap>
                  <ogc:Filter xmlns:ogc="http://www.opengis.net/ogc">
                    <ogc:Literal>3</ogc:Literal>
                  </ogc:Filter>
                </se:Gap>
              </se:GraphicStroke>
            </se:Stroke>
          </se:LineSymbolizer>
        </se:Rule>
      </se:FeatureTypeStyle>
    </UserStyle>
  </NamedLayer>
</StyledLayerDescriptor>

this doesn't work for any complex multi symbols like lines with graphical hashes and so. even when i try using Qgis plugin " geoserver explorer " to publish it the published SLD file ignore that multi symbols or font symbol at all however when i go to SLD cookbook and take some code samples for some complex symbols but in v 1.0.0 it kinda works, but i need the advanced options in GUI that QGIS has.

i read in that thread that QGIS produce 1.1.0 while geoserver validate 1.0.0 SLDs but i'm looking for an update for that issue or is it still ?thread link


texts: Managing non spatial tables with geoserver

texts: I am developing an application involving PostGis+Geoserver+QGIS, the application is basically a QGIS plugin, as client of Geoserver. We need to manage forms related to some spatial data (more than one form per spatial item). I would like to know if there is any way using Geoserver to manage the editing of these forms stored in a non spatial table. In summary, Is Geoserver capable  of managing non spatial tables (including read and write)


texts: Create primary and foreign key on Spatialite

texts:
I try to create a database on spatialite (for various reason, I can't use PostGis). I need to link a table to another one using a foreign key.
The problem with spatialite is that it is not possible to define a primary key afterwards, so I create a new table (tMAJ_TBE_REB_15₂) with the Primary key « C_SGE_15 » (it is field with no null value, and unique values) and then I try to copy the content from the other one (which have got the default Primary Key PUIK).

When I start this command, I get this error message.

Une erreur est apparue lors de l'exécution du code Python :
Traceback (most recent call last):
  File "C:/Documents and Settings/XXXXXXX/.qgis2/python/plugins\QspatiaLite\QspatiaLiteApp.py", line 178, in runQuery
    header,data=self.db.executeQuery(query)
  File "C:/Documents and Settings/XXXXXXX/.qgis2/python/plugins\QspatiaLite\Classes\Database.py", line 51, in executeQuery
    cursor.execute(query,params)
IntegrityError: column C_SGE_15 is not unique

Version de Python :
2.7.4 (default, Apr  6 2013, 19:54:46) [MSC v.1500 32 bit (Intel)]
Version de QGIS :
2.8.2-Wien Wien, 1b929ef
Chemin vers Python : ['C:/Documents and Settings/XXXXXXX/.qgis2/python/plugins\\processing', 'C:/PROGRA~1/QGISWI~1/apps/qgis/./python', u'C:/Documents and Settings/XXXXXXX/.qgis2/python', u'C:/Documents and Settings/XXXXXXX/.qgis2/python/plugins', 'C:/PROGRA~1/QGISWI~1/apps/qgis/./python/plugins', 'C:\\PROGRA~1\\QGISWI~1\\bin\\python27.zip', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\DLLs', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\plat-win', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\lib-tk', 'C:\\PROGRA~1\\QGISWI~1\\bin', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\site-packages', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\site-packages\\PIL', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\site-packages\\jinja2-2.7.2-py2.7.egg', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\site-packages\\markupsafe-0.23-py2.7-win32.egg', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\site-packages\\python_dateutil-2.2-py2.7.egg', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\site-packages\\pytz-2014.2-py2.7.egg', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\site-packages\\win32', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\site-packages\\win32\\lib', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\site-packages\\Pythonwin', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\site-packages\\Shapely-1.2.18-py2.7-win32.egg', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\site-packages\\six-1.6.1-py2.7.egg', 'C:\\PROGRA~1\\QGISWI~1\\apps\\Python27\\lib\\site-packages\\wx-2.8-msw-unicode', u'C:/Documents and Settings/XXXXXXX/.qgis2//python', 'C:\\Documents and Settings\\XXXXXXX\\.qgis2\\python\\plugins\\DigitizingTools\\tools', 'C:\\Documents and Settings\\XXXXXXX\\.qgis2\\python\\plugins\\mmqgis/forms', 'C:\\PROGRA~1\\QGISWI~1\\apps\\qgis\\python\\plugins\\fTools\\tools', 'C:/XXX/MAJ_TBE_REB']

To resume my idea, I would like to have a geographic table link using the field C_SGE_15 to another table  tMAJ_TBE_REB_15₂ which contain description, etc. The goal after is to create a form easy to use.
I do this on Qspatialite and Qgis 2.8.2.

CREATE TABLE 'tMAJ_TBE_REB_15_2' (
'C_SGE_09' TEXT,
'MODIF_15' INTEGER,
'TypModif15' TEXT,
'C_SGE_15' TEXT NOT NULL PRIMARY KEY,
'NOM_09' TEXT,
'NOM_15' TEXT,
'TBE_09'  INTEGER,
'TBE_15'  INTEGER,
'REB_09' INTEGER,
'REB_15' INTEGER,
'Espec_pres' TEXT,
'Obs_2009' TEXT,
'Obs_2015' TEXT)

INSERT INTO   "tMAJ_TBE_REB_15_1" (
'C_SGE_09',
'MODIF_15',
'TypModif15',
'C_SGE_15',
'NOM_09',
'NOM_15',
'TBE_09',
'TBE_15',
'REB_09',
'REB_15',
'Espec_pres',
'Obs_2009',
'Obs_2015')
SELECT 'C_SGE_09',
'MODIF_15',
'TypModif15',
'C_SGE_15',
'NOM_09',
'NOM_15',
'TBE_09',
'TBE_15',
'REB_09',
'REB_15',
'Espec_pres',
'Obs_2009',
'Obs_2015',
FROM  "Modifs_MAJ_TBE_REB_15"

Does anybody see from where the error come from?


texts:
I would like to construct a map tool application that allows only to add some point features.

How do I complete the following code to activate snap mode while editing a layer with PyQGIS 3?

from qgis.core import *
from qgis.gui import *
from PyQt5.QtGui import *
from PyQt5.QtCore import *

class AddPointTool(QgsMapTool):
   def __init__(self, mapCanvas, layer):
       QgsMapTool.__init__(self, mapCanvas)
       self.layer=layer
       self.setCursor(Qt.CrossCursor)

   def canvasReleaseEvent(self, event):
       point=self.toLayerCoordinates(self.layer, event.pos())
       fields=self.layer.dataProvider().fields()
       feature=QgsFeature()
       feature.setFields(fields)
       feature.setGeometry(QgsGeometry.fromPointXY(point))

       self.layer.addFeature(feature)
       self.layer.updateExtents()
       self.layer.triggerRepaint()

canvas=iface.mapCanvas()
layer=iface.activeLayer()

layer.startEditing()
layer.triggerRepaint()
tool=AddPointTool(canvas,layer)
canvas.setMapTool(tool)


texts: QGIS Select Single Feature not correct in attribute table

texts:
For some reason, QGIS is not displaying the correct row of my data when I click on a feature with "Select Single Feature" or even the other select tools.

I have a point database linked from mysql.  There is a column named ASR with a unique number.  When I plot all of the points, they are correct - everything is where it should be.  When I label the points they are all correct - ASRs correspond to points in the correct locations.  When I open the attribute table, select a line, and then zoom to that feature, it works fine.  Correct point, correct labeling, all good.

HOWEVER, if I just display the map of all the features and try clicking on one with the 'Show selected only' checked, it doesn't actually select the right point, and moreover it shows the wrong row.  In most cases if I select a rectangle of 20 or 30 points, it actually highlights 20 or 30 OTHER points, and displays various results in the attribute table.

I'm totally confused.  Almost everything is working right and its clear that the record data and point data are correct and in agreement, but selecting features by mouse is totally hosed.

Any ideas?


texts: How to use PostgreSQL/PostGIS with ArcGIS for Desktop?

texts:
Has anyone used PostgreSQL with PostGIS successfully?

I am a semi-experienced ArcGIS user and want to build out a database and connect it to ArcCatalog.

I've scoured the internet for resources but nobody has really given a good walkthrough.


texts:
When I tried to add a wms server in ArcCatalog, the following message came up: "The URL you passed in uses a protocol that is not supported." . The WMS Server is the following gis.ktimanet.gr/wms/wmsopen/wmsserver.aspx?. It's the Greek Cadastre WMS Server. I use ArcGis 10.2. Exactly the same message came up in QGIS.

Any ideas?


texts: Using PyQGIS processing runalg clipvectorsbypolygon

texts:
I try to use clip vector by polygon in pyqgis.

My code (minimal working example)

import sys
from qgis.core import *
import qgis.utils

QgsApplication.setPrefixPath("/usr",True)
qgs = QgsApplication([],True)
qgs.initQgis()

# load base layer
blayer = QgsVectorLayer("lokale_orte.shp", "lokale_orte", "ogr")
if not blayer.isValid():
        print "Layer failed to load!"
else:
        print "success:  load layer"

# load clip layer
clayer = QgsVectorLayer("kleine_kreise.shp", "kleine_kreise", "ogr")
if not clayer.isValid():
        print "Layer failed to load!"
else:
        print "success:  load layer"

# Prepare processing framework
sys.path.append('/usr/share/qgis/python/plugins') # Dir where Processing is located
from processing.core.Processing import Processing
Processing.initialize()
from processing.tools import *

# run clip
general.runalg('gdalogr:clipvectorsbypolygon', 'lokale_orte.shp',                 'kleine_kreise.shp', None, None)

The complete output and error message is:

success:  load layer
success:  load layer

Error in sys.excepthook:

Traceback (most recent call last):
  File "/usr/lib/python2.7/dist-packages/qgis/utils.py", line 196, in qgis_excepthook
    showException(type, value, tb, None, messagebar=True)

  File "/usr/lib/python2.7/dist-packages/qgis/utils.py", line 107, in showException
    open_stack_dialog(type, value, tb, msg)

  File "/usr/lib/python2.7/dist-packages/qgis/utils.py", line 142, in open_stack_dialog
    iface.messageBar().popWidget()

AttributeError: 'NoneType' object has no attribute 'messageBar'

Original exception was:

Traceback (most recent call last):
  File "a.py", line 34, in <module>
    general.runalg('gdalogr:clipvectorsbypolygon', 'lokale_orte.shp',
'kleine_kreise.shp', None, None)

  File "/usr/share/qgis/python/plugins/processing/tools/general.py", line 75, in runalg
    alg = Processing.runAlgorithm(algOrName, None, *args, **kwargs)

  File "/usr/share/qgis/python/plugins/processing
/core/Processing.py", line 300, in runAlgorithm
    ret = runalg(alg, progress)

  File "/usr/share/qgis/python/plugins/processing
/gui/AlgorithmExecutor.py", line 52, in runalg
    progress.error(e.msg)

AttributeError: 'NoneType' object has no attribute 'error'

I guess my problem is the runalg call, according to http://docs.qgis.org/2.6/de/docs/user_manual/processing_algs/gdalogr/ogr_geoprocessing/clipvectorsbypolygon.html

I think I can write None for output layer to force QGIS to write the output in an memory layer. But what about options? What to write? Where is a proper documentation?


texts:
I have seen some examples of this but they are not working for this situation.

What I need to be able to do is when I click on a marker get the index of that marker to write to the console. The reason I need this info is so that I can then go and later update the popup for this marker.

I have a JSBIN of this at http://jsbin.com/hepeli/4/edit?js,console,output

When you create a marker it makes a popup that displays an index, but this is actually just the length of an array that the markers are being pushed into.

When you click a marker with the delete tool it does write to console the index but I can't get it to do this when I just click on the marker.


texts: Change or remove default style in GeoServer tiles layers

texts:
I have publish vector tiles layers in GeoServer. It's displaying with default style as shown in image

and

.

I want to change style for vector tiles in GeoServer.


texts:
Using ArcGIS Desktop v10 i am editing a "tank" layer.

It is selected as the only editable layer

It has three features all of which have circle as the diagrammatic property.

I use the circle tool and construct a new circle.

It draws nicely.

Then I save edit

Then I stop editing the layer.

Visually it turns into an oval!

Thoughts?


texts:
I'm totally a newbie with geo tools. I'm using Sharpmap as I'm working with C#.

I have a map where I can put an item on a UTM Location with E,N Zone and number. So the item appears on the desired position.

The problem comes when I want to move the item to a new position. So, for example, I want to increase easting by 2000 meters.

How can I calculate if my item has crossed to a new zone? The easting ranges differ from equator to the poles, am I right? Is there a way to calculate easting ranges of a certain zone?

What is the best way to work with location points in different zones -- for example, to calculate distance between them?


texts: Arcpy replaceDataSource for different cases

texts:
I have a layer in my gdb file called lines.
I also have a *.lyr file lines-N with symbology for the lines layer.
I need lines-N to have the same data source as lines, using arcpy.replaceDataSource.

My code:

kyl=arcpy.GetParameterAsText(0)
desc=arcpy.Describe ("lines")
path = desc.path
for lyr in arcpy.mapping.ListLayers(mxd, "",df):
      if lyr.name == "lines-N":
          if lyr.supports ("DATASOURCE"):
              lyr.replaceDataSource (path,"FILEGDB_WORKSPACE","",True)

produces this error:

return convertArcObjectToPythonObject(self._arc_object.replaceDataSource(*gp_fixargs((workspace_path, workspace_type, dataset_name, validate), True)))
ValueError: Layer: Unexpected error

What am I doing wrong? I'm trying to avoid manually adding the file path like so:

lyr.replaceDataSource (r "D:\Folder\file.gdb","FILEGDB_WORKSPACE","",True)

UPDATE:
I'm running it from the ArcMap Python window to test it first, and then adding it to a script that performs similar operations on other layers.
ArcGIS version is 10.2.2


texts: How to select featureMember by gml:id in WFS 2.0

texts:
In WFS 1.1.0 I could filter like this:

<ogc:Filter>
<ogc:GmlObjectId gml:id="5"/>
</ogc:Filter>

Is there an equivalent in WFS 2.0?

I've tried:

<fes:Filter>
    <fes:PropertyIsEqualTo>
        <fes:ValueReference>gml:id</fes:ValueReference>
        <fes:Literal>B6_Aardgastransportleidingen.2</fes:Literal>
    </fes:PropertyIsEqualTo>
</fes:Filter>

But this returns an exception

<ows:ExceptionText>Illegal property name: gml:id for feature type B64_nutsvoorzieningen:B6_Aardgastransportleidingen</ows:ExceptionText>


texts: I want to download Portal for ArcGIS, but I am not finding the link to download it. Where I can download Portal for ARCGIS?


texts:
I need a query that return the area of two polygons that have been split by a line.

This is what I did so far:

SELECT ST_Area(ST_Split(surface.geom, limite_a.geom)) from surface,
limite_a   where surface.id = 945 and limite_a.id=123


texts:
Trying to center the map at a certain coordinates but map only centers at 0,0 instead. Coordinates are in decimal degrees. Code below.

function foo(data) {

        var lon = data.address.longitude;
        var lat = data.address.latitude;
        var coord = [lat,lon];

        var point = new esri.geometry.Point(coord)

        console.log(point)

        map.centerAt(point);
    };


texts: Temporally reducing image bands using Google Earth Engine

texts:
I am trying to temporally reduce an image from 365 bands (daily) to 12 bands (monthly).
Context: I imported an image collection with Soil Moisture (SM) that had 8 values per day (every three hours). I first reduced it to average daily values to calculate SM anomalies with this formula: SManom = (SM-SMmean)/ SMstdev
Then I created a new variable (irrigation events) that is 1 when there is an irrigation event and 0 if not. So now I have an image with one band for each day, and I want to reduce it to accumulated monthly values. This is the code that I have for now, below it I will post what I tried.
(This is actually not what I want to achieve (accumulate irrigation events) but it is a working example, otherwise the entire code would be too long. )
var dataset = ee.ImageCollection("NASA/GLDAS/V021/NOAH/G025/T3H"),
    geometry = ee.Geometry.Polygon(
        [[[21.742985018165424, 52.321076580783036],
          [21.742985018165424, 44.02344004522736],
          [40.28790689316542, 44.02344004522736],
          [40.28790689316542, 52.321076580783036]]], null, false);

//Import precipitation and soil moisture for a year.
var pcp = dataset.select('Rainf_tavg').filterDate('2019-01-01', '2020-01-01');
var sm = dataset.select('SoilMoi0_10cm_inst').filterDate('2019-01-01', '2020-01-01');

//Reduce SM to daily values and compute mean, std dev, and anomalies.
var sm_daily = ee.ImageCollection.fromImages(days.map(function(d) {
      var filtered = sm.filter(ee.Filter.calendarRange({
        start: d,
        field: 'day_of_year'
      }));
      return filtered.mean().set('day', d);
    }));

var sm_mean = sm_daily.mean();
var sm_stdev = sm_daily.reduce(ee.Reducer.stdDev());
var sm_img = sm_daily.toBands();
var sm_anom = (sm_img.subtract(sm_mean)).divide(sm_stdev);

//Reducing precipitation to daily values.
var pcp_daily = ee.ImageCollection.fromImages(days.map(function(d) {
   var filtered = pcp.filter(ee.Filter.calendarRange({
     start: d,
     field: 'day_of_year'
   }));
   return filtered.mean().multiply(86400).set('day', d);
}));

var pcp_img = pcp_daily.toBands();

//Computing irrigation events.
var irrigation_event = pcp_img.lte(1).and(sm_anom.gte(1));

I tried to adapt the part of the code that I used to reduce to daily values but the fact that now it is an image rather than an image collection is problematic. The best solution I have found is to use slice on the image in order to return images from specifically selected bands and creating an image collection from those images. I think this could work but I'm having trouble selecting bands in an iterative manner. I also tried creating an array with lists that contain sequences of the days of each month but that didn't work since the lists have to be the same length.
var months = ee.List.sequence(1, 12);
//Here I want to map this function to the variable months but I can't seem to figure out how to select
//the right bands in every iteration.
var irr_jan_acc = ee.ImageCollection.fromImages(months.map(function(g) {
  var filtered = irrigation_event.slice(g);
  return filtered.reduce(ee.Reducer.sum());
}));

I am quite new to GEE and especially to iterations.


texts: How to preprocess Bing Singleband Satellite image for eliminating different exposure tiles problem

texts:
I have used around 200 sqkm high resolution single band colour Bing image for object based image analysis on ecognition developer. But due to different tile based exposure made it worse to classify. the segmentation goes different for two areas. Can anyone suggest how to pre process the image before segmentation. I have added an example of exposure problem through image. []

I have used ecognition developer provided equalizing tools like Gamma Correction, Histogram, linear equalizer and standard deviation. But nothing works.


texts: Create a GEOSGeometry from a FeatureCollection in GeoDango

texts:
Is it possible create a GEOSGeometry from a FeatureCollection(GeoJSON)?
For instance: GEOSGeometry(json.dumps(fc)) where:

fc = { "type": "FeatureCollection",
"features": [
  { "type": "Feature",
    "geometry": {"type": "Point", "coordinates": [102.0, 0.5]},
    "properties": {"prop0": "value0"}
    },
  { "type": "Feature",
    "geometry": {
      "type": "LineString",
      "coordinates": [
        [102.0, 0.0], [103.0, 1.0], [104.0, 0.0], [105.0, 1.0]
        ]
      },
    "properties": {
      "prop0": "value0",
      "prop1": 0.0
      }
    },
  { "type": "Feature",
     "geometry": {
       "type": "Polygon",
       "coordinates": [
         [ [100.0, 0.0], [101.0, 0.0], [101.0, 1.0],
           [100.0, 1.0], [100.0, 0.0] ]
         ]
     },
     "properties": {
       "prop0": "value0",
       "prop1": {"this": "that"}
       }
     }
   ]
 }


texts:
I am using Mapbox GL JS and querying a source on a common ID.

var tempId= 55;

var relatedFeatures=map.querySourceFeatures('featureSource', {
    sourceLayer: "sourceLayerId",
    filter: ['in', "FID",tempId]
 });

This returns an array of geojson features. If I pass this directly to turf union like below it throws an error

turf.union(relatedFeatures)

I can iterate over the array like below and union the features together but this seems unnecessary. Any suggestion for doing this more efficiently?

var featuresJoined=relatedFeatures[1]
  if(relatedFeatures.length>1){
    $.each(relatedFeatures,function(i,item){
      if(i>1){
        featuresJoined=turf.union(featuresJoined,item)
      }
    })
  }


texts: Modifying width and precision of database using QGIS?

texts:
I have a database in txt format made of several columns of real numbers. Two column are integer, 17 are real (rounded to the first decimal) and 17 are real rounded to the second decimal.

If I import it through a txt format QGIS identified the integer columns properly, but also assign a width of 24 and a field precision of 15.

How can I limit it?


texts: How do I get a polyline Displayed in Mapbox Visual for Power BI?

texts: All of the Mapbox Examples involve point layers where you add a csv and drop lat/lngs into power bi.  I have also tried creating a custom map in studio and adding and styling my polyline.  However when I load the custom map in power bi it shows the correct basemap but my polyline is not there. Are there extra steps in enabling this layer in my custom map?


texts: How to edit layers which are using via OSM Bright in Tilemill

texts:
How could i edit a map which is edited via OSM Bright in TileMill? Layers are in TM added via PostGis. OSM Bright it's simply built according to informations here. I need delete some roads also delete labels. I have deleted and edited .shp file in ArcGisfrom project folder of TileMill but after saving of the .shp file and starting TM there are no any changes. Thanks for help.

I think it could be resolve by editing of .shp but it's still without changes. Do you think that changes could be shown after "re-postGis" to pgsql.


texts:
Within CartoDB my CartoCSS looks something like

#table_copy_1 [ lease_rate <= 55.15] {
   marker-width: 35.0; marker-fill: #081B47
}

Can I add another variable (size) so that the display filters by both lease_rate and size? My idea is that it would look like this:

#table_copy_1 [ lease_rate <= 55.15, size <= 40000] {
   marker-width: 35.0; marker-fill: #081B47
}


texts: I have some MODIS imagery that I want to display in true colour, however when I add the data it gives me 13 subdatsets to choose from. I've never worked with this before and I'm not sure which subdatasets get which results. It doesn't work with the normal band combinations that I'm used to. My goal is to look at the change in fall foliage over a five year period. Any information on how to work with this data?


texts: Labelling the GeoTIFF canvas symbology in QGIS3

texts:
I used the Interpolation (Nearest neighbor) feature.
After watching this tutorial:
https://www.youtube.com/watch?v=By-K5C7jDQc&ab_channel=HansvanderKwast
I could set the color rendering for my values.
I also put labels next to them, because I want them visible next to the color band, as we used to see by contour lines and so on.
Unfortunately, they don't come out, as you can see below:

I though that the thread below will be helpful:
How to add text on top of symbols for vector layers in QGIS
But I can't see the label option in the Properties panel at all.
I want something like you to see below:

Is it possible with GeoTIFF band rendering?


texts: Is there a way to control the legends of cartographic representations with Python?

texts: I have a ArcGIS feature class with 3 cartographic representations.  Two of them represent the features with line offsets.  These representations are determined based on multiple fields.  One field determines the color of the line, the other field determines the offset.  In the legend this means that each color is repeated once for each offset.  This can be controlled using Layer Properties>Symbology>Legend Options to make it so that each line color is only shown once in the legend.  Is there a way to control this setting with Python?


texts: Draw a polygon in OpenLayers with coordinates from a variable?

texts:
I have the polygon coordinates defined as a variable. This draws nothing when I use the variable, when the same coordinates are used as such in the polygonFeatureT definition, the polygon is drawn correctly.

With point coordinates as variable this works, provided that the quotes are not used.

What am I doing wrong here?

The code is below:

var polyCoords = '[3143090.603086447, 9928281.393790578], [3283734.7351311715, 9928892.890016861], [3181003.3691158947, 9849398.380600277], [3143090.603086447, 9928281.393790578]';

// this Point var works only when no quotes are used in definiton!
var pointCoords = [3229617.319105267, 9916160.39109719];

var polygonFeatureT = new ol.Feature(new ol.geom.Polygon([[polyCoords]]));
var pointFeatureK = new ol.Feature(new ol.geom.Point(pointCoords));

var map = new ol.Map({

        target: 'map',
        layers: [
          new ol.layer.Tile({
              source: new ol.source.OSM()
          }),
               new ol.layer.Vector({
                   source: new ol.source.Vector({
                       features: [polygonFeatureT, pointFeatureK]
                   }),
...


texts:
Some of the ones I'm familiar with (in no particular order):

 Google Maps (http://code.google.com/apis/maps/index.html)
 Bing Maps (http://www.microsoft.com/maps/developers/)
 Yahoo Maps (http://developer.yahoo.com/maps/)
 MapQuest (http://developer.mapquest.com/)
 OpenLayers (http://openlayers.org/)
 WorldWind Java (http://worldwind.arc.nasa.gov/java/)
 OpenMap (http://www.openmap.org/)

Are there others that I'm missing?


texts: Measure length in specific units like 'aln' arcmap

texts: I have an old map in aln, 1 aln=0.5934 meters. I would like to draw lines using directly this measure unit without having to convert in meters each time to draw the lines. Is there any way to customize the unit type in arcmap 10.1?


texts: Why the quality of image is changed as it is added to a mosaic

texts:
I have added an image to a mosaic dataset but I have noticed that the quality of the image has been changed automatically.
How to keep the original quality of the image after adding it to the mosaic dataset?


texts: Failure to Add Layer from a WFS Provider QGIS 2.18.10

texts:
I just installed QGIS 2.18.10. When I tried to add a layer from a WFS provider, it failed citing Bad Request as the reason. Exact message:

Falló la solicitud de red GML Getfeature con el error: Error downloading http://ide.igvsb.gob.ve/geoserver/escala250000/wfs/?SERVICE=WFS&VERSION=1.0.0&REQUEST=GetFeature&TYPENAME=escala250000:ba_4189_orla_de_venezuela_250000_20161116&SRSNAME=EPSG:4189 - server replied: Forbidden


texts: Best practices for multiple admins on ArcGIS Server and IIS

texts:
What are best practices for having multiple Administrators of ArcGIS Server?  I ask in terms of one admin starting and stopping services while another admin is deleting/creating/stopping different services and other similar situations?

If best, line out scenarios to avoid and best practices.  Thanks.


texts: How to analyse point data as covariate for areal data

texts:
I have two data sets:

a) a set of polygons with areal data variables.
b) point measures of various co-variates of interest (note I do not necessarily have a point for each polygon)

How do I correctly analyse this data set for correlations between an areal variable (dependent variable) and point measures (independent variable)?

Unsure what techniques I need to use. Working in R.

Edit on 7 July 2014:

What I'm really wondering is - if you have areal data i.e. some variable that depends on an area and not points - so something related to population density for example - and you want to analyze that versus a point measurement - rainfall at specific measurement points for example. What is the correct/best way to do this ?
I've thought about it since and I see two approaches:
1) Krige the point data and get an average value for each polygon (I presume by picking multiple sample points within the polygon and getting an average) - then build polygon based models using the kriged estimates as covariates.
2) Treat the areal data as point estimates of a field taking the centroid of each polygon as the location for each point estimate, then build a geostatistical model using all data as point estimates.
Both these approaches strike me as "incorrect"  - what I'm asking really is which is least "incorrect" ?


texts:
I am attempting to learn ArcObjects using a tutorial that uses VB.NET and I am trying to convert it to C#.  I have Visual C# 2008 Express Edition.  This is an ArcMap Add-In project that creates a button that displays a message box with the name of the top most layer in the data frame.  Here is the VB.NET code:

Dim pMxDoc As IMxDocument
pMxDoc = My.ArcMap.Application.Document
Dim pMap As IMap
pMap = pMxDoc.FocusMap
Dim pFLayer As IFeatureLayer2
pFLayer = pMap.Layer(0)

Dim strDisplayField As String
strDisplayField = pFLayer.DisplayField

MsgBox(strDisplayField, vbOKOnly, "Test")

My conversion in C# so far:

using ESRI.ArcGIS.ArcMapUI:
using ESRI.ArcGIS.Carto;

IMxDocument pMxDoc = ArcMap.Application.Document as IMxDocument;
IMap pMap = pMxDoc.FocusMap;
IFeatureLayer2 pFLayer = pMap.

There is no 'Layer' property.  Only 'LayerCount' to choose from.
I even tried this:

IFeatureLayer2 pFLayer = ArcMap.Application.Document as IFeatureLayer2;
pFLayer = pMap.

this didn't work either.
I have looked all over the ArcMap Object Diagram and I can't find the FeatureLayer2 class anywhere.  According to ArcObjects SDK for .NET there is a Layer member for the IMap interface, but I don't have access to that.  Any help would be grealty appreciated.


texts: QGIS geometry generator - chainage nodes and labels along line at ellipsoidical distance

texts:
All points along line solution proposed in other question or plugin or processing toolbox use projected distance and not ellipsoidal.

I work with very long lines and is quite important that the points along line are based on ellipsoidal measures and not projected/plane.

Any idea?


texts: OpenLayers - Can't remove zoom control buttons, other controls simply lay on top

texts:
I have a simple OpenLayers script with a Google map as the base layer.  By default the map has a + and - zoom control button, which I can't find how to get rid of.

I've added an empty controls array to the map options, but it still appears.  If I add any other controls, for example PanZoomBar, it just sits ontop of the + and - buttons, semi-obscuring them. The effect is the same whether I add PanZoomBar in the controls array or with map.addControl.

Here is what it looks like by default:

And with added PanZoomBar:

And here's my code:

<!DOCTYPE html>

    <title>Map test</title>

Map Test

Can anybody help me get rid of the + and - buttons?


texts: Measuring area of raster within polygons using QGIS

texts:
I am trying to use QGIS3.10 to analyze the area of forest within regions. It looks like this

where the polygons belongs to WFS layer and Rasters belongs to WMS layer.
The legends fo Raster suggests that the raster contains different catagories, looks like this

I would like to find the percentage of green areas in polygons based on percentage coverage info provided in Raster. I have tried converting raster to polygons (with no polygons outcome), raster zonal analysis (with error message) and Extract layer content (no desired outcome).

Can anyone suggest me a way to extract the data from it?

I have knowledge in Python, but I have no idea on how to handle this kind of data.

Dataset used in this case:

http://geo.stat.fi/geoserver/postialue/wfs
-> Paavo-postinumeroalueet 2019
https://image.discomap.eea.europa.eu/arcgis/services/GioLandPublic/HRL_TreeCoverDensity_2015/MapServer/WMSServer
-> Tree Cover Density 2015


texts: Excluding self-intersections when performing a linestring difference in shapely

texts:
I have a (complex) LineString (with 1 self-intersection):
ls_self_intersection = LineString([(1, 2), (4, 2), (4, 3), (3, 3), (3, 1)])

and I have a (simple) LineString (with no self-intersection):
ls = LineString([(2, 1), (2, 3)])

when I do a difference:
ls_diff = ls_self_intersection.difference(ls)

the result is a MultiLineString
mls = mapping(ls_diff)

with the following content
{
    'type': 'MultiLineString',
    'coordinates': (
        ((1.0, 2.0), (2.0, 2.0)),
        ((2.0, 2.0), (3.0, 2.0)),
        ((3.0, 2.0), (4.0, 2.0), (4.0, 3.0), (3.0, 3.0), (3.0, 2.0)),
        ((3.0, 2.0), (3.0, 1.0))
    )
}

I was expecting this result (excluding the self-intersection):
{
    'type': 'MultiLineString',
    'coordinates': (
        ((1.0, 2.0), (2.0, 2.0)),
        ((2.0, 2.0), (4.0, 2.0), (4.0, 3.0), (3.0, 3.0), (3.0, 1.0)))
    )
}

How can I exclude the self-intersection when performing the difference operation?


texts: Are features of new labelling tool available to control via Python API?

texts:
I am writing a plugin using Python and currently using methods that set the label similar to those seen under Layer Properties. Some of the features under the new labelling engine available on the Toolbar would be useful, such as (not) "Label every part of a multi-part feature" and avoiding label collisions. Also the look of the labels generated using the latter method are usually better. Is there any way of setting these programatically in Python?

I understand the Toolbar labelling system will supersede the version under Layer Properties, so if they are not available in the API, will they be in the future?


texts: Where to find QGIS geoprocessing results?

texts:
Is there a way to access in QGIS 2.6.1 something like the Geoprocessing Results window found in ArcGIS?

The History and Log window in QGIS appears to show the geoprocessing results from SAGA, OTB and other algorithms but nothing from the Raster Calculator, which is what I am really after. There does not appear to be anything useful in the metadata, or in the layer properties.

Is there a way to find out what raster magic I performed, or should I resign myself playing around with the raster calculator and eventually figure out what I did in the past?


texts: Where should I start advanced GIS learning

texts:
I am relatively familiar with most of the GIS stuff. I develop GIS applications as well with open source offerings. Still I find some lapses in my theoretical knowledge in GIS. What resources can you propose which might be good for an advanced GIS learner with a touch of software development?

Edited

One area I think I am weak in is coordinate systems and the theory behind coordinate systems.

Update:
"I will mostly appreciate if you could point out one or two books as well"


texts: Clip raster function leaves overlap when clipping to polygon using ArcObjects?

texts:
I'm testing a portion of a Desktop Add-In that uses a raster clip function to quickly clip a topo raster to a quadrangle polygon (Data management was far too slow as I will be clipping many rasters). Everything works fine except that the output raster has a bit of an overlap from the polygon extent. I have noticed in a PCS the overlap is uneven and in a GCS it is uniform.
When I use the Image Analysis window to clip the raster by selecting the quad polygon the output is perfectly flush and I'm wondering what my program is doing differently.
ArcGIS 10.1 SP1, C#.NET/4.0
    protected override void OnClick()
    {
        mxdoc = ArcMap.Application.Document as IMxDocument;

        IWorkspaceFactory rwsf = new RasterWorkspaceFactoryClass();
        IWorkspaceFactory fwsf = new FileGDBWorkspaceFactoryClass();
        IRasterWorkspace rws = (IRasterWorkspace)rwsf.OpenFromFile(topoFolderPath, ArcMap.Application.hWnd);
        IRasterDataset rds = rws.OpenRasterDataset(topoBaseName);

        IFeatureLayer2 quadLayer2 = new FeatureLayerClass();
        IFeatureWorkspace fws = (IFeatureWorkspace)fwsf.OpenFromFile(quadFolderPath, ArcMap.Application.hWnd);
        quadLayer2.FeatureClass = fws.OpenFeatureClass(quadBaseName);
        IFeatureLayer quadLayer = (IFeatureLayer)quadLayer2;

        IClipFunctionArguments rasterFuncArgs = new ClipFunctionArguments() as IClipFunctionArguments;
        rasterFuncArgs.Raster = rds;
        rasterFuncArgs.ClippingType = esriRasterClippingType.esriRasterClippingOutside;
        rasterFuncArgs.ClippingGeometry = quadLayer.AreaOfInterest;

        IRasterFunction clipFunction = new ClipFunctionClass();
        IFunctionRasterDataset funcRasterDataset = new FunctionRasterDataset();
        IFunctionRasterDatasetName funcRasterDatasetName = new FunctionRasterDatasetName() as IFunctionRasterDatasetName;
        funcRasterDatasetName.FullName = @"C:\Users\dcanady\Desktop\clip.afr";
        funcRasterDataset.FullName = (IName)funcRasterDatasetName;
        funcRasterDataset.Init(clipFunction, rasterFuncArgs);

        rds = funcRasterDataset as IRasterDataset;
        IRasterLayer rLayer = new RasterLayerClass();
        rLayer.CreateFromDataset(rds);
        mxdoc.AddLayer(rLayer as ILayer);
    }

Above is the clip result behind the quad polygon I'd like it to be flush with (UTM 14) and a topo basemap layer.


texts: Aggregating Virtual Fields in the QGIS Layout

texts:
I have an atlas creating maps of parks. These parks contain different areas with different types of surfaces, which need to be managed differently.
In my map-layout I'd like to summarise the area of all different kinds of surfaces per park and I'm using the following function to display that (example for all Lawns in a certain Park):
[%aggregate(layer:='surfaces',aggregate:='sum',expression:="srea",  filter:="parknumber"=attribute( @atlas_feature, 'parknumber' ) AND "Type"='Law')%]
Needless to say, that the surface layer and the park layer are linked in a 1:n relation via the parknumber.
What seems not to work, is making the "Area" field a virtual field, so when geometries change, the area calculation changes automatically and the information in my map-layout changes automatically, too.
Or phrased differently, the aggregate function does not find/accept the virtual fields in the table.
Is there anything I could do, to get this solved?


texts: Can not find Points in Polygons and Spatial Query in QGIS 2.18.9

texts: I am using some old tutorial with Points in Polygon and Spatial Query but can not find these options in QGIS. Do I need to turn on some options to see them?


texts: Selecting layers of a raster stack using a sliding windows as a step in R in for loop

texts:
I would like to create a loop that every time selects the first ten layers of a raster stack, calculates their mean and saves it to a layer of a new raster stack. Also, I would like to have a step that works like a sliding window. So, in the first loop, it will select layers from 1:10, in the second loop it will select layers from 11:20, and so on.
So far, I have written the following code:

file_list=stack(myfile.nc)
mn = stack()
#
for(i in seq(1,length(file_list@layers), by=10)){
  mn[[i]] <- mean[[i]] # supposed to calculate the mean of the first ten layers
}


texts:
I want to use array_first(@myArray) to get the first value of myArray in a data defined box.
myArray is defined in the project variables, but I can't find the correct syntax. I tried :

(1, 2, 'abc')
[1, 2, 'abc']
array_first( array(1,2, 'abc')) returns 1 but as soon as I define myArray as array(1,2, 'abc'), it says "can't convert into list" (translation from french).


texts: Calculating distance in feet between points in a Pandas Dataframe

texts:
I am pulling my data from an AGOL feature service using the arcgis python api. I pull the data into a spatial dataframe then merge it in various ways in pandas, do various column calculations, move columns around, etc. I just can't figure out how to apply arcgis.geometry.distance or point.distance_to to whole columns in pandas to get a distance between the two columns in feet. The relevant columns of my data are in the following format:

SHAPE_1: {'x': -118.78663214499994, 'y': 46.70276817100006, 'spatialReference': {'wkt': 'GEOGCS["NAD83(NSRS2007)",DATUM["D_",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]]'}}

SHAPE_2: {'x': -118.78642560599997, 'y': 46.70243596700004, 'spatialReference': {'wkt': 'GEOGCS["NAD83(NSRS2007)",DATUM["D_",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]]'}}

Distance: ??


texts: no Get More in Plugin Manager in qgis 2 + ubuntu 13.04

texts: After upgrade of qgis 1.8 to qgis 2.0 I couldn't find the choice "Get more" in Plugin Manager and I can't install plugins.


texts: PostGIS: Find outers and inners inside MultiPolygon geometries

texts:
I'm trying to make a request to find outers and inners inside a multipolygon geometries. With this request I need to return a column with outers as GeoJSON and a column with inners as GeoJSON.
An image to explain my problem:

Multipolygon geom contains 5 polygons (A, B, C, D, E). A, B, E are outers and C, D are inners of A.
I can find number of rings with ST_NRings(geom), I can split multipolygon into polygons with ST_Dump(geom) but no luck with ST_Within.
Here is my query:
SELECT ST_AsGeoJSON(???) as outers, ST_AsGeoJSON(???) as inners
FROM(
    SELECT (ST_Dump(geom)).geom as d_geom1, (ST_Dump(geom)).geom as d_geom2
    FROM(SELECT * FROM schema.table WHERE ST_NRings(geom) > 1) AS Q1
) AS Q2

If I add ST_Within(d_geom1, d_geom2), result is true.
I don't understand how to loop through each polygon and determine if polygon[x] is inside polygon[y].
And how to return geom after result of ST_Within!


texts: TIleMill multiple design feature of the label road

texts:
Does TileMill have possibility to design features/labels in my case based on more than one command? I have this carto CSS code

#mainroad_label[type='primary'][zoom>13], {
  text-name:'[name]';
  text-face-name:@sans;
  text-placement:line;
  text-fill:@road_text;
  text-halo-fill:@road_halo;
  text-halo-radius:1;
  text-min-distance:60;
  text-size:9.5;

But i need one type of "primary" road at the zoom level from 12! SO i have written this code:

 #mainroad_label[type='primary'][zoom>13], {
          text-name:'[name]';
          text-face-name:@sans;
          text-placement:line;
          text-fill:@road_text;
          text-halo-fill:@road_halo;
          text-halo-radius:1;
          text-min-distance:60;
          text-size:9.5;
    [name='Leave Road'] [zoom>=12] {
      text-name: '[name]';
      text-face-name:@sans;
      text-placement:line;
      text-fill:@road_text;
      text-halo-fill:@road_halo;
      text-halo-radius:1;
      text-size:10; }
}

It did not show me this road at Z12 because "Leave Road" is "primary" road and obviously all "primary" roads are shown from zoom>13. When i have changed it to zoom>=12 it showed me all "primary" roads but i need show only "Leave Road" primary road at this zoom level. Is there any way how can i do it? Simply i would write something like [type='primary'], [name='Leave Road'] {.....} The result will be only this one "primary" road with "Leave Road" name shown at my target zoom level. Thanks


texts: Overpass API query error: query ran out of memory

texts: I am using a simple query the one that is in the examples which returns all data in the bounding box. Currently using overpass-api.de/api/. I have been trying for a couple of days thinking that server may be overloaded but i still get this error


texts: How to georeference a map - which is in projected coordinate system - using geographic lat long?

texts:
We have some hard copy maps which were projected in Gauss Kruger coordinate system (in the legend it has been mentioned). But, the coordinates were given in Geographic lat long (Beijing 1954). In Global mapper, if we enter the geographic lat long and change the projection, it will change the coordinates accordingly & will detect the Zone too. Is there any such thing in QGIS or even in ArcGIS.

We have tried by assigning geographic lat long (Beijing 1954) to the data frame (Default CRS in the case of QGIS) and entered the geographic lat long and tried with the Target CRS as Gauss Kruger projection. But still, it doesn't convert the coordinates correctly.

All i can see is, first we have to georeference the image in Geographic lat long and then we have to reproject (Warp in case of QGIS) the image. Is there any other options available ?


texts: How to determine the size of the cropped tiles we are dividing in using GDAL?

texts:
I tried the following code to divide my tiff file into smaller chunks of tiff files. I want my tile-size to be 30m*30m. However I am not sure if it is possible to do that. Does the tilesize parameter signify the length of the tile here?
import os, sys
from osgeo import gdal

dset = gdal.Open(sys.argv[1])

width = dset.RasterXSize
height = dset.RasterYSize

print width, 'x', height

tilesize = 5000

for i in range(0, width, tilesize):
    for j in range(0, height, tilesize):
        w = min(i+tilesize, width) - i
        h = min(j+tilesize, height) - j
        gdaltranString = "gdal_translate -of GTIFF -srcwin "+str(i)+", "+str(j)+", "+str(w)+", " \
            +str(h)+" " + sys.argv[1] + " " + sys.argv[2] + "_"+str(i)+"_"+str(j)+".tif"
        os.system(gdaltranString)


texts: Exporting the Google Street View coverage layer?

texts: How can I export the coverage layer of Google Street View (the blueish layer that appears when dragging the pegman on Google Maps) to (for example) GeoJSON format to work with it in different places?


texts: Geoprocessing Ungroups Layers in ArcGIS 10.1

texts:
Is there a fix for this, am I only imagining this, or is it a bug:

ArcMAP 10.1 for Desktop Build 3035 on a 64 bit Windows 7 system

When I run a geoprocessing task that modifies the input data (such as Snap), if the input layer is in a group, the group is removed once the process is complete! I have lost groups with detailed definition queries to this problem.

I suppose removing the layer from the group before processing will help alleviate this. However, is there a more definite fix? I can easily see myself forgetting to separate that layer and subsequently losing my group.

Plus, there is no 'save map every x minutes' feature (that I know of), so closing without saving will bring back the groups, but will negate any changes I made while working.


texts: ArcGIS 10.1: Subtype/Domain Functionality

texts:
Is there any way to create more complex conditional domains in a table than a single subtype and associated domains? I would like to get the same functionality that more than one subtype in a single table would provide but ArcGIS only allows one. I want to be able to use multiple domains in one field that are part of a subtype, where each domain might only allow other domains used in the next field.

Any ideas?


texts:
I am tracing the Postgis topological network with Postgresql via:

WITH RECURSIVE path AS (
 SELECT edge_id,start_node, end_node
   FROM simulation.edge_data
   WHERE start_node = 1 OR end_node = 1 OR start_node = 8 OR end_node = 8
 UNION
 SELECT e.edge_id, e.start_node, e.end_node
   FROM simulation.edge_data e
   INNER JOIN path s ON
     (s.start_node = e.start_node
     OR s.start_node = e.end_node
     OR s.end_node = e.end_node
     OR s.end_node = e.start_node) )
SELECT * FROM path;

I want to improve the query in only visiting every node once. Also this network is full of gas pipelines. If the node is a valve that is closed i want to not further track in that direction. This returns me all nodes that have access to gas but not all edges. What am I doing wrong?

I do this with:

WITH RECURSIVE path AS (
 SELECT edge_id,start_node, end_node
   FROM simulation.edge_data
   WHERE start_node = 1 OR end_node = 1 OR start_node = 8 OR end_node = 8
 UNION
 SELECT e.edge_id, e.start_node, e.end_node
   FROM simulation.edge_data e
   INNER JOIN path s ON
     (SELECT * FROM topology.validateAttribute('topology','valve','status','closed',e.start_node,1)) = false
     AND (s.start_node = e.start_node
          OR s.start_node = e.end_node
          OR s.end_node = e.end_node
          OR s.end_node = e.start_node) )
SELECT * FROM path;

This is my graph, where purple represents the gas sources. The yellow marked are closed valves. The result of the above query is missing the two edges 2 and 3, but includes the correct nodes:


texts: Removing overlapping parts of polygons from the same layer in QGIS

texts:
I'm working on a favourability mapping based on geologic inputs.
One of those inputs is composed of fault zones which in my case are represented as one unique polygon layer with one attribute that represents the favorability of a certain fault to be open or not defined in my case by 5 "scores.
The problem I have is that the different polygons composing the layer overlap, meaning that I have two or more fault polygons with different favorability in some areas. You can see an example here:

Therefore, what I would achieve is to produce a new layer where the overlaps disappear and where the priority is given to the highest favorability.
I'm doing that by exporting each favorability class to a new layer and then sequentially use the "difference" tool to eventually produce what I need, but it is a tedious process and I'm sure that QGIS has some tools that could help but for the moment I did find nothing that really works ("SAGA Vector polygon tools", "v.clean", "Eliminate sliver polygons", "Check Geometries" do not provide results that can even partially help me).


texts:
Using two different time of images in the same places. Can i develop an automated change detection for forest area?

What / which are the suitable software to used?


texts: Remove values below a given value with GDAL for DEM geo tif files

texts:
I have a DEM in GeoTiff format. The processing have made areas that have really low values that I would like to remove.

Are there any GDAL tools that could set all values below -30 to NODATA or -30.


texts: Finding gaps and overlaps in same feature class with Data Reviewer?

texts:
I am trying to find gaps and overlaps in the same feature class (e.g. political boundary) with Data Reviewer.

The check I am using is Polygon Overlap/Gap Is Sliver (Feature on Feature Checks), and it requires me input 2 feature classes in the check window, so I input the political boundary in the Feature Class 1 and input it again in Feature Class 2.

Is this the right way or procedure?


texts: For some reason my sub isn't 'linking' in MapBasic, I'm pretty new to this, but why is the sub assuming it begins on line 110 rather than where it actually says sub?


texts:
I am new to GIS .
I would like to create a geodatabase in ArcGIS that stores information about different regions in a country. In each region there are few villages that are selected to monitor their  land cover change by taking images (jpeg) of those sites every 2 years. The database should have all the images taken each two years and different features found in that region such as services ( roads, rivers and villages). It should be able to add new information every two years.


texts: Carto: Can I use some “if/then/else” structure in popup custom HTML?

texts:
I believe I am asking a question similar to this question that was asked a few years ago about popup functionality in the old CartoDB. I just tried my best to replicate the solution to that question in the newer Carto Builder, but I'm failing. (I tried using the # and ^ symbols in what I thought were similar places to where they were used in the solution to the question above; doing so didn't break the HTML, but it didn't work either).

My algorithm is:

If <field_A> is not null, show fields A, B, and C in the popup
If <field_D> is not null, show fields D, E, and F in the popup

Or, an alternate algorithm could be:

If <field_A> is not null, show fields A, B, and C in the popup
If <field_A> is null, show fields D, E, and F in the popup instead

This is the default HTML that Carto created for me to simply show fields A, B, and C in the popup. I'd like to edit this HTML to match my algorithm above. My arrows show where I think I need to edit the HTML, but I'm at a loss for how to actually edit it:

field_A
    <-- I think this is where I need to specify "if field_A is not null..."

{{field_A}}
     <-- And this is where I need to specify "...then show field_A"
                                                                  <-- And then here is where I'd need to specify "...else show field_D"

field_B
    <-- Repeat the above for the next set of fields

{{field_B}}

field_C
    <-- Repeat the above for the last set of fields

{{field_C}}


texts:
I would like to download MODIS npp images for Mexico but I cannot figure out which files correspond to the country. How can I associate each tile to its location without looking at the million of browsing images. I am an experienced python and R user if that helps.

I am manually downloading the images from this link: http://e4ftl01.cr.usgs.gov/MOLT/MOD17A3.055/2000.01.01/


texts: Raster of missing values inside polygon in R

texts:
I have a raster rrepresenting values from 0 to 50 for the world. There are large areas with missing values NA values, including:

All oceans
Large inland regions such as deserts

The problem is that when I plot r all NAs appear with the same colour regardless of whether they are oceans or deserts. This does not look nice, so I want to display inland NAs in a different colour than oceans.

I was thinking perhaps I could grab a map representing continents (e.g. map("world")), and create a layer of all NAs in r that fall inside continents. Then, I can just overlap this new raster on top of r in a different colour.

How can I create a new raster with all the inland missing values in a raster?


texts: Adding shapefile layer in PostGis creates two tables

texts:
When I try to upload a shapefile to my PostGIS database, it creates two tables (I can see them in QGIS). I import them by drag and dropping them on the PostGIS connection.
If I delete one of them, both disappear. While it doesn't hurt, it seems worrying.
For a reason I don't understand, it happens especially with Lambert-93 tables.
If I am doing it wrong, how can I fix that?


texts: Checking if selected features are adjacent using OpenLayers

texts:
I need to be able to merge features on the condition that the selected features are adjacent to each other.

How can I achieve this using Open Layers version 5?


texts: Moving feature classes between FGDB's creates unknown topology errors

texts: I am attempting to merge together multiple large feature classes from different FDGB's into a single feature class and am getting unknown topology errors when checking for gap and overlap rules.  I have run topology checks on these feature classes before and the layers are clean, however for some reason when I merge the different datasets together the resolution of the features seems to diminish and create topology errors that weren't there prior to the merge.  I am thinking this could be an issue caused by geoprocessing the feature classes from different FDGB's but can't confirm this is true.  Any ideas on issues that may be causing this would be extremely helpful.  I am running ArcMap 10.2 with 64 bit processing.


texts: New QGIS versions and nightly builds for mac?

texts:
When I head over at Kyngchaos QGIS page I notice the latest version isn't available. Also when it comes to nightly builds it seems like the ones over at Dakota Cartography hasn't been updated since september 2016.

Are there any other sources of recent versions for Mac without compiling it all from scratch?


texts:
I'm working on a python script for ArcGIS 10.2/10.3. What I want is to automatically fill in a field of a point layer based on the location of that point in relation to a polygon layer.

Say I have a two layers: cities and regions in a database, and the cities layer has a field called "region_name". Then all I want is that I can fill in this field "region_name" based on the location of the city points within the regions.
I know this question has been asked before but I'm stuck at a specific point:
I know I can either use Select by location or a Spatial Join but what bothers me is that this means that I get a new layer with the result while what I need is the same "cities" layer (without additional fields) just having the field "region" filled-in.

It sounds so simple so I'm wondering whether I have just missed something.
Also since I have a lot of cities and regions it really needs to be done programmatically (not manually) but if you tell me it's easier if I use something else than Python (I'm working with an enterprise geodatabase) I'm happy to change.


texts: Coordinate input does not seem to be correctly processed when Georeferencing in QGIS?

texts:
I am trying to georeference an old scanned map. I divided the scan in several parts. The original map uses an uncommon coordinate system (Equidistant Conic Projection according to De L'isle with Parallel lines at 50° and 53°), which I tried to manually create by using the following proj4 code:

Proj4: +proj=eqdc +lat_0=30 +lon_0=10 +lat_1=50 +lat_2=53 +x_0=0 +y_0=0 +ellps=intl +units:m +no_defs

To georeference the map, I used the qgis georeference plugin. However, something seems to go wrong in processing the coordinates. I entered the coordinates in degrees and used the following transformation settings (sorry in German, problem still occurs when changing to other coordinate systems):  .

Below in the picture you can see the coordinate format (does this seem correctly specified?).

When I now run the georeferencing, and load the result into qgis, it does not display the coordinates correctly. In particular (compare mouse position/coordinate change), it seems that the picture is shrinked (no change in coordinates when mousing over the picture). Coordinates are displayed in Decimal Degrees. The pictures display two different mouse positions on the map.

What am I doing wrong?

I am a beginner in this field.


texts:
I have a dump from a version of PostGIS that is no longer supported.  It dumps the value of a (lat, long) as 01010000003D2CD49AE6A549409AB1683A3B19F1BF

How can I decode that number (without using PostGIS) to determine the original lat long?


texts: How to get raster layer pixels contained in a shape defined in a vector layer

texts:
I am new in QGis python plugins dev and I have a "maybe stupid" question (in the title :-) ).

I try to retrieve pixels (from a raster layer) which are "contained" in a polygon defined in a vector layer.
How can I achieve that via a python plugin ? The user has to select the vector layer and the raster layer using the GUI of the plugin (which is already coded).

Thanks in advance !


texts: Creating an SQL expression and writing/reading it to/from a file

texts:
I'm working on a way to take the number for projects that are marked as active in a spreadsheet and put that into a SQL expression that I can use to select from. So far I've got this

dirname = os.path.dirname(__file__)
gdbname = 'G:\\GIS\\General Maps\\AltumProjects\\AltumProjects.gdb'
env.workspace = gdbname
env.overwriteOutput = True
Riv = gdbname + '\\RivParcels'
df = pd.read_csv(dirname + '/ActiveProjects.csv')
newdf = df.loc[df['Status'] == 'Active', 'Project'].reset_index()
print(newdf)
Proj = np.array(newdf['Project'])
with open("Output.txt", "w") as text_file:
    print(f"\"Project IN ({Proj})\"", file=text_file)
with fileinput.FileInput("Output.txt", inplace=True, backup='.bak') as file:
    for line in file:
        print(line.replace("' ", "', "), end='')
with fileinput.FileInput("Output.txt", inplace=True, backup='.bak') as file:
    for line in file:
        print(line.replace("[", ""), end='')
with fileinput.FileInput("Output.txt", inplace=True, backup='.bak') as file:
    for line in file:
        print(line.replace("]", ""), end='')
expr = open('./Output.txt').readline()

Which creates a text file that is:

"Project IN ('C1335', 'C1410', 'C1237')"

This all works fine but when I add it to SelectLayerByAttribute:

Projects = arcpy.MakeFeatureLayer_management('Projects', 'Projects.lyr')
Env = arcpy.SelectLayerByAttribute_management(Projects, "NEW_SELECTION", expr)
arcpy.CopyFeatures_management(Env, gdbname + '\\ActiveEnvProjects')
EnvProjects = gdbname + '\\ActiveEnvProjects'
SelectEnv = arcpy.SelectLayerByLocation_management(Riv, 'COMPLETELY_CONTAINS', EnvProjects, 0, 'NEW_SELECTION')
arcpy.CopyFeatures_management(SelectEnv, gdbname + '\\EnvProjectParcels')
print("Environmental Projects Done")

I get:

  ERROR 000358: Invalid expression

If I add this in manually as a string it works fine but I can't get it working like this. Is there something with how I'm writing the expression wrong?


texts: Impossible to sum rasters because of nodata pixels

texts:
I have a set of rasters (.vrt) with daily soil moisture data. I want to sum the pixels of all rasters in order to have a measure by month. However, the whole world is not covered each day, which results in nodata values at places where measures exist for the other days.

What I want to do, is sum the values of each raster. However, it seems that each time a nodata pixel is in the sum, the resulting pixel is directly categorized as nodata. I would like to have the opposite: ignoring all nodata values and summing the rest.

I thought of 2 ways of solving the problem:

summing rasters ignoring nodata values
converting nodata pixels to value 0, then sum all the rasters

Unfortunately, I can't find any tool to do this.

Can anyone help me?


texts:
I'm trying to make a matrix that shows the cost from a node to the rest of the nodes in a network.

The final result I'm looking for is not the cost, but only know if a car can or can not reach each vertex.

I've made a test in a small area of my network, and it seems to work perfectly. I get the cost from a node to the rest, and if can't reach it, only leaves the value as NULL, then I can set this value to false.

The problem is where I try to extend this methodology to the rest of the network, I get an error showing that the target vertex was not found.

Why is pgRouting saying this, if I'm doing the same analysis?

The steps I've followed are:

--creating topology
ALTER TABLE "3prueba" ADD COLUMN "source" integer;
ALTER TABLE "3prueba" ADD COLUMN "target" integer;
SELECT pgr_createTopology('3prueba', 0.00001, 'geom', 'id');
--(I already have the cost field)
--indexes
CREATE INDEX vials3_source_indx ON calles("source");
CREATE INDEX vials3_target_indx ON calles("target");
--and
create table "3catchment82" as
select
    id,
    the_geom,
    (select sum(cost) from (
       SELECT * FROM pgr_dijkstra('
       SELECT id AS id,
          source::int4 AS source,
          target::int4 AS target,
          coste::float8 AS cost
       FROM "3prueba" where mostrar = 1
       ', --<-- show only roads
       82, --<-- specify node
       id::int4,
       false,
       false)) as foo ) as cost
from "3prueba_vertices_pgr"

This are the steps I've followed.

The main idea is get the costs only for the roads with the field "mostrar=1", as I show in the image.  If the road is closed (mostrar=0) it can not reach the node. So, the red road is closed, and can not reach nodes 58, 59, 5...
After this,I'll change the roads. So, I can get a matrix when some roads are closed, other matrix when other roads are closed, etc.

When I apply this analysis to the rest of the network, I get the error.

Any Idea of what am I doing wrong, or what would be the correct way to do this?

After @underdark advice, it seems to work.

I've run the model applying a high cost to the 'closed' roads, and seems to be ok.

I've seen that there are 3 types of output:

In green: nodes that the source reaches.
In red: nodes that source do not reach.
In orange: nodes disconnected. The values was NULL.

So, I think the problem was with this type of nodes.


texts: How do I select all points from a range of dates?

texts: Is there a way to select points from a range of dates to create a new layer? I want to clean up my data by making all points from previous years a separate layer but doing a selection by attributes and choosing each date individually is taking way too long.It would be perfect if there was a way to query out all points from 2013, for example.


texts: Emitting map click using ArcGIS API for JavaScript?

texts:
I am trying to emit a map click from the MapView object using its emit method. I have successfully done this in the 3x API from the Map object.  This same method, however, is not working for me in the 4x API.  The Map object no longer has an emit method, and I can't get the MapView object's emit method to work. It is not generating any errors for me to debug.  This is the code I am currently using:

var mapPt = evt.result.feature.geometry;
app.mapView.emit("click", {mapPoint:mapPt});

Based on the API documentation, it is not clear to me why this isn't working.  I have also tried specifying the screenpoint X and Y.  This had no effect.  I tried simulating a mouse click using vanilla javascript and with jquery, but neither method interacts with the app or MapView objects.  The click event triggered outside of the API appears to be superficial to the DOM.

I am trying to emit a click event after the search widget returns a result, as I want the popup behavior/content to be identical to what is created using a map click.  I am trying to avoid having to re-code all the popup templates for the various layers specifically for the search widget when I have already specified the popup templates when I defined the feature layers.


texts: How to config geoserver to read style of vector layer from postgis/SQL store?

texts:
I've imported vector layers from map info pro to both SQL server and post gis using easy loader.
In both cases I've included style column in the exported table.
When I create new store in geoserver and see the list of layers , each layer has the mi_style column in the list of features at the bottom of the screen . But when I publish the layer it still have the default style of geoserver . How can I set it to read the style from the column ?


texts: Is there any readily available implementation of Mapserver (or Geoserver)- Postgres-Openlayers combination?

texts:
I am developing an web GIS for map data presentation for public.

Is there any readily available implementation of Mapserver (or Geoserver)- Postgres-Openlayers combination or similar?

More Details:

Overall system requirements: Uploading, managing and publishing map
   data

Back-end :

Upload geospatial data into the administration interface.  create new
layers and organize.

Admin should be able to:

Organize the mapping layers in a navigation menu
Edit the name, legend, and
Choose the order of the mapping layers
Change the thematic colors of layers
Customizable icons for map layers in the table of contents ( should
reflect in the published map also)
Should be able to add few attributes from the front end.

Front-end should allow the public to:

Add google satellite layer.
Zoom directly to counties from the navigation  pane.
Turn ON/OFF layer visibility
Turn ON/OFF layer legend
Turn ON/OFF layer labels
Change the layer order(up-down) in the navigation
Export the map to PDF
Print layout option(with or without legend, labels etc.)

Please suggest.


texts: Is there a QGIS Python binding for generating a grid programmatically without using the GUI?


texts: Classification of NDVI using Google Earth Engine

texts:
I want to classify my NDVI raster using training data in Google Earth Engine. For that i am using the following code. But it is returning some error like ndvi is not a usable band.
Please help in this

// Cloud Masking function
function maskS2clouds(image) {
  var qa = image.select('QA60');

  // Bits 10 and 11 are clouds and cirrus, respectively.
  var cloudBitMask = ee.Number(2).pow(10).int();
  var cirrusBitMask = ee.Number(2).pow(11).int();

  // Both flags should be set to zero, indicating clear conditions
  var mask = qa.bitwiseAnd(cloudBitMask).eq(0).and(
             qa.bitwiseAnd(cirrusBitMask).eq(0));

  // Return the masked and scaled data.
  return image.updateMask(mask).divide(10000);
}

// AOI of Study Area
var boundary = ee.FeatureCollection('ft:1rbhNtC1TqDBvY9Rt2BZR-DjhpIPuC3nU5kmz49WW');//Jorhat Boundary
//var boundary = ee.FeatureCollection('ft:1ABffZYEE4XhMTOXsoSfWENKXBK2fQfOrdMplEaVo');//Midnapur Boundary

// Import of Images (Sentinel 2 multispectral)
var image = ee.ImageCollection(sent2img
  .filterDate("2017-12-01","2018-01-30")
  .filterBounds(boundary)
  .map(maskS2clouds)
  .sort("CLOUD_COVERAGE_ASSESSMENT")
  .median()
);

// Preprocessing
var mosaic = image.mosaic()
var clip = mosaic.clip(boundary);
print(clip);

// FCC creation and visualisation of AOI
Map.addLayer(clip, {bands: ['B8','B4','B3'], min: 0, max: 0.3},'clip');

// NDVI  calculation
var ndvi = clip.expression(
    ' ((NIR - RED) / (NIR + RED))', {
      'NIR': clip.select('B8'),
      'RED': clip.select('B4'),
}).rename('nd');

print(ndvi);

// Colour Palette
var palette = ['FFFFFF', 'CE7E45', 'DF923D', 'F1B555', 'FCD163', '99B718',
               '74A901', '66A000', '529400', '3E8601', '207401', '056201',
               '004C00', '023B01', '012E01', '011D01', '011301'];

// Map display
Map.addLayer(ndvi, {min:0, max:1, palette: palette},"NDVI");

// Training Classes
var newfc = wb.merge(plantation)
              .merge(agriland1)
              .merge(agriland2)
              .merge(habitation)
              .merge(sandyarea)
              .merge(deepplantation);
var bands = ['nd'];
var training = clip.select(bands).sampleRegions({
  collection: newfc,
  properties: ['class'],
  scale: 20,
  geometries:true
});

var classifier = ee.Classifier.cart().train({
  features: training,
  classProperty: 'class',
  inputProperties: bands
});

var classified = clip.select(bands).classify(classifier);

Map.addLayer(
  classified,
  {min: 1, max: 7, palette: ['#0d1898', '#ff0841', '#138b11','#fff81c','#154c17','#529400','#F1B555']},
  'classification');


texts: Map to KML tool issue: legend is cut off in Google Earth

texts:
I am facing this issue with MAP to KML tool in ArcGIS Desktop (10.2.2): I have created the legend in Map Layout as indicated here then exported to KML. But the legend is partially cut off in Google Earth. I have already tried increasing spacing and halo, but it does not solve the problem.
Does anyone have experienced the same issue and/or can suggest solutions?


texts: Algorithm "Merge Layers" from SAGA not showing in QGIS Modeler

texts:
I'm trying to merge multiple vector-layers via the modeler in QGIS 3.14.1 as a part of a bigger model.
As a matter of fact the QGIS-native "merge layers" algorithm doesn't work here, for the layers to be merged have some different data types, due to the layer creation by a different software, which should not be the topic here.
So I need to use the SAGA-algorithm "merge vector layers", which used to work out fine in previous versions of QGIS (up to 3.12), but now it does not show natively in the toolbox or the modeler.
After I check the Box "show algorithms with issues" in the processing preferences, the mentioned SAGA-algorithm shows up in the toolbox, marked with a special warning (see picture below), but it works as usual without any errors or problems.

But: In the modeler the algorithm is still not available, so I still can't use it for my model.
Anyone knows a solution for this?
I tested it both on a Windows and on a Linux System, it's exactly the same issue.
The use of the QGIS-native "merge layers" algorithm is not an alternative here, for I cannot influence the creation of the input layers.


texts: How to use node.js with ArcGIS Web AppBuilder to list files in a widget

texts:
I have built an ArcGIS Web AppBuilder
and a custom popup (widget.js). My goal is to list every file in a specific directory (in 'wwwroot' folder) to the popup user.

I didn't find any built-in methods, so I'm attempting to use node.js.
I installed node glob lib into my app:
    npm install --prefix C:\inetpub\wwwroot\WABCustomApp\libs glob
(first I tried to install only glob module in the dir, that I copied all nodejs dependencies for another try.)
I configured init.js to load this libary:
    window.path + 'libs/nodejs/node_modules/npm/node_modules/glob/glob.js'

And received the following error:
Uncaught ReferenceError: module is not defined at glob.js:41

How would you make Web AppBuilder work with nodejs?

If you wouldn't use this approach, what would you suggest for listing dirs in popups?


texts: Selecting all the pixels in a raster coverage that have a specific pixel value? (in PostGIS raster)

texts:
I have a land use map in the format of TIFF. I imported this file as a table in PostGIS using raster2pgsql. The land use file has 7 different classes.

How can I select all those pixels that have the value (class) of 5? The table has only two columns (rid, rast).


texts: Is it possible to withhold attribute data from a shapefile?

texts: Basically, what I'm trying to figure out is if its possible to withhold or lock attribute data from a shapefile if you are sending the shapefile to someone?


texts: Why when I use the ReadAsArray() Shell restart in this code?

texts:
I want to read the .pix file in the 3 different  Array in python and I use this program but when I run that  nothing do And shell be restart and program closed .

from scipy import ndimage
import os,csv,math,subprocess
from osgeo import gdal
from gdalconst import *
os.chdir(r'"D:\RADARREXAMPLE\outputsac ')
def importpix (a):
    # open the image
    fn = str(a) + ".pix"
    ds = gdal.Open(fn, GA_ReadOnly)
    print ds
    # get image size
    rows = ds.RasterYSize
    cols = ds.RasterXSize
    band = ds.RasterCount
    # Set NoData Value
    R = ds.GetRasterBand(1)
    G = ds.GetRasterBand(2)
    X = ds.GetRasterBand(3)
    B = ds.GetRasterBand(4)
    c11 = R.ReadAsArray(0, 0, cols, rows)
    c22 = G.ReadAsArray(0, 0, cols, rows)
    c33 = B.ReadAsArray(0, 0, cols, rows)
    c11=pow(c11,2)
    c22=pow(c22,2)
    c33=pow(c33,2)
    # Delete the band and dataset objects
    ds =None
    R = None
    B = None
    X = None
    G = None
    # Write the files to ENVI format
    driver = gdal.GetDriverByName('ENVI')
    ds = driver.Create('c11.bin', cols, rows, 1, GDT_Float32,['SUFFIX=ADD', 'INTERLEAVE=BSQ'])
    ds.GetRasterBand(1).WriteArray(c11, 0, 0)
    ds = None
    ds = driver.Create('c22.bin', cols, rows, 1, GDT_Float32,['SUFFIX=ADD', 'INTERLEAVE=BSQ'])
    ds.GetRasterBand(1).WriteArray(c22, 0, 0)
    ds = None
    ds = driver.Create('c33.bin', cols, rows, 1, GDT_Float32,['SUFFIX=ADD', 'INTERLEAVE=BSQ'])
    ds.GetRasterBand(1).WriteArray(c33, 0, 0)
    ds = None
    # Write the accompanying config.txt file
    file = open('config.txt', 'w')
    text = "Nrow\n" + str(rows) + "\n ---------\nNcol\n" +str(cols) + "\n ---------\nPolarCase\nmonostatic\n ---------\nPolarType\nfull\n"
    file.write(text)
    file = None
>>>
>>> importpix ('test11')
<osgeo.gdal.Dataset; proxy of <Swig Object of type 'GDALDatasetShadow *' at 0x00000000032A6150> >
13977 7985 4

>>> ================================ RESTART ================================


texts: i have an unhandled exception in asp.net

texts:
i am going to create a website and load a gis map to it.

i am using arcgis 10 and i created a web application in arcgis server 10

i created a new WebSit in visual studio 2010 and add a map resource manager.

add a MapResourceItem for this.

after that add a map control and set MapResourceManager property.

i config ArcGISIdentity too.

but unfortunately when i run application occurs this error

  Error processing server request

  [SoapException: Error processing server request]
     System.Web.Services.Protocols.SoapHttpClientProtocol.ReadResponse(SoapClientMessage message, WebResponse response, Stream responseStream, Boolean asyncCall) +4159
     System.Web.Services.Protocols.SoapHttpClientProtocol.Invoke(String methodName, Object[] parameters) +328
     ESRI.ArcGIS.ADF.ArcGISServer.MapServerProxy.ComputeScale(MapDescription MapDescription, ImageDisplay MapImageDisplay) +117
     ESRI.ArcGIS.ADF.Web.DataSources.ArcGISServer.MapFunctionality.GetScale(Envelope extent, Int32 mapWidth, Int32 mapHeight) +564
     ESRI.ArcGIS.ADF.Web.UI.WebControls.Map.get_Scale() +794
     ESRI.ArcGIS.ADF.Web.UI.WebControls.Map.setExtent(Envelope value, Boolean addSetExtentCallbackResult) +2556
     ESRI.ArcGIS.ADF.Web.UI.WebControls.Map.set_Extent(Envelope value) +39
     ESRI.ArcGIS.ADF.Web.UI.WebControls.Map.validatePrimaryMapResource() +2586
     ESRI.ArcGIS.ADF.Web.UI.WebControls.Map._OnLoad(EventArgs e) +115
     ESRI.ArcGIS.ADF.Web.UI.WebControls.Map.OnLoad(EventArgs e) +37
     System.Web.UI.Control.LoadRecursive() +134
     System.Web.UI.Control.LoadRecursive() +232
     System.Web.UI.Control.LoadRecursive() +232
     System.Web.UI.Page.ProcessRequestMain(Boolean includeStagesBeforeAsyncPoint, Boolean includeStagesAfterAsyncPoint) +3412


texts: Huge processing time for gdalwarp command extracting ROI by extents

texts:
I'm running the below gdalwarp command to extract a 1024 x 1024 ROI from a large raster which has Dimensions X: 34396 Y: 36707 Bands: 4

Input raster has previously been pan-sharpened and orthorectified from original World View source, admitedly it's a fairly large file at 19GB, but initial steps did not take too long, however the below is taking nearly 36 hours to complete within an OSGeo4W shell?

gdalwarp -t_srs EPSG:32630 -r cubic -wo SOURCE_EXTRA=1000 -tr 0.5 0.5 -r cubic -of GTiff -te 419175.301303 5577359.141256 421223.301303 5579407.141256 F:/GIS/orthorectified.tif F:/GIS/extract_sat_WV_5.tif

The output is a 16mb tif

This is longer than I've ever waited to process a raster, I tested a similar process within Global Mapper (GIS application) and it took 2 minutes to cut using a shapefile and export to tif within the interface.


texts: Why does QGIS not detect the CSV/DBF column headers?

texts: I am trying to import a CSV or DBF into QGIS so I can perform a join on a shapefile.  However, QGIS is not detecting the first row as header names.  I know I've done similar joins before and there's never been a problem detecting column headers.


texts:
I'm using a PostGIS View as a layer to join multiple lookup tables.
I published this view through Geoserver to edit using WFS-T.

My first question is that this is right strategy using a VIEW as a layer for Geoserver's WFS-T or not?
Is there any way to serve normal table as a layer using some alternatives to lookup tables?

Second question, if I have to use view as a Geoserver's layer in this case, is it possible to edit this view using WFS-T?


texts: QgsSpatialIndex.nearestNeighbor() always returns 0 id

texts:
I'm creating a set of features and adding them to a QgsSpatialIndex to find the ones closest to a certain location.

Why does nearestNeighbor() always return [0L] instead of the correct feature id?

import random

index = QgsSpatialIndex()

for i in range(0,10):
    feat = QgsFeature()
    pt = QgsPoint(
        random.randint(0,10),
        random.randint(0,10))
    feat.setGeometry(QgsGeometry.fromPoint(pt))
    index.insertFeature(feat)

nearest = index.nearestNeighbor(QgsPoint(5,5),1)
print nearest


texts: How can I set the output feature type name dynamically?

texts:
How can I set the output feature type name so that it consists of a constant part and a variable part, such as "Circle_205", where the number part could be different each time?

When using a fanout, this is easy to set, but without a fanout, I do not know how to do this, unfortunately.


texts: ArcGIS - Help Using Oracle Function in Definition Query - No ArcSDE

texts:
I have set up a Query Layer pointing at an Oracle table with a SDO.Geometry column, and that is working and behaving fine.
However, when I am attempting to set up a definition query that calls a Oracle function, I get an error message saying "There was an error with the expression.  Failed to parse the where clause".
I am fairly confident the SQL expression I am using is valid:  it works fine inside of SqlDeveloper.  Here is what the sql looks like:
select * from gis_meters where gis_getmeterstatus(METERNAME,'2013/08/08 10:25:00','2013/08/09 10:25:00') = 0

In the Query Builder of ArcMap, I am only adding, as the query layer alread has the select *:
gis_getmeterstatus(METERNAME,'2013/08/08 10:25:00','2013/08/09 10:25:00') = 0

After some research, I found some ESRI documentation stating "Coverages, shapefiles, and other nongeodatabase file-based data sources do not support subqueries".  I am guessing that the call to the Oracle function (gis_getmeterstatus) is a subquery - and thats why its not working for me?  Since I don't have SDE installed in the database, its not considered a geodatabase.
Can anyone help me find a way to use this function in my definition query?


texts: Google maps on iOS - query string to achieve what I need

texts:
I would like Google maps to show indoor maps on iOS. I use Google API on Android, but accordingly to their site, Google maps SDK for iOS does not support indoor mapping.

Thus, I think the closest thing I can come is showing URLs like this in embedded browser:

  https://maps.google.com/maps?q= %string% &output=embed'

I would like the query string to achieve the following:

Multiple markers (I think one is max)
Enable position location (so user can see himself as arrow/dot)
Hybrid view


texts: Python: How to export Shapely objects to DXF file

texts: Using shapely package for Python how to export the resulting shapely objects such as buffer to a DXF file?


texts: How to remove decimal places from currency formatting in ArcGIS?

texts: For years I have put up with the 2 decimal places that ArcGIS adds to the "currency" formatting, but most of the time I do not want them and end up removing them in my exported map images. I cannot see a way to edit the number of decimals - does anyone have any ideas how to get rid of these decimals?