texts: Buffering in WGS 84 / Pseudo-Mercator 3857 using QGIS?

texts:
I have read all the buffer issues topics but haven't found solution to my problem.

I have two layers:

Border layer of Tokyo ward - WGS 84 / Pseudo-Mercator 3857
A CSV text layer with x,y coordinates of buildings which I saved in the same CRS.

Now, I need to buffer one (central) point from the point layer on the scale of the whole Tokyo to make some calculations. But when I add buffer and set the distance (cca 14 km) buffer actually has the distance of cca 11km.

I tried re-saving and opening the point layer, re-doing the layer couple of time. Followed all the instructions on other similar questions. But the same problem occurs. I have done the same on couple of cities more and there was no problem in making the buffer correctly.

Any other inputs?


texts:
I have a QGIS project with several shapefiles like the one in following image. I want to know how many times the shapefiles overlap and where is produced a gradient.
I think the best option is to transform the polygons in a grid line 1km*1km with data presence/absence to produce the gradient. But I am not pretty sure if this is possible and how to do it.


texts:
I am trying to make a stacked area graph in ggplot2:

Zooplankton<-read.table("clipboard".header=TRUE)

TransZoop   Zoop    DensZoop
0   Copepoda    0.066
0   Nauplii 0.425
0   Cladocera   0.199
0   Rotifera    0.269
0   others  0.041
217 Copepoda    0.045
217 Nauplii 0.567
217 Cladocera   0.110
217 Rotifera    0.267
217 others  0.011
462 Copepoda    0.046
462 Nauplii 0.325
462 Cladocera   0.522
462 Rotifera    0.090
462 others  0.017
848 Copepoda    0.043
848 Nauplii 0.020
848 Cladocera   0.928
848 Rotifera    0.009
848 others  0.000
1248    Copepoda    0.019
1248    Nauplii 0.013
1248    Cladocera   0.964
1248    Rotifera    0.003
1248    others  0.000
1765    Copepoda    0.017
1765    Nauplii 0.003
1765    Cladocera   0.980
1765    Rotifera    0.000
1765    others  0.000
3216    Copepoda    0.030
3216    Nauplii 0.003
3216    Cladocera   0.962
3216    Rotifera    0.004
3216    others  0.000

When I use the command:

Z<-ggplot(Zooplankton,aes(x=TransZoop,y=DensZoop,group=Zoop,fill=Zoop))+geom_area(position="fill")+scale_fill_brewer(palette="Dark2")

It turns out like this:

But when I try to change some graphical parameters, the plot disappears:

Z<-ggplot(Zooplankton,aes(x=TransZoop,y=DensZoop,group=Zoop,fill=Zoop))+geom_area(position="fill")+scale_fill_brewer(palette="Dark2") + labs(list(x="Distance from margin (m)",y="Zooplankton"))+theme(axis.text.x=element_text(size=30),axis.text.y=element_text(size=30),axis.title.x=element_text(size=30),axis.title.y=element_text(size=30,angle=90),legend.text=element_text(size=30),legend.title=element_blank(),legend.key.size=unit(2,"cm"),legend.key=element_rect(colour="black",linetype="solid",size=1.5),panel.border=element_rect(colour="black"))

What am I doing wrong?


texts:
I would like to use the new ST_TileEnvelope() function from PostGIS 3.0. My current version is 2.4.4.

I would like to update to PostGIS 3.0 in pgAdmin. Is this possible?

I have tried
ALTER EXTENSION postgis UPDATE TO '3.0.0';

This returns the error...

ERROR:  extension "postgis" has no update path from version "2.4.4" to version "3.0.0"

I am on Windows 10.

I initially installed our older version of postGIS (2.4.4) with no postGIS file on my computer. Nevertheless, I downloaded postgis-3.0.0 and I can now find this folder in the Windows File Explorer.


texts: How can I create a circle feature that will have a constant radius in pixels even on zoom in and zoom out?

texts:
Creating SVG icons have also solve it with much greater flexibility:

var svgXmlCircle = '<svg version="1.1" id="Layer_1" xmlns="http://www.w3.org/2000/svg" ' +
        'xmlns:xlink="http://www.w3.org/1999/xlink" width="55" height="55" ' +
        ' xml:space="preserve">'+
        '<circle fill="#FF0000" cx="27.5" cy="14" r="14" stroke="black" stroke-width="1"/>' +
        '</svg>';

    var svgIcon = new Image();
    svgIcon.src = 'data:image/svg+xml,' + escape(svgXmlCircle);

    return [
        new ol.style.Style({
            image: new ol.style.Icon(({
                img: svgIcon, // for new implem
                imgSize:[55,55], // for new implem
            }))
        }),
        new ol.style.Style({
            image: new ol.style.Icon(({
                src: 'images/map-' + itemType + '.svg',
                scale: .1,
                //center bottom of the marker
                anchor: [0.5, 1],
                anchorXUnits: 'fraction',
                anchorYUnits: 'fraction',
                opacity: 1
            }))
        })
    ];


texts: Why won't ArcMap load after borrowing a license?

texts:
We have an ArcGIS license manager on a server on our network. We have to occasionally take the server down for maintenance, so we have the production team borrow their licenses so that they can keep working.

ArcMap will not launch using the borrowed license while still connected to the network (from which the license server has been removed). The program just hangs on the splash screen. However, if we unplug the machines from the network ArcMap launches without issue. If we then plug the machine back into the network with ArcMap launched, several features of the software operate more slowly (e.g. accessing feature properties). We need to maintain a network connection in order to work with data located on another server.

Is ArcMap ignoring the checked-out license because it thinks it can still find the license manager on the network?

Is there a way in which we can check-out the licenses and use them without having to disconnect from the network?

NOTE: The solution provided by @PolyGeo worked. I did end up having to have ESRI reauthorize our licences as we had to reinstall the License Manager while all licenses were borrowed, which ended up causing a conflict.


texts: Exporting layout file from ArcGIS Pro project to *.pagx file using ArcPy?

texts:
I'm trying to implement a tool to import *.pagx files exported from ArcGIS Pro  to a  project, modify text elements and export the altered layouts back to the *.pagx files. It is pretty easy to export these files via ArcGIS Pro GUI, but I'm struggling to find a way to export these files via ArcPy.

So far I used the following code to modify layout from the individual *.pagx file:

import arcpy

layout = arcpy.mp.ConvertLayoutFileToLayout(r"C:\temp\templates\8x11 Landscape.pagx")
for element in layout.listElements('TEXT_ELEMENT'):
    if element.name == 'Text 4':
        element.text = "New text"

but there is no way to save the modified layout object as a file.
Also I tried to loop through the layouts, imported to the Pro project

import arcpy
aprx = arcpy.mp.ArcGISProject(r"C:\Users\username\Documents\ArcGIS\Projects\Layouts\Layouts.aprx")

for layout in aprx.listLayouts():
    for element in layout.listElements('TEXT_ELEMENT'):
        if element.name == 'Text 4':
            element.text = "New text"
aprx.save()

Using the save() method I can save the project, but again there is no way to export the *.pagx files. I heard somewhere that I could export these *.pagx files via the .Net API only, but I'm still not 100% sure if such export feature is not implemented in arcpy, or I overlooked something.

Update:
Tried to employ jsontocim module, located in "C:\Program Files\ArcGIS\Pro\Resources\ArcPy\arcpy\cim\cimloader\jsontocim.py"

from arcpy.mp import ConvertLayoutFileToLayout
from arcpy.cim.cimloader.cimtojson import CimJsonEncoder
from json import dumps
layout = ConvertLayoutFileToLayout(r"C:\temp\templates\8x11 Landscape.pagx")
cim = layout.getDefinition('V2')
outJSON = dumps(cim, cls=CimJsonEncoder, indent=2)
print(outJSON)

It produced the output, resembling the structure of the *.pagx file, exported via GUI ,but some parts are missing or their structure differs from the original file.


texts:
Examples of GIS Visualisation

Example:
Dr. John Snow (1813-1858), a legendary figure in the history of public health, epidemiology and anesthesiology created one of the first uses of geographical analysis on a static map to solve a problem of stopping the spread of cholera through Soho, London, England.

http://upload.wikimedia.org/wikipedia/commons/thumb/2/27/Snow-cholera-map-1.jpg/643px-Snow-cholera-map-1.jpg

The first two cholera cases in 1848 occurred in the same room but a week apart at No. 8, New Lane (in red) near Gainsford Street in the community of Horsleydown... water pumps were causing the spread and over a period of time it was found after collating deaths and visualising them.

Good Examples can be ancient, historical or cutting edge.

source:
http://www.ph.ucla.edu/epi/snow/1859map/snowindex.html


texts:
Looking for a way to add an AutoCad-type circular inset in ArcGIS 10.1. Have tried many different methods but the best we have been able to do is shown here

using "insert picture" after various image/raster processes.


texts: How to display features one by one in Map Book?

texts:
I'm using "Data Driven Pages" Tool to create a map book in ArcMap. I want to display buffer circles one by one in each page, but the "Data Driven Pages" Tool displays all features around it.

How can I display one buffer circle without others in one page?

Can ArcMap do that? Or I need to learn ArcPy?

What I get:

What I want:


texts: I have created my own customized tile. Now I am facing some problem with its indexing or naming convention. I am giving an clear example. I have taken an area where 6 tiles were created. In Google Maps or OSM in the same extent its giving 6 tiles and the same tiles. So there is no referencing or projection problem. Even the X or the folder name is also matching. I am getting trouble with the Y column or Image name. Its not matching. I have taken a close look at it. Interestingly I found in Google or OSM the Y column is incrementing from top to bottom where as tilemill created tiles number those incrementing from bottom to top. That is why the naming of tiles are not matching. What should I do?


texts: Resize rasters to add empty data around edges

texts:
I have a set of DEMs that have different extents due to different data being available at each of the sample points these DEMs were generated from.  I need to use these to create a set of uniformly sized rasters with the same spatial extents.  This would mean creating a new raster that has a null or 0 value everywhere except when it can use a value from the original (which is smaller).  I've played with con but it doesn't seem to do what I want.

Another way I've tried to approach this is to create the larger raster at the resolution I want, then use the raster calculator to go through and set values when it can extract them from the smaller original raster, leaving them blank/0/null otherwise.  But I haven't been able to figure out how to implement that because I don't know how to refer to whether you're outside the extent of a raster or not.  Please let me know if I need to reword part of this to be more clear.


texts: QGIS Authentication System for plugin manager

texts:
Does anyone have any details/documentation or information yet on the new Authentication System for Plugin Manager added to QGIS v2.14?

We would love some more information, and in particular if it resolves previously known issues with regard to the storage of login details on a secure network?


texts: Extracting basic statistics (raster data) from different latitudes using R?

texts:
I want to extract the mean, min and max stats from latitude -50 to 0 and from 0 to 50 etc. No idea how to do this.

I would prefer using R for this as I have a lot of these to do.


texts: MapServer propertyname=fieldname only works when using gml as output

texts:
Is it possible to make the WFS propertyname work for MapServer GeoJSON outputformat too?

Example:

curl "http://localhost/cgi-bin/mapserv?map=$MAP&SERVICE=WFS&version=2.0.0&request=GetFeature&typenames=lines&count=2&outputformat=geojson&propertyname=myid"


texts: Import error: no module named _gdal_array

texts:
I am using python 2.7 and have gdal installed. I want to read an entire .jpg image at once and so I used statement
import gdal
from gdalconst import *
driver = gdal.GetDriverByName('JPEG')
driver.Register()
dataset = gdal.Open('1.jpg', GA_ReadOnly)

but I am getting an import error for gdal_array.
ImportError                               Traceback (most recent call last)
<ipython-input-7-6a16047a395b> in <module>()
----> 1 data = band.ReadAsArray(0, 0, cols, rows)

C:\Users\DELL\Anaconda\lib\site-packages\osgeo\gdal.pyc in ReadAsArray(self, xoff,        yoff, win_xsize, win_ysize, buf_xsize, buf_ysize, buf_obj)
1174     def ReadAsArray(self, xoff=0, yoff=0, win_xsize=None, win_ysize=None,
1175                     buf_xsize=None, buf_ysize=None, buf_obj=None):
-> 1176         import gdalnumeric
1177
1178         return gdalnumeric.BandReadAsArray( self, xoff, yoff,

C:\Users\DELL\Anaconda\lib\site-packages\osgeo\gdalnumeric.py in <module>()
----> 1 from gdal_array import *
2 from numpy import *

C:\Users\DELL\Anaconda\lib\site-packages\osgeo\gdal_array.py in <module>()
 24                 fp.close()
 25             return _mod
 ---> 26     _gdal_array = swig_import_helper()
 27     del swig_import_helper
 28 else:

 C:\Users\DELL\Anaconda\lib\site-packages\osgeo\gdal_array.py in swig_import_helper()
 16             fp, pathname, description = imp.find_module('_gdal_array', [dirname(__file__)])
 17         except ImportError:
 ---> 18             import _gdal_array
 19             return _gdal_array
 20         if fp is not None:

 ImportError: No module named _gdal_array

Please tell me what could possible be going wrong. I am new to using GDAL.


texts: Adding a well discernible symbol on top of a heat map in QGIS without duplicating a point layer

texts:
I have a point layer in QGIS 3.14 on which I applied a heat map symbology:

I would like to add a well defined / discernible shape (circles here) on top of the heatmap for the underlying points to actually make the points discernible for the user. This would help him to select the points on a web app based on this QGIS project file, because when I select to "Weight points by" an intensity attribute, the heatmap blends with the basemap layer and some spots may be really hard to see where the intensity values are low as shown for 4 dots on this second picture:

Is this possible without having to actually duplicate the layer?
For the moment, the only workaround I found, as you can see in the layer window, is to duplicate the layer with simple circles symbols on top of the heatmap.


texts: Activating graphic creation from button using ArcGIS API for JavaScript?

texts:
I need to have 2 graphic layers. One is used to set a starting location on a network, and the other is used to mark a stopping location.

There should only be a maximum of 1 feature in each layer. So, a single starting point, and a single stopping point.

I'd like the user to be able to click a button called, for example, Set Starting Location which would activate that graphic layer and allow the user to put a point on the map. Then there'd be another button that they could click called, for example, Set Stopping Location which would work the same way.

How do I make it so one button activates one graphic layer, and another button actives the other graphic layer? I'd also want to allow only 1 feature per layer, which I assume I'd use some type of count function.

The geometry of these graphics would eventually be used in a Python geoprocessing tool.


texts:
I'm using the US Census WMS service describe here.  Using code similar to that found in this example, but specifying a data format of "image/svg+xml", I've been able to get the TIGERweb WMS service to give me SVG that is rendered in my map.  However, each shape in the SVG just the actual path, with no identifying data.

Is there any way to get TIGERweb to annotate the shapes it serves so I can tell what they are (and process them accordingly with JS, CSS, etc.)?  Like, I would like to get an SVG where each path is the boundary of a census block, tract, etc., and where the each path element in the SVG has an attribute giving the block number, FIPS code, or whatever piece of info will allow me to identify what the shape represents.

As it is, I can't see how one could use the TIGERweb SVG data for anything except a plain rendering of the boundaries of counties, tracts, or what have you.  Without any identifying information, there's no way to spiff up the shapes with special styling to represent data associated with the census areas.


texts: Identify polygons greater than n-th percentile and perform further analyses using arcpy

texts:
My goal is

to identify tracts that are above the 99 percentile of a variable (xvar),
to dissolve contiguous tracts by summing a variable (yvar), and
to get the coordinates of the centroids (longitude and latitude) of the merged polygons.

At the very first stage for picking the above-99th-percentile tracts, I have got some errors. The following is my current work:

import arcpy
import numpy as np

input = r'D:\7362.shp'
myarray = arcpy.da.FeatureClassToNumPyArray(input, ('xvar'))

print myarray
# Printed output:
# [(5.204770523,) (8.6471839055,) (2.1095837756,) ..., (1.4788840302,
# (6.24183499,) (3.6710648163,)]

p99 = np.percentile(myarray, 99)

# use cursor to create the indicator field
# Do I have to create the abovep99 field in the first place???
with arcpy.da.UpdateCursor(input, ['xvar','abovep99']) as cursor:
    for row in cursor:
        if row[0] > p99:
            row[1] = 1
        else:
            row[1] = 0
        cursor.updateRow(row)

The error messages in the percentile calculation p99 = np.percentile(myarray, 99) are following:

Traceback (most recent call last):
File "<pyshell#11>", line 1, in <module>
p99 = np.percentile(myarray, 99)
File "C:\Python27\ArcGIS10.3\lib\site-packages\numpy\lib\function_base.py", line 3096, in percentile
return _compute_qth_percentile(sorted, q, axis, out)
File "C:\Python27\ArcGIS10.3\lib\site-packages\numpy\lib\function_base.py", line 3132, in _compute_qth_percentile
return add.reduce(sorted[indexer]*weights, axis=axis, out=out)/sumval
TypeError: unsupported operand type(s) for *: 'numpy.ndarray' and 'numpy.ndarray'

How do I fix this problem for the percentile calculation?


texts: Adding Oracle Spatial View as Layer in QGIS?

texts:
Is it possible to add an Oracle Spatial View on QGIS MAP?

My view is registered on MDSYS metadata but, when I connect with Oracle DB, it is disabled on the list.

QGIS Desktop 2.2.0

Oracle Database 11g Enterprise Edition Release 11.2.0.4.0 - 64bit Production
PL/SQL Release 11.2.0.4.0 - Production
CORE    11.2.0.4.0      Production
TNS for 64-bit Windows: Version 11.2.0.4.0 - Production
NLSRTL Version 11.2.0.4.0 - Production


texts: Output of QGIS raster clipper not as expected?

texts:
My goal is to clip a raster with a polygon in QGIS 1.8. Based on an article I read, this should be a pretty straightforward task. However, the output appears to be empty. The raster produced is the same size as the original. The metadata for the output raster shows that it is the same extent as the original, with the same pixel size, same number of bands, but there is no information about the content of the bands.

I've attached a screenshot of the parameters I am using for the clip. I have done this two other ways, all with similar results. The first is with the NoData value set to 32767, which is the same as the NoData value of the input raster. The second is with 'Create an output alpha band' checked, which creates a 2-band raster that is otherwise identical to the output described above.

I would very much appreciate it if someone could critique my parameters, and suggest what I might be doing wrong.

My OS, if that is important, is Ubuntu 12.04, and I'm using the KDE desktop.


texts: GeoTIFF - multiple tiles with different georeferencing

texts:
Can I have GeoTIFF file with multiple tiles where each tile has different bounds? The projection is the same for all tiles.
If it is possible, how can I create such file?
I have tried:
When I create tiled tif with:  gdalbuildvrt out.vrt 0_0.tif 2_2.tif and gdal_translate -of GTiff -co "TILED=YES" out.vrt mosaic.tif the resulting TIFF when listed with gdalinfo acts as a single gereferenced image. There are only one corner coordinates, but I would expect two sets, one for each tile. If I open this image in QGis, I get one large image with tiles but places where no tile is are "black". I expected that places with no tile will be "empty".
I have also tried gdal_merge.py -o output.tif -of gtiff 0_0.tif 2_2.tif", but with the same result.


texts: Identify shape of the polygons: elongation, roundness, etc

texts:
I am generalizing geological map. In the initial step I am dealing with small polygons which should be either deleted, enlarged, aggregated and etc.

Before applying operations mentioned above I would like to start with the shape analysis, such as elongated polygons or more round polygons and so on, as this may be key factor to deal with.

Geological maps have quite complex shapes which should be considered and maintained as much as possible.

Is there any way that identifies the shapes of the polygons based on the geometry?

I am using ArcGIS, Python basically.


texts: Read GPX file from desktop in Openlayers 3

texts:
I'm aware that Openlayers 3 has the ability to read and parse gpx files using the drag/drop event and I've incorporated this functionality into my map. What I'd also like to do is give the user the ability to 'upload' the file in a more conventional way by browsing through their filesystem.

I've hooked up the basic functionality and the resulting file can be output to a div as text but how do I get the final step of instructing Openlayers to read this information and subsequently display it on a map?

There seems to be very little information regarding this and reading the ol-debug.js file hasn't shed any light. Can anyone give me any pointers please?


texts: Exporting map to image with annotations created using PyQGIS

texts:
Using QGIS 3.16, I have been able to create annotations following this solution. The problem is when saving the map as an image, the annotation will not render. Here are the portions of my code that create the annotation and save the image:
from qgis.utils import *
from qgis.core import *
from qgis.gui import *
from PyQt5.QtGui import *
from PyQt5.QtCore import *

Layer1 = QgsProject.instance().mapLayersByName('Layer1')

# MAP SETTINGS
img = QImage(QSize(715, 850), QImage.Format_ARGB32_Premultiplied)
color = QColor(255, 255, 255, 255)
img.fill(color.rgba())
p = QPainter()
p.begin(img)
p.setRenderHint(QPainter.Antialiasing)
ms = QgsMapSettings()
ms.setBackgroundColor(color)
ext = Layer1[0].extent()
ms.setExtent(ext)

# TITLE ANNOTATION
title = QgsTextAnnotation()
title.setFrameSize(QSizeF(100, 25))
title.setFrameOffsetFromReferencePoint(QPoint(0, 0))
title.setDocument(QTextDocument('TEST'))
title.setMapLayer(Layer1[0])
title.setMapPositionCrs(QgsCoordinateReferenceSystem(Layer1[0].crs()))
title.setMapPosition(QgsPointXY(ext.xMinimum(), ext.yMaximum()))
tsymb = QgsMarkerSymbol()
tsymb.setSize(0);
title.setMarkerSymbol(tsymb)
QgsProject.instance().annotationManager().addAnnotation(title)

# OUTPUT
ms.setLayers([Layer1[0]])
ms.setOutputSize(img.size())
render = QgsMapRendererCustomPainterJob(ms, p)
render.start()
render.waitForFinished()
p.end()
img.save('C:/Users/user1/Desktop/test.jpeg')

What am I missing? It is odd that the annotation does appear in the map (canvas) when created via pyqgis and even when comparing it to one created using the GUI 'Text Annotation' they look exactly the same, but when exporting the map as a jpeg (using the GUI 'export map to image' option), only the one created with the GUI is rendered.


texts: Fitting Time-series Values with a Gaussian Function

texts:
I'm trying to fit a stack of NDVI values to a Gaussian model to allow for determining dates of certain NDVI values using Python and NumPy/SciPy. I've attempted to do this with scipy.optimize.curve_fit.

The following is the function I'm using when applying curve_fit to the stack

def func(x, *p):
    A, mu, sigma = p
    return A * np.exp(-(x-mu)**2/(2.*sigma**2))

The parameters that I used for the curve fitting are the following

p0 = [1., 0., 1.]  # initial guess for the fitting coefficients
newX = np.linspace(date[0], date[-1], (date[-1] - date[0]) + 1)  # range of days of the year between first image in stack and last

To get the stack of values for each pixel and apply the curve fit, I used the following nested loop

for x in range(imageHeight):
    for y in range(imageWidth):
        lai = imgStack[x, y]
        popt, pcov = curve_fit(func, date, lai, p0)
        yFit = func(newX, *popt)

When I run the script, the output of yFit is always just an array of 0.0. The values of the popt variable are always the same as p0 which leads me to believe that the curve was not fit properly. Is the model not suitable for a stack of 20 NDVI values?

Credit to this post for most of this code.


texts: In R: Convert from WGS84 lat/long to UTM and back again to WGS84 lat/long

texts:
I am looking for a way in R to project WGS84 latitudes and longitudes to a UTM grid and then convert the UTM coordinates back to the WGS84 latitudes and longitudes. I have found a number of ways to project to a UTM grid (mapproj::mapproject, rgdal::spTransform), but going from the grid back to latitudes and longitudes seems to be more complicated. Is there a projection function that offers the inverse?

I would like to use find the reverse of this:

newdata <-mapproject(dataset$lon, dataset$lat, projection = "lambert", parameters = c(mean(dataset$lon), mean(dataset$lat)))


texts: Inconsistent display of HTML table font size in QGIS layouts

texts:
I'm using QGIS 3.10.1, and am having trouble with getting HTML tables to render correctly in a layout.  I generated the HTML table using the "TABLEIZER!" online tool, and just pasted the results into the HTML box.  The CSS looks like something this:

<style type="text/css">
    table.tableizer-table {
        font-size: 12px;
        border: 1px solid #CCC;
        font-family: Arial, Helvetica, sans-serif;
    }
</style>

I first created a new HTML object on the layout, and pasted the HTML content from Tabelizer into the code box.

Then I copy/paste the same HTML table object (including the same HTML values) into different layouts.  In the different layouts, the font sizes are being displayed in dramatically different in sizes.  Finally, I copied one of the tables with "appropriately" sized fonts and pasted it back on the layout with the fonts that are being displayed too large.

The image below shows how on the same layout, 2 tables with the exact same HTML code, is being displayed differently.  Both the HTML and CSS are located in the "HTML Source" box.  Any ideas what the heck is going on, or how to get the "big font" HTML table to refresh to look correct (smaller text)?


texts:
Does anyone have any idea how I could download all QGIS plugins without having to do that from the interface one by one (possible but tedious).

I am doing this for a friend who uses an offline computer, I'd like to send him a USB drive with the new version of QGIS and the plugins.


texts: Entering digitizing mode makes QGIS freeze

texts:
When I toggle the editing mode QGIS freezes. The layer is stored in PostGIS.
I'm using QGIS 3.10.


texts:
I'm having problems installing the gdal library for python. I use a Windows 7 and have WinPython 2.7.9.2 which runs Python 2.7.9. I'd like to use gdal + Python for raster manipulation.

I've read that PIP is the preferred way for installing packages on windows and that its packed with Python 2.7.9. I've edited my PATH environment variable to include pip and python:

;c:\~\WinPython-32bit-2.7.9.2\python-2.7.9\Scripts\;c:\~WinPython-32bit-2.7.9.2\python-2.7.9\

so I can now run both python and pip directly from command prompt.

The error appears when I type

c:\>pip install gdal

which gives:

C:\~\AppData\Local\Programs\Common\Microsoft\Visual C++ for Python\9.0\VC\Include\xlocale(342) : warning C4530: C++ exception handler used, but unwind semantics are not enabled. Specify /EHsc
extensions/gdal_wrap.cpp(2855) : fatal error C1083: Cannot open include file: 'cpl_port.h': No such file or directory
error: command 'C:\\~\\AppData\\Local\\Programs\\Common\\Microsoft\\Visual C++ for Python\\9.0\\VC\\Bin\\cl.exe' failed with exit status 2

So apparently something called 'cpl_port.h' is missing. This is also described in this question where it's suggested that header files cannot be found. The solution provided in that thread is however LINUX based and incompatible with my Windows.

I'm not quite sure how to proceed from here so should I forego PIP and try another route?

I've verified that PIP works by installing some other packages, so the problem seems to be with gdal.


texts:
I'm trying to import a KML file into an existing GeoDjango app (using a Postgres 9.3 database). This is my Django model:

class PCT(models.Model):
    code = models.CharField(max_length=3, primary_key=True)
    boundary = models.GeometryField(null=True, blank=True)
    objects = models.GeoManager()

And this is my code [updated to add the relevant part!]:

    filename = 'myfile.KML'
    ds = DataSource(filename)
    lyr = ds[0]
    for feat in lyr:
        code = feat.get('NAME')
        geom = feat.geom
        pct = PCT.objects.get(code=code)
        pct.boundary = geom.json
        pct.save()

But when I run it I get this error:

TypeError: Cannot set PCT GeometryProxy (GEOMETRY) with value of type: <class 'django.contrib.gis.gdal.geometries.Polygon'>

What am I doing wrong?

Note that some of the KML features are Polygon and some are MultiPolygon (which I think makes it difficult to use LayerMapping to import the data, which I had ben considering as an alternative approach).

UPDATE:

This seems to work (keeping the boundary field as a GeometryField):

    filename = 'myfile.KML'
    ds = DataSource(filename)
    layer_mapping = {
        'code': 'Name',
        'boundary': 'Unknown'
    }
    lm = LayerMapping(PCT, filename, layer_mapping, transform=False)
    lm.save(strict=True)


texts:
I uploaded the mmqgis plugin to merge layers together, but when I tried to open mmqgis> combine> merge layers there is a python error message. Details read, "AttributeError: 'QListWidget' object has no attribute 'findText'"

Traceback (most recent call last):
  File "C:/Users/VolAero/.qgis2/python/plugins\mmqgis\mmqgis_menu.py", line 316,
   in merge dialog = mmqgis_merge_dialog(self.iface)

  File "C:/Users/VolAero/.qgis2/python/plugins\mmqgis\mmqgis_dialogs.py", line 1714,
   in __init__mmqgis_load_combo_box_with_vector_layers(self.iface, self.sourcelayers, True)

  File "C:/Users/VolAero/.qgis2/python/plugins\mmqgis\mmqgis_dialogs.py", line 2610,
   in mmqgis_load_combo_box_with_vector_layers combo_index = combo_box.findText(layer.name())

AttributeError: 'QListWidget' object has no attribute 'findText'

Does anyone have the solution?


texts:
I wonder if there is a "unified" source for map color conventions, that GIS map builders are encouraged to follow?

EDIT

To clear up my question a bit. I'm more interested in the color with respect to the object it represents.  For example, when you look at a map, and see "blue", the assumption could be water; however, it gets a bit shady when talking about let's say property maps: "is it pink or grey for industrial?". It'll be nice if there were such conventions whether loose or rigid about associating colors to map objects


texts: Writing GeoPandas data frame to PostGIS?

texts:
Similar question has to this has been asked before Adding GeoPandas Dataframe to PostGIS table?
I have followed the example but I am running into an issue I cannot solve. There is some geometry problem. What am I doing wrong?

I have defined a simple geometry in PostgreSQL using PostGIS:

alter table helloworld
ADD COLUMN geom geometry(Point,4326)

This is the error I get running the script

ProgrammingError: (psycopg2.ProgrammingError) can't adapt type 'Point'
[SQL: INSERT INTO helloworld (hello, world, geometry) VALUES (%(hello)s, %(world)s, ST_GeomFromEWKT(%(geometry)s))]
[parameters: ({'hello': 'hello1', 'world': 'world1', 'geometry': <shapely.geometry.point.Point object at 0x000001F7C5B21BE0>}, {'hello': 'hello2', 'world': 'world2', 'geometry': <shapely.geometry.point.Point object at 0x000001F7C5B21518>})]

Here is my script.

# Creating SQLAlchemy's engine to use
engine = create_engine('postgresql://postgres:postgres@127.0.0.1:5432/geospatial')

# ====== Writing table ====== #
# Creating a simple pandas DataFrame with two columns
liste_hello = ['hello1','hello2']
liste_world = ['world1','world2']
geom=['POINT(7.14111328125 51.957807388715516)','POINT(9.426269531249998 52.24125614966341)']

# Create DataFrame
df1= pd.DataFrame(data = {'hello' : liste_hello, 'world': liste_world, 'geom':geom })
df1['geometry'] = df1['geom'].apply(wkt.loads)

# Create GeoDataFrame with geom from dataframe
gdf1 = geopandas.GeoDataFrame(df1, geometry='geometry')

#https://gis.stackexchange.com/questions/239198/geopandas-dataframe-to-postgis-table-help
gdf1['geom'] = gdf1['geometry'].apply(lambda x: WKTElement(x.wkt, srid=4326))
gdf1

#drop the geometry column as it is now duplicative
gdf1.drop('geom', 1, inplace=True)
gdf1

# Write to db
gdf1.to_sql('helloworld', engine, if_exists='append', index=False, dtype={'geometry': Geometry(geometry_type='POINT', srid= 4326)})


texts:
I'm trying to make a topographic map of an area in Iceland, using SRTM data from here. There are other data sources, but this seems fine for what I'm trying to do.
The unedited, raw HGT file looks like this when imported into ArcGIS Pro:

I'm trying to find out how to convert this into a nice looking shaded relief map, similar to something like this.
After searching through the pain that is the ArcGIS documentation, I concluded that I had to do something with either the HillShade or Slope functions in the 3D Analyst toolbox. However, when I try both of these tools, the produce some very weird results.
For example, the HillShade function on default settings:

There is also this weird triangular prism produced:

And then the Slope function on default settings:

Even when changing the colours used it looks wrong.
How do I make a Topo map in ArcGIS Pro?


texts: QGIS - Attachment widget: copy file, use camera?

texts:
With a QGIS Project in QField, you can configure the Attachment Widget to store the picture you take (or select) in a folder together with the project.

Is it possible with QGIS to use a built-in camera with the widget? (rather than taking the picture, browsing to the directory)
Is it possible to copy a selected file to a folder stored with the project?


texts: Transformation not available in QGIS despite adding it to proj.db and adding GSB file to proj folder

texts:
In QGIS v3.10, I am trying to add a transformation from NAD83 (EPSG:2958) and NAD83 CSRS (EPSG:26917) for Ontario, Canada, using the GSB file from here.

I have added my grid transformation file here: C:\Program Files\QGIS 3.10\share\proj\ON83CSv1.gsb
(after reading this post)
I addded an entry to C:\Program Files\QGIS 3.10\share\proj\proj.db:grid_transformation.
When I go to Settings>Options>Default Datum Transformations, and select From EPSG:2958 to EPSG:26917, I do not see the transformation that I added.
I believe I am adding the entry to the correct table, because when I edit another entry in this table for an existing transformation, I see my edit within the QGIS GUI.

Does anyone know:

Why I am not seeing my transformation in QGIS despite it being in the database (do I need to edit other tables?)
Is there an easier way to add a transformation with a GSB file in QGIS?


texts: Calculating distance between latitude and longitude Python and CSV

texts:
I wanting to compare the distance between two sets of lat/long and
I need help inputting my lat/long fields into the Haversine formula.

I am able to print my two sets of lat/longs (print floats1, floats2, floats3, and floats4). I am unsure how to take these variables and use them as my inputs for the Haversine definition.

I was able to find similar code on GIS Stack Exchange, but I am fairly new to python and having some trouble.

import csv
import math

def haversine(lon1, lat1, lon2, lat2):
    """
    Calculate the great circle distance between two points
    on the earth (specified in decimal degrees).
    Source: http://gis.stackexchange.com/a/56589/15183
    """
    # convert decimal degrees to radians
    lon1, lat1, lon2, lat2 = map(math.radians, [lon1, lat1, lon2, lat2])
    # haversine formula
    dlon = lon2 - lon1
    dlat = lat2 - lat1
    a = math.sin(dlat/2)**2 + math.cos(lat1) * math.cos(lat2) * math.sin(dlon/2)**2
    c = 2 * math.asin(math.sqrt(a))
    km = 6367 * c
    return km

with open(r'PathtoFile.csv', 'rb') as f1:
    reader1 = csv.reader(f1, delimiter=',', quoting=csv.QUOTE_NONE)
    header1 = reader1.next()
    for row1 in reader1:
        floats1 = map(float, row1[16:]) #lon1
        floats2 = map(float, row1[17:]) #lat1
        floats3 = map(float, row1[18:]) #lon2
        floats4 = map(float, row1[19:]) #lat2
        print(floats1)
        print(floats2)
        print(floats3)
        print(floats4)
        print haversine(floats1[1],floats1[0],floats3[1],floats4[0])


texts: Using row.setValue() function with arcgisscripting?

texts:
I'm trying to set column value NULL. I tried different ways to do that, I used CalculateField_management and UpdateCursor but i did not get success. when i used CalculateField_management. It says success no error but no data updated. Then i used UpdateCursor and it always gave me below error. Please help me in this.

'2016-01-06 11:20:37.563000
 Updating field values of  C:/Bhaskar/CO_LOAD/CO JAN/SIR 13996/SIR 9999/TEST_STAGE_ELEC_22.sde
Cursor created.Row extracted.Curent field values is  TESTING. Updating field values of  ABC_TEST
Unable to add field:  (<type 'exceptions.RuntimeError'>, RuntimeError('ERROR
99999: Error executing function.',), <traceback object at 0x3DB01148>)>>> '

I am using bellow code in PYTHON 2.5.1 and ArcGIS 9.3.1.

import xlrd
import os.path
import DomainAssignmnetScriptForSde
import Tkinter
import tkFileDialog
import ctypes
import re
import sys, string, os, arcgisscripting
from datetime import datetime

Tkinter.Tk().withdraw() # Close the root window
#open the Dailog box to browse the database
databasePath =tkFileDialog.askopenfilename(title='Select Sde File')

file=open("ValidationReport.txt","a")
currentDateTime = datetime.now()
print str(currentDateTime)
file.write(str(currentDateTime) + '\n')
file.write("..............Validation Strated..........\n")

# Create the Geoprocessor object
gp  = arcgisscripting.create()

# Load required toolboxes...
gp.AddToolbox("C:\Program Files (x86)\ArcGIS\ArcToolbox\Toolboxes\Data Management Tools.tbx")

def UpdateField():
try:

    fclass= databasePath + "/SEU_GIS_ELEC.Cameras"
    gp.workspace=fclass
    #gp.workspace=databasePath
    #gp.CalculateField_management("SEU_GIS_ELEC.SecUGConductor", "OPERATINGVOLTAGE", None, "PYTHON")

    rows = gp.UpdateCursor(fclass)
    print "Cursor created."
    row = rows.Next()
    print "Row extracted."
    while row:
        print "Curent field values is  " + str(row.getValue("LOCATIONLINK"))
        #row.OPERATINGVOLTAGE = None
        #myVal='ABC_TEST'
        #row.setValue("LOCATIONLINK", myVal)
        row.setValue("LOCATIONLINK", None)
        #print row.OPERATINGVOLTAGE
        print "Updating field values of  " + str(row.getValue("LOCATIONLINK"))
        rows.updateRow(row)
        print "Updating field values of  " + str(row.getValue("LOCATIONLINK"))

        row = rows.Next()
    del row
    del rows
    file.write("........Processing finished.........."  +'\n')
except:
    print "Unable to add field: ", sys.exc_info()

UpdateField()
file.write("..............Validation finished..........\n")

file.close()


texts:
I have a line shapefile tc_line with 166 segments and a point shapefile tc_point with 8 points, and I loaded them into PostGIS database with DB manager:

and the following were their attribute table:

table tc_line (part)

table tc_point

Because these points are just nearby the line, not right on it, I try to use ST_ClosestPoint first, then use ST_Split to split the line.

@amball, I tested your code in the edited answer and the result showed:

SELECT ST_Split(ln.geom, ST_Collect(newpt_geom))
FROM tc_line AS ln,
    (SELECT ST_ClosestPoint(ST_Collect(ln.geom), pt.geom) AS newpt_geom
     FROM tc_line AS ln, tc_point AS pt GROUP BY pt.geom) t1 GROUP BY ln.geom;

The result was strange because it had 166 rows just the same as the rows original line had, which means that ST_split didn't work. So I loaded the result as new layer to check it.

But In QGIS I could not see it even I set its CRS, and as I wanted to export it as ESRI shapefile, error message showed up:

Export to vector file failed.
Error: creation of layer failed (OGR error:Geometry type of `Geometry Collection' not supported in shapefiles.
Type can be overridden with a layer creation option
of SHPT=POINT/ARC/POLYGON/MULTIPOINT/POINTZ/ARCZ/POLYGONZ/MULTIPOINTZ.
)

So I tried to split the line with 2 stages:

ST_ClosestPoint
ST_split

In the first stage I ran the code:

    SELECT
    tc_point.id as pointid,
    St_ClosestPoint(ST_collect(tc_line.geom), tc_point.geom) as closest,
    st_astext(St_ClosestPoint(st_collect(tc_line.geom), tc_point.geom)) as text
    FROM tc_point, tc_line GROUP BY tc_point.id;

The result were fine: (8 closest points)

and I successfully loaded as new layer tc_proj, and its attribute table:

Then I ran stage 2:

SELECT ST_Split(ln.geom, ST_Collect(tc_proj.geom))
from tc_line as ln, tc_proj GROUP BY ln.geom;

error showed up:

Splitting a Line by a GeometryCollection is unsupported

I am newbie to SQL and SQL open source such as PostGIS. How to work the above task successfully? I am working with PostGIS 2.2.1

Update#1

@amball, I tested your following code,

SELECT ST_Split(ln.geom, newpt_geom)
FROM tc_line AS ln,
    (SELECT ST_ClosestPoint(ST_Collect(ln.geom), pt.geom) AS newpt_geom
     FROM tc_line AS ln, tc_point AS pt
     GROUP BY pt.geom) t1
WHERE ST_Intersects(ln.geom, newpt_geom);

and its result,

I don't know why there were 4 rows in the result, and I could see nothing using Create a view, also, Load as new layer outputted nothing. Thus, I still can't determine whether result is correct.


texts: Designing small letters under words in QGIS

texts:
I want to design small letters under words in QGIS when I designed the labels (attach picture).

Is there any way to do it?


texts: Geoserver ImageMosaic: error connecting to WorldImage

texts:
I would like to use the ImageMosaic plugin to set up a Geoserver layer which refers to a WorldImage, just like in the official tutorial.

As a simple test case, I put a custom folder into my Geoserver data directory with one World Image (image + WorldFile + projection). When I specify the folder location to a ImageMosaic Store in Geoserver, it simply displays the error message:

Could not list layers for this store, an error occurred retrieving them: Argument "value" should not be null.

On the other hand, when I put my image into a default folder that comes with the openGeo suite (like file:data/dem/ with the Medform Digital Elevation Data), the thing works as it should.

So I suspect there is something going wrong with user rights, can anybody enlighten me?

PS: My Geoserver error log

    INFO [geoserver.web] - Getting list of coverages for saved store file:mosaic
    java.lang.RuntimeException: Could not list layers for this store, an error occurred retrieving them: Argument "value" should not be null.
    at org.geoserver.web.data.layer.NewLayerPageProvider.getItemsInternal(NewLayerPageProvider.java:123)
    ...
    Caused by: java.lang.IllegalArgumentException: Argument "value" should not be null.
    at org.geotools.util.SoftValueHashMap.ensureNotNull(SoftValueHashMap.java:140)
    ...


texts: How to identify coordinate system of a point in Guinea?

texts:
I need help to identify the coordinate system of the below point in Guinea

{"Name":"Simandou",
 "Latitude":0.29482446877975,
 "Longitude": 0.08190694324693,
 "k1":0.22528363047002,
 "k2":1.4976525821596}

I am not able to identify the unit of the Lat and Long but I looked for the Name on Google Earth and I found the equivalent in decimal degrees which is 8.656198, -8.887153?

How do I convert the Lat/Long of the point to decimal degrees?


texts: Extracting features from WMS layer in QGIS

texts:
I have a WMS layer in QGIS and I want to select and extract certain features in order to vectorize them and perform clipping with other vector layers.

Is there any way of doing this?

I do not have access to the WMS source files.


texts: Iterate zoom to layer then export as PDF

texts:
I'm looking to automate a process in ArcGIS that I believe needs Python. However, I have zero experience with Python only with ModelBuilder.

Here is the process:

I have 400 individual layers. I have a standard page setup in the layout view. I need to zoom to each layer in the layout data frame and then export each one as its own PDF. How would I automate this?


texts:
I am using infowindow to display the results of identifyTask that has two layers from my map service so it returns multiple features.
The infowindow has more button in each feature it identifies so that when I click the more button I get the data of current feature. How can I get the current selected feature? I am using this example.
 I want to get incase the popup identifies more than one feature to get current feature How can this be done.


texts: Create Address Locator with HERE/Navstreets Data

texts:
I just create Address Locator : Dual Range by using Streets feature class from HERE data.

But, in this HERE data, Streets Feature Class not have field for City and State. If you notice, Left City or Place. Right City or Place, Left State and Right State is blank.

In my finding, to get the full complete address in searching, we have to have City and State field.

Here data also have feature class Adminbndy1, Adminbndy2, Adminbndy3, Adminbndy4 and Adminbndy5.

But how to make use these feature classes to have relationship with Street.

How can I get information on City and State to get Full Address in searching box and complete my Address Locator?


texts:
I have a map of Ireland and I want to create a 10km x 10km grid over the whole country. I used the 'create grid' feature from the toolbox and put 10km as my horizontal and vertical spacing. However, when I insert a scale bar it shows me that the grid boxes are not 10km apart but are closer to 6.25km apart. (see pic 1)

The CRS I used for creating the grid is EPSG:3857 - WGS 84/ Psuedo-Mercator.

Does anyone know why this is happening and how I can make a grid that is actually 10km x 10km?

I have also noticed that when I change the scale bar units to 'map units' suddenly 10 grid squares do fit within 100000 units (see pic 2)


texts: Class Extension Problem - RelatedObjectChanged not firing

texts:
I'm trying to set up class extensions to have 3 feature classes communicate. I have feature classes A, B, and C. A has a simple, forward relationship with B, B has a simple, forward relationship with C. The goal is to have attribute updates in A force attribute updates in B, which in turn force attribute updates in C.

I have one class extension on feature class C that implements IRelatedObjectClassEvents2 RelatedObjectChanged, which listens for changes in feature class B. When I manually make an update to feature class B, as expected, the correct updates are made in feature class C. So that part is working.

Likewise, I have a class extension on feature class B that implements IRelatedObjectClassEvents2 RelatedObjectChanged, which listens for changes in feature class A. Also as expected, when I manually make updates to feature class A, the expected changes are made on feature class B. So this part is also working.

The problem I run into is that when I manually update Feature Class A, Feature Class B gets updated properly, but the extension on Feature Class C doesn't fire the RelatedObjectChanged event. As a result, my updates to feature class A aren't cascading all the way to feature class C.

What am I missing that is causing the event to not fire on feature class C?


texts: Is there any easy way to convert GRIB1 files to CSV in Linux? I've seen that supposedly GRIB2 can readily be converted to CSV, and I've tried to convert the Grib1 files to Grib2 using grib1to2 but its a nightmare as the documentation and software is no good. Is there any other way to get a CSV or any other easily readable format out of GRIB1?


texts: Performing land cover classification over whole country

texts:
I am doing land cover classification over the whole country of Djibouti.
I am downloading Sentinel 2 images for this purpose. I have realised doing land cover classification over a large area needs downloading several images as I was used to smaller regions prior to this.
My question is, is there a simpler way to do land cover classification over a whole country without downloading several images, doing atmospheric correction (and other preprocessing procedures) where such kind of procedures may be faced with time constraints.
I hear coding may assist but I am not good in this area. Anyone who has had such a similar workload can advise me. I am using QGIS 3.14


texts: In QGIS add new attribute table features from updated csv

texts: In QGIS I created a point layer from a csv using "add delimited text layer" (the comma icon). I set style and viewing scales the way I want for the layer. Since then the csv has been expanded with many new points. Is there a procedure for updating the attribute file with the new data?


texts: Where are imported Mosaic Dataset images stored on ArcSDE server?

texts:
I am trying to learn more about Mosaic Datasets and specifically how they interact with ArcSDE (SQL Server 2008). I read in How raster data is stored in a geodatabase that "...In the case of the mosaic dataset, this table remains empty, since the raster data for a mosaic dataset is not stored in the DBMS; instead, the raster data is referenced from an image file..."

Presumably this means that if I perform the Mosaic Dataset -> Add Rasters operation using files on the same machine as SDE then a file location is inserted into the murky depths of SDE. However, what if I connect to SDE remotely and Add Rasters from a different machine? I can only assume that SDE has a default directory to which the files are copied but I can't find anything about it, even in Mosaic datasets in a geodatabase in SQL Server.

I did find a table named <schema>.AMD_<mosaic name>_CAT which includes a field Uri ("A serialized Uniform Resource Identifier (URI) that is used to re-create the raster dataset") but it means nothing to me.

The test server in this case just has a single partition, but what if I wanted SDE to store mosaic dataset images on a separate drive for performance? How would I go about changing the default location, or even finding out what it is?

Any input very much appreciated.


texts: Rasters added to mosaic dataset only display at large acale

texts:
I have a set of 357 rasters in .jp2 that I've added to a mosaic dataset. All but 7 of them display properly, but the 7 only display at larger scales, around 1:6500 and under. I've tried adding them several times, deleting duplicates, overwriting, building overviews, etc., but can't figure out how to get the 7 to display normally like all the rest. Any clues?

Using ArcGIS 10.0, sp5, Windows 7 Enterprise.


texts: Add layers to a feature group with a function/loop

texts:
I have a bunch of markers that I would like to add to a feature group. Currently I'm doing it manually:

1) First adding all the markers (around 50)

var shelter1 = L.marker([54.962725, 12.548215], {icon: shelterIcon}).on('click', clickZoom);
var shelter2 = L.marker([54.944250, 12.515722], {icon: shelterIcon}).on('click', clickZoom);
var shelter3 = L.marker([55.008244, 12.299891], {icon: shelterIcon}).on('click', clickZoom);
...

2) Then I define a feature group:

var shelterMarkers = new L.FeatureGroup();

3) And finally I add the markers in my feature group:

shelterMarkers.addLayer(shelter1);
shelterMarkers.addLayer(shelter2);
shelterMarkers.addLayer(shelter3);

Instead of doing 3) manually, I would like to create a function that does it automatically. Something like this:

var i = 1;
text = "";

while (i < 52) {
    text += "shelterMarkers.addLayer(shelter" + i + ")" + ";" + " ";
    i++;
    console.log(text)
}

When I look at the output, it looks pretty much like 3) (which successfully adds markers to the map) but for some reason the markers are not added to my map.


texts: Offline editing solution for small GIS project with multiple editors?

texts:
I have a specific use-case that I cannot seem to find any solutions for online, so any and all advice would be most welcome. Some basic details:

I am part of a small academic organization (a research lab), and our research leans heavily on GIS. We are all field scientists who specialize in our particular discipline, NOT computer scientists, so while we are all competent at the analytical side of desktop GIS, as soon as we get into the weeds of database programs, web maps, online servers, etc. we are lost. Because we are a bunch of broke academics, we mostly use QGIS (currently 2.18 but switching to 3.x as soon as there's an LTR) and I'd like to keep it that way.

Our research requires us to be offline (or with only tenuous and expensive internet access) for MONTHS at a time. Most GIS multi-user solutions seem to presume that people are coming back to the office and syncing their data every night. For us, this is not realistic.

As such, currently I spend about a week after every field campaign manually integrating 3 or 4 updated/edited versions of the same shapefile to reflect everyone's changes. We also have an absurd proliferation of shapefiles that I need to get under control.

It looks like Spatialite and GeoPackage both allow multiple geometries to live in the same file/database, so that is very appealing for dealing with the second problem (shapefile proliferation). However, how do these file formats handle data syncronization?

Is there a way for us to work with a single, simple database of features (e.g. Spatialite) while in our home GIS lab, then take 3 or 4 local copies of that database to our field sites and make a few hundred edits per user, and then sync all of those edited versions back together at the end of the field season?


texts:
I recently use QGIS 3.0.0, in Windows 8. Today, I work with symbol. When I use SVG symbol, the symbol selector window become so long, even under the taskbar, I cannot resize or drag up the window (it won't allow me). No, it's not just covered by taskbar, it's under the monitor boundary.

Previously I used QGIS 2.18 .12 and never met this problem.

Maybe this is just silly question, but I need help.


texts:
In my application I have a gis file to render an (java BufferedImage) image of countries with geotools which is sent directly to the browser. The colours of the countries (plus name and additional icons) are derived from the database. Depending on the database values, different rules and visualizers are applied to the Features from the gis file.

public BufferedImage renderOverviewImage(FileDataStore dataStore) {

    FeatureSource featureSource = dataStore.getFeatureSource();
    ContentFeatureCollection featureCollection = (ContentFeatureCollection) featureSource.getFeatures();
    SimpleFeatureIterator sfi = featureCollection.features();
    Map<FeatureId, PojoType> featureIdMap = new HashMap<FeatureId, PojoType>();

    BoundingBox boundingBox = null;

    try {
        while (sfi.hasNext()) {
            SimpleFeature simpleFeature = sfi.next();

            // this part here is actually a bit more complicated, the map will be enriched
            // with pojos which determine the eventual colour
            PojoType pojo = gatherInformationFromDatabase(simpleFeature);
            featureIdMap.put(simpleFeature.getIdentifier(), pojo);
        }
    } finally {
        sfi.close();
    }

    // Create a map context and add our shapefile to it
    MapContent map = new MapContent();

    // Create a basic Style to render the features
    Style mapStyle = createMapStyle(featureIdMap);

    // Add the features and the associated Style object to
    // the MapContext as a new MapLayer
    map.layers().clear();
    Layer someLayer = new FeatureLayer(featureSource, mapStyle);
    map.addLayer(someLayer);

    boundingBox = map.getMaxBounds();

    BufferedImage image = renderImage(map, 600, boundingBox);
    map.dispose();

    return image;

}

private Style createMapStyle(Map<FeatureId, PojoType> featureIdMap) {

    List<Rule> ruleList = new ArrayList<Rule>();

    for (FeatureId featureId : featureIdMap.keySet()) {
        Rule rule;
        PojoType pojo = featureIdMap.get(featureId);

        // create a partially opaque outline stroke
        Stroke stroke = styleFactory.createStroke(filterFactory.literal(outlineColor), filterFactory.literal(1), filterFactory.literal(0.5));

        // create a partial opaque fill
        Fill fillstyleFactory.createFill(filterFactory.literal(fillColor), filterFactory.literal(0.5));

        PolygonSymbolizer polygonSym = styleFactory.createPolygonSymbolizer(stroke, fill, null);
        Rule rule = styleFactory.createRule();
        rule.symbolizers().add(polygonSym);

        // add other symbolizers here...

        // tie this rule to the featureId
        rule.setFilter(ff2.id(featureId));
        ruleList.add(rule);

    }

    FeatureTypeStyle fts = styleFactory.createFeatureTypeStyle(ruleList.toArray(new Rule[] {}));
    Style style = styleFactory.createStyle();
    style.featureTypeStyles().add(fts);

    return style;
}

Now I need to add geo-referenced icons to the image map. The coordinates are also stored in the database.

Setting the AnchorPoint does not help because the anchorpoint is only a point within the rendered image.

The comments in org.geotools.styling.Graphic class say:
"The 'hot spot' to use for positioning the rendering at a point must either be inherent from the external format or be defined to be the  "central point" of the graphic."

Is it possible to add a PointVisualizer with Graphic and then set the location coordinates for it? Or would I have to create artificial Features with fixed coordinates in order to render an image graphic at specific location?


texts: Converting Shapefile points to Lat/Long?

texts:
Quick Overview of the Project: I have 9000 records of individuals, each with a coordinate pair in a csv. I also have a shapefile/folder(dbf, prj, sbn, sbx, shp, shx) of municipalities within the range of coordinates. I'd like to know which municipality(shape) in the shapefile corresponds to each record in the csv

Problem: When I look at the shapefile in python, the grid isn't lat/long. I subscribed to a trial version of arcGIS to be able to load up the coordinates and look at them, and they map correctly, so I must just be missing something. Attached is my code and the map printed by Python.

How can I transpose my csv coordinates to the map or vice versa?

import shapefile as shp
import matplotlib.pyplot as plt
from shapely.geometry import Point

import csv

n = 0
m = 1

coorDict = {}
coor = ''
fileName = 'TestLatLong'

sf = shp.Reader("wcmun/wcmun")
sfRec = sf.records()
plt.figure()
for shape in sf.shapeRecords():
    x = [i[0] for i in shape.shape.points[:]]
    y = [i[1] for i in shape.shape.points[:]]
    plt.plot(x,y)
plt.show()
plt.figure()

n = 0
#Then create a dictionary of lat/longs from 'TestLatLong.csv'

The output of all shape records and data in the first shape record above
For comparison, all the data should be roughly in a rectangle from
(40.881321, -73.926644) to (41.442483, -73.542874), and as you can see this map isn't.


texts: Geoserver proxying custom tiles through a WMS

texts:
I have a set of tiles that I was given via disk (So no way to recreate them, they are already generated/cached).  There is a clear pattern/style to the naming of the tiles (zoomLevel-x-y.gif).  The bounds are actually configured the same way that "Google" organizes their tiles(from top left, not bottom-left. ugh).  More information on this specific format and how I request them in OpenLayers can be found http://trac.osgeo.org/openlayers/wiki/UsingCustomTiles

My question comes into play as I want to provide these tiles via Geoserver.  I have a client who is only able to access data and layers via wms calls (don't ask). The file system fetching that openlayers uses won't work for him (as well as being able to use a top-left origin instead of bottom-right).

So, from my limited knowledge of geoserver, it seems like it should be able to front these tiles to the users via a WMS server, where it takes the request, translates it to a file path/url, and goes from there.  Sort of a proxy functionality.  Am I completly off base here?  And if not, how can I even start to go about this?

edited for clarity, and a bit more background information

Thanks!


texts: Value NAN in netcdf4 file for lat&lon using python

texts:
I have two file netcdf4 ( latlon and other file without latlon)

i open the 2 file in the same code and i choose one variable and the two variable lon and lat
when i write lon i have the value nan i don't know why

the 2 files have same dimensions ( size) and same dimensions for variables all i need is to show the value of lat&lon when i open the 2 file how can i do it ?

In [54]: import numpy as np
    ...: import netCDF4
    ...:
    ...:
    ...: fic1='/data/latlon_+000.0_globe.nc'
    ...:
    ...: nc1 = netCDF4.Dataset(fic1,'r')
    ...:
    ...: lon = nc1.variables['lon'][:]
    ...: lat = nc1.variables['lat'][:]
    ...:
    ...: fic='/data/Mmultic3kmNC4_msg04_201905080200
    ...: .nc'
    ...:
    ...:
    ...: nc = netCDF4.Dataset(fic,'r')
    ...: IR=nc.variables['IR_108'][:]
    ...:
    ...:

In [55]:

In [55]: lon
Out[55]:
masked_array(
  data=[[nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan],
        ...,
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan]],
  mask=False,
  fill_value=1e+20,
  dtype=float32)

In [56]: lat
Out[56]:
masked_array(
  data=[[nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan],
        ...,
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan]],
  mask=False,
  fill_value=1e+20,
  dtype=float32)

In [57]:

how i can get value of lat & lon ?


texts:
I am not being able to extract or interact with the svg tiles data on mouse hover/click using leaflet.
I have a tile server that serves svg tiles eg:https://sometileserver.com/roads/{z}/{x}/{y}.svg, the svgs returned in this have some attributes embedded in each path drawn on them. I am using leaflet in my vue app to show the tiles.
example of how I am loading my tiles:
const tileUrl = "http://{s}.tile.openstreetmap.org/{z}/{x}/{y}.png";
const tiles = L.tileLayer(tileUrl, {
        attribution
      });
      var myMap = L.map("mapContainer", {
        layers: tiles
      }).setView([0,0], 12);

const roads = L.tileLayer('https://sometileserver.com/roads/{z}/{x}/{y}.svg', {
        attribution
      }).addTo(map)

At this point leaflet shows the svg tiles on top of the OSM base map which is great. however these svg tiles has many roads drawn on them and each such road/path has some attributes about it embedded into.
an example of one such svg tile:
<ns0:svg xmlns:ns0="http://www.w3.org/2000/svg" height="256pt" version="1.1" viewBox="0 0 256 256" width="256pt">
<ns0:g id="surface4396">
<ns0:path d="M 254 13 L 252 13 " road_status="[1.04,4.09]" road_id="abc" road_name="unknown" style="fill:none;stroke-width:3.0;stroke-linecap:butt;stroke-linejoin:miter;stroke:rgb(0%,100%,0%);stroke-opacity:1;stroke-miterlimit:10;" />
</ns0:g>
</ns0:svg>

My goal is to show those attributes/values(eg:road_id,road_name) when I hover over each road/path drawn on this tile.


texts:
I'm new to qgis and in here I want to find a path between two selected points on the map(Roads-vector layer). The points are selected by the user, using mouse clicks.

So here I used the astar algorithm to find path between two points.

*********************************astar.py************************************
import heapq

class AStar(object):
      def init(self, graphAstar):
                self.graphAstar = graphAstar

def heuristic(self, node, start, end):
    raise NotImplementedError

def search(self, start, end):
openset = set()
    closedset = set()
    current = start
    openHeap = []
    openset.add(current)
    openHeap.append((0,current))
    while openset:
    temp = heapq.heappop(openHeap)
    current = temp[1]
        if current == end:
            path = []
            while current.parent:
                path.append(current)
                current = current.parent
            path.append(current)
            return path[::-1]
        openset.remove(current)
        closedset.add(current)
        for node in self.graphAstar[current]:
            if node in closedset:
                continue
            if node in openset:
                new_g = current.gg + current.move_cost(node)
                if node.gg > new_g:
                    node.gg = new_g
                    node.parent = current
            else:
                node.gg = current.gg + current.move_cost(node)
                node.H = self.heuristic(node, start, end)
                node.parent = current
                openset.add(node)
        heapq.heappush(openHeap, (node.H,node))
    return None

class AStarNode(object):
    def __init__(self):
        self.gg = 0
        self.H = 0
        self.parent = None

    def move_cost(self, other):
        raise NotImplementedError

*****************************astar_grid.py*******************************

from astar import AStar, AStarNode
from math import sqrt

class AStarGrid(AStar):
    def heuristic(self, node, start, end):
        return sqrt((end.x - node.x)**2 + (end.y - node.y)**2)

class AStarGridNode(AStarNode):
    def __init__(self, x, y):
        self.x, self.y = x, y
        super(AStarGridNode, self).__init__()

    def move_cost(self, other):
        diagonal = abs(self.x - other.x) == 1 and abs(self.y - other.y) == 1
        return 14 if diagonal else 10

and in the main code, the following method is used to create graph from vector layer.

****************************plugin.py************************************

def make_graph(self, mapinfo):
    nodes = [[AStarGridNode(x, y) for y in range(mapinfo['height'])] for x in range(mapinfo['width'])]
    graphAstar = {}
    for x, y in product(range(mapinfo['width']), range(mapinfo['height'])):
        node = nodes[x][y]
        graphAstar[node] = []
        for i, j in product([-1, 0, 1], [-1, 0, 1]):
                if not (0 <= x + i < mapinfo['width']): continue
                if not (0 <= y + j < mapinfo['height']): continue
                graphAstar[nodes[x][y]].append(nodes[x+i][y+j])
    return graphAstar, nodes

And I called that method in FindRoutes method..

def findRoutes(self):
# QMessageBox.information( self.iface.mainWindow(),"Info", "in findRoutes function" )
# vl = qgis.utils.iface.mapCanvas().currentLayer()
vl=self.canvas.currentLayer()
director = QgsLineVectorLayerDirector( vl, -1, '', '', '', 3 )
properter = QgsDistanceArcProperter()
director.addProperter( properter )
# crs = qgis.utils.iface.mapCanvas().mapRenderer().destinationCrs()
crs = self.canvas.mapRenderer().destinationCrs()
builder = QgsGraphBuilder( crs )

global x1
global y1
global x2
global y2
pStart = QgsPoint( x1, y1 )
pStop = QgsPoint( x2, y2 )
graphAstar, nodes = self.make_graph({ "width": 8, "height": 8 })
paths = AStarGrid(graphAstar)

start, end = ??
path = paths.search(start, end)

How do I pass the start and end coordinates to the function above?

Passing them just as coordinates (start, end = pStart, pStop) does not work.

How do add them to the graph created as nodes?


texts: How to crop a vector shapefile according to raster?

texts: I am trying to clip a vector shapefile- points according to a DEM in QGIS. I used GDAL function and polygonized the DEM. Then I ran the Vector-Geoprocessing tools-clip in order to crop the vector data points. The whole QGIS crashed couple of times while doing this and I have never achieved any results. Could you please suggest any alternative steps for me? Many thanks!


texts:
I would like to get a center point approximation of a set of polygons so I can do a NN lookup to use for assignment instead of a point in polygon query.

Lets say I have this set of polygons:

I would like to calculate a set of centerpoints that approximates the polygon when used in a NN lookup:

Does anybody have a clue to achieve this (Python/postgis preferred)? And does anybody have a pointer to any algorithms that can achieve this?


texts: Compare fields of two tables and calculate field

texts:
Can you help me to compare fields in two tables?
I have two tables, which contain ID and Measure fields.
I need to take stroke by stroke values from field1(table1) and find all values in field1(table2) which are less than sought-for value. And then take largest of them and put it on stroke field2(table1).
Example showed on picture below.

It's will be great, if it would be python script.


texts: QGIS 2.10 crashes while working with Qgis2web plugin

texts: Currently I am working on Qgis2Web for editing and generating web maps from my vector data. But after opening the plug in window. QGIS got crashed.


texts:
I'm trying to change the name of the layers in the tree.
This is part of my code:

var layer_group1 = new OpenLayers.Layer.WMS("layer_group1",
                    url,
                    {
                        LAYERS: ['namespace:layer_a','namespace:layer_b','namespace:layer_c'],
                        format: "image/png",
                        transparent: "true",
                        projection: 'EPSG:4326'
                    },
                    {

                        buffer: 0,
                        displayOutsideMaxExtent: true,
                        isBaseLayer: false,
                        displayInLayerSwitcher: false,
                        yx : {'EPSG:4326' : true}

                    }
                );

var layer_group2 = new OpenLayers.Layer.WMS("layer_group2",
                    url,
                    {
                        LAYERS: 'namespace:layer_d',
                        format: "image/png",
                        transparent: "true",
                        projection: 'EPSG:4326'
                    },
                    {
                        opacity:0.5,
                        buffer: 3,
                        displayOutsideMaxExtent: true,
                        displayInLayerSwitcher: false,
                        yx : {'EPSG:4326' : true}

                    }
                );

    mapPanel = new GeoExt.MapPanel({
        border: true,
        region: "center",
        // we do not want all overlays, to try the OverlayLayerContainer
        map: mapa,
        center:[-6613943,-3793944],
        zoom: 7,
        layers: [ capaBase, layer_group1, layer_group2, gmap,gsat,vector]
    }
)

I get the tree:

layer_group1
   |_ namespace:layer_a
   |_ namespace:layer_b
   |_ namespace:layer_c

layer_group2
   |_namespace:layer_d

But I want to expect:

layer_group1
   |_ Layer a
   |_ Layer b
   |_ Layer c

layer_group2
   |_ Layer d

How can I change the name of the layer in the tree?

PS: don't hesitate if you need more code.


texts: Why can't I run QGIS android on my Samsung Galaxy Tab S2?

texts:
I wanted to be able to run QGIS on my Galaxy Tab S2.  Play store says that the device is not compatible.  But the support page lists the Tab S (first one) as supported, but doesn't mention the Tab S2.

What is it about my device that prevents the download via play store, and is there a workaround?

This was one of the reasons I bought the device, but I was mistaken in the support wiki between a Samsung Galaxy S2 (a phone) and a Tab S2?


texts: Archiving activation in ArcGIS Web AppBuilder

texts:
When I allow archiving, ArcGIS also keeps old records of the data. When I publish this application, can I deactivate the current record of any detail and activate the old record?
For example;
I have bus stops on the routes. In the past, the route code of the A station was x101. With the update made, the current route code of the A station is x102. Can I activate the code of this stop as x101 in the web environment after a certain day?


texts: How to convert Common Alerting Protocol (CAP)  to WKT

texts: I've come across some data being sent in the CAP (http://docs.oasis-open.org/emergency/cap/v1.2/CAP-v1.2.html) format, which in turn contains an <area> element with either a polygon or circle definition.  I want to ultimately get this into a SQL Server geometry, but I think parsing it into WKT would facilitate that.  Are there any common algorithms already out there that deals with this CAP format geometries, or am I on my own for parsing it out?


texts: best setting for exporting a map ( jpeg ) that will be used in an A4 report

texts:
Often i find some of the maps i make will make there way into reports. This can be either taken up a whole page of A4. Or being 50 percent of the Page or sometimes. 25.

Do people have any page settings they would use


texts: How to calculate a missing (one or more) lat/long value of a bounding rectangle?

texts:
Is it possible to calculate an unknown (one ore more) latitude/longitude value of a bounding rectangle. The log file that I need to process provides a pair of lat/long values for Upper left and Bottom right of a bounding rectangle. However sometimes (don't know why) the logging source corrupts the data resulting in one or more of coordinate values logged as 0000.00. A snapshot of what data might look like is given below.

My questions:

Is it possible to to figure out the missing coord point if only of the values is missing.
Is it possible to to figure out the missing coord point if more than one values is missing.

PS: I may have access to GPS heading and speed. Any help will be greatly appreciated.


texts:
I used the following code to label the line features(roads) in a layer.It works in the console, but not by python code(application) It works for the points too.

What did I do wrong?

 def labell(self):
    layer = self.iface.mapCanvas().currentLayer()
    palyrr = QgsPalLayerSettings()
    palyrr.readFromLayer(layer)
    palyrr.enabled = True
    palyrr.fieldName = 'name'
    palyrr.setDataDefinedProperty(QgsPalLayerSettings.Size,True,True,'08','')
    palyrr.writeToLayer(layer)
    layer.commitChanges()
    self.iface.mapCanvas().refresh()


texts: QGIS projected coords vs lat/longs in coordinate capture

texts:
My project CRS is 23031 (ED50,UTM 31N).
Lat long coords displayed in Coordinate Capture are WGS84. How can I get ED50 lat longs when using coordinate capture? Is there a GUI coordinate transform utility that will go back and forwards between projected coords and lat longs in my choice of datum or spheroid?


texts:
I'm working on OpenStreetMap data on Postgis. I have in a table "roads_polygons_label" containing all the polygons and roads in a given area and for each of them a label saying what it is (e.g. building, residential, grassland, according to OSM database). The problem is that some polygons are nested. For instance polygon A labeled as "residential" contains polygons B and C both labeled as "building". I want to separate A from B and C  to get three records:

polygon A\(BuC) labeled as "residential"
polygon B labeled as "building"
polygon C labeled as "building"

I managed to get a table "overlapping" where I join "roads_polygons_label" with itself on an ST_Intersects query to have for each larger polygon a row for each polygon it contains. So in my basic example it would give two records for polygon A:

A - B
A - C

Then I wanted to use ST_Union to make a big multipolygon made of all polygons nested in one larger polygon, so here a polygon made of BuC. But I would need to loop on each larger polygon and I don't know how to do that. I can create a table "large_polys" listing all polygons containing smaller polygons, but then I don't know how to join "overlapping" with "large-polys" in the loop.
I have read the following posts but they did not provide a solution for me problem:

Separate polygons based on intersection using PostGIS

SQL Server - Select all polygons nested within a larger polygon

Edit : illustration of what I mean:

Maybe "nested" is not the right term.
A, B and C have different osm ids. I also have several polygons with the same id and different area, but that's another issue.


texts: Could someone suggest a good GeoExt Tutorial?

texts:
In my research I've found these two, which are only explaining the basics.

Developing OGC Compliant Web Applications with GeoExt
GeoExt Workshop (Pdf)

And of course the native (but short) tutorials and examples on the GeoExt website.

Can anyone know about other resources?


texts: How to connect Postgis database with/in Remote Servers...?

texts:
I am preparing a 'Web based map' application using 'Geoserver-Openlayers-Postgis'.

I am trying to connect my Postgis database(located at 10.50.110.11 location) with the following command from 10.40.48.22 location....

$host = 'localhost/10.50.110.11';
$port = '5432';
$database = 'template_postgis_20';
$user = 'postgres';
$password = 'user';

$connectString = 'host=' . $host . ' port=' . $port . ' dbname=' . $database .
    ' user=' . $user . ' password=' . $password;

$link = pg_connect($connectString);
if (!$link)
{
    die('Error: Could not connect: ' . pg_last_error());
    echo "Could NOT connect POSTGRES...";
}

This is yielding a blank response.


texts: Cell size changing when cropping a raster

texts:
I am cropping a large raster with the extent of a smaller one.

The printout for large:

class       : RasterLayer
dimensions  : 6327, 12697, 80333919  (nrow, ncol, ncell)
resolution  : 463.3127, 463.3127  (x, y)
extent      : -9632735, -3750053, 4869417, 7800796  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=sinu +lon_0=0 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs
names       : large
values      : 0, 4976998  (min, max)

the printout for small:

class       : RasterLayer
dimensions  : 2400, 2400, 5760000  (nrow, ncol, ncell)
resolution  : 463.3127, 463.3127  (x, y)
extent      : -4447802, -3335852, 5559753, 6671703  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=sinu +lon_0=0 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs
names       : small

Therefore, I want the cropped raster to have the same dimensions, resolution and extent of small.

I am attempting this like so:

library(raster)

large <- raster ("C:/large.tif")
small <- raster ("C:/small.tif")

clip_extent = extent(ras)

#clip
clipped <- crop(large, small)
extent(clipped) <- extent(small)

but the output changes the resolution and dimensions:

class       : RasterLayer
dimensions  : 2400, 1506, 3614400  (nrow, ncol, ncell)
resolution  : 738.347, 463.3127  (x, y)
extent      : -4447802, -3335852, 5559753, 6671703  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=sinu +lon_0=0 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs
values      : 187927, 2240075  (min, max)


texts: Updated version of QGIS coding and compilation guide?

texts:
is there an updated version of QGIS coding and compilation guide? The current version is clearly outdated. I have been trying to compile the first example, 1_hello_world_qgis_style, though to no avail regrettably. I then adjusted my code as kindly suggested by Tiggey here. I suppose something must be wrong with the include files as I'm still unable to compile.

My platform is Windows/cygwin; I'm using the OSGeo4W64 packages. This is the code I'm using:

//
// QGIS Includes
//
#include <qgsapplication.h>
#include <qgsproviderregistry.h>
#include <qgssinglesymbolrendererv2.h>
#include <qgsmaplayerregistry.h>
#include <qgsvectorlayer.h>
#include <qgsmapcanvas.h>
//
// Qt Includes
//
#include <QString>
#include <QApplication>
#include <QWidget>

int main(int argc, char ** argv)
{
  // Start the Application
  QgsApplication app(argc, argv, true);

  QString myPluginsDir        = "/cygdrive/c/OSGeo4W64/apps/qgis";
  QString myLayerPath         = "/cygdrive/c/Users/defloriv/Documents/qgis/qgis_sample_data/shapefiles";
  QString myLayerBaseName     = "regions";
  QString myProviderName      = "ogr";

// Instantiate Provider Registry
QgsProviderRegistry * preg = QgsProviderRegistry::instance(myPluginsDir);

QgsVectorLayer * mypLayer =  new QgsVectorLayer(myLayerPath, myLayerBaseName, myProviderName);

QgsSingleSymbolRendererV2 *mypRenderer = new QgsSingleSymbolRendererV2(QgsSymbolV2::defaultSymbol(mypLayer->geometryType()));

QList <QgsMapCanvasLayer> myLayerSet;

mypLayer->setRendererV2(mypRenderer);
if (mypLayer->isValid())
{
    qDebug("Layer is valid");
}
else
{
    qDebug("Layer is NOT valid");
}

// Add the Vector Layer to the Layer Registry
QgsMapLayerRegistry::instance()->addMapLayer(mypLayer, TRUE);
// Add the Layer to the Layer Set
myLayerSet.append(QgsMapCanvasLayer(mypLayer, TRUE));

// Create the Map Canvas
QgsMapCanvas * mypMapCanvas = new QgsMapCanvas(0, 0);
mypMapCanvas->setExtent(mypLayer->extent());
mypMapCanvas->enableAntiAliasing(true);
mypMapCanvas->setCanvasColor(QColor(255, 255, 255));
mypMapCanvas->freeze(false);
// Set the Map Canvas Layer Set
mypMapCanvas->setLayerSet(myLayerSet);
mypMapCanvas->setVisible(true);
mypMapCanvas->refresh();

return app.exec();
}


texts: How to point out a place in iframe using latitude and longitude

texts:
I want point out a place using latitude and longitude in following frame,

<iframe width="75%"
        height="550"
        frameborder="0"
        scrolling="no"
        marginheight="0"
        marginwidth="0"
        src="http://www.openstreetmap.org/export/embed.html?bbox=79.84942674636841%2C6.899900944350287%2C79.85650777816772%2C6.90436374803887&amp;layers=ND"
        style="border: 1px solid black">
</iframe>

How should I point out?


texts: Using ArcPy to retrieve maps posted on an ArcGIS portal?

texts:
I'm using ArcPy and python 3.5 to connect to a few different instances of ArcGIS portal.

The goal is to produce a report listing the maps, layers, and datasources that each portal instance has to compare differences between them.

My problem is, I can't figure out how to retrieve or list the maps from python.  I can connect to the portals and output the data from the GetPortalDescription function but I'm not sure where to go next.

import arcpy
import arcgis
import json
from urllib.parse import urlparse
from arcgis.gis import GIS

def parse_url(url):
    parsed_uri = urlparse( url)
    hostname=parsed_uri.hostname
    return hostname

def output_portal(portalUrl, portalAdminUser, portalAdminUserPassword):
    #gis = GIS(portalUrl, portalAdminUser, portalAdminUserPassword)
    gis=GIS()
    url = arcpy.GetActivePortalURL()
    print("Current Portal: "+url)
    domain=parse_url(url)
    portalDesc =arcpy.GetPortalDescription()
    outfile=open(domain+".json","w")
    outfile.write(json.dumps(portalDesc, indent=4, sort_keys=True))
    outfile.close()
    # get all maps posted to the portal
    # for each map print all the layers in the portal


texts:
I have a script that shows certain features that are within a user-set distance of a point. I would like, however, to have those features that are "in the buffer zone" to be listed at the end when the layer is added in the form of a table. I know MakeTableView_management, but that only adds a table to the Table of Contents and doesn't open it so that doesn't really help me. Listing out rows also doesn't help since it doesn't appear clear in the message window. Therefore, is there some way to make a Table View, and have all the rows and columns of that table be displayed as a message? Or even a way to add the table and program it to automatically open?

This is what I have tried, but all it did was list the name of my table and not show anything as a message:

SourceCount = arcpy.MakeTableView_management("Selected_Sources", "Source_Table")
if arcpy.Exists("Selected_Sources"):
    arcpy.AddMessage("Creating a source layer to add to the map...")
    SourceLyr = arcpy.mapping.Layer("Selected_Sources")
    arcpy.mapping.AddLayer(df, SourceLyr, "TOP")
    arcpy.ApplySymbologyFromLayer_management("Selected_Sources", "E:/GIS/File_Organization/Data/Sources/Source_Colors.lyr")
    arcpy.AddMessage("Sources surrounding your point have been added to the map and are listed below.")
    arcpy.AddMessage("{0}".format(SourceCount))
else:
    pass


texts: Seeking a gvSIG plugin to go from shp to kml to cityGML

texts:
I'm looking for a gvSIG plugin to convert .shp to .kml and then to .cityGML.
I'm working with gvSIG + 3D extension and I have .shp format data but I need it in .cityGML format.


texts: QGis polygon to points different CRS problems

texts:
I'm trying to join the values from polygons onto a grid of points. I had this problem solved here by @Underdark , however a few layers weren't working. It seemed they were all British National Grid CRS [EPSG:27700] whereas my project & points CRS is WGS84. The error was:

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

So I tried re-saving the BNG layers as WGS84 as per @Underdark's guide here, the result of which looks fine, however after then adding polygons to points and saving the output as CSV (as I did successfully for the rest), there are no data in the field column.

Can anyone think why this might be? Or what I could do to solve it?

Thanks in advance!


texts: Making join with SQL Server table using ArcPy?

texts:
I'd like to join .shp to table by using AddJoin_management but I'd like to make join to an SQL Server table (not to *.dbf).

Is it possible to do?

It should look something like this:

AddJoin_management ("my_layer.shp", "ID", SQL-SERVER-TABLE, "ID")


texts: Joining table data with expression in QGIS

texts:
I use QGIS 3.4 as a beginner.

I have a vector layer 'Address' with the fields:

"Id"(int)
"postalcode" (str)
"straatname"(str)
"housenr"(int)

and a CSV-file with the fields

"section" (int)
"postalcode" (str)
"housenr_low" (int)
"housenr_high" (int)

I want to add the field "section" from the CSV -ile to the 'Address' layer where Address.postalcode = CSV.postcalcode and Address.housenr >= CSV.housenr_low and Address.housenr <= CSV.housenr_high.

I tried "join with the attribute table" and "Join attributes by field value" without success.

How can I do this?


texts: How to place two data driven pages on one map layout?

texts:
I need to place two data frames (DF) with different data driven pages (DDP) on one map layout. i have two different Strip Map Index Features (SMIF) for each data frame.
When i turn next page in DDP series for DF1 this work only for DF1. If i want next DDP in DF2 i need to setup DDP and choose SMIF2.
Is it possible on one map?
Dynamic text work only current DDP this is problem too.
P.S. I don`t have Map Atlas Extension.


texts: Use of File Geodatabase uploaded to ArcGIS Online

texts:
I have a few feature classes in a file geodatabase. Instead of uploading the feature classes to ArcGIS Online individually as shapefile, I was trying to upload the entire zipped geodatabase to ArcGIS Online.

MyFGD.gdb
  .FeatCls1
  .FeatCls2
  .FeatCls3
  .FeatCls4

The File Geodatabase was uploaded successfully and it was honored File Geodatabase "MyFGD", and only one AGOL Feature Layer with the same name as the FGB was also created and hosted "MyFGD".

However, I could not Open the hosted Feature Layer to Map Viewer. It failed with a "Invalid URL" error msg.

How can we add the feature classes in a file geodatabase uploaded to ArcGIS Online to Map Viewer?

Are those feature classes published and hosted on ArcGIS Online individually as AGOL Feature Layer?


texts: Find a geographic concentration of values?

texts:
I have a database table (PostGIS/Postgres) of points (dozens to thousands), each of which have a Value and a Position (latitude/longitude in a geography datatype).

I am trying to figure out whether I have too much Value concentrated in any one place.

More specifically, I am trying to find an efficient approach to figure out if any circle of radius r, encloses a total summed Value over a threshold MAX_ACCEPTABLE_VALUE_IN_CIRCLE.

This can be done in PostGIS SQL or the data can be pulled and an algorithm run locally in Python, it doesn't really matter.

Are there any preexisting approaches or optimizations that might help?


texts: QGIS Equivalent of MapInfo's 'Combine Objects using Column'

texts:
MapInfo has a function that takes a table, and creates a new table populated combines them by a specified value.

So, combining on the name field, and specifying that the Area field should be summed...

Name                  Area
Ben's Field           100
Ben's Field           50
Tom's Field           25

Becomes

Name                  Area
Ben's Field           150
Tom's Field           25

How can one achieve this in QGIS?


texts:
I know little about geographic data.  I am trying to get a file of elevation data for the United Kingdom so that I can import it into a CAD program and convert it to an *.stl or *.dxf format.  I intend to use a 3D CNC router to make a physical model of this data.

Can anyone point me anywhere or explain how I would acquire data in this format of data I can convert?


texts: How to set high value clustered feature to the front?

texts:
Here is the clustered vector layer I just added to the map.

var defaultStyle = new OpenLayers.Style({
        strokeColor     : "#FFFFFF",
        strokeOpacity   : 1,
        strokeWidth     : 2,
        fillColor       : settings.color,
        fillOpacity     : 0.9,
        cursor          : "pointer",
        pointRadius     : "${point_radius}",
        pointerEvents   : "visiblePainted",
        label           : "${label}",
        fontColor       : "#000000",
        fontSize        : "12px",
        fontFamily      : "verdana",
        fontWeight      : "bold",
        labelAlign      : "center"
    }, {
        context: {
            label: function(feature) {
                var sumValue = 0;
                if(feature.attributes.count>1) {
                    for (var i=0; i<feature.cluster.length; i++) {
                        sumValue++;
                    }
                    return sumValue;
                }else{
                    return 1;
                }
            },
            point_radius: function(feature) {
                var pix = 10;
                if(feature.cluster) {
                    pix =Math.min(feature.attributes.count, 30) * 2.5 + 8;
                }
                return pix;
            }
        }
    });

    var onSelectStyle = new OpenLayers.Style({
        strokeColor     : "#000000",
        fontColor       : "#FFFFFF",
        pointerEvents   : "visiblePainted",
        fontSize        : "15px"
    });

    var strategy = new OpenLayers.Strategy.Cluster({distance: 30, threshold: 2});

    var ly=new OpenLayers.Layer.Vector(GSMapLayerName,{
            strategies: [strategy],
            styleMap:new OpenLayers.StyleMap({
                        "default":defaultStyle,
                        "select": onSelectStyle
                    }),
            displayInLayerSwitcher:false
        });
    map.addLayer(ly);

I am getting the below output

Look at the image the feature which have 9 is behind the 3 and 5. This is not a good to view.
 So, How can I set this feature to the front of all other features?

Note: I am using Openlayers 2.11

Thanks for any reply.


texts: Extract HDF5 group name and apply it as a variable name using R

texts:
I'm new to using HDF5 files and need to extract information from several thousand such files, each with 6 groups of interest (e.g., ('/gt1l', '/gt1r', '/gt2l', '/gt2r', '/gt3l', '/gt3r').  I would like to collate lat/lon information from each of the groups and place them in a dataframe. Within the dataframe I need one unique_id column representing the name of the HDF5 file and another column representing the group id.
So for each row I would want to know what HDF5 file the data is from and which HDF5 group the data was extracted from.
Is there a straightforward way to append the group name in the row column?
I've looked around for a solution to this but have not yet been able to identify one.
Attached is the code I've been using to access the data.
fp <- "test_data.h5"

hdf = H5Fopen(fp)

HDF5 FILE
        name /
    filename

                name     otype dclass dim
0 METADATA           H5I_GROUP
1 ancillary_data     H5I_GROUP
2 gt1l               H5I_GROUP
3 gt1r               H5I_GROUP
4 gt2l               H5I_GROUP
5 gt2r               H5I_GROUP
6 gt3l               H5I_GROUP
7 gt3r               H5I_GROUP
8 orbit_info         H5I_GROUP
9 quality_assessment H5I_GROUP

#Extract the variables from the HDF5 file
lat <- hdf$gt1l$sea_ice_segments$latitude
lon <- hdf$gt1l$sea_ice_segments$longitude
group <- paste group name here

#Create a new dataframe containing variables of interest
df <- as.data.frame(lat, lon, group_name)


texts: How to convert vector format (geopackage, shapefile) to raster format (geotiff) using rasterio without template raster

texts: I'm looking to create a raster given just a geopackage (or shapefile, or other vector format). I don't have a template raster, and all I wish to do is burn the shapes into a raster using rasterio. Without the template raster, I'm uncertain how to generate the required inputs (transform, mainly). I'd rather not construct it by hand.


texts: How could I assign point values to a raster stream segment?

texts:
I need to assign a point value in raster format to a raster stream segment, which contains the point. See bellow:

The expected result is a segment, whose value is the point value. How could I?


texts: Including snapped distance in Network Analyst route calculation?

texts:
In designing a routing tool, I'm attempting to minimize the amount of surface streets (non-motorways) needed to reach the closest facilities. The only options available in Network Analyst are snap to closest feature or snap to first feature (within some tolerance radius) when using the AddLocations tool.

As you can see from the image below, the code will snap the input location via part 1 every time (either by closest or by first). The favorable direction is to snap via path 2 to avoid that southern loop.

I believe I can minimize this issue by running multiple solves on the results of a Generate Near table for the 5 closest features then checking the cost for the cheapest route. It does leave the issue of a single continuous road feature that might loop back like that.

Does anyone have any ideas of how to route using part 2?

(Note: The Network dataset shown below is from OpenStreetMap and prohibits travel on motorways and motorway links. All the remaining roadways are to be treated equally in this scenario.)


texts: Extracting desired area from group of shapefiles, area remains same as original file

texts:
I have an output of a classified image as a polygon shapefile (call it File A, the larger, irregular shape in the screenshot below). I used the clip (Analysis) tool to extract a rectangular area in the center (inner green polygon in the screenshot).

My output is shown below. When I checked the area (square meters) of the purple layer I found that the area of output is almost equal to the area of file A. I would expect it to be much smaller. Am I doing anything wrong?


texts: Gdal.open crashes python for some .tif images

texts:
I just upgraded to GDAL 2.2.0 and now I'm having issues with my python kernel crashing when I run gdal.Open() on certain files that it used to work on when I had GDAL 2.1.2. It seems to only be crashing on images I have created with gdal_translate or gdalwarp, but not all of them. This is the error I keep receiving in my terminal window after a crash:

x = gdal.Open("tiffpath.tif")
TIFFReadDirectory: Warning, Unknown field with tag 42113 (0xa481) encountered.

These are the commands I'm using to create the GTiffs that crash my kernel:

subprocess.call(["gdalwarp", "-of", "GTiff", "-r", "bilinear", "-srcnodata", "0", "infile1path", "infile1path", "tiffoutpath"])

gdalTranslate = r'/usr/local/Cellar/gdal2/2.2.0/bin/gdal_translate'
subset = (520926.0, 7339961.0, 686430.0, 7223749.0,"path_in","path_out")
trans = gdalTranslate + ' -of GTiff -projwin %s %s %s %s %s %s' %(subset)
os.system(trans)

I can open these images using PIL Image.open(). However, I still need gdal.Open to work so I can retrieve projection information.

I found that tag 42113 is associated with this: 'a special pixel value to mark geospatial areas for which no info is available'. Is there a glitch in the way the new gdal_translate and gdalwarp are creating geotiffs, an error in the way gdal.Open() functions with some .tif images, or something else I might be doing wrong?


texts: QGIS export to PDF layer is disappearing

texts:
I would like to reopen/look again at this question. I have a aerial raster file with a vector file over the top indicating specific points on the ground. I created my map, set the rotation to zero, not using WMS-layers, toggled the transparency, exported multiple times; the output PDF file continues to loose written data. The raster looks good, but the shapefile only shows up on the top half.
Has anyone else ran into this problem?
QGIS raster background with red points.

Exported Geo-referenced PDF