texts: Merging overlapping polygons with SpatiaLite

texts:
The green points, buffered by one of the point attributes, should result in the three areas outlined in red (a table with three rows):

How can I do this with SpatiaLite?
If it's not possible, is there any other way to do this in QGIS?


texts: Extracting raster values from polygon perimeters using ArcGIS Spatial Analyst?

texts:
I am currently working on a project to identify local ponds or water bodies that have high probability of species occurrence. I'm using BioClim climate data, SSURGO soil data, and a DEM layer for a MaxEnt analysis. On top of that I have a global water presence layer from the EU Global Surface Water Explorer. My problem at this point is that the soil layers do not have data where water is present, even if they did I would think they would be moderately different from the outlying perimeter of the water bodies. At this point, the water bodies have been converted to a polygon layer, but I need to extract probability values from the MaxEnt output around the perimeter of these water bodies. I tried using a buffer style extraction, but I got an error that I didn't have enough memory. I could deal with that if necessary, but is there another way around this issue?

I do have the Spatial Analyst toolbox.


texts:
I am working on modeling surface run-off due to a rain event (e.g., 5 year rain for a duration of 30 mins.) over an approx. area of 30 km2 for which I have a DEM and necessary soil as well as land-cover data (Manning's n values). I am going to implement Green-Ampt infiltration model in ArcGIS 10.2 (with Python scripting) and have already computed rasters like flow direction, slope, flow accumulation, stream network, stream length, etc. with hydrology tools.

At each time-step, I have an excess water depth raster and velocity raster (with slope raster from slope tool) at each cell with which I can calculate the discharge volume that is flowing out from a cell to its neighbor. Now, my problem is:

For the next time-step I have to add those discharge volumes to the
corresponding cells to ensure mass balance but how can I determine
to which neighboring cell is the discharge flowing? I only know
that ArcGIS is calculating the steepest gradient for each cell but
how can I ascertain which out of 8 neighbors does that gradient
correspond to?
For a time step and a cell, how can I define a run-off water depth?
Is it simply (outflow volume - inflow volume)/area of cell?


texts: When trying to save a vector in QGIS a blank dropdown window appears

texts:
I am learning and using QGIS Brighton and have hit a obstacle. Everything has been going fine until I started using the experimental generalizer to smooth the contours that I had extracted from some DEM's but it may be something else.

When I create a vector, it all works well until I right click to save. I get a blank dropdown screen. I uninstalled generalizer and then QGIS did a reinstall to no avail. Several times. This blank drop down is starting to apply to other menus (contour extraction, attribute table,   qgis-bin etc.) it is a little bit annoying as you can imagine. Any ideas in how to resolve?


texts: Exporting list using Google Earth Engine

texts:
I want to extract lat/lon values in GEE and I found the following code on stack exchange which does the job. However, when I try to print the lat/ lon values in GEE, it starts to crash. I am going to take a guess here and say it most likely because there is a large number of values ( I could be wrong). Therefore I am trying to export the lists (var point_list in this instance) into a csv, but not quite sure to go on doing that.

Any ideas?

var dataset = ee.ImageCollection('USDA/NASS/CDL')
              .filter(ee.Filter.date('2017-01-01', '2018-12-31'))
              .first()
var ROI = ee.Geometry.Polygon(
[[[-98.32711839622601,43.77454409032633],
  [-98.32711839622601,43.296678036299085],
  [-97.40426683372601,43.3126679077047],
  [-97.44821214622601,43.80626716084386],
  [-98.32711839622601,43.77454409032633]]]

var image_clipped = dataset.clip(ROI)

var proj = image_clipped.select('cropland').projection()

var latlon = ee.Image.pixelLonLat().reproject(proj)

var coords = latlon.select(['longitude','latitude']).reduceRegion({
             reducer:ee.Reducer.toList(),
             geometry:ROI
             })

var lat = ee.List(coords.get('latitude'))
var lon = ee.List(coords.get('longitude'))

var point_list = lon.zip(lat)

var coords = cropLandcover.pixelLonLat()


texts: Extract JPEG out of raster based on rectangle location in QGIS

texts:
I have 10 overlaying rasters (both one- and multi-band) - same area, different date. I also have shapefile layer with 12 rectangle polygons.
What I want is to extract proportion of raster inside each of rectangle and save it as JPEG image. Sure, it's easy to do it step by step. I select one polygon, clip the raster ten-times for each raster and then move to another polygon. The procedure is clear but please, is there any way how to make it automatically? I mean always select one polygon, extract it 10 times and then select another one - together make 120 pictures ideally in JPEG and save it to folder.
Maybe I could use Graphical Modeler but I'm not sure how to include the condition i+1 for selecting each polygon separately (something like for cycle). Probably it will be necessary to use Python Console but I'm not familiar with scripting in QGIS so far, so I'll need some hint.


texts: How can I zoom to 55% using arcpy in the layout view, instead of doing it manually with "Zoom to percent" button? I have 55 maps and I prefer doing it automatically with arcpy code.


texts: Difference between azimuth degrees when i use st_project in postgis

texts:
I have two points p1 and p2.

 SELECT ST_GeomFromText('POINT(92.78488163381019 55.9757874622428)',  4326) as p1,
        ST_GeomFromText('POINT(92.78482645521103 55.97564239375992)', 4326) as p2

degree between points equal 200.824
and distance between points equal 16.492

SELECT degrees(ST_Azimuth(p1, p2)) as degree,
       ST_Distance_Sphere(p1, p2)  as distance

But i have different degree when i use st_Project with this arguments:

SELECT degrees(ST_Azimuth(p1, _p2)),
       ST_Distance_Sphere(p1, p2)
FROM (
  SELECT (ST_Project(p1, distance, radians(degree))::geometry as _p2 FROM ...
)

It returns:

  st_distance_sphere  ||       degrees
======================||====================
  16.4661410575154    ||  214.150581571718

Distance difference between (p1,p2) and (p1,_p2) changing only by 0.02 meter. But degree changing by 13 degrees.

Why degrees are not equal between other?


texts: How to calculate raster area or number of cells for multiple classes in a raster?

texts:
I have a raster with many classes and I want to know the area or the number of cells for each class.
I found quite a few posts on this topic (see links below). None of the suggestions matched my wishes. I do not want to convert to vector since my raster is fairly big. I am looking for a tool similar to the summarize functions available for vector data, with a table as output.

The raster is a 32bit float tif with 5 unique values (1.0, 2.0, 3.0, 4.0 and 5.0) as classes as a result of classification with QGIS SAGA-Toolbox.

QGIS solutions are preferred. Related posts:

Calculating raster area?

Calculating area of vegetation in raster file?

Calculating area of rasters in QGIS?


texts:
I am having trouble initializing the GRASS plugin.I've tried un-installing and re-installing QGIS already. But the following message keeps appearing.

 GRASS init error : GRASS was not found in '/usr/local/opt/osgeo-grass/grass-base' (GISBASE), provider and plugin will not work.


texts: Is there a way using ArcGIS 10.3, when you have points for Surface Well and Bottom Hole... with API numbers, match up the API numbers and have a "LATERAL LINE" be drawn without hand drawing in each one?


texts: Calculating shape perimeter/area per zone in ArcGis 10

texts: How can I calculate geometrical features of polygon (perimeter and area) per zone (similar to zonal stats for rasters)? I have a square fishnet (approx 20x35, 1 square kilometer each) and I need to calculate forest area and boundary length per each zone.


texts: Image size posted from OpenLayers 3 to MapServer

texts:
I have the following map file

MAP
  PROJECTION
    "init=epsg:4326"
  END # PROJECTION

  FONTSET "H:/mapserver/data/fonts/fonts.list"
  EXTENT 19 33 30 42
  IMAGETYPE "png24"
  NAME "airways"
  SIZE 1280 1024

  WEB
    METADATA
      "wms_title" "Airways"
      "wms_onlineresource" "http://localhost/cgi-bin/mapserv.exe/wfs?"
      "wms_enable_request" "*"
      "wms_srs" "EPSG:4326"
      "wms_format" "image/png"
    END
  END #WEB

  LAYER
    NAME "airways"
    DATA "H:/mapserver/data/airways.shp"
    TYPE LINE
    LABELITEM "NAME"
      PROJECTION
        "init=epsg:4326"
      END # PROJECTION

    CLASS
      NAME "airways"
      LABEL
        COLOR 255 255 255
        TYPE TRUETYPE
        FONT Vera
        SIZE 7
        ANTIALIAS TRUE
        POSITION CC
        ANGLE AUTO

        STYLE
          GEOMTRANSFORM 'labelpoly'
          COLOR 0 0 0
        END # STYLE
      END # end of label

      STYLE
        COLOR 0 0 0
        WIDTH 0.25
      END # STYLE
    END # CLASS

  END # LAYER

END # MAP

and I'm using Openlayers 3 with following code

function Initialise() {

  var sourceMQ = new ol.source.MapQuest({layer: 'osm'});
  var layerMQ = new ol.layer.Tile({source: sourceMQ});

  var wmsSource = new ol.source.ImageWMS({
    url: 'http://localhost/cgi-bin/mapserv.exe?',
    params: {
      'SERVICE':'WMS',
      'map': 'H:/mapserver/data/airways.map',
      'LAYERS': 'airways',
      'VERSION':'1.3.0',
      'REQUEST':'GetMap',
      'FORMAT':"image/png",
      'TRANSPARENT':'true',
      'CRS':'EPSG:4326',
      'STYLES':'',
      'WIDTH':'1920',
      'HEIGHT':'1200',
      'BBOX':'19, 34, 29, 42'
    },
    serverType: 'mapserver'
  });
  var wmsLayer = new ol.layer.Image({
    source: wmsSource
  });

  var view = new ol.View({
    projection: 'EPSG:4326',
    center: [23.48, 37.55],
    zoom: 7
  });

  var map = new ol.Map({
    target: 'map',
    layers: [layerMQ, wmsLayer],
    view: view,
    controls: ol.control.defaults().extend([
      new ol.control.ScaleLine({units: 'nautical'})
    ])
  });

  map.on('pointermove', getCoordinates);
}

function getCoordinates(evt) {
  var out = ol.coordinate.toStringHDMS(evt.coordinate);
  $('#coordinates').text(out);
}

The html code is

Open layers test

The css file that I use is

.map
  {
     height: 800px;
     width: 1300px;
     border: red thin solid;
  }

  #coordinates {
    width: 100%;
    height: 25px;
  }

The request send to mapserver captured from Chrome developers tools is the following

http://localhost/cgi-bin/mapserv.exe?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&FORMAT=image/png&TRANSPARENT=true&map=H:/mapserver/data/airways.map&LAYERS=airways&CRS=EPSG:4326&STYLES&WIDTH=1950&HEIGHT=1200&BBOX=29.970703125,14.23828125,43.154296875,35.66162109375

The image size that is requested is 1950 by 1200 but the div that the map lies in is 1300 by 800. If the div width is increased then the requested image width exceeds 2048 which is the maximum image dimension for mapserver.
Also WIDTH and HEIGHT that are set in params in javascript are ignored.

So my question is how openlayers create the width and height request that is send to mapserver and how can I make openlayers send a request regarding width and height in correct limits.


texts: Visualising relationships in ArcGIS Desktop?

texts:
I have two polygon feature classes related through a key field. I can see the relationship through the attribute table but I'd like to visualise the relationship.

The ideal process would be if click on a polygon, I could run a tool that would highlight all other polygons related. I can workout that way to set it up though.

I'm using ArcMap 10.5.


texts: Import Aerial Images from Google into QGIS?

texts:
QGIS 2.6.1, Windows 7

I am looking to import an image from Google Earth into an existing QGIS project. Is it possible to download geo-referenced images from google maps (or google earth) and import into QGIS?
What CRS are these likely to be in?


texts: Python- Make use table from a MDB with name changes every week (name+date)

texts:
I have been trying to make a script to automate some steps in a process that uses tables from different databases.  One of them changes its name every week - name+date eg: "MegaDataBase_16062019.mdb"

How can I make use of the table in it?

I've tried arcpy.ListWorkspaces, something like:

 datain = arcpy.ListWorkspaces("\\hugeserver\mainwork\databases\MegaDataBase_*","Access")
 tableat = datain+".mdb\WorthyTable"

and that got an error, because the result in "datain" was:

  u"["\hugeserver\mainwork\databases\MegaDataBase_16062019.mdb"]

I was linking to

  u"["\hugeserver\mainwork\databases\MegaDataBase_16062019.mdb"]\WorthyTable

So, the wildcard* works, and I find the MDB I want, but can't make use of it to build the path I need

What I need from that table, is to export it with some fieldmappings to a .GDB on my machine

I used * or ? not quite sure now...tomorrow I can confirm it, and share more info if needed, as well as test your suggestions


texts: How to display label in symbol in ArcGIS?

texts:
As you see on the map below, there are labels in circle coloured symbols. How can I show them in this format in ArcGIS?

If you look at the legend of this map it shows the different categories with different size and color, so my request is a little bit complicated. And there are 3 different labels displaying on the map. How can I do this kind of map?


texts:
When I try to calculate a new field. I receive the following error:

idx = self.layer.fieldNameIndex(attributeName)
AttributeError: 'QgsVectorLayer' object has no attribute 'fieldNameIndex'

What´s the problem? I'm using QGIS 3.0.1 64x in windows.


texts: ECW to KML from ArcGIS Desktop failing to open in Google Earth?

texts:
i have some .ECW file and i want to convert to KML to used in another program (AutoCAD)..
so i convert this file using conversion tools which is ArcGIS Desktop 10.2 provide..
but i have problem in here..
when i just convert directly (without set the output image setting) it was successful the KML output size file is arround 500KB but the file is became blur because low resolution in this file..
so decide to set the resolution in the layer to KML tool..
here is my setting.

image size: 100,000
DPI: 100,000

that file is successfully converted the file size is around 300 mb.. but when tried to open it in google earth it just show square box and X red..

How do I do it?


texts:
I have a Leaflet map that is customized by an end-user by panning and zooming to an area of interest.
I'd like to then take all of the relevant parameters (projection, zoom level, bounds, etc) and submit them to a backend service that needs to project lat/long coordinates onto the screen coordinates (X, Y) of the Leaflet map.
(I do understand I could just ask Leaflet to project these points on the client-side, but this is not my goal. I do want to get the same answer as if I had done the projection on the client.)
Probably on the backend I will use something like pyproj. How can I generate a PROJ string or similar from the Leaflet map that conveys what the backend needs to do the projection?


texts:
I'm trying to view some Environment Agency (EA) WFS data in ArcMap 10.2.2. I've installed the Data Interoperability extension and enabled it. I've added a new connection with the url below. It should contain the EA's main river dataset. When I load the WFS, the layer appears in the table of contents, but contains no features.

EA Main River WFS URL: http://environment.data.gov.uk/ds/wfs?SERVICE=WFS&INTERFACE=ENVIRONMENTWFS--25dde009-ba7d-40de-8380-c5c3bb32ccdc&VERSION=1.0.0

(WFS available from http://environment.data.gov.uk/ds/catalogue/#/catalogue)


texts: Cannot change an SRTM raster rendering with single band pseudocolour

texts:
I'm trying to colorize a SRTM raster in QGIS 2.10 following a tutorial found on the Ireland QGIS User Group Blog. My problem is that when I try to set the render type as single band pseudocolour in the Style Preferences, nothing happens and the raster rendering remains unchanged. Then when I open again the Style Preferences, everything is default-set.

What am I doing wrong ?


texts: Calculate distance between polygons on two separate layers using ArcGIS for Desktop?

texts:
For my purposes, I need to find for each polygon in one feature dataset (red polygons) the outemost (from edge to edge) polygon in the second feature datasets (green polygons). Paradoxically, seeking the opposite function of the Near function in ArcGIS.

How would this be done?

om/en/arcmap/10.3/tools/analysis-toolbox/near.htm


texts: Grid-like lines appearing after using Curvature tool

texts: I have been attempting to use the Spatial Analyst Curvature tool with mixed results.  I originally tried to use it with a raster that was projected using a 1983 Teale Albers projection and the results were black, indicating the lowest values.  I then transformed the raster into UTM projection and ended up with grid-like lines on the raster which have values that are obviously not accurate.  Does anyone know why these lines have appeared or what I'm doing wrong.  Thanks.


texts: Raster to vector, matching another vector in QGIS?

texts:
I have taken a pdf of a map and managed to import it into my map using the QGIS geolocator. Now I want to convert the resulting raster layer to a set of polygons that match up with the existing vector layers. When I polygonize the raster layer, it's a mess.

Can someone get me started in the right direction?

Map

Map with raster (w/ transparency)

Map with polygonized raster

Or is it easier to use the raster overlay to note by hand which of the small areas fall within the the respective raster's areas?


texts: Mapfish printing application from Openlayers : disableScaleLocking

texts:
Mapfish printing application from Openlayers: I get the following error when I use a different scale than the scales defined in the config.yaml even though I am using  disableScaleLocking: true

  Error while generating PDF:
  org.mapfish.print.InvalidJsonValueException: spec.pages[0].scale has an invalid value: 65000

I thought that disableScaleLocking: true would enable me to print a map at any scale and not only the ones defined in config.yaml . Did I get it wrong?

Here are the first lines of my config.yaml :

dpis: [75, 150, 300]

scales:
   ... (list of predefined scales)

hosts:
    - !localMatch
        dummy: true
    - !ipMatch
        ip: 0.0.0.0
        mask: 0.0.0.0

disableScaleLocking: true

layouts:
      ...


texts:
I set up a batch file to launch pycharm on Windows in accordance with https://qtibia.ro/2016/06/09/how-to-setup-pycharm-for-qgis-development-under-linux-and-windows/

@echo off
SET OSGEO4W_ROOT=C:\OSGeo4W64
SET PYCHARM="C:\Program Files\JetBrains\PyCharm Community Edition 2018.2\bin\pycharm.exe"

call "%OSGEO4W_ROOT%"\bin\o4w_env.bat
call "%OSGEO4W_ROOT%"\apps\grass\grass-7.4.0\etc\env.bat
@echo off
path %PATH%;%OSGEO4W_ROOT%\apps\qgis\bin
path %PATH%;%OSGEO4W_ROOT%\apps\grass\grass-7.4.0\lib

set PYTHONPATH=%PYTHONPATH%;%OSGEO4W_ROOT%\apps\qgis\python
set PYTHONPATH=%PYTHONPATH%;%OSGEO4W_ROOT%\apps\Python27\Lib\site-packages
set QGIS_PREFIX_PATH=%OSGEO4W_ROOT%\apps\qgis

start "PyCharm aware of QGIS" /B %PYCHARM% %*

Pycharm recognizes the OSGeo4W64 python interpreter but qgis.core is still not importable.
These are the paths that are recognized by the interpreter:

I tried to change the batch file in accordance with: http://spatialgalaxy.net/2014/10/09/a-quick-guide-to-getting-started-with-pyqgis-on-windows/ and https://nathanw.net/2014/05/10/setting-up-pycharm-for-pyqgis-and-qt/ but then the OSGeo4W64 python interpreter was not recognized. I also tried adding additional paths and installing available packages through pycharm but so far nothing has worked. The error message that I'm receiving when running my code reads as follows:


texts: How to access and modify properties of one single layer with ArcPy

texts:
In General
I have to access and modify a specific property, e.g. data source, of one specific layer with arcpy. I know the name of the mxd-file, the dataframe and the layer. How can I access, for modification, this one single layer directly without the common for layer in arcpy.mapping.ListLayers(...) loop as suggested in other answers (e.g. in the manual or in case 18433)? I also don't manage to use UpdateLayer: this command allocates properties by fetching them from another layer, however, I don't have another layer which could deliver the desired properties.
Example
import arcpy, os
workspace_mxd = r"L:\project\home"
workspace_old = r"L:\project\scenario1"
workspace_new = r"L:\project\scenario2"
mxdfile = os.path.join(workspace_mxd,'mxd','map.mxd')
datasource = os.path.join(workspace_old,'data','points.shp')
layername = 'village'
arcpy.env.workspace = workspace_mxd
arcpy.env.overwriteOutput=True
mxd = arcpy.mapping.MapDocument(mxdfile)
df = arcpy.mapping.ListDataFrames(mxd)
for dataframe in df:
    for lyr in arcpy.mapping.ListLayers(mxd, layername, dataframe):
        lyrDatasource = lyr.dataSource if lyr.supports('dataSource') else 'N/A'
        if lyrDatasource == datasource:
            lyr.findAndReplaceWorkspacePath(workspace_old,workspace_new,False)
mxd.save()
del mxd

By this approach, the property data source is modified for one specific layer. The approach looks inflated for this purpose. E.g. the two for loops are, in the specific case, for one dataframe and one layer only. Isn't there a shortcut to do what is finally done with lyr.findAndReplaceWorkspacePath(workspace_old,workspace_new,False)?


texts: Set up the center of a map with Google Layer as basemap in OpenLayers and GeoExt?

texts: I looked on similar questions but in vain. The map is working fine but has a center of (0,0) as returned in the Firebug console not (-8, 32) that I want. I did the projection transformation.. Here is my code :


texts:
I have 4 raster layers contained in a raster stack within program R. The layers represent occurrence surfaces for two species in two seasons (i.e. species one summer and winter and species 2 in summer and winter). Values for all raster layers range from 0 to 10 and represent a relative ranking of habitat, worst to best, respectively. For each layer I want to make polygon/s that encompass raster values greater than a specified number, for example 8.

Within a given layer there is a patchy distribution of 'good' habitat (i.e. > 8) which will results in multiple polygons for each layer.

This process is similar to creating animal home ranges from kernel densities and I was surprised to not see previous posts on the topic, although I may be searching the wrong terms...

My hope is to then calculate the overlap of seasonal polygons among species (i.e. overlap of species 1 and 2 in summer), likley using intersect from the raster package, although I welcome any suggestions.


texts: Creating Minimum Convex Polygon - Home Range from Points in QGIS

texts:
I just downloaded QGIS (3.0.2) and trying to create minimum convex polygon around points using minimum bounding geometry (convex hull). The videos and help documents I see online appear to be from older editions of QGIS. I cannot find comparable menu options in this edition.

Can this be done in this version of QGIS?

If so, I seek simple instructions (I'm new to GIS).


texts: Network Quantity Penalized of Distance in sDNA

texts:
I'm following up on a question asked on the calculation of the Network Quantity Penalized of Distance (NQPD) with sDNA.

I can calculate the value for the link asked in Why my calculated NQPD value is different than sDNA calucated?; e.g. the links where the value is 0.031481.

However, for the adjacent values I obtain a different value when calculating the NQPD by hand.

For instance, for the links with a NQPD of 0.024074 I obtain:

1/(180+90) * 4 (four links with 90 degree turn)
1/(180+0) * 1 (one link with 0 degree 'turn')
1/(180+0/3) * 1 (self-contribution)
 = 0.025925926

To obtain the value calculated by sDNA I need to adjust the self-contribution of the link to

1/(180+0)*2/3

Is the self-contribution of these links calculated differently?


texts: Python, GDAL and building raster attribute tables

texts:
I have an integer raster for which I would like to build a raster attribute table using Python and GDAL. I can create a GDAL raster attribute table in Python as follows:

>>> rat = gdal.RasterAttributeTable()

This works fine, as we can see:

>>> rat
<osgeo.gdal.RasterAttributeTable; proxy of <Swig Object of type 'GDALRasterAttributeTableShadow *' at 0x0000000002A53D50> >

The table thus created has no rows or columns:

>>> rat.GetRowCount()
0
>>> rat.GetColumnCount()
0

I create a column called "Value" to store the unique values in the raster:

>>> rat.CreateColumn("Value", gdalconst.GFT_Integer, gdalconst.GFU_MinMax)
0

This is fine, and the column count is updated:

>>> rat.GetColumnCount()
1

Now I have to add values (records) to the column for it to be of any use. I can get a list of unique values from the raster's band like so:

>>> data = band.ReadAsArray(0, 0, dataset.RasterXSize, dataset.RasterYSize)
>>> vals = list(numpy.unique(data))
>>> vals
[3, 7, 8, 10, 11, 12, 13, 14, 17, 18, 20, 22, 23, 25, 27, 28, 41, 45, 52, 56]

What I would like to do is create a for loop to loop through vals and populate the column in the attribute table. I thought I could do something like this:

>>> for i in range(len(vals)):
        rat.SetValueAsInt(i, 0, vals[i])

...where i is the row (record), 0 is the field index and vals[i] is the integer value I want to insert. But it causes an error:

Traceback (most recent call last):
  File "<pyshell#32>", line 2, in <module>
    rat.SetValueAsInt(i, 0, vals[i])
  File "C:\Python27\lib\site-packages\osgeo\gdal.py", line 1139, in SetValueAsInt
    return _gdal.RasterAttributeTable_SetValueAsInt(self, *args)
TypeError: in method 'RasterAttributeTable_SetValueAsInt', argument 4 of type 'int'

The error is caused because I use vals[i] in the call to SetValueAsInt() rather than using an integer directly. For example, rat.SetValueAsInt(0, 0, 0) works fine, but is useless for populating the column if I just want to loop over the list of unique values.

Is this a known issue? Google has so far not been very useful. What can I do to get around this problem?


texts: Can I install GIS Software on RSS (informix server) running 11.7.fc6

texts:
I currently use ArcGIS Server 10.0 and consume data through a db2 database. I recently installed an RSS informix server and was hoping to pull data from this since it is updated real-time versus the nighly update of the db2.

RSS Informix Server Version: 11.7.fc6 OS: AIX 6.1
GIS Server: 10.0

If I stay at Server 10.0, I would have to install ArcSDE (I assume) on the Informix Server. Can SDE be installed on an Informix Server running 11.7.fc6 OS: AIX 6.1?

**
Edit: I see ArcSDE 10 isn't compatible with 11.7.FC6.
But is it possible to install ArcSDE 10.1 on Informix (assuming Informix version is supported by 10.1) while still remaining at 10.0 Server with capability of publishing the Informix data as a service on Server 10.0?
**

If I upgraded to 10.1 Server, the Informix capability is already there without installing SDE, BUT this link shows support for only FC8 and FC4. How would FC6 react or would I need an upgrade to FC8?

http://resources.arcgis.com/en/help/system-requirements/10.1/index.html#//015100000071000000


texts: Extract polygon pixels using lazy reading in rasterio

texts:
I have several hundred tif files of ~50mb each where I need to read only a few hundred pixels within a polygon. I currently do this by loading each tif into memory and using rasterio.features.geometry_mask to get the pixels of interest. But reading in the large tifs takes a very long time.

import fiona
import rasterio

with fiona.open('polygon.geojson') as f:
    geom=f[0]['geometry']

all_data=[]
for tif_filename in tif_file_list:
    with rasterio.open(tif_filename) as raster_obj:
        polygon_mask = rasterio.features.geometry_mask(geometries=[geom],
                                               out_shape=(raster_obj.height, raster_obj.width),
                                               transform=raster_obj.transform,
                                               all_touched=False,
                                               invert=True)

        all_data.append(raster_obj.read()[polygon_mask])

Is there any way I can use the Window option to only read the portion of the raster I want? For example I can read just the extent of the area of interest:

with fiona.open('polygon.geojson') as f:
    b = f.bounds

polygon_window = rasterio.windows.from_bounds(left=b[0],
                                              bottom=b[1],
                                              right=b[2],
                                              top=b[3],
                                              transform = raster_obj.transform)

extent_data = raster_obj.read(1, window=polygon_window)

This takes a fraction of the time, but I can't figure out how to then get the pixel data from only within the polygon.


texts: Categorising raster data based on shapefile boundaries

texts:
I've tried searching for an answer but honestly I'm not even sure how to word my search.
I'm using Arc GIS Pro. I have two layers that I want to analyse.
The first layer is a shapefile and it covers three countries, split up into regions. So, it records the boundaries between those regions. The attribute table has fields for Country and Province.
The second layer is an ESRI grid raster. It contains population info across the study region. I have already clipped it to match the area of the first layer.
I'd like to find the population of each Province from the first layer. So I guess I want a new layer, prov_pop, where the attribute table has: Country, Province, Population. Ultimately, I want to create a map that shows the province boundaries and is colour-graded to population size at a provincial level.
I'm new to GIS and really not sure where to start to produce this.


texts: Query the distance to the nearest location through a join table

texts:
I'm looking to select the distance from a given location to the nearest location in a joining table, and order by it. In my use case I have a listings table which is joined to locations table through listing_location. The locations table has two columns for coordinate - latitude and longitude (both doubles).
SELECT
  listings.*
  ST_Distance(ST_Point(?, ?), ST_ClosestPoint(:unsure:, ST_Point(?, ?)), FALSE) AS distance)
JOIN listing_location
  ON listing_location.listing_id = listings.id
JOIN locations
  ON locations.id = listing_location.location_id
ORDER BY distance ASC

The question marks would be bound to the longitude and latitude of the given location, but I'm unsure what to pass into :unsure: in order to select the joined locations. Am I missing something, or thinking about this problem in the wrong way? Is this problem solvable when my location data is stored as latitude/longitude double columns?


texts: R giving error when attempting map algebra

texts:
I'm trying to manipulate some downloaded rasters using R, and keep getting the same error when I try to do any type of map algebra.

The example below uses a raster file downloaded from http://www.csc.noaa.gov/digitalcoast/data/ccapregional/ geotiff and .grd files, regardless of whether I downloaded them or generated them previously.

require(raster)
test1=raster("nh_2006_ccap_land_cover.img")
test2=test1+1

Gives the error:

Error in matrix(unlist(ini), ncol = 2, byrow = TRUE) :
  'data' must be of a vector type, was 'NULL'
In addition: There were 11 warnings (use warnings() to see them)

On the other hand, I'm not encountering problems with map algebra for rasters generated on the fly.  So I'm not having any problems with this example from the documentation:

r1 <- raster(ncols=10, nrows=10)
r1[] <- runif(ncell(r1))
r2 <- setValues(r1, 1:ncell(r1) / ncell(r1) )
r3 <- r1 + r2
r2 <- r1 / 10
r3 <- r1 * (r2 - 1 + r1^2 / r2)
plot(r3)

What's going on here?


texts: Printing a degree symbol ("º") yields "unknown symbol" (question mark in romb) via Python

texts:
I am working with QGIS 2.14.3, I just create a JSON file for detecting text, numbers, bearing and distance with OCR process using python 3 with pytorch, successfully I completed for creating a JSON file. When I load my JSON file
the degree symbol is not displaying instead, it produces like in below snap:

Even I used # -*- coding: utf-8 -*- in my code still getting same error


texts: MapInfo Styles in Oracle Spatial Database

texts:
I have geographic data in an Oracle database. The data was imported by EasyLoader, which is a MapInfo-Tool. I want to publish this data in a Web Map Service. I have tried to use MapServer, but I have problems handling the styles.

I know that the in MapServer, it is possible to get the correct styles by OGR, if they are in a TAB-File.

CONNECTIONTYPE OGR
CONNECTION "Cables.TAB"
STYLEITEM "AUTO"

The same data is in the database and I created an new layer section:

CONNECTIONTYPE PLUGIN
PLUGIN "msplugin_oracle.dll"
STYLEITEM "AUTO"

An error occurs

General error message. 'STYLEITEM AUTO' not supported for this data source.

If I define a manual style with

CLASS
    STYLE
        COLOR 0 0 0
        WIDTH 2.5
    END
END

the lines are shown.

So here is the question: Is there any possibility to use the AUTO StyleItem with an oracle spatial database to decode the mapinfo style?

Or do you have any suggestions, how I could use this information for styling my items?


texts: Error message received when clipping raster

texts:
QGIS 2.10.1 Pisa

Message -> Error: Computed -srcwin 8208 142185 0 0 has negative width and/or height.

Tried with several different files and types of raster files [GeoTIFF, MrSis, ECW], the same msg with different numbers.


texts:
My layers are present in postgis db which are being served to openlayers JS through geoserver. I have my layers corresponding to polygons,lines etc. whose attributes are stored in the same table itself. Like Halls table contain columns (hid, hall_geom,hall_name) etc. I want hall_name of each hall_geom(polygon) to be placed as a label on it in the map. I have tried SLD-styling form geoserver but couldn't make it successful(I think the styles given are for shape files and not tables in db).

Please show me a right way


texts: Arc ToolBox Display paramter inputs after user stops script

texts: I am using arcpy 10.0, python and arctoolbox. I have a question regarding the parameters for the arc toolbox. In my script I have a function that does a check on an input file. The user is charged with selecting the appropriate input file. Say, if they select a shapefile instead of a text file, they get an error. A box pops up saying they chose the wrong input file. At that point, the script just stops and the input parameters box has alos gone away. Is there a way to make the input parameters box appear again on erroneous user input? So instead of running the program again on the users error, the input parameters box automatically reappears so the user can enter in a correct input file. Any knowledge on the subject would be greatly appreciated. Thanks.


texts: Create adjacent buffers around grid points

texts:
I've a file of grid points (ERA5 meteorological info grids) and I'd like to convert the points to adjacent polygons so I can create nice maps from them.
The resolution of my grids is around 30 km x 30 km, however unfortunately it's not unique, neither in km, neither in degrees.
I could calculate distance of grid points in another program and add desired height & widths as new column to apply them as buffers' dimension, however I'd also be interested in a solution completely done in QGIS.
Sample data:

location
lat
lon

38.6°N/27.8°E204.9m
38.64167
27.77344

38.6°N/28.1°E357.3m
38.64167
28.125

38.6°N/28.5°E462.3m
38.64167
28.47656

38.6°N/28.8°E619.0m
38.64167
28.82812

38.4°N/27.8°E630.4m
38.36064
27.77344

38.4°N/28.1°E795.8m
38.36064
28.125

38.4°N/28.5°E651.6m
38.36064
28.47656

38.4°N/28.8°E663.0m
38.36064
28.82812

38.1°N/27.7°E328.0m
38.07961
27.66667

38.1°N/28.0°E605.1m
38.07961
28

38.1°N/28.3°E747.2m
38.07961
28.33333

38.1°N/28.7°E741.9m
38.07961
28.66667

37.8°N/27.7°E176.4m
37.79858
27.66667

37.8°N/28.0°E188.5m
37.79858
28

37.8°N/28.3°E437.6m
37.79858
28.33333

37.8°N/28.7°E655.0m
37.79858
28.66667

Desired output (Example from another file with 0.1 x 0.1 degree resolution):


texts:
My task is to calculate the area available for specific plant species in a natural reserve that has been reconnected to the tidal influence.

First step is to create an inundation map for a natural reserve.
I have a height map in .dxf format (annotation,point, polyline,polygone and multipatch layers). I can see the height of the points in the attributes layer.
I have available the data (for one year) for the tide levels from a connected tide gauge for the fluctuating water level. I'm now combining this with the plant preferences for inundation time. I expect to get a result like for species 1 the height between 1.5 an 1.7 meter is suitable, for species 2 between 1.7 and 1.8 and so on.

How do I implement this in the map and find the total suitable area for each species? I am running Arcgis 9.3 with the Spatial Analyst and 3D Analyst extensions.


texts: How to route on OSM without using any already released routing engine?

texts:
I'm a computer science student and a complete beginner in GIS and stuff. I'm trying create my own web-based multi-modal journey planner as my thesis. And as it is a thesis I can't use any ready-made routing engines like graphserver or opentripplanner, etc. though I can use third party applications for data conversions, like map data into SQL or KML, etc. I've done some research and this what I've learned so far.

that I can download map data from OpenStreetMap (OSM) and convert them into SQL like SQLlite or PostgreSQL for querying. (in RTree index structures)
that I can also convert GTFS data into SQL, again for querying.
that I can link those two data sets through stops (those loading and unloading stations for each mode of transit).
that I should convert those two data sets into a graph then link them in order to perform multi-modal routing. (I don't know if there are other solutions)
that I can convert the resulting route into KML so I can overlay/display them on Google Maps or other API's, for presentation. (forgot where I read this, but I know I read about this somewhere)

My question is how can I convert an OSM database into a graph data... How can I perform routing on those graph data... this is where I got lost, I still don't know how can I perform routing algorithms from those given data sets as input/source. I tried reading the source code of graphserver, but it did no good to me... besides it needs to be in Linux system to work, and I'm currently working on Windows... I have few understanding on Java... but if I have to, I can learn it, if it is what it takes so I can grasp the steps in performing routing. (as I see.. most of the routing engines are written in Java)


texts: How to skip null images in Google Earth Engine?

texts:
I've called sentinel2 image collection but I don't know how can I skip null images instead of filter.neq():
code link: https://code.earthengine.google.com/b762152c33e6eddca96c2006e62eed54
var s22 = ee.ImageCollection("COPERNICUS/S2")
.filterDate(start_date,end_date)
.filterBounds(table)
.filter(ee.Filter.lessThan('CLOUDY_PIXEL_PERCENTAGE',45))
.filter(ee.Filter.neq('PRODUCT_ID','S2A_MSIL1C_20180518T071621_N0206_R006_T39SUR_20180518T093201'));


texts:
I am interested in how to create a base map for Leaflet offline with XAMPP.

Is this just a picture of what you draw from the Internet and then specifies the path with L.tileLayer?

I'm still a bit new about Leaflet and JavaScript


texts:
For a project I'm working on, I need to represent network flows over 1,000 zones with directionality.

Any suggestions on displaying the information in a nice and non overwhelming way?
Any examples on displaying network flows in heavy urban areas?
Color suggestions and thickness schemes?
any "not to dos" when representing network flows?


texts:
I created a yaml config file and ran it with mapproxy-util serve-develop my-config.yaml
The file contains this:

services:
  demo:
  tms:
  kml:
  wmts:
     capabilities: http://localhost:8080/wmts/1.0.0/WMTSCapabilities.xml
  wms:

layers:
- name: photographiesaeriennes
  title: Photographies aériennes
  sources:
  - ign_cache
caches:
  ign_cache:
    grids:
    - ign_grid
    sources:
    - ign_source
sources:
  ign_source:
    type: tile
    grid: ign_grid
    url: https://wxs.ign.fr/pratique/geoportail/wmts?SERVICE=WMTS&REQUEST=GetTile&VERSION=1.0.0&LAYER=ORTHOIMAGERY.ORTHOPHOTOS&STYLE=normal&FORMAT=image/jpeg&TILEMATRIXSET=PM&TILEMATRIX=%(z)s&TILEROW=%(y)s&TILECOL=%(x)s
grids:
  ign_grid:
    base: GLOBAL_MERCATOR
    srs: EPSG:3857
    origin: nw
    res:
    - 156543.0339280410
    - 78271.5169640205
    - 39135.7584820102
    - 19567.8792410051
    - 9783.9396205026
    - 4891.9698102513
    - 2445.9849051256
    - 1222.9924525628
    - 611.4962262814
    - 305.7481131407
    - 152.8740565704
    - 76.4370282852
    - 38.2185141426
    - 19.1092570713
    - 9.5546285356
    - 4.7773142678
    - 2.3886571339
    - 1.1943285670
    - 0.5971642835
    - 0.2985821417

and when I run it I get this error :

mapproxy.config - WARNING - unknown 'capabilities' in services.wmts

I can still use this proxy but it loads very slow compared to the source and
while it runs in QGIS I get many [Errno 32] Broken pipe messages.

Did I write my config file wrong ?


texts:
I am having a little trouble debugging an ex-coworkers code.

I am receiving this error:

RuntimeError: FieldMappings: Error in adding table to field mappings

The code:

CIP_FDFC_Tuple =  [('Cadastre', 'Installation', 'Cadastre\\Installation'),
              ('Cadastre', 'LandParcel', 'Cadastre\\LandParcel'),
              ('General', 'Wall', 'General\\Wall')]

env.workspace =r"C:\\Users\\1468394150E\\AppData\\Roaming\\ESRI\\Desktop10.1\\ArcCatalog\\CIP-GIO FUSED.sde" #sys.argv[1]
inputgdb = env.workspace

#Loop through CIP datasets and feature classes
#env.workspace = inputgdb
for (a,b,c) in CIP_FDFC_Tuple:
    inFeature = inputgdb + "\\" + c
    arcpy.AddMessage( "inFeature: "+ inFeature)

    #Remove GlobalID within SDE feature tables
    fieldmappings = arcpy.FieldMappings()
    fieldmappings.addTable(c)
    fieldmappings.removeFieldMap(fieldmappings.findFieldMapIndex("GLOBALID"))
    #Copying Feature class as a layer into fused gdb
    arcpy.AddMessage("Copying " + c + " feature layer records")
    arcpy.FeatureClassToFeatureClass_conversion( inFeature, tempWorkspace, c,  whereClauseAll)


texts:
I recently did an sudo apt dist-upgrade on two Ubuntu 16.04 computers. I didn't pay especially close attention to the messages and now QGIS is gone. So I do
sudo apt install qgis and I get qgis : Depends: python-qgis (= 2.14.17+dfsg-1~xenial1) but it is not going to be installed

Any idea where my problem might be?


texts: Recommendations on an MVVM framework for ArcGIS Desktop add-ins and/or ArcGIS WPF Runtime

texts:
Is anyone using an MVVM framework with ArcGIS Desktop add-ins or the ArcGIS Runtime for WPF?

I am leaning towards Caliburn.Micro simply because it seems to be the most popular "lightweight" MVVM framework, but I was wondering if anyone had any experience with it or others (or if they are using none at all) before diving in.

My main requirements would be that it work with Visual C# 2008 Express, allow for targeting .NET 3.5 (and easily transition to .NET 4.0), allow for WPF Windows to be created, displayed as a well-behaved modal child dialog, and destroyed inside the ArcGIS Desktop applications, and provide logging, file open/save dialog, and message box services.

Looking for light-weight, minimal dependencies, promotes best design practices, etc. Everything else would be pretty much icing on the cake.

Or does everyone use Windows Forms still? :)


texts: Is there a "SnapToLevels" property in ArcGIS Runtime SDK for .NET?

texts:
In ArcGIS Runtime SDK for WPF, I can set SnapToLevels as True, so when I zoom the map, it will snap to scale levels if the map contains cached layers.

However, I can't find anything like SnapToLevels in ArcGIS Runtime SDK for .NET. Is there a property like that?

If so, will it work on touch device such as Surface Pro 3?  (SnapToLevels in ArcGIS Runtime SDK for WPF doesn't work when pinching with fingers.)

I've posted this question to GeoNet, but no one has answered until now.


texts:
I have a combination method to combine two raster( raster x and y) .
the method is : x+y / 1- min(|x|,|y|)
my problem is about "," character .whats the "," means in the method and  how can i use this method in raster calculator ?


texts: Opening files with ".int" extension using ArcMap

texts:
I'm trying to perform some experiments on datasets from sites of Coleambally Irrigation Area (CIA: http://dx.doi.org/10.4225/08/5111AC0BF1229) and Lower Gwydir Catchment (LGC: http://dx.doi.org/10.4225/08/5111AD2B7FEE6). However, the downloaded Landsat and MODIS data have extensions of '.bil' or '.int' without '.hdr' files. Therefore, I tried to create '.hdr' file according to the provided data description as (for LGC landsat):

NCOLS 3200
NROWS 2720
NBANDS 6
NBITS 8
LAYOUT BSQ
BYTEORDER I
SKIPBYTES 0
MAPUNITS DECIMAL DEGREES
ULXMAP 148.8815
ULYMAP -28.7455
XDIM 0.00025
YDIM 0.00025

After I created the '.hdr' file, data can be opened in ArcMap, however, they appeared to be problematic as:

Does anyone have an idea of what I did wrong? How can I correctly open such file?


texts:
I'm looking for some functionality beyond what ESRI already provides in the JSAPI.  For example, the ability to display the LOD in a tooltip on the zoom control.  It seems like there should be a dojo+jsapi answer to openlayers+geoext but I have not found it.

Can anyone point me to Any Open Source Extensions to ESRI JSAPI?


texts: Create a column in one attribute table based on a column in another attribute

texts: In QGIS, I have a river shapefile that does not include the basin that is in. I have a separate shapefile that does. Is there a way to get that river shapefile to have a new column of basin data in the attribute table based on the separate shapefile's basin? The data in the separate shapefile may not have the same length as the shapefile in the river shapefile and the only way to link the two shapefiles is based on geographic location.


texts:
I have a raster file and I would like to save it to be able to open it with a attribute table in ArcMap.
When I use the command writeRaster(raster, "test_output11", format = "GTiff"), I don´t have the aux files such as .tfw, .xml, .ovr, .cpg, .dbf.
Should I produce all the files manually by changing format = "GTiff"? I can ´t find a general command for it.


texts: Using Raster>Analysis>DEM to display hillshade produces the message "gdaldem.exe had stopped working" so I asked the community. I'm new to QGIS, and I now realise after a little experimentation that I should have used Raster>Terrain Analysis>Hillshade which now works OK - just in case another newbie has encountered the same problem. Why the former action should crash gdaldem.exe remains a slight concern though.


texts: Monthly Average Solar Maps (R.Sun) in Grass

texts: I'm trying to create monthly average maps as described in R.Sun in Grass. I have my DSM. In reading the page on OSGEO it mentions that it is possible without showing how. Ideas? I assume I will need to write a script.


texts: Showing progress of raster geoprocessing tools with ArcPy?

texts:
Background: I am doing DEM resampling for a pretty big raster (341020x341020 cells, ~ 10 GB file size, ~ 433 GB uncompressed size) with arcpy. My Python script (uses 64-bit Geoprocessing) has been running for more than half a day with the output file size is ~ 20 GB so far and I have no idea how long I will have to wait for it to finish.

import arcpy
arcpy.env.workspace = r"C:/GIS"

inRaster = 'dem_2m.tif'
resampleRaster = 'dem_1m.tif'
newCellSize = '1'
resamplingMethod = 'BILINEAR'

arcpy.Resample_management(in_raster=inRaster, out_raster=resampleRaster,
                          cell_size=newCellSize, resampling_type=resamplingMethod)

Question: Is there a way to show the percentage of completion when the script is running? Something similar to this when we use ArcToolbox within ArcGIS interface

System info: Windows 10 Enterprise, ArcGIS 10.6.1 with 64-bit Geoprocessing tools, Core i7-8850H (12 cores), 32 GB RAM, 300 GB free SSD disk space


texts: How to calculate mean of GPS data assigned to 30 m segment in QGIS?

texts:
I have collected some temperature data from a Google Car dataset with each data having Lat/Lon and temperature values. The data is stored in CSV file. I already fixed the GPS lat/lon by using "project" and "interpolate" using this shapely Snap points.

I also have a shapefile of the road the car took divided into 30 m segments.
What I want to know, how can I snap/group GPS data to these 30 m segments in QGIS?
I am using QGIS 3.16 and all data is projected in ESPG 4326.


texts: Free API or database for Continent, country, city, region/state?

texts:
Is there a database or API that can provide: Continent, country, city and region/state? It should also have the country codes in ISO format. It would also be great if if could all be in one file for easy access. Ideally, the API/database shouldn't have too many fields as I would prefer concise data. It would be a good bonus if it included support for different languages.

I found that GeoNames contains too many places and landmarks. Same with Yahoo GeoPlanet. GeoDataSource is too scattered between different files. I prefer something that is not too big and concise and possibly with an hierarchical structure.

Actually I found one here!!!!

http://dev.maxmind.com/geoip/geoip2/geolite2/


texts: Model Builder not completely honoring lyr file

texts:
I have created a series of map templates for one of our departments that each have an associated model.

As part of the model I also created a lyr file for symbology and labeling (the subject parcel layer portion).

This is referenced in the output feature parameters of my model.

The tool runs correctly and the Subject Parcel_layer gets added to the map template. The overall symbology from the layer file is correct, however, the labeling preferences I've set up are not (font size, type, etc.)

As a side note, the Subject Parcel Layer does not exist in the template, it gets added through the model.

Any ideas on why the labeling portion is not being honored?


texts: Semtech LR1110 GPS coordinates from GNSS

texts:
I am thinking of using Semtech LR1110 for designing tracking device in field with LoRa. Since this one has GNSS & LoRa on SoC, so going with LR1110.

But I have a concerns about LR1110.

How can I use GNSS data to get accurate GPS location ?
Is there any open source (free) way to get GPS coordinates from GNSS data ?


texts: Label layer as symbol colors with QGIS 3

texts:
Is it possible to label a layer using same colors as the categorized symbology of that layer?

It is possible using the Rule-based labeling, but it can be difficult for a layer with many categories and maybe it can be done using the Data-defined override, but I don't know how to get the color used for every category.


texts: QGis through python script : no map exported (map will be printed here)

texts:
I am very new to GIS in general and even more to QGis. I have a pretty simple task which is to import all shapefiles in one folder, style them and finally export them.

I have wrote python code that does it one by one just well from the console but when I try to launch it from a pyton class, I do not obtain the shapefile layer printed as I did from the console. Only the elements such as the legend, label appear and a big "Map will be printed here" !

I guess it is because I'm not dealing explicitly with the shapefile I just imported but I didn't find examples doing so on the documentation...

Here is my class so far :

from glob import glob
from os import path
from qgis.core import *
from qgis.gui import *
from PyQt4.QtCore import *
from PyQt4.QtGui import *
from PyQt4.QtXml import *
import PyQt4

    class Loader:
        def __init__(self, iface):
            """Initialize using the qgis.utils.iface
            object passed from the console.

            """
            self.iface = iface

        def load_shapefiles(self, shp_path):
            """Load all shapefiles found in shp_path"""
            print "Loading shapes from %s" % path.join(shp_path, "*.shp")
            shps = glob(path.join(shp_path, "*.shp"))
            for shp in shps:
                (shpdir, shpfile) = path.split(shp)
                layer = self.iface.addVectorLayer(shp, shpfile, 'ogr' )
                layer.loadNamedStyle('C:\FOLDER\USER\shape\\shapestyle.qml')

                #Create composition
                mapRenderer = self.iface.mapCanvas().mapRenderer()
                c = QgsComposition(mapRenderer)
                c.add
                x, y = 0, 0
                w, h = c.paperWidth(), c.paperHeight()

                composerMap = QgsComposerMap(c, x,y,w,h)
                c.addItem(composerMap)
                composerLabel = QgsComposerLabel(c)
                composerLabel.setText("Hello world")
                composerLabel.adjustSizeToText()
                c.addItem(composerLabel)

                legend = QgsComposerLegend(c)
                legend.model().setLayerSet(mapRenderer.layerSet())
                c.addItem(legend)

                item = QgsComposerScaleBar(c)
                item.setStyle('Numeric') # optionally modify the style
                item.setComposerMap(composerMap)
                item.applyDefaultSize()
                c.addItem(item)

                composerLabel.setItemPosition(20,10, 100, 30)
                dpi = c.printResolution()
                dpmm = dpi / 25.4
                width = int(dpmm * c.paperWidth())
                height = int(dpmm * c.paperHeight())

                # create output image and initialize it
                image = PyQt4.QtGui.QImage(PyQt4.QtCore.QSize(width, height), PyQt4.QtGui.QImage.Format_ARGB32)
                image.setDotsPerMeterX(dpmm * 1000)
                image.setDotsPerMeterY(dpmm * 1000)
                color = PyQt4.QtGui.QColor(155, 155, 155)
                image.fill(color.rgb())

                # render the composition
                imagePainter = QPainter(image)
                sourceArea = QRectF(0, 0, c.paperWidth(), c.paperHeight())
                targetArea = QRectF(0, 0, width, height)
                c.render(imagePainter, targetArea, sourceArea)
                imagePainter.end()

                image.save("C:\FOLDER\USER\\outfromscript.png", "png")


texts: Updating missing values in one table based on values in another

texts:
I am trying to accomplish something that seems very straight forward but I can't seem to figure it out.  I have two stand alone tables.  The first is an output report containing a field of land use codes.  The other is a lookup table containing definitions of those codes.  I need to update a field in the lookup table with any codes that are present in the report but missing in the lookup table.  See Code below.  I am getting stuck on trying to retrieve the values from the search cursor.

import arcpy
import pythonaddins
import csv, arcpy, os
from arcpy import env

table = "S:\\Zone_GP_Update\\Zoning.gdb\\CrystalReportViewer"
field = "Gen_Plan_Primary"
table2 = "S:\\Zone_GP_Update\\Zoning.gdb\\GenPlanLU"
field2 = "LANDUSECOD"
companyList = []

with arcpy.da.SearchCursor(table,field)as cursor:
    companyList = [row for row in cursor]

with arcpy.da.InsertCursor(table2,field2) as tablecursor:
    for tablerow in tablecursor:
        if tablerow[0] in companyList:
            continue
        else:
            tablerow[1] = row.getValue(field)
        cursor.insertRow(tablerow)


texts: How to set layers, scale, center automatically and default zoom based on JSON file config?

texts:
I want to know how I can init my OpenLayers 3 Configuring the Viewer by JSON file (layers, scale, center and default zoom) at the start of my application from a JSON configuration file.

My app.js:

var ol = require('openlayers');
var jQuery = require('jquery');

map = new ol.Map({
  target: 'map',
  controls: ol.control.defaults({
    zoom: false,
    attribution: false
  }),
  layers: [
    new ol.layer.Tile({
      source: new ol.source.OSM()
    })
  ],
  view: new ol.View({
    center: ol.proj.transform(
      [3.6656, 45.9192],
      'EPSG:4326',
      'EPSG:3857'
    ),
    zoom: 6
  })
});


texts: QGIS loads without its python module, postgreSQL, and plPython3u

texts:
In the (C:\Program Files\PostgreSQL\9.6) the "installation-notes.htm" for postgreSQL 9.6 show python 3.3 as being the compatible version.  I installed it, put it into my windows 7 path/variables, and can run python 3.3 at the command line.  After issuing a CREATE LANGUAGE plpython3u, I was able to write plpython3u functions.  But, QGIS uses python 2.7.5.  When I load QGIS, it will no longer load its python console because it explains that python33.dll is not compatible.  I suppose I could try reverting to python 2.7, in windows 7, but is that supposed to be compatible with postgreSQL 9.6?

Combining the advice of the two responses, I edited my environmental variables to remove any mention of python3.3.  QGis loads, and with the python console! I am not sure whether the python 3.3 environmental variables may be needed in order for plpython3u functions to work in postgreSQL.


texts: SQL 2012 - How can I convert MultiLineString to Geography?

texts:
I have a Geometry LineString that I wish to convert to Geography.
Since the Geometry isn't valid , When I use MakeValid on it, SQL convert it to MultiLineString and than the Geography creation is failing with error :

The specified input does not represent a valid geography instance. Use MakeValid to convert the instance to a valid instance. Note that MakeValid may cause the points of a spatial instance to shift slightly.

What can I do in order to convert the Geometry to Geography ?

Geometry text : http://pastebin.com/Ly8ex41U

Geometry text after MakeValid : http://pastebin.com/HGNFKnS3


texts: One popup for multiple layers in ArcGIS JS API

texts:
I am trying to have multiple feature layers that all share one popup in an ArcGIS JS map.  As far as I know, the only way to define a popup is on a specific layer, when you do something like:
const myFeatureLayer = new FeatureLayer({
  url: '...',
  popupTemplate: {
    ...popupTemplateStuffHere
  }
})

When a map has multiple layers, each with their own popups, an ArcGIS map will allow you to toggle through them.  But I can't seem to find a way to have all the content from all the different popups end up on a single popup element.
Take this codepen example.  There are 2 layers, each describing some data according to county in the US.  Each layer has a popup, so when you click a county, you can see the two popups, but only by tabbing through them with the arrows in the bottom right (< 1 of 2 >, < 2 of 2 >, etc.)  And each popup has data that comes from fields specific to that layer.  How can I conglomerate them?  I do know that similar things can be done with a Feature widget, but I want this in a popup that opens when the map is clicked.


texts: GeoTools & Shapefile: How to rename column (attribute) name

texts:
I'm accessing exiting shapefile with:

File file = new File("file_name.shp");
Map<String, Object> map = new HashMap<String, Object>();
map.put("url", file.toURI().toURL());
DataStore dataStore = DataStoreFinder.getDataStore(map);

or directly the associated dbf:

FileDataStore myData = FileDataStoreFinder.getDataStore(file);
SimpleFeatureSource source = myData.getFeatureSource();
SimpleFeatureType schema = source.getSchema();

but in both ways I don't know how to proceed in order to rename the column (attribute) names. I haven't found solution in the documentation.


texts:
Following this question, I am trying to use QgsRasterHistogram to obtain a histogram of my raster image. However, I cannot figure out how to change bin count.

I looked at the docs here, and it would seem to me that when calling

l = iface.activeLayer()
p = l.dataProvider()
p.initHistogram(QgsRasterHistogram(),1,100)
h = p.histogram(1)

the second parameter in initHistogram is the band number, and the third is bin count. However, no matter what number I put as third parameter,

h.binCount

always returns 1000. What am I missing?


texts: How to restrict or block certain paths in NetworkX graphs?

texts:
I am trying to calculate shortest path between 2 points using Dijkstra and A Star algorithms (in a directed NetworkX graph).

At the moment it works fine and I can see the calculated path but I would like to find a way of restricting certain paths.

For example if we have following nodes:

nodes = [1,2,3,4]

With these edges:

edges = ( (1,2),(2,3),(3,4) )

Is there a way of blocking/restricting  1 -> 2 -> 3 but still allow 2 -> 3 & 1 -> 2.

This would mean that:

can travel from 1 to 2
can travel from 2 to 3
cannot travel from 1 to 3 .. directly or indirectly (i.e. restrict 1->2->3 path).

Can this be achieved in NetworkX.. if not is there another graph library in Python that would allow this?


texts: Showing specific kernel density contours using QGIS Animove?

texts:
This is more a graphics question (I calculated all the important stuff in R), but is it possible to edit the Kernel Density Estimation (KDE) output layers given by the Animove plugin to show only the 50% and 95% contours?

Currently it gives me somthing like a heatmap that I cannot change.

@csk -  see attachment. Also I just realised that the created layer isn't saved with the project, so it seems like Animove calculates/creates temporary data-layers...?


texts:
The result from this query

Select pointlayer.geometry from pointlayer

return Point(40 30)

I need to add quotation marks to the result of my select query so I can use ST_GeomFromText. Something like this
I am hoping a query similar to this..

Select (' + pointlayer.geometry + ') from pointlayer

Would produce 'Point(40 30)'

So I can use the geometry like so

Select ST_GeomFromText((' + pointlayer.geometry + ')) AS geom from pointlayer

Would it be possible to concatonate in the select query to use the ST function?


texts:
I have done a bunch of internet research on this topic with no such luck. The closest I got to a solution is Removing extra string from text column in Attribute Table to make numeric column from number part using ArcGIS Field Calculator?.

I have an attribute field "UniqueID" with information formatted as "NE_MH143". other information in this field start with "SE_", "NW_" and "SW_". I would like to reduce this information to remove the first three characters of the string.

Can this be done with ArcGIS field calculator?


texts: Determining bare earth DEM from unclassified LAS file?

texts:
I have data in LAS format with RGB values created from aerial photogrammetry using a UAV. I am  trying to find a solution to extract the bare earth DEM from the point cloud.

I have tried SAGA, Fusion, MCC-LIDAR, but it is seems they need the LAS file to be already classified (which it naturally isn't). Can anyone point me in the right direction with a brief explanation of the process?

Generally, I would need to process about 100 mill points at a time (can tile them if needed).


texts: Getting a list of all neighbors within a specific distance for every node in QGIS?

texts:
I have a small task I would like solve in QGIS. I have a single layer with ~ 30000 nodes.
How can I create a list for every node in this layer which names all of its neighbors within a specific distance?

The only helpful tutorial I could find was about how to find neighbor polygons in a layer:

http://www.qgistutorials.com/de/docs/find_neighbor_polygons.html

Unfortunately I am just a beginner in python.

How do I adapt the mentioned script for my task?

from qgis.utils import iface
from PyQt4.QtCore import QVariant

# Replace the values below with values from your layer.
# For example, if your identifier field is called 'XYZ', then change the line
# below to _NAME_FIELD = 'XYZ'
_NAME_FIELD = 'NAME'
# Replace the value below with the field name that you want to sum up.
# For example, if the # field that you want to sum up is called 'VALUES', then
# change the line below to _SUM_FIELD = 'VALUES'
_SUM_FIELD = 'POP_EST'

# Names of the new fields to be added to the layer
_NEW_NEIGHBORS_FIELD = 'NEIGHBORS'
_NEW_SUM_FIELD = 'SUM'

layer = iface.activeLayer()

# Create 2 new fields in the layer that will hold the list of neighbors and sum
# of the chosen field.
layer.startEditing()
layer.dataProvider().addAttributes(
        [QgsField(_NEW_NEIGHBORS_FIELD, QVariant.String),
         QgsField(_NEW_SUM_FIELD, QVariant.Int)])
layer.updateFields()
# Create a dictionary of all features
feature_dict = {f.id(): f for f in layer.getFeatures()}

# Build a spatial index
index = QgsSpatialIndex()
for f in feature_dict.values():
    index.insertFeature(f)

# Loop through all features and find features that touch each feature
for f in feature_dict.values():
    print 'Working on %s' % f[_NAME_FIELD]
    geom = f.geometry()
    # Find all features that intersect the bounding box of the current feature.
    # We use spatial index to find the features intersecting the bounding box
    # of the current feature. This will narrow down the features that we need
    # to check neighboring features.
    intersecting_ids = index.intersects(geom.boundingBox())
    # Initalize neighbors list and sum
    neighbors = []
    neighbors_sum = 0
    for intersecting_id in intersecting_ids:
        # Look up the feature from the dictionary
        intersecting_f = feature_dict[intersecting_id]

        # For our purpose we consider a feature as 'neighbor' if it touches or
        # intersects a feature. We use the 'disjoint' predicate to satisfy
        # these conditions. So if a feature is not disjoint, it is a neighbor.
        if (f != intersecting_f and
            not intersecting_f.geometry().disjoint(geom)):
            neighbors.append(intersecting_f[_NAME_FIELD])
            neighbors_sum += intersecting_f[_SUM_FIELD]
    f[_NEW_NEIGHBORS_FIELD] = ','.join(neighbors)
    f[_NEW_SUM_FIELD] = neighbors_sum
    # Update the layer with new attribute values.
    layer.updateFeature(f)

layer.commitChanges()
print 'Processing complete.'


texts: qgis:JoinAttributeTable is giving a Syntax Error in QGIS 3.0

texts:
I am new to QGIS 3.0 and using "qgis:Joinattributetable" for joining .csv and .shp file. Below is my code, but I am getting Syntax Error.

AttributeError: module Processing has no attribute 'runandLoadResults'

I don't know how to fix this. Any suggestions?

vectorLyr=QgsVectorLayer("C:/Test/demand/Demand.shp","BuildingLayer","ogr")
vectorLyr.isValid()
uri='file:///C:/Test/data.csv?delimiter=,'
infoLyr=QgsVectorLayer(uri,'Population','delimitedtext')
infoLyr.isValid()
QgsProject.instance().addMapLayer(vectorLyr)
QgsProject.instance().addMapLayer(infoLyr)
shpField='FID_1'
csvField='FID'
result=processing.runandLoadResults('qgis:joinattributetable',vectorLyr,infoLyr,shpField,csvField,None)


texts: Adding attribute columns together QGIS 2.4

texts: I am in QGIS trying to combine a Zip code and City (separate columns) together and my methods are not working. Here is a screenshot of the issue.


texts: Using base layer of two different projections (Spherical Mercator and WGS84) in OpenLayers

texts:
How do I add baselayer with two different projections?

I have added overlayer with a different projection(EPSG:4326), but I need to add base layer with various projection

Google Maps
My vector layer with EPSG:4326 projection

I can see everything is fine with Google Maps but when I switch base layer to EPSG:4326 layer everything is white.

window.onload = function init(){
                var lon = 84;
                var lat = 28;
                var zoom = 8;
                var url = "http://10.6.7.6/geoserver/wms";

                var mapControls = [
                    new OpenLayers.Control.Navigation(),
                    new OpenLayers.Control.PanZoomBar(),
                    new OpenLayers.Control.LayerSwitcher({'ascending':false}),
                    new OpenLayers.Control.ScaleLine(),
                    new OpenLayers.Control.MousePosition(),
                    new OpenLayers.Control.OverviewMap(),
                    new OpenLayers.Control.KeyboardDefaults()
                ];

                var map = new OpenLayers.Map(
                'map',
                {
                    numZoomLevels: 21,
                    projection: new OpenLayers.Projection("EPSG:900913"),
                    displayProjection: new OpenLayers.Projection("EPSG:4326"),
                    controls: mapControls
                });

                var dsb = new OpenLayers.Layer.WMS(
                "Nepal Border",
                url,
                {layers: 'GisData:IR'}
            );

                var gmap_gs = new OpenLayers.Layer.Google(
                "Google Streets",// the default
                {'sphericalMercator': true}

            );

                var dG = new OpenLayers.Layer.WMS(
                "2G Down Sites",
                url,
                {layers: 'GisData:dG', transparent: "true",format: "image/png"},
                {
                    isBaseLayer: false,
                    buffer: 0
                });

                map.addLayers([dsb, dG,gmap_gs]);
                var proj = new OpenLayers.Projection("EPSG:4326");
                var point = new OpenLayers.LonLat(lon,lat);

                point.transform(proj, map.getProjectionObject());

                map.setCenter(point, zoom);

            }


texts:
I have a DGN file made by Workshop Paint Pro and would like to view it in ArcMap 10 so that I can create a risk map.

What are the procedures I need to follow?


texts: How to add a halo to a text in the composer?

texts:
I would like to add texts with halo, in the composer (not in the map). Do you know a way to do it?

I'm using QGIS 2.0.1


texts: How to substract 2 dem.tiff rasters through GDAL?

texts:
I'm a beginner in GDAL and don't understand the logic of it work completly.

I have a folder with 2 DEM-tif files and I need to do this:

gdal_calc.py -A D:WORK\gdal\dem05.tif -B D:WORK\gdal\dem31_size05.tif --outfile=demdiff_0531

The DEM have different sizes and different pixel sizes (so I'm getting errors), but the same location and the same coordinate systems.

"Error! Dimensions of file D:WORK\gdal\dem31.tif (25521, 21308) are different from other files (28986, 13493). Cannot proceed"

And when I did "gdal_translate -outsize" made them the same size I got this:

– Traceback (most recent call last):
File "C:\Program Files (x86)\GDAL\gdal_calc.py", line 335, in <module> main()
File "C:\Program Files (x86)\GDAL\gdal_calc.py", line 332, in main doit(opts, args)
File "C:\Program Files (x86)\GDAL\gdal_calc.py", line 168, in doit myOut.SetGeoTransform(myFiles[0].GetGeoTransform())
AttributeError: 'NoneType' object has no attribute 'SetGeoTransform'

I will have different DEMs in the future so I need to know how to do it in general.


texts:
I have a layer (currently the source is a csv file) with 130,000 points.  Many of these points are co-located in a building room.

How can I identify points that are within a small distance of another point on the same layer?  Basically I want to add a buffer around each points and then only display points where the buffer is touching another.


texts: How to sum or even perform more complex mathematical operations (multiply, add a constant, ...) two GDAL grids?

texts:
I have two grids with same cell size, same starting position (top left) and same grid size. I want to sum the grids and create a sum grid as an result:

e.g. grid A + grid B = grid C

How can I achieve that with GDAL?

What happens, if my grids are not identical (1) in grid size (2) in starting position and (3) in cell size.

Thank you very much for the answer.


texts: Define a legend for a MapBox Studio Project

texts:
Given that it is possible to request the TileJSON for a MapBox project at https://api.mapbox.com/v4/{mapId}.json?access_token={token}, I was wondering if it was possible to define a legend for MapBox projects created using either MapBox Studio or MapBox Studio Classic.

Normally, I just define legends in the front-end of whichever site, but for this use case I was hoping to expose all tilesets through a public service accessed by any number of clients.

I know this was possible using TileMill.  Perhaps by editing, in the case of Studio Classic, the project.yml file?


texts: How to convert a polygon into a set of points within its border using MapInfo?

texts: I need to convert a polygon into a set of points in MapInfo with a specific inter distance, I can't insert them manually as I need to have a large number of points (up to 10,000) at each run.


texts: ArcGIS GeoProcessing Model - Store List of Altered Features in Database

texts:
I am trying to build a GeoProcessing model which will find any parcels that have had their geometry values changed and store the corresponding parcel numbers in a separate non-spatial table.  I am using ArcGIS 10.

Does anyone know of any python scripts that would accomplish these tasks?  One potential solution I found is the Compare Layer to Snapshot tool.  Even though this is listed in the ArcGIS Resource Center, I can not find it in my version of ArcMap.

Any advice would be greatly appreciated.


texts: How to iterate through attribute visibility with PyQGIS?

texts: I have a point shapefile with a timestamp that I would like to turn off and turn back on 1 by 1 with a new element of the attribute table appearing every second


texts: File Geodatabase to map package or map document

texts:
I have a file geodatabase which is output of a certain model builder. I need to make it as input of another tool that takes only map documents or map packages.

Is there any tool in ArcMap or ArcGIS Pro that creates map packages or MXD from geodatabse.gdb ?


texts: Heat map / Density map from dynamic points table in MapServer (GeoServer)?

texts:
How can I create a heat map (some would call a density map) from a points layer in MapServer (or GeoServer)?
(preferably able to choose Interpolation method, NN, IDW, etc and color map with a transparent color)

The points are stored in PostGIS table and the table is dynamic, meaning I need to create the heat map on the fly.


texts: Print composer: Draw coordinates in right format and language

texts:
When drawing the coordinates of a map in WGS84 in decimal format, they appear in "German", which means they look like 81.3°O (O for "Ost" which means "East" in German). I am using a German version, but need them in English, indicating "E" and "N".

Is there a way to change this?


texts: Setting right coordinate system for northern Australia?

texts:
I have shapefile (points) which I got from adding xy data in ArcMap from Excel table. I had Latitude/Longitude data. It is area of northern Australia.

For example (5 points):

  Longitude
  136.882
  142.356
  143.049
  142.984
  145.086

  Latitude
  -12.249
  -11.151
  -12.659
  -12.881
  -15.135

What coordinate system should I define here and how?


texts: Input value into new field based on value in OID field

texts:
I need to input PofC in the PNT_TYP field (created in the script below) in my feature class (assigned to the variable verts) the odd OID values will have PofC input in PNT_TYP and the even with PofT. I am not sure how to do this. I have tried a number of things but this code does not seem to be working. If I have to run it in ArcMap then this will not work as I am trying to make it part of a larger script? if that is the case how do I do this so that I can run it from a script? I am not getting an error message, it just prints out script completed at the end like the script says to but the PNT_TYP field is blank in the feature class (as this script converts it to a lyr file) and I do not see any lyr file. I have never used one so I do not know if it would show in the catalogue or not. I have updated this to include the whole script but it is just the end where I need that one column to be updated with the two values based on the OID value that I am not sure if it is working as I see nothing to indicate it is.
# imports
import arcpy, os
arcpy.env.overwriteOutput = True
gdb = raw_input("please input the path of the GDB you are working with. >> ")
arcpy.env.workspace = gdb
in_fc = raw_input("Please input the fc that requires utm coordinates. >> ").upper().strip()
out_fc = raw_input("Please enter the name of the new start and end point feature class. eg) RoadNo_amdNo_pnts >> ").upper().strip()

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

arcpy.FeatureVerticesToPoints_management(in_fc, out_fc, "BOTH_ENDS") #creating point fc from line fc
# coord = arcpy.AddField_management ("verts", "PNT_TYP","TEXT", "","",25)
# code using modulus possibly to determine PofC and PofT.

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

arcpy.AddXY_management(out_fc) #creating points for each vertex

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

verts = raw_input("What do you want to call your final feature class? Eg. R18564_amd5_UTM >> ")
y_field = 'POINT_Y' #column of y values
x_field = 'POINT_X' #column of x values
input_format = 'SHAPE'
output_format = 'UTM'

arcpy.ConvertCoordinateNotation_management(out_fc, verts, x_field, y_field, input_format, output_format)#takes out_fc (from vertices tool) and writes out to verts

coord = arcpy.AddField_management (verts, "Coords","TEXT", "","",25)
expression = '!UTM!.split(" ")[1] + ", " + !UTM!.split(" ")[2]'
arcpy.CalculateField_management(verts, "Coords", expression, "PYTHON_9.3")

in_fld = arcpy.AddField_management (verts, "PNT_TYP", "TEXT", "","",25)

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

in_feature = os.path.join(gdb, verts)

feat_lyr = raw_input("What would you like to call the feature layer? >> "

arcpy.MakeFeatureLayer_management (verts, feat_lyr)
arcpy.SelectLayerByAttribute_management (feat_lyr, "NEW_SELECTION", 'MOD("OID", 2) = 0')
arcpy.CalculateField_management(feat_lyr, "PNT_TYP", '!PNT_TYP! = "PofT"', "PYTHON_9.3")

arcpy.SelectLayerByAttribute_management (feat_lyr, "NEW_SELECTION", 'MOD("OID", 2) > 0')
arcpy.CalculateField_management(feat_lyr, "PNT_TYP", '!PNT_TYP! = "PofC"', "PYTHON_9.3")

print ("script complete.")


texts: Making a 3D map with shapefiles, lidar and aerial image data

texts:
I want to try to make a simple 3D map for a small area (10 square miles) with some housing, some roads and some vegetation on it. My geodata is as follows:

shapes with houses/roads/vegetation.
LIDAR point cloud data over the area (can easily be built to a dem).
georeferenced aerial image.

Any idea where to start to bundle these geodata together and make a neat looking 3D map?


texts: Fill attribute-table on feature-creation based on intersection with other layer

texts:
I am digitising line-features based on WGS-84 point-coordinates from a PDF which cannot be copied from.

The final line layer should contain lines which each connect two points and have the following attributes:

Name of start point - name of end point
value in [ft]
value in [m]
remark (same for all features)
date of creation

Currently I have a csv-based point layer which simply contains the name and coordinates of each point.

Creating the point layer is no problem, but creating the lines is a bit of a hassle. Either I manually create a .csv containing all existing connections between points like the following table and run connect by lines:

Or I digitise the lines manually but fill in the values (except for the value in [ft]) automatically.

Personally I would prefer option 2, since option 1 involves more manual work and the process has to be repeated regularly (and I figure I could transfer the solution to other tasks), but I can't get the expression right which should fill in the name of the line. All other values are created flawlessly (though they are not taken from another layer based on geometry).

I added the WKT of each point to the point-layer. Currently I am running the following, which returns NULL, despite it seems to be correct to me:

attribute(get_feature('point-layer','geom',geom_to_wkt(start_point($geometry))),'Waypoint')||' - '||attribute(get_feature('point-layer','geom',geom_to_wkt(end_point($geometry))),'Waypoint')

Any ideas where I went wrong?


texts:
I use openlayers in my web app.I want to set map extent on specific country to show cities of that country but i don't know which option of map can show specific country in map!.
How can i show specific country on map and can show cities of that country?
To setting map on specific bound when i use below code i can set extent to desire extent!

var bounds = new OpenLayers.Bounds();
    bounds.extend(new OpenLayers.LonLat(8.91,55.18).transform(new OpenLayers.Projection("EPSG:4326"),
    new OpenLayers.Projection('EPSG:900913')));
    bounds.extend(new OpenLayers.LonLat(16.19,56.22).transform(new OpenLayers.Projection("EPSG:4326"),
    new OpenLayers.Projection('EPSG:900913')));
    map = new OpenLayers.Map('map_element', {
        projection: new OpenLayers.Projection('EPSG:900913'),
        zoomToExtent:bounds
  });
map.setCenter(new OpenLayers.LonLat(-1.788, 53.571).transform(
    new OpenLayers.Projection("EPSG:4326"),
    new OpenLayers.Projection('EPSG:900913')
), zoom);

But when i use zoomToExtent and zoomTo function map can't show my desire bound!

var bounds = new OpenLayers.Bounds();
    bounds.extend(new OpenLayers.LonLat(8.91,55.18).transform(new OpenLayers.Projection("EPSG:4326"),
    new OpenLayers.Projection('EPSG:900913')));
    bounds.extend(new OpenLayers.LonLat(16.19,56.22).transform(new OpenLayers.Projection("EPSG:4326"),
    new OpenLayers.Projection('EPSG:900913')));
    map = new OpenLayers.Map('map_element', {
        projection: new OpenLayers.Projection('EPSG:900913')
  });
map.zoomTo(10);
map.zoomToExtent(bounds);

why in the above code i can't zoom to desire bound?


texts: Visual Studio 2010 templates for ArcGIS missing

texts:
It seems that I do not have ArcGIS templates in VS 2010 templates folder -

  \Program Files (x86)\Microsoft Visual Studio
  10.0\Common7\IDE\ProjectTemplates\CSharp\

and as a result ArcGIS templates do not show up as installed templates for New Project. To get these templates will I need to reinstall VS or can I get the zipped ArcGIS templates from somewhere?
Thank you,
-Dee


texts: How to modify the coordinates of a shapefile geometry using GeoTools?

texts:
I have a shapefile with some polygons. I would like to change the height for this elements using Java. The rest of the attributes should be the same. I've thought using GeoTools but I haven't found how to implement this.

Any help?


texts: Trouble applying bfastmonitor function in R

texts:
I am trying to apply a BFAST monitor to my time series but keep getting the following error. What does it mean and how can I fix this?
    Error in if (!foundBreak & any(abs(obj$process) > boundary)) { :
       missing value where TRUE/FALSE needed

Below is my code. The values in sales_ex are the yearly mean NDVI values for a region between 2000 and 2020.
sales_ex <- c(df[2])
sales_ex

# [1] 0.7765 0.7989 0.7764 0.7952 0.7935 0.7987 0.8048 0.7930 0.7949 0.8053 0.8100 0.8055 0.7864
# [14] 0.8040 0.8001 0.8109 0.8126 0.8063 0.7980 0.8188 0.8123

tsales_ex <- ts(as.numeric(unlist(sales_ex)), start=c(2000), end=c(2020))
plot(tsales_ex)

bfm1_ex <- bfastmonitor(tsales_ex, response ~ trend + harmon, order = 3, start = c(2020))

# Error in if (!foundBreak & any(abs(obj$process) > boundary)) { :
#   missing value where TRUE/FALSE needed

plot(bfm1_ex)


texts: Leaflet HTTP.Get bad request; not loading tiles

texts:
I am trying to walk through the Leaflet quick start tutorial and am having a really simple issue: None of my tiles are loading.

I get this error at the image url and an HTTP Bad Request in the inspector:

  {"message":"id must be a string with the format account.handle"}

My code is copied and pasted straight from the tutorial:

<link rel="stylesheet" href="//cdn.leafletjs.com/leaflet-0.7.3/leaflet.css" />


texts: list of coordinates with unknown projection

texts:
I've a list of coordinates which were recorded in a last years aerial survey. Besides this I've the flightlog. The list of coordinates was created "manually" by writing down the position into a spreadsheet.

The projection of the flightlog data is: WGS84, EPSG4326

This files works fine and the locations seem to be correct.

Now I'd like to create a point-shapefile from my excel-spreadsheet. The problem right now, is that I don't have any information of the used projection. When I try to use EPSG4326 I get a gap (see screenshot).

The coordinates are as follows:

S074233 E0341051

S074109 E0340115

S073729 E0340115

S071631 E0342037

S071628 E0342837

S071409 E0343446

S071525 E0343546

I'm a bit perplex. Maybe important to know that they maybe used Skydemon or Garmin as GPS recorder. The area of interest is in Tanzania. I've access to ArcGIS and QGIS.


texts: Changing pixel values using QGIS raster calculator?

texts: How could I convert values of 32767 to 0 using raster calculator in QGIS?


texts: Use Java to automate service creation and display up to date data on map?

texts:
Our team is using ArcGIS and our requirement is

To create a web application in JAVA with a map displaying various data from DB.
The data gets updated periodically in the DB  and the map should display it as and when the data is ready in DB.
Numerous fields of the DB should be selectable on the map on a demand basis.

Right now we are creating a mxd file and publishing it as service and showing it on the ArcGIS for Dexktop.

So I am pretty unsure if I should use the ArcGIS for Java SDK for all of my above requirements or go with ArcGIS for Javascript API.
Please let me know the best possible way to achieve this.


texts: Generating a txt file from attribute table, multiple shapefiles in a folder

texts:
Edited for clarification. Thisis the process I currently use:
In ArcMap 10.6 I perform an intersect between 2 shapefiles, one is a drivetime area and one is a set of points in a geographical area. The intersect identifies which points lie within the drivetime area. I am not saving the result in the gdb but in a folder as a shapefile. ArcMap automatically adds this shapefile to the table of contents. I right click on the generated shapefile and select open attribute table. In the window that opens, I navigate to the column that I want "pat_id" then highlight it, then right click and select summarize. I change the output from a dbf to a text file then point it to a folder on my local machine, giving it a unique name. The pat_id field is a list of unique integers, it is the only column (field) that I need. I need to do this about 1400 times. I am a novice with python and model building, my attempts (which are many) have not resulted in an output that matches the process above. There is some baseline knowledge for model building that I have not acquired yet, and all the snippets I have tried normally lead to no output being generated. I am editing this question again to see if I can get help with a solution. Ideally I would like to build a scrip so that I can export the single field from each shapefile attribute table to a text file. I my be butchering the terms, but I do not know how to explain it another way.


texts: Computing distance between two points following network and doing it for 100.000 points?

texts:
My starting data :

a collection of points which are buildings
a collection of points which are electrical transformation plant
a network of vectors which is the public electric network from an operator
buildings are not connected to the electric network

I use QGIS to find which are the buildings which are located less than 250 meters
from an electrical transformation plant by following the network. I need this information more 100.000 buildings.

To do it by hand :

trace a line between the building and the nearest electrical network segment,
compute the distance by the network between this entry point and the nearest tranformation plant.
I have found the RoadGraph plugin (http://gis-lab.info/qa/road-graph-eng.html) which works for that
add the two distances and check if < 250m.

Now, my question is: how to automatize it for the 99.999 other buildings ?
I know Python a lot but I am a beginner with QGIS. I did not find a way to automatically interact
with the RoadGraph plugin and source code seems not to be available anymore.

Maybe I am doing it wrong ?


texts: QgsVectorFileWriter creating empty shapefile and invalid GEOJson

texts:
A script I wrote recently as an answer here was written and was working in QGIS 2.14.?, but now that I've updated to 2.18.3, the script is not working.

I'm not getting any error messages saying anything isn't working correctly, it is just creating an output point layer with 0 features. I'm saving the shapefile to file, and the shape files are not empty. The .dbf file is about 500Mb and the .shp file is about 10Mb, but when loaded in QGIS there are no features. EDIT more details on saving as a shapefile. Similar to saving as a GEOJson, when the script initially creates the shapefile it loads in with 0 features. If I close and reopen QGIS, the shapefile changes and if I load it again, then all of my expected features are there.

I tried outputting to different formats like GEOJson, and with GEOJson it created a perfectly valid .geojson file, but it left out the closing bracket for the features object and it left out the closing curly brace for the whole geojson file. I manually added in the closing bracket and curly brace and was able to load the geojson file into QGIS and everything looks correct from that point on. EDIT More details on saving as GEOJson. When I run the script the file is missing the closing bracket and curly brace, so it fails to load the layer, but if I close QGIS, then the file changes and if I open it again, then the closing bracket and closing curly brace have been added. It's as if QGIS is somehow preventing QgsVectorFileWriter from finalizing the file for some reason.

I've also tried completely reinstalling QGIS from the OSGeo4W 64bit installer (I'm using this installer instead of the standalone QGIS installer because I need the FileGDB libraries)

I tried to reproduce the error with a simpler script, so I took the example scripts here and tried to run those. Unfortunately (or fortunately depending on your perspective) those very simple scripts worked to generate valid shapefiles and geojson files.

So, why is the QgsVectorFileWriter failing to write out valid shape files and geojson files in my script?

The script from the answer linked above is copied below for convenience.

##measure_layer=vector
##output_layer=output vector

from qgis.core import *
from qgis.utils import *
from PyQt4.QtCore import *
from PyQt4.QtGui import *

def cleanUp(writer):
    del writer
    iface.messageBar().clearWidgets()

def addVertices( linestring, writer, inputFeature ):
    type = linestring.geometryType()
    if type == "MultiLineString" and linestring.isMeasure():
        geomCount = linestring.numGeometries()
        for i in range(0, geomCount):
            line = linestring.geometryN(i)
            addVertices(line, writer, inputFeature)
    elif type == "LineString" and linestring.isMeasure():
        vertexCount = linestring.vertexCount()
        for i in range(0,vertexCount):
            v = linestring.pointN(i)
            feature = QgsFeature(fields)
            m_value = v.m()
            feature.setGeometry( QgsGeometry.fromPoint( QgsPoint(v.x(), v.y()) ) )
            attributes = inputFeature.attributes()
            attributes.append(m_value)
            feature.setAttributes(attributes)
            writer.addFeature(feature)
    else:
        cleanUp(writer)
        raise ValueError("Input layer must be a linestring or multilinestring and contain M-values")
    return

iface.messageBar().clearWidgets()
progressMessageBar = iface.messageBar()

progress = QProgressBar()
progress.setMaximum(100)
progressMessageBar.pushWidget(progress)

inlayer = processing.getObject(measure_layer)
provider = inlayer.dataProvider()
fields = provider.fields()
fields.append( QgsField("MValue", QVariant.Double) )
writer = QgsVectorFileWriter(output_layer, 'UTF-8', fields,  QGis.WKBPoint , inlayer.crs())

features = inlayer.getFeatures()
featureCount = inlayer.featureCount()
featureIndex = 0

for f in features:
    percent = (featureIndex/float(featureCount)) * 100
    progress.setValue(percent)
    g = f.geometry()
    lines = g.geometry()
    addVertices(lines, writer, f)
    featureIndex +=1

cleanUp(writer)


texts:
I have a DTM raster which has too high a spatial resolution for my intended use, and I need it more generalized. Can't put my finger on how much more, but considerably more.

My first thought was to smooth its peaks and valleys with a gaussian filter, but after a day of repeated passes (it's quite a large area, so each pass takes a lot of time) it was nowhere near how much generalized I needed it to be. Given that I'll later need to do the same with other similar areas, this approach is unfeasible.

I thought that, instead of gathering information from surrounding pixels to change the middle pixel (which requires me to process almost all the pixels) I'd just make a median filter and apply the value to all cells in the kernel, thus making the process vastly faster.

But I have no idea how methodologically sound this is, and how much I'd be losing in terms of precision. This is purely for aesthetic purposes, so the resulting height values do not need to be overtly accurate, but I do need the resulting DTM to largely resemble the area's topography, as if it had been created from a smaller resolution sensor (or close enough).

If this approach is gibberish, which would yield good results in this scenario, without being too computationally heavy?


texts:
I am displaying a bunch of points on a map that have a field with an image name to popup when the point is clicked. The images are in a subfolder.

Is there an easy way to make the popup go to a new window with the target="_blank" tag?

I am using qgis2web, so I have this code snip to work with, it seems like this would be the right place to put a new window request, but instead I spawn a new tab for every feature in the map*

code from qgis2web

    var popupContent = '
\

\
                ' + (feature.properties['Image'] !== null ? '' : '') + '	\

\

';
    layer.bindPopup(popupContent, {maxHeight: 100});
}

end of qgis2web code

*because this gets called  by the GeoJSON loader...


texts: ogr2ogr: define own false easting value?

texts:
I have a shapefile with something similar to EPSG Projection 31256 - MGI / Austria GK East. However, the EPSG code states a false easting value of 0, while my prj file contains 750000 (that's probably also the reason why ArcGIS doesn't suggest any transformation parameters and thus won't reproject my shapefile?).

How can I define my false easting value when I want to re-project my shapefile from epsg:31256 to WGS 84 / UTM zone 33N (epsg:32633)?

Do I miss something?


texts: Converting GeoJSON to WKT with OpenLayers2


texts: I have a problem with my date format. how can I convert day of year values to mm-dd format? For instance if my doy have a value of 161 this will be converted to 06-10 (June 10 for instance). And this can work with leap and non-leap years.


texts:
I would like to be able to take my SQL queries for PostGIS and make each row for the geom attribute part of an sp object so I can plot it with R's leaflet bindings. I can't really share what my database but I will share how I learned how to connect and what I am trying to do to bind each row with a geom into a single sp object.  When I use readWKT for a single row the object class is SpatialPolygons.

#make my query to return multiple rows with just the geom attribute
query <- "SELECT ST_AsText(geom) FROM assembly"

#connect to db, send query, recieve table, and clear the result
drv <- dbDriver("PostgreSQL")
con <- dbConnect(drv, host = "localhost", port = "5432",
                 dbname = "anc", user = "postgres")
res <- dbSendQuery(con, query)
df <- dbFetch(res) ; dbClearResult(res)

Here is where I give it a shot to bring the rows together into a list object similar to what I've found with this answer.

#preallocate vector for combining spatial objects into a list
n <- dim(df)[1]
map_obj <- vector(mode = "list", length = n)

for(i in seq(n)) {
  map_obj[i] <- readWKT(df[i,1])
}

#Put the spatial objects together?
Polygons(map_obj,1:n)

An output of str(readWKT(df[1,1])):

Formal class 'SpatialPolygons' [package "sp"] with 4 slots
  ..@ polygons   :List of 1
  .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
  .. .. .. ..@ Polygons :List of 1
  .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
  .. .. .. .. .. .. ..@ labpt  : num [1:2] 1684530 2627594
  .. .. .. .. .. .. ..@ area   : num 3.54e+08
  .. .. .. .. .. .. ..@ hole   : logi FALSE
  .. .. .. .. .. .. ..@ ringDir: int 1
  .. .. .. .. .. .. ..@ coords : num [1:820, 1:2] 1690044 1690045 1690046 1690046 1690047 ...
  .. .. .. .. .. .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. .. .. .. .. .. ..$ : NULL
  .. .. .. .. .. .. .. .. ..$ : chr [1:2] "x" "y"
  .. .. .. ..@ plotOrder: int 1
  .. .. .. ..@ labpt    : num [1:2] 1684530 2627594
  .. .. .. ..@ ID       : chr "1"
  .. .. .. ..@ area     : num 3.54e+08
  ..@ plotOrder  : int 1
  ..@ bbox       : num [1:2, 1:2] 1670993 2612705 1698389 2642795
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "x" "y"
  .. .. ..$ : chr [1:2] "min" "max"
  ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slot
  .. .. ..@ projargs: chr NA


texts: Converting full road names to abbreviation using ArcGIS field calculator?

texts:
I have a data table with a column 'names' that have street names with street type (Road, Lane, Street etc.). I need to get the names in a format like in column 'new_names'using field calculator in ArcGIS.
How can I do this?