texts: QGIS: Making circles of specific diameter from attribute in point layer

texts:
I have seen a similar question at Creating different diameter circles round points using QGIS but it is very old am seeking up-to-date and possibly more conclusive guidance.
I would like to create a new layer with circles of a specific diameter based on an attribute within a point layer. This relates to trees, specifically BS5837 In the UK.
For example with a DBH, (Diameter at Breast Height of tree stem), reading of 250mm I want to see a circle with a radius of 3m. This circle is typically calculated and then drawn manually in aCAD. I'm looking to improve this workflow.
The DBH/RPA calculation is simple. RPA=DBH*12 normally I work this into meters for CAD (/1000).
Answer notes: The calculation Ive used in QGIS is slightly different to compensate for units (DBH*0.012) would work either way though. The projections were key in getting the correctly sized circles.
Many of the suggested methods, including those received here utilise buffer in Vector Tools this is likely the most robust solution however I have an issue that my projection (EPSG3857 - WGS84 Pseudo) is in degrees. In the answer provided below to arrive at a correctly formatted output I used EPSG:23030 - ED50 / UTM zone 30N.
My ultimate goal is to have some form of template style/formula/symbology/buffer that runs through the point data drawing all RPAs.
Answer notes: My next goal will be to condensed the newly learned workflow down into some form of template style or project i can import data into.
Many thanks to Håvard Tveite and others for input.


texts: Create points along line using CCTV data (irregular interval!)

texts:
I am trying to create points along a line using CCTV data for our wastewater network using ArcMap 10.3.

There is no field within our pipe to specify upstream/downstream so I have determined the flow of the pipe using the Utility Network Analyst tool.

All measurements within the CCTV Inspections are taken from the upstream manhole.  I would like to have a point created at each of these inspection points.
(A future challenge would be to then include continuity information to create a line rather than point.)

Below is a simplified database diagram.

I am unsure how to progress with this ... tracing or linking using the linear referencing tool??

Any suggestions on the next steps I should take?


texts: Add a new polygon to an existing shapefile from WKT in QGIS

texts:
I have a layer of polygons in QGIS;
I want to add programmatically (PyQGIS) a new polygon to this layer using WKT.
Example:

   mywkt = POLYGON ((522762 376413 0,530635 368671 0,522237 353713 0,522762 376413 0))

I want to add this polygon to my layer.

I try this :

layer.startEditing()
geom = QgsGeometry()
geom = QgsGeometry.fromWkt(mywkt)
feat = QgsFeature()
feat.setGeometry(geom)
layer.dataProvider().addFeatures([feat])
layer.updateExtents()
layer.reload()


texts: How to mosaic images with different CRS with rasterio

texts:
I am doing mosaicing of 14 raster images (.img format) to completely cover a province.see original spatial location of separately opened tiles(). Mosiac operation is completed but 4 tiles are mosiaced on wrong spatial location. whereas other 10 tiles are mosiaced accurately (). when i am opening separate files all files are getting opened at their origina spatial location. please help me to resolve this issue.
I am using this code for mosaicing of 14 raster images.

import rasterio
from rasterio.merge import merge
from rasterio.plot import show
import glob
import os

dirpath = r"path"
outpath = r"outpath\Mosiac.tif"

search_criteria = "*.img"

q = os.path.join(dirpath, search_criteria)
print (q)
img_fps = glob.glob(q)
img_fps

src_files_to_mosiac = []

for sim in img_fps:
src = rasterio.open(sim)
src_files_to_mosiac.append(src)

src_files_to_mosiac

mos, out_trans = merge(src_files_to_mosiac)

show(mos, cmap='terrain')

with rasterio.open(outpath,"w", driver ='Gtiff',count=3,
            height= mos.shape[1],
            width= mos.shape[2],
            transform= out_trans,
            crs= src.crs,
            dtype= src.dtypes[0]) as dest:
dest.write(mos)


texts: QGIS - Raster Calculator VS Raster Difference (SAGA)

texts: I am trying to calculate the difference between two rasters.  If i use Raster Calculator in QGIS I get the correct result (when manually checked using the information tool and the pixel values for the new layer and input layers).  When I use Raster Difference tool in the SAGA tools (within qGIS) I get a very similar but not exact result when it comes down to the decimal place. I'm wondering what the difference is between these two tools?  I would use the raster calculator but the method I have been directed to is the rater difference tool so I was just wondering why they are outputting a slightly different result with the same inputs.


texts: How to associate a point with a file on Google Earth Engine

texts:
Trying to plot CTD data from an oceanography cruise on Google Earth Engine so that it can be compared to MODIS ocean color.

The trouble is, each measurement is a separate CSV, and I would like to associate the CSV with each point where the measurement took place.

I was thinking of making points (213 of them), each named after the file (ex: var poi = ee.Geometry.Point(long, lat), {name: 'measurement_1.csv'}; and then calling the files—is there anyone who has done this?


texts: Is  only available in GeoServer from version 2.9?

texts:
I need to visualize PostGIS' polygons as GeoServer WMS with specific style such that the polygons' stroke have different color between inner and outerside of the polygon.
According to this reference, a <PerpendicularOffset> tag can be used. Below is the illustration from the reference.

I tried to add that tag in my GeoServer 2.7.0 but does not work. The stroke was rendered with usual symmetric stroke, (inner and outer polygon with similar color and width). Below is the result (i am expecting the gray stroke in only inner side) but fails:

The reference apparently is for "2.9.x", and when I check the same reference spot (LineSymbolizer) of lower version (2.8.x and 2.7.2), the explanation of <PerpendicularOffset> tag in LineSymbolizer does not exist.

Is <PerpendicularOffset> only available in GeoServer from version 2.9 ?
If <PerpendicularOffset> is only available in GeoServer from version 2.9, then what are the possible options for me using GeoServer 2.7.x ?
What's the difference between using <PerpendicularOffset> tag in GeoServer 2.9.0 and using Fluxomajic ? (Fluxomajic was also discussed here).

I'm using:

GeoServer 2.7.0 on Apache Tomcat 7.0.61 on Windows Server 2008 R2
Datacenter 64-bit on VMWare
JVM: Oracle Corporation: 1.7.0_80 (Java HotSpot(TM) 64-Bit Server VM)
PostGIS - PostgreSQL 9.4.1 (on Windows Server 2008 R2 Datacenter 64-bit on VMWare)

I tried the blending modes, with the following success result:

The above uses 4 FeatureTypeStyle, more or less based on the example you refer in the answer by @AndreAime.

I also tried the Fluxomajic but blending modes are far better.


texts: Overriding updateParameters()-Values in Python toolbox

texts:
I'm just getting started with Python toolboxes in ArcGIS and I'm pretty new to Python in general. I've created a python toolbox which

takes an raster input file
derives some values based on the name of the input file (via a custom function)
suggests several parameter values
processes the files in some way (not relevant to my question)

My issue is, that I cant seem to override the suggested parameter values when I run the tool. The values always switch back to the original values as suggested via the custom function. My guess is that whenever I change the suggested parameter, the function updateParameters() is run which in turn overrides my manual input. Is there a way to solve this dilemma?

I've reduced my toolbox (.pyt-File) to a minimal example (excluding the geoprocessing part) in the following code. In this example, the function suggestName() returns a new filename with the ending ".shp" and a path to the scratch folder. Imagine the function suggestName() to be slightly more complex and returning several values, not just one.

import arcpy

def suggestName(inputname):
    "Extract the Sheetnumber and Year from a given map Name"
    return inputname+".shp"

class Toolbox(object):
    def __init__(self):
        """Toolbox"""
        self.label = "Toolbox"
        self.alias = "Toolbox"

        # List of tool classes associated with this toolbox
        self.tools = [Tool]

class Tool(object):
    def __init__(self):
        """Some Info"""
        self.label = "Tool"
        self.description = "Tool"
        self.canRunInBackground = False

    def getParameterInfo(self):
        """"""

        ## Raster map to extract features from
        param0 = arcpy.Parameter(
            name = "inputRaster",
            displayName = "Inpurt Raster File (TA25)",
            parameterType = "Required",
            direction = "Input",
            datatype = "DERasterDataset"
            )

        ## Output feature class
        param1 = arcpy.Parameter(
            displayName="Output Features",
            name="out_features",
            datatype="GPFeatureLayer",
            parameterType="Requred",
            direction="Output")

        return [param0,param1]

    def updateParameters(self, parameters):
        """"""

        if parameters[0].value:
            parameters[1].enabled = True

            newname = suggestName(arcpy.Describe(parameters[0].valueAsText).basename)

            parameters[1].value = newname
        return

    def execute(self, parameters, messages):
        """"""
        return

Edit:

It seems like I wasn't very clear in clarifying my issue, I'll try to do so now. As you can see in the screenshot below, the tool suggests the name "test.shp" as the output file. But I cannot override this parameter. If I try to do so (e.g. change the name to "test_2.shp") the parameter switches back to the original name ("test.shp") once I click anywhere outside the box.


texts: Using surface - slope analysis tool box. All values for the slope calculation are in between 80 and 89. There are no values below 80. What am I doing wrong?


texts:
I would like to compute shortest path between 2,5 millions pairs.

I created table with id, origin, the nearest node of origin, destination, the nearest node of destination.
I also created index on source, target, id.
I have table with roads network of one city - also with indexes. All geoms have gist index.

When I tried to use dijkstra, the execute of query is too slowly. It would take over a week.
I have idea, I will add where condition (f.e. id <10000) and I will make insert to new table. Therefore I can do a few queries in the same time.
I noticed that using a larger index (few columns) speeds up the query.

Below my code:

SELECT a."id" as id
, SUM(a.cost) AS "Aggregated Cost"
FROM ( SELECT "id"
, (pgr_dijkstra( 'select id, source, target, shape_leng as cost
from roads', "mz_node", "tz_node", false )).*
FROM pair_with_near where id < 1000) AS a
JOIN roads AS b ON a.edge = b.id GROUP BY a."id"

I don't have more ideas, how can I make it faster.
I read about buffer creating but I can't do it.
The same with THIS

Do you have any idea?


texts: Creating multi-image single band single file geotiff using GDAL with NumPy?

texts:
I have a geotiff file that is single band but contains multiple subdatasets.  I use gdal with numpy to open a specific subdataset and mask some of its values.

Example:

$ gdalinfo sample.tif
...
Subdatasets:
SUBDATASET_1_NAME=GTIFF_DIR:1:sample.tif
SUBDATASET_1_DESC=...
SUBDATASET_2_NAME=GTIFF_DIR:2:sample.tif
SUBDATASET_2_DESC=...

# pseudocode
sds = gdal.Open('GTIFF_DIR:2:sample.tif')
sds_mask = apply mask to sds

I want to put sds_mask back into the original geotiff or a copy of it.  Is this possible with gdal?

I tried writing out all the subdatasets, modifying one of the subdatasets, and then merging them back together with gdal_translate and gdal_merge.py with the -separate option but it does not create a geotiff with subdatasets.


texts: Can I use QGIS plugins in PyQGIS as a stand alone script without iface?

texts:
Can I use a qgis plugin especially (openlayers plugin) in stand alone pyqgis script?

I don't find examples, just code using iface but that cannot used in stand alone script.


texts: Couldn't change ArcGIS Web App Builder application spatial reference (wkid 102100)

texts:
I couldn't change ArcGIS Web App Builder application spatial reference (wkid 102100) to wkid 4326 even though I changed it in the configuration.
Is there anywhere else where I have to make changes to a config file?
The one below is in WebAppBuilderForArcGIS2.6\server\apps\2\config.json

  "mapOptions": {
    "extent": {
      "xmin": -73.3563082698,
      "ymin": -33.485687785,
      "xmax": 61.373452187,
      "ymax": 29.4752271893,
      "spatialReference": {
        "wkid": 4326
      }
    }
  },

What I am using:

Web App Builder 2.6
Javascript API 3.22


texts: Attribute stats query with Execute SQL tool in QGIS Modeler

texts:
I'm trying to get a simple temporary table with the surface of different categories of a polygon shapefile, classified by category. My final aim is to hide an SQL query in a tool with optional variables to give it back as a user-friendly model tool.

Problem : i get some error messages and cannot do the task to the end.
Is there a syntax error in my parameters ?

Below details of what I do :

This is my test data where I want to get in 2 different row the sum of S and C items areas.

In the database menu, in thr SQL window with virtual layers, i do the request on my layers to check my formula - it works :

SELECT  layer.ID as CAT, SUM(area(layer.geometry)) as SURF FROM layer GROUP BY CAT;

I get back the following table :

CAT SURF
C   139327
S   185247

I would now like to use the modeler to be able to ask the user what layer he wants and what column to use for my GROUP BY. I'm trying to use the "Execute SQL" tool but with no success. Below the model I use (it's really basic)

Setup of Execute SQL

I get the following error :

  Préparation de l'algorithme : qgis:executesql_1 Exécution de Exécuter
  SQL [1/1] Paramètres en entrée: { INPUT_DATASOURCES:
  ['//.../surf.shp'], INPUT_GEOMETRY_CRS: None,
  INPUT_GEOMETRY_FIELD: 'geometry', INPUT_GEOMETRY_TYPE: 4, INPUT_QUERY:
  'SELECT layer.ID as CAT, SUM(area(layer.geometry)) as SURF FROM layer
  GROUP BY CAT;', INPUT_UID_FIELD: 'TYPE', OUTPUT: 'memory:' } Traceback
  (most recent call last): File
  "C:/OSGEO4~1/apps/qgis-ltr/./python/plugins\processing\algs\qgis\ExecuteSQL.py",
  line 169, in processAlgorithm raise
  QgsProcessingException(vLayer.dataProvider().error().message())
  _core.QgsProcessingException: virtual: Referenced table layer in query not found!


texts: Iterate Dataset and Workspace to Add Rasters to Mosaic Dataset?

texts:
My raster mosaic datasets are in a file gdb ordered by block. For example Block_470 and Block_480.

My workspace has images ordered by block with subfolders with in each block. For example workspace Block_470 contains subfolders BW and RGB. This is the same for workspace Block_480. Meaning workspace Block_480 will also have subfolders BW and RGB.

I need to add the tiff files in workspace Block_470/BW to raster mosaic dataset Block_470 and the tiff files in workspace Block_480/BW will need to be added to raster mosaic dataset Block_480.

How can I iterate through so that the right tiff files in the subfolders go into the corresponding mosaic datasets?

I was thinking something like this:

import arcpy

workspaces = arcpy.ListWorkspaces("//nyx-galli/eGIS/Data/Raster/Imagery/Aerial/SGIC/2008_Mapcon/Block4*", "Folder")
datasets = arcpy.ListDatasets("D:/TestRasterMosaic/BW.gdb/Block_4*", "Mosaic")

for i in workspaces:
    workspaces2 = arcpy.ListWorkspaces("", "Folder")

After building this model did not work for me:

Model Builder Image


texts: Whats steps can I suggest to achieve the best Geolocation Result

texts:
We are using Geolocation (getCurrentPosition()) in a website to determine a users position when using our site from a mobile device.

I want to write an article explaining how the user can obtain the best results.

Am I correct in assuming:

Enabling GPS will yield the best result when in rural areas (less buildings to obscure line of sight to the satelites)
Enabling Wi-Fi will yield the best results when in urban areas (generally more Wi-Fi hotspots available)
Is it true that Android phones have better results from silently harvesting Wi-Fi hotspot details?

Any links to reference material on this are appreciated


texts: Running parallel search cursors in Python?

texts:
I have a multipoint feature class and a polyline feature class (transects) filled with the same number of objects. Basically, each set of points corresponds to a single transect.

I want to run the intersect analysis tool for each individual set of rows, so each transect and set of points should have an output. The output would be a set points that is "snapped" to the transect polyline. Then at the end of the code, I have as many intersect outputs as I do transects.

How can I get two search cursors to run simultaneously? I need a search cursor that can take one row from each file, run the intersect tool, and then move to the next set of rows.

A nested loop does not work (unsurprisingly). Is there a workaround for this?

EDIT:
I don't want to compare two datasets, I actually want to analyze corresponding features of both. The formatting here didn't copy well... apologies.

# Import system modules
  import arcpy
  import numpy as np
  import sys, os, string, arcgisscripting
  arcpy.env.overwriteOutput=True

  #Define local variables
  workspace=
  transect=
  point_file=
  output_file="clusters.shp"

  sr=arcpy.Describe(point_file).spatialReference
  arcpy.env.overwriteOutput=True

  gp=arcgisscripting.create()
  output_fc=gp.CreateFeatureClass_management(workspace,output_file,     "MULTIPOINT", "#", "DISABLED", "ENABLED",sr)

 #Perform fits and intersections
 with arcpy.da.SearchCursor(transect,["SHAPE@"]) as cursor1:
     for row in cursor1:
         profile=row[0]
 with arcpy.da.SearchCursor(point_file,["SHAPE@"]) as cursor2:
     for row in cursor2:
         point=row[0]

        #Create Swath Profile
        clusterTolerance=0.5
          swath_points=arcpy.Intersect_analysis([point,profile],output_fc,"",clusterTolerance,"point")


texts: Near analysis for nearest feature in a distinct angle

texts:
I want to calculate the nearest distance from a point feature to a line feature in a distinct angle.

When using the Near analysis tool in ArcGIS I get the nearest distance in a angle where the nearest feature lies. In the picture below for -170°.

But I want to be able to calculate the distance for different angles for example for 90° or 180°.
Someone an idea how I could solve this? I'm also capable of creating my own script with python.


texts: How do I convert Eastings and Northings projection to WSG84 in geopandas?

texts:
I'm very new to projections and struggling to convert a map dataset from OSGB36 (Eastings and Northings) to WGS84 (lat/longs).

I've used geopandas to import a shapefile. Now I want to export it as GeoJSON so I can then use Folium to map it on leaflet.js.

If I check the crs variable with .crs I get:

{u'datum': u'OSGB36',
u'k': 0.999601272,
u'lat_0': 49,
u'lon_0': -2,
u'no_defs': True,
u'proj': u'tmerc',
u'units': u'm',
u'x_0': 400000,
u'y_0': -100000}

Trouble is if I'm a bit lost on the attributes to set to convert it to WGS84:

staffs.to_crs({'datum':"wgs84"})

RuntimeError: projection not named

Can you help?


texts:
I calculated a series of Kernel Density Estimations in CrimeStat (because more parameters can be set than in the built-in tool in ArcGIS) and then imported the results into ArcMap. I wanted to combine all the Kernel Density Estimations to get an overall result by using the Raster Calculator so I can assign different weights to each Kernel Density Estimation. But the Raster Calculator only accepts actual rasters, while crimestat can only output vector grids.
Is one of these solutions possible:

convert the vector grids to raster grids without loosing the original values? (when I use Polygon or Feature to Raster the grid looks totally different and I don't know why, it should look exactly the same only with a different format)
use a different tool instead of the Raster Calculator that accepts vector grids? I thought of a Spatial Join, but I need to weigh every layer differently and also the Spatial Join takes forever because of all the little grid cells.


texts: Python script for constructing a WHERE clause from user input

texts:
I'm trying to construct a where clause from user input that will eventually be passed through as a SQL statement in the select tool or the select by attributes tool. I've started a script thats been adapted from some previous posts, however, I still don't have it correct. the way its currently written the resulting SQL statement will look something like:

"123" = 'ABC' OR 'DEF'

I need the resulting statement to look more like

"\123\" = 'ABC' OR "\123\"='DEF'

Important: the values (i.e. ABC DEF) will be user inputs and the number of values entered as arguments may vary from user to user. My current script is below

 def whereClause(table, field, values):
    """Takes a semicolon-delimited list of values and constructs a SQL WHERE
    clause to select those values within a given field and table."""

    # Add field delimiters
    fieldDelimited = arcpy.AddFieldDelimiters(arcpy.Describe(table).path, field)

    # Split multivalue at semicolons and strip quotes
    valueList = [value[1:-1]
        if (value.startswith("'") and value.endswith("'"))
        else value for value in values.split(';')]

    # Determine field type
    fieldType = arcpy.ListFields(table, field)[0].type

    # Add single-quotes for string field values
    if str(fieldType) == 'String':
        valueList = ["'%s'" % value for value in valueList]

    # Format WHERE clause in the form of an IN statement
    whereClause = "%s = %s" % (fieldDelimited, 'OR' .join(valueList))
    return whereClause

def outName(input,post="Out",fileExtension="shp"):
    """Returns output name."""
    outName=os.path.basename(input).split(".")[0]+post+"."+fileExtension
    return outName

# Create SQL Statement
InputFC = arcpy.GetParametersAsText(0)
outputFC=outName(InputFC,"_Select")
Field ="CO_NAME"
Values = arcpy.GetParameterAsText(1) #user input
SQL=whereClause(InputFC,Field,Values)
print SQL


texts: Converting/using ArcGIS *.style file for/in QGIS (QML or SLD)?

texts:
I have downloaded a dataset from a government organisation. it comes with a .style file.

I use QGIS 1.8 on windows 7.

Is there a way to use the .style file as a qml or sld? I have tried simply changing the extension but to no avail.


texts: How to calculate corner(edge) coordinates for geo image

texts:
I want to add image layer to Azure map. To do that I need coordinates for the corners. With help of this site http://viglino.github.io/Map-georeferencer  I am able to get center, rotation/angle and scale but having trouble to calculate coordiates for the corners. To calclate coordinates I tried:

TOP RIGHT :
Top_Right.x = center.x + ((width / 2) *(scale.x) * cos(angle)) - ((height / 2) *(scale.y) * sin(angle))
Top_Right.y = center.y + ((width / 2) *(scale.x) * sin(angle)) + ((height / 2) *(scale.y) * cos(angle))

TOP LEFT :
Top_Left.x = center.x - ((width / 2) *(scale.x) * cos(angle)) - ((height / 2) *(scale.y) * sin(angle))
Top_Left.y = center.y - ((width / 2) *(scale.x) * sin(angle)) + ((height / 2) *(scale.y) * cos(angle))

BOTTOM LEFT :
Bot_Left.x = center.x - ((width / 2) *(scale.x) * cos(angle)) + ((height / 2) *(scale.y) * sin(angle))
Bot_Left.y = center.y - ((width / 2) *(scale.x) * sin(angle)) - ((height / 2) *(scale.y) * cos(angle))

BOTTOM RIGHT:
Bot_Right.x = center.x + ((width / 2) *(scale.x) * cos(angle)) + ((height / 2) *(scale.y) * sin(angle))
Bot_Right.y = center.y + ((width / 2) *(scale.x) * sin(angle)) - ((height / 2) *(scale.y) * cos(angle))

But output is incorrect.


texts: create POV-Ray image from DEM-file in grass

texts:
I am working on r.out.pov Tool of grass 6.4 creating a POV-ray(beautiful 3D) image from an elevation raster.

You can see this link for more information about 3D analysis of GRASS.
r.out.pov manual
Can anyone help me(step by step for windows) creating such an image of my DEM file? and some help about targa file (.tga).


texts: How to calculate RMSE, skewness, and kurtosis of raster image?

texts:
I'm doing a project to analyze the impact different interpolation techniques plays on DEM accuracy.

In order to do so, I need to calculate the RMSE, skewness, and kurtosis.

I have:

an original source DEM,
a DEM generated from 10 different interpolation techniques,
a DEM generated by the subtraction of above two DEMs.

How do I calculate the three parameters?


texts: Calculating multiple polygon areas coverage within multiple polygons in QGIS?

texts:
I have shapefiles for school attendance zones (SAZs) and Census Tracts (CTs), where the red lines below are CTs and the lighter lines are SAZs. I would like to calculate the proportion (%) of each CT covered by each SAZ. For example, if one CT was covered by two SAZs (SAZ1 and SAZ2), I want to be able to estimate that 28% of the CT is covered by SAZ1 and 72% of it is covered by SAZ2.

I have not tried anything yet, since I am fairly new at this and cannot think of a solution. I think I would start out by calculating the areas of each polygon in sq. miles, ft, meters, etc., then do a join of some sort (union? intersect?), but I am stuck.


texts: Adding new actions to identify feature context menu

texts: How can I add new actions to identify feature context menu? I mean when selecting a feature with identify tool and right click on the selected feature, here I want to add new actions to the popup menu.


texts: Isolation and interpretation of vineyard multispectral data

texts: I have 3 band multispectral raster imagery over a vineyard from which i have created a number of indices e.g. NDVI. The data is fairly high resolution (5cm/pixel) so shows a lot of along vine and across vine variation.  For visualisation (and actual use and actioning of the data) i would like to be able to isolate the actual vines from the background (they obviously stand out really well as they are irrigated and are easily separated with the correct thresholds set) and then interpolate the average between the vine rows to show broad scale variation across the vineyard.  Does anyone have any ideas how i might do this in the QGIS GUI environment?


texts: How to disable automatic cropping of WMS layers in QGIS?

texts:
When importing a WMS layer, QGIS automatically limits the extents of the image to the CRS area of use.

Is there any (GDAL) setting/variable to prevent this annoying behavior?


texts: Working with Map Series in Map View to assist Editing in ArcGIS Pro

texts:
In ArcMap I sometimes used Data Driven Pages to actually compare features and make edits.  For example, I have two feature classes consisting of the same list of addresses geocoded via different locator services.  If feature class A is more accurate over all but feature class B has a few points that are more accurate, I'd want to move those points in A to match those of B.  I could do this easily in ArcMap; just turn on Data Driven Pages and click the arrow on the toolbar to cycle through all the features and do my manual inspection of the two features.

In ArcGIS Pro it seems like the replacement for Data Driven Pages, the Map Series, is only available in Layout View.  If I want to do the same thing in ArcGIS Pro, I have to go to Layout View, add the map series, click the first index, go back over to map view to perform any edits I may need to make, click back over to Layout view, click the next index in the series, repeat.  Lots of extra clicking for the same functionality.

Is there a way to enable Map Series in Map View?


texts: Manually digitalize points in QGIS with automatic ID

texts:
I would like to create a few thousand points by clicking on specific points in my background map.
I have created a shapefile-layer for this task. But when I click on the map, QGIS always asks for a point ID. Is it possible to give these points an ID automatically so I don't have to type one? Do the points need an ID? If I just click enter, the point is saved without an ID. Is that a problem?


texts: Calculating shortest distance between two tables using PostGIS

texts:
I am trying to calculate shortest distance from each point to polygon where idfield match and add the distance to a field in the point table:

For the two "A" points the distance is the same, which is incorrect. I understand why but not how to fix it. If I had only one A point, it would work.
update public.tbl2
SET dist = sub.dist
from (

select  tbl2.idfield,
        st_shortestline(tbl1.geom, tbl2.geom) geom,
        round(st_length(st_shortestline(tbl1.geom, tbl2.geom))::numeric, 0) dist
from public.tbl1 , public.tbl2
where tbl1.idfield = tbl2.idfield
) sub
where sub.idfield = tbl2.idfield


texts: How to find measures of points that are close to route in ArcGIS 10.2?

texts: What I have is a route layer of a river and points layer with associated temperature data (temperature of the river) that were collected using an airplane. The path of the points does not match the path of the river, but the points are in the vicinity of the river. My goal is to calculate river km locations of those points along that river layer. I could do this manually, but there are over 600 points. Is there a way to find locations for all of those points on the route layer that is river?


texts:
I want to process a 16 bit satellite image(.geotiff) file in python using GDAL and other imaging libraries like PIL. But when I display image, it is not displayed properly as it is 16 bit. Is there any way to display it properly in python? Please help. Code is given below.

Algorithm is as below:

Open geotiff raster using GDAL library.
Get it in the form of array and process as per requirement.
Write processed image into output raster.

This code works and output raster (image) can be seen in other remote sensing software like Arcmap, QGIS, ENVI. But during processing, I need to see resultant image in python. Due to 16 bit nature of image, it shows only black and white (no other grey shade) color in image.

Is there any way that I can see output image in python only instead of going to remote sensing software?

Code is given below. (write raster part is removed from code)

from PIL import Image,ImageOps
import numpy as np
import scipy as sp
import pywt
from osgeo import gdal
from osgeo import gdal_array
from osgeo import osr
import cv2
#from osgeo.gdalconst import *
#import matplotlib.pylab as plt

gdal.AllRegister()

################## Read Raster #################
inRaster='C:\python\Results\Database\Temporary data base\clipped_26July2017\HV_clip9Nov.tif'

inDS=gdal.Open(inRaster,1)
geoTransform = inDS.GetGeoTransform()
band=inDS.GetRasterBand(1)
datatype=band.DataType
proj = inDS.GetProjection()
rows = inDS.RasterYSize
cols=inDS.RasterXSize
data=band.ReadAsArray(0,0,cols,rows)#extraction of data to be processed#
im=Image.fromarray(data)#-----creation of image from raster data----#

ca1,(ch1,cv1,cd1)=pywt.dwt2(im,'db1');
ca2,(ch2,cv2,cd2)=pywt.dwt2(ca1,'db1')
ca3,(ch3,cv3,cd3)=pywt.dwt2(ca2,'db1')
coeff3=ca3,(ch3,cv3,cd3)
c_op3=pywt.idwt2(coeff3,'db1');
c_op3=Image.fromarray(c_op3)
c_op3.show()
coeff2=ca2,(ch2,cv2,cd2)
c_op2=pywt.idwt2(coeff2,'db1');
c_op2=Image.fromarray(c_op2)
c_op2.show()
coeff1=ca1,(ch1,cv1,cd1)
c_op1=pywt.idwt2(coeff1,'db1');
[u,v]=np.shape(c_op1)
c_op1=np.delete(c_op1,u-1,0)#deletion of last row, appended by pywt due to odd size of image
c_op1=np.delete(c_op1,v-1,1)#deletion of last column, appended by pywt due to odd size of image


texts: How to match the colors of different RGB mosaic rasters?

texts:
I am beginner in QGIS. I want to classify many  RGB rasters (see picture) but their colors don´t match. As you can see, the green polygon shows the same land cover class (meadow) but in each tile the color for this class is totally different.

I would like to "re-paint" them or match their histograms somehow. They overlap sligthly (like a mosaic), maybe not enough for matching overlapping features. I would like to take a sample from one image, tell the program that the sample is a meadow, and re-paint all the other meadows in the other images. Is it possible to do this?

I tried to do it in arcGIS, but always got an error. Answers in QGIS or R are welcomed.

I am running QGIS 2.12.3 in a PC.


texts:
I am attempting to draw a WMS layer consisting of multiple polygons of different colours and an opacity of 0.4.

The relevant part of my SLD is

      <PolygonSymbolizer>
        <Fill>
          <CssParameter name="fill"><ogc:PropertyName>colour</ogc:PropertyName></CssParameter>
          <CssParameter name="fill-opacity">0.4</CssParameter>
        </Fill>
        <Stroke>
          <CssParameter name="stroke">#000000</CssParameter>
          <CssParameter name="stroke-width">1</CssParameter>
        </Stroke>
      </PolygonSymbolizer>

Which draws the layer using the correct colours, but completely opaque.

The layer is being created with OpenLayers 2 as follows:

self.createBoundingAreaLayer = function (typeCode, boundingAreaIdentifier) {
    return new OpenLayers.Layer.WMS(
        typeCode,
        esd.pathUtility.geoServerUrl + "/ESD/wms",
        {
            LAYERS: 'ESD:areas',
            viewparams: 'parent:' + boundingAreaIdentifier + ';typecode:' + typeCode,
            format_options: 'antialias:none',
            format: 'image/gif',
            transparent: true
        },
        {
            singleTile: true,
            //opacity: 0.4,
            ratio: 1,
            isBaseLayer: false,
            transitionEffect: 'resize',
            yx: { 'EPSG:4326': false }
        }
    );
};

and if I uncomment the opacity: 0.4 line from there then it works. However, that isn't exactly what I want because it also changes the opacity of the stroke making it hard to see.

Is my SLD syntax correct for setting the opacity of the fill only? If so then why is this not working? Everything is currently running locally and I've rebooted my PC, so I don't think it's a caching issue.

Possibly relevant information:
The style already existed and I added the fill-opacity parameter.
I select the SLD by setting the default style for the layer in GeoServer. I never specify it in code.


texts: Displaying multiple symbols/patches for the one layer in an ArcMap legend

texts:
I'm mapping incidents of an invasive species on a map, where the observations during surveying were recorded as points and also a line where there was a long continuous strip of the species along a bank. These two layers are displayed in the legend as this:

1

The two different symbols (point and line) ideally will sit beside each other in the legend or at least it will be communicated that both the point and the line represent the species, in one legend item.

I've achieved this by converting to graphics and ungrouping before (see knotweed example below).

2

However, as I am using data driven pages to output several maps at once, I don't want to break my dynamic legend or for the converted graphics to appear across all the pages (this species is only relevant to this map), plus I expect the map to have to be updated in the future.

Can something like in the second examples above be achieved without converting to graphics?


texts:
How can an intersection of 3D object with 3D object be done in R, for example a 3D line with 3D polygon ? Saw that "sf" provides some possibilities but didn't get the desired result. What was accomplished:

library(dplyr)
library(sf)

# create 3D polygon
pol1 <- data.frame(lon = c(6,6,6,6,6),
                   lat = c(3,1,1,3,3),
                   alt = c(2,2,4,4,2)) %>%
st_as_sf(coords = c("lon", "lat", "alt"), dim = "XYZ") %>%
summarise(do_union = FALSE) %>%
st_cast("POLYGON")

#create 3D line
lns1 <- matrix(c(5,8,2,2,2,4), ,  3) %>% st_linestring(dim = "XYZ")

# trying to find the point where they intersect
>st_intersection(lns1, pol1)
GEOMETRYCOLLECTION EMPTY

Is there a better way to do this ? Thankful for any suggestions


texts: Getting URL for specific tile on XYZ layer?

texts:
I am using custom XYZ tiles layer served from a local server. It works great most of the time, but some of the tiles appear to be missing or damaged (show up as white squares).

Is there a way to figure out what are X/Y/Z numbers for a specific tile visible on the screen so I can go to that URL and investigate what is wrong with the file?

I am using QGIS 3.6.2 on Windows.


texts: Using Spatial Adjustment to georeference shapefile?

texts:
I am having trouble using spatial adjustment to georeference a shapefile that I created from a CAD drawing to a basemap of satellite imagery.

At first, the shapefile had an unknown coordinate system, and I created about 10 displacement links and the adjust option was still greyed out in the drop-down menu.

I tried defining the projection of the shapefile to NAD 83 UTM Zone 17N, because I would need it in that projection eventually, but the adjust option is still greyed out.

It seems that I can't get the option no matter how many displacement links I add as well.

Am I missing something?


texts:
If I use the standard buffer tool in ArcMap, I can set the units for the buffers. I'm trying to execute a python script to do the same thing, but the resulting buffers are huge and incorrect, indicating that the units are wrong. Is there a good way to set the units in the script, or is this just the result of the projection I'm using? My data is all being projected into a state plane system using feet, but the feature I'm trying to buffer is in WGS 1984 (degrees).

I tried to research a bit but didn't find much.

#assign variable bus routes
routes = r"C:\Users\bobby\Desktop\gis data\Bus redesign.gdb\BusRoutes"
buffDists = [300, 600, 900]
#create arcpybuffer

for buff in buffDists:
    arcpy.Buffer_analysis(routes, routes + str(buff), buff)
    print (str(buff) + " is processing")

print("script complete")


texts: What might be the reasons for polyline not splitting at the intersection even after using the command 'Split Lines with lines" in vector overlay tools


texts:
I'm using r.patch (within QGIS, GRASS 7) to patch two Float32 rasters both of which have a respective filesize of under 200MB.

When I save the memory layer raster created by the r.patch function, the new GTIFF is float64 and has a filesize in excess of 7GB.

How can I prevent this, or how can I reduce the created file size?  I've tried using GDAL-Translate to compress to JPEG but this just created a file with all nan values.


texts: Calculating best memory allocation for GeoServer

texts:
I am having lots of OutOfMemory exceptions with GeoServer 2.6.0, and now I'm trying to increase the memory allocated to it and Tomcat. How can I figure out the best allocation of memory to set in the configuration?

This is a virtual machine with 16Gb of RAM.


texts: Keeping certain scale of OS layers while zooming

texts:
On QGIS I have a WMS OS Premium Raster layer supplied by an external company. As I zoom in and out (i.e. change scale), the map changes to reflect the different OS map scales.

How can I stop this happening and ensure that, for example, the 1:10,000 scale OS map is displayed regardless of zoom/scale level in QGIS?

I've searched on the internet for a while and found no answers.


texts: st_clip and gdal_translate: "Error browsing database for PostGIS Raster properties" - Why?

texts:
I'm using shell-script, PostGIS 2.0 and GDAL for trying to clip raster apart, using an overlapping raster and a the_geom-object. At first I create a view to clip the overlapping part of the raster. It gives me several raster-tiles as results, which I bundle with st_union as one raster. It works all fine and I can see results in the tables.

psql -h localhost -p 5432 -d testlauf -U postgres <<EOF
CREATE OR REPLACE VIEW profilschnitt.cliptest AS
SELECT st_union(st_clip(a.rast, b.the_geom))
FROM rastertable a
JOIN tiletable b
ON (st_envelope(a.rast) && b.the_geom);
EOF

Now I want to see how the unified raster looks in GIS and try to export the raster via gdal_translate into SAGA-rasterformat:

gdal_translate -ot Float32 -of SAGA PG:"host=localhost port=5432 dbname=testlauf user=postgres password=postgresdb schema=profilschnitt table=cliptest mode=2" ~/Desktop/cliptest.sdat

Problem is, that GDAL tells me the following:

  ERROR 1: Error browsing database for PostGIS Raster properties
  GDALOpen failed - 1

My db-connection is set up correctly and gdal-translate already worked with other raster-views/-tables too. So I assume that the problem must be somewhere within the query-view...


texts: Use map projections to find closest point over dateline

texts:
I want to identify closest points to a polygon over the dateline in an OGR/Python script. I found this workaround for PostGIS but I'm not 100% certain about the applicability in my case. This may be because I'm not very good at reading SQL but I'm also not convinced about the fixed longitude values to detect dateline intersections (-178.75 lon) in the script. I experimented with map projections because I want to process data automatically, regardless of the input data and the distance to the dateline.

Problems don't occur when the longitutes are outside the -180 / +180 degree range (as long as coordinates are adjacent), but when, for instance, the polygon is at -179° longitude and a close point is at +179° longitude. The first image shows a case, where the closest points are correctly identified. Imagine the orange line to be the dateline and the green lines the connectors between points and closest point on the polygon outline. The second case is what I would like to consider.

This is OK:

This is not OK:

I experimented with map projections. I defined a central meridian for the longitude of the polygon's centroid. However, this didn't solve the issue as in the 'not OK' case, the point which is causing problems simply has a longitude of +350° while the polygon has a central longitude of 0°. For calculations, points and polygon seem far off even though they are very close on a spheroid.
Can you help me find a workaround for this? Are map projections suitable for the consideration of dateline intersections?


texts: How can I install the greenbrown package to calculate a trend of a NDVI time series?

texts:
I have a NDVI series and I would like to calculate the trend and then test the significance in R.
I thought I could use the greenbrown-package, but I struggle with installing it.
This is the code I used:
setwd("xxx")

install.packages("strucchange")
install.packages("bfast")
install.packages("phenopix")
install.packages("quantreg")
install.packages("ggplot2")
install.packages("raster")

library(raster)
library(ggplot2)
library(strucchange)
library(bfast)
library(phenopix)
library(quantreg)

install.packages("Kendall")
install.packages("fields")
install.packages("ncdf4")

library(Kendall)
library(fields)
library(ncdf4)

#install.packages("greenbrown", repos="http://R-Forge.R-project.org")

# Rtools

PATH="${RTOOLS40_HOME}\usr\bin;${PATH}"
writeLines('PATH="${RTOOLS40_HOME}\\usr\\bin;${PATH}"', con = "~/.Renviron")
Sys.which("make")
install.packages("jsonlite", type = "source") #Beispiel

install.packages("greenbrown_2.4.3.tar.gz", repos = NULL, type="source")

library(greenbrown)

### Laden der NDVI-Zeitreihe

img <- stack("versuchsgebiet_ndvi_serie.tif")
print(img)

plotRGB(img, r = 1, g = 15, b = 36, stretch = "lin")

I couldn't install the greenbrown package just like that, so I tried to download it but I keep getting the error message:
install.packages("greenbrown_2.4.3.tar.gz", repos = NULL, type="source")
WARNING: Rtools is required to build R packages but is not currently installed. Please download and install the appropriate version of Rtools before proceeding:

https://cran.rstudio.com/bin/windows/Rtools/
Installing package into ‘C:/Users/ronja/Documents/R/win-library/3.5’
(as ‘lib’ is unspecified)
* installing *source* package 'greenbrown' ...
** R
** data
** inst
** byte-compile and prepare package for lazy loading
Error in loadNamespace(j <- i[[1L]], c(lib.loc, .libPaths()), versionCheck = vI[[j]]) :
  there is no package called 'ggplot2'
ERROR: lazy loading failed for package 'greenbrown'
* removing 'C:/Users/ronja/Documents/R/win-library/3.5/greenbrown'
In R CMD INSTALL
Warning in install.packages :
  installation of package ‘greenbrown_2.4.3.tar.gz’ had non-zero exit status
> library(greenbrown)
Fehler in library(greenbrown) : es gibt kein Paket namens ‘greenbrown’


texts: How to get pre-change Object state from a Class Extension change event in ArcObjects?

texts:
I have an Esri based system with a class extension registered for the feature classes. This extension implements IObjectClassEvents, and handles OnChange. I get passed the modified object, but after the change was made (e.g. an object was moved from within ArcMap).

However, I also need to know the state of the IObject before the change was made. How do I go about that?


texts: Installing and using ArcGIS REST API Connector plug-in?

texts:
I downloaded this plug-in  http://plugins.qgis.org/plugins/connector/ from the plugins link at QGIS.org and it unzipped into a folder entitled "connector" in my "downloads" folder (running Mac Yosemite 10.10.5) but as with everything to do with QGIS, it comes with no read-me instructions. I have no idea what to do with the files that were downloaded.

Is there some folder I am supposed to put them in? or do I go to the QGIS app and point it at the plug-in somehow?

Meanwhile I tried signing up for a public account at ArcGIS and that is a hopeless affair. Is that necessary to use this plugin or does the plugin have a path to public maps there without a sign-in?

ArcGIS allowed me to enter information, it knows my username so I can't register it again but it won't allow me to sign in. It never sent an acknowledgement to my email but when I click forgot password (which I haven't, but what else can I do) it sends me a note that my account has been blocked so it has the right email.

Is there something about GIS in virtually every aspect that makes it the most intransparent uninelligble unintuitive platform to get running or am I just denser than the average user?


texts: Copying the standard HTML design of function descriptions in the Function editor

texts:
I am creating my own function in the function editor and want to copy the HTML design of the standard function description, but I can't find it anywhere.
Where can I find it?
UPDATE:
I found a description of the area function in the QGIS repository, and it contains a description that corresponds exactly to the screenshot below. Apparently this description is somehow handled separately to give the description a standard style.
UPDATE 2:
It looks like I found the right piece of code that handles the JSON data for the HTML style I want.


texts:
I am new to GIS. As a part of my research, I want to figure out among OSSIM, Saga API, Terralib and Mapnik which one will best satisfy the following functionalities and why (all for desktop GIS application).

Loading/Saving major GIS files (library OGR/GDAL)
Database compatibility.
Vectorial element linked to a field of the database
Georeferencing.
Zooming, panning
Multiple vector Layers
Performances (loading of large shapefile, quick zoom, quick panning, tile support, etc)
User defined symbol
Zone selection (polygon, circle, square, etc.)
Animations
OGC compatibility

Also which all can be used with Qt.


texts: Is there any way to create a legend from label classes or annotations in ArcGIS 10.2.2?

texts:
I have a single "geonames" layer with 56 label classes and another few layers with around 10 label classes. Each class has a different font and colour or graphic. I would like to add most of these labels along a regular symbol in my legend.  I am building my legends by converting them to graphics but it would still be nice to not to have to create a representative label for each label class manually.  In case of symbology I simply make a legend, tweak the font and size of the symbols, export and work with the static graphics. I am wondering if anyone has a script or a procedure that can make a dump of all the label symbols, each a representative of each class.  (I believe this could be obtained programatically since ArcMap is already doing it - Labels Tab -> Expression -> Verify -> resulting Expression Verification shows a sample label graphic.) I would like to shomehow spit out this same graphic along with the label class name onto the map layout.

Programatic VB.NET / ArcObjects solution showing how retrieve the label class name and it's symbol and would be acceptable.

I found a similar post here asking a very similar question about version 9.3 and although marked as answered it doesn't really provide a satisfactory solution. I am using 10.2.2

Figure 1: ArcMap Label Expression Verification showing that it already generates the graphic:

Figure 2:  A example I made manually showing what I would like to get from a layer with 4 label classes:


texts: Generate cache tiles for EPSG:27700 for use with OpenSpace

texts:
I am using GeoServer 2.1.0 (Windows Jetty local install for dev), OpenLayers and Ordnance Survey openspace API and have a set of shape data that I am using as an overlay.

Everything is working as expected and drawing correctly but I need to create cached tiles to increase drawing speed.

In understand that the default cache tile generation is either EPSG:4326 or EPSG:900913 and to add EPSG:27700. I understand that this is achieved by creating a geowebcache.xml file containing the definition.

Unfortunately I haven't found an example of a geowebcache.xml file and so if any one could help me I would be most grateful.

BTW: I'm not a cartographer just a developer!


texts: Why do ArcPy SQL expressions use three double quotes?

texts:
I am using ArcGIS Pro 2.4 and I am trying to understand why ArcPy SQL expressions use three double quotes. ESRI's "Specifying a query in Python" page, https://pro.arcgis.com/en/pro-app/arcpy/get-started/specifying-a-query.htm,
gives two examples of using SQL statements (both of which seemingly have where clauses of "roadclass" = 2 wrapped in three double quotes):

Example without field delimiters:

import arcpy
fc = "D:/St_Johns/data.gdb/roads"

# Create a search cursor using an SQL expression
cursor = arcpy.da.SearchCursor(fc, ("roadclass", "name"), """"roadclass" = 2""")
for row in cursor:
    print(row[1])

Example with field field delimiters:

import arcpy
fc = "D:/St_Johns/data.gdb/roads"
fieldname = "roadclass"

# Create field name with the proper delimiters
whereclause = """{} = 2""".format(arcpy.AddFieldDelimiters(fc, fieldname))

# Create a search cursor using an SQL expression
cursor = arcpy.da.SearchCursor(fc, ("roadclass", "name"), whereclause)
for row in cursor:
    print(row[1])

It is my impression that SQL statements in a SearchCursor constructor in ArcGIS Pro should have string values enclosed in single quotes and that the overall expression should be in double quotes as in " PARCEL = 'A2003KSW' ". That syntax make sense, however I am not understanding the jump from one to three sets of double quotes.


texts: Trouble with adding WFS features in ArcMap

texts:
I add the following WFS address to ArcCatalog just fine and the different feature classes show up in ArcCatalog.  However, when I attempt to drag or add a feature into the TOC of ArcMap, it gives me a "This dataset appears to be invalid" error and the feature gets added to the TOC but is empty with no feature records.

Could I be missing a step here?  I'm using ArcGIS 10.1.  Here is the WFS address:

http://ogi.state.ok.us/geoserver/wfs?service=WFS&version=1.0.0&request=GetFeature&


texts: Not Able to Enable Snapping Management on top of a feature service

texts:
Can you please take a look at This Demo and let me know why I am not able to set up and enable Snapping on top of following feature Service? What I need is snapping to all features in mouse hover/over of each?

var map;
require(["esri/map",
  "esri/SnappingManager",
  "esri/layers/FeatureLayer",
  "dojo/keys",
  "dojo/domReady!"
], function(Map, SnappingManager, FeatureLayer, keys) {
  map = new Map("map", {
    basemap: "topo",
    center: [-88.158805, 41.786075],
    zoom: 18
  });

  var propertyLayer = new FeatureLayer("https://sampleserver6.arcgisonline.com/arcgis/rest/services/CommunityAddressing/FeatureServer/0", {
    mode: FeatureLayer.MODE_ONDEMAND,
    outFields: ["*"]
  });
  map.addLayers([propertyLayer]);

  var snapManager = map.enableSnapping({
    snapKey: keys.CTRL
  });
  var layerInfos = [{
    layer: propertyLayer,
    snapToVertex: true
  }];
  snapManager.setLayerInfos(layerInfos);

});


texts: Deleting text from Map Items in QGIS print composer?

texts:
I have been using QGIS for some time, and what I am experiencing is not something that I have had issues with before. It appears to be a gitch.

I created a map in QGIS 3.6.2. Then, I exported it as a new print layout in the print composer. At some point during the creation of a print layout, I was working with text labels to create a north arrow/final touches, and didn't notice the random symbology that was created in the map item (see photo). Now, it appears I cannot select these symbols to delete them and they appear when I close and reopen print composer and even when I restart QGIS and export a new print layout from my original map. They appear to auto generate, but I cannot seem to be able to select them. Again, they are not in the original QGIS map, only in the print composer/map layout.

How do I remove these text label errors that are somehow embedded in my print layouts?


texts: Reversing coordinates of MultiPoint using Shapely?


texts: Can I join a CSV field to one of any matching values in a single cell within a shapefile?

texts:
I'm working with the Global Administrative Boundaries shapefile, and am have a CSV of values I'd like to join to it.

My issue, however, is that there are multiple ways to spell the region names. Look at the Catalonia example below: there are 6 different ways to spell the region name, depending on the language.

My question is this: If the join field value is "Catalonia" in my CSV, how would I go about joining my data to the shapefile on any of the values in the VARNAME_1 column? Or should I split VARNAME_1 into multiple columns, and then somehow join on "one of any, if any" values?

Note 1: There isn't a consistent number of languages in the VARNAME_1 column.
Note 2: Using QGIS 2.18.2, MacOS


texts:
In composer I can select multiple items but if I change their attributes the change is only applied to the first item selected.

Example:

Select multiple text items using shift select (Mac OS X), change font colour or size. Change is only applied to the first item selected all others remain unchanged.

Is there away to apply changes to all, or the equivalent of paste attributes to text or set a character style globally?


texts: 3D interactive terrain visualization with GRASS 7

texts:
GRASS 6 had an interesting and interactive terrain visualization tool which was based on Tcl/Tk and called NViz. After my Ubuntu 15.10 updated itself for GRASS 7.0, the NViz interactive tool is missing.

What can I do to install NViz / interactively visualize terrains?

I am using Ubuntu 15.10 / 64 bit / GRASS 7.0.


texts:
I am trying to add a the draw interaction to my openlayers 4.1.0 VectorTile layer. I have been able to add the draw interaction to a Vector layer, but not a VectorTile layer. Below is my js code:

Please note:

When cur_page = "page1" the draw interaction is working.

When cur_page = "page2" the draw interaction will draw the polygon, and on completion of drawing the polygon (either by double clicking or clicking on the first polygon point) an error will be generated and the polygon will disappear.

var defaultFeature = new ol.style.Style({
    stroke: new ol.style.Stroke({
        color: '#0000ff',
        width: 2
    })
});

var selectedFeature = new ol.style.Style({
    stroke: new ol.style.Stroke({
        color: '#ff0000',
        width: 2
    })
});

if(cur_page == "page1"){
    var latitude = -99.624;
    var longitude = 38.6855;
    var zoom = 3;
    var raster = new ol.layer.Tile({
        source: new ol.source.OSM()
    });

    var source = new ol.source.Vector({wrapX: false});

    var vector = new ol.layer.Vector({
        id: 'drawlayer',
        source: source,
        style: defaultFeature
    });
    var layers = [raster, vector];
}

if(cur_page == "page2"){
    var latitude = -97.587453;
    var longitude = 35.471398;
    var zoom = 11;
    var layer = 'test%3Alayer';
    var projection_epsg_no = '4326';

    var source = new ol.source.VectorTile({
        tilePixelRatio: 1,
        tileGrid: ol.tilegrid.createXYZ({maxZoom: 19}),
        format: new ol.format.MVT(),
        url: 'http://192.168.10.15:8080/geoserver/gwc/service/tms/1.0.0/' + layer + '@EPSG%3A'+projection_epsg_no+'@pbf/{z}/{x}/{-y}.pbf'
    });

    var layers = [
        new ol.layer.Tile({
            source: new ol.source.OSM()
        }),
        new ol.layer.VectorTile({
            style:defaultFeature,
            source: source
        })
    ];
}

var projection = "EPSG:4326";
var latLong = ol.proj.fromLonLat([latitude, longitude]);

var map = new ol.Map({
    layers: layers,
    target: 'map',
    view: new ol.View({
        center: latLong,
        zoom: zoom
    })
});

$(document).on("click", ".map-side-panel", function() {
    var inter_type = $(this).data("interaction-type");
    var inter_descr = $(this).data("interaction-description");

    addInteraction(inter_type, inter_descr);
});

var draw;
var modify;
var selectInteraction;
var writer = new ol.format.GeoJSON();
function addInteraction(inter_type, inter_descr) {

    if(inter_type == 'Draw'){
        clearInteraction();
        draw = new ol.interaction.Draw({
            source: source,
            type: (inter_descr)
          });
          map.addInteraction(draw);
    } else if(inter_type == 'Modify'){
        ...
    } else if(inter_type == 'Source'){
        ...
    } else if(inter_type == 'Data'){
        ...
    } else{
        clearInteraction();
    }
}

function clearInteraction() {
    map.removeInteraction(draw);
    map.removeInteraction(selectInteraction);
    map.removeInteraction(modify);
}

This is the error I receiving in Chrome developer tools console:

ol-debug.js:64479 Uncaught TypeError: this.source_.addFeature is not a function
    at ol.interaction.Draw.finishDrawing (ol-debug.js:64479)
    at ol.interaction.Draw.handleUpEvent_ (ol-debug.js:64180)
    at ol.interaction.Draw.ol.interaction.Pointer.handleEvent (ol-debug.js:18743)
    at ol.interaction.Draw.handleEvent (ol-debug.js:64129)
    at ol.Map.handleMapBrowserEvent (ol-debug.js:31688)
    at ol.MapBrowserEventHandler.boundListener (ol-debug.js:6973)
    at ol.MapBrowserEventHandler.ol.events.EventTarget.dispatchEvent (ol-debug.js:7393)
    at ol.MapBrowserEventHandler.handlePointerUp_ (ol-debug.js:12691)
    at ol.pointer.PointerEventHandler.boundListener (ol-debug.js:6973)
    at ol.pointer.PointerEventHandler.ol.events.EventTarget.dispatchEvent (ol-debug.js:7393)

How do I add a draw interaction on my VectorTile layer?


texts: Scale national economic data down to a watershed

texts: Economic Input Output databases such as e.g. the World Input-Output Database by the European Commission contain economic (and often environmental) statistics. The spatial resolution of data provided in these databases is usually given on a country level e.g. "the turnover of the agricultural sector in the Netherlands per year is X milion €". I would like to link national trade data to individual watersheds. However obviosuly, the boundaries of a watershed rarely coincide with the political boundaries of a country. Does anyone have an idea/a suggestion, how I can relate national economic data to individual watersheds within the respective country?


texts:
I'm looking for a dataset or an API that provides information about toll costs of european roads.

Something similar to ViaMichelin, that is what we are using now, but preferably free of charge.

Does somebody know some free dataset with this information?

I was thinking in OSM, but I do not know if offers this information.


texts:
I am currently creating an online GeoPortal for open data with the idea being that a user sources open data from a source such as the following: http://www.opendataphilly.org/opendata/resource/60/aerial-photography-2008/

The user then logs into the portal and using one of the menus can import this data as XML metadata using this example and the portal take this file and convert it into WMS ready to then display to the user.

Is this possible? And if so how would I go about achieving this?


texts: Uploading table (shape) to asset using Google Earth Engine?

texts:
I am trying to upload a table (shape) to my assets and it has been uploading for 7 days.

Is it normal?

I do not understand why spends so many time.


texts: How can I find start and end points from a road layer?

texts:
I have a shapefile with geometries MultiLineString and Point. I have loaded the shapefile into PostGIS. From that table I cannot find any (lon,lat) except geometrical values with some SRID,.... columns.

How about knowing of nodes, edges, start point and end point for a particular layer that has above geometries; and to create a start point and end point from a shape file

I have tried few shapefiles, finally stopped with

Error:Start vertex was not found.


texts:
I keep receiving this error whenever I try to use this API call:

https://{username}.carto.com/api/v2/sql?q=SELECT+*+FROM+OBS_GetDemographicSnapshot(CDB_LatLng(40.760410,-73.964242))

I do use my username at the front so that is correct, however I don't know how to change the url so that it will work. Is there further documentation I could use to get a more in depth understanding? The basic documentation for this paticular api isn't that great with saying how I should be using functions.

Here is the complete error:
{"error":["function obs_getdemographicsnapshot(geometry) does not exist"],"hint":"No function matches the given name and argument types. You might need to add explicit type casts."}


texts: How to Reference a specific layer in a WMS Feed

texts:
I am looking to access HRDPS.CONTINENTAL - Air temperature (°C) to use in an ArcGISOnline Map from the WMS Feed http://geo.weather.gc.ca/geomet/?lang=E&service=WMS&request=GetCapabilities

Does anyone know how to pull in the specific layer that is available in the WMS.

The data is the Government of Canada GeoMet Geospatial Web Services


texts: QGIS - Rendering of two overlapping attributes of one line feature

texts: I have a network of routes and some polylines belong to 2 different routes at the same time (see picture). Is it possible to display these lines in a way that I do not need to setup a rule-based rendering for the exception of the depicted case? I would like to achieve that both lines lie next to each other like in subway system maps...


texts:
Can anyone help me with the distance matrix tool in QGIS? Is there a tutorial about it? It sometimes crashes when I use it and the values that it produces seem strange. I am a beginner in QGIS, so I extensely searched for a good tutorial, paper or document about the distance matrix plugin but I could not find something about it, which really helped me.

So I still have several questions:

Do my coordinates have to be all in the same CRS format?
How can I produce the distance from a csv where I have many points and I want to get the distances from each point to each point in kilometers?

Are the results in degrees? How is that possible that a distance is measured in degrees?

I hope you can help me and am always happy for any constructive ideas or suggestions!


texts:
Completely new to GIS so this may be the most basic thing in the world however I've followed several other posts on the same question and it does not seem to be working.

My end goal(hopefully) is to take the files(.shp, .dbf, etc) that I've downloaded from the NHESP and turn certain habitats into lat and long values in a .csv.

Basically, I'm trying to determine if bee-hive distance to core habitats etc makes an impact on health, ability to survive, etc(if anyone cares).

I've tried using the field calculator to generate an X & Y coordinate, as well as the "save vector layer as" and changing the geometry to "as_xy."

Is there something I'm missing, or is it possible the files I have simply cannot generate or do not contain the information needed for lat/long?


texts: Automatically aggregate German ZIP codes by the first two digits in a shapefile

texts:
I have a shapefile which contains all German Zip-Codes in a a five-digit format (e.g. 80637). For what I intend to do, this is way to much detail. I would like to aggregate all Zip-Code areas by the first two digits (e.g. 80).

My approach is ordering ZIP codes in the QGIS attributes table numerically, then selecting them manually, then joining them. This is somewhat cumbersome I'll have to do this a hundert time.

Is there a way of doing this programmatically? I tried automate this process following this instruction, but I failed at importing sextante.


texts: 'Result' object has no attribute 'GetMessage', ArcGIS Example

texts:
i'm using the following ArcGIS example to calculate the surface and volume of a raster. I added the script to my toolbox:

'''****************************************************************************
Name: Surface Volume Example
Description: This script demonstrates how to use the
             Surface Volume tool.

****************************************************************************'''
# Import system modules
import arcpy
from arcpy import env
import exceptions, sys, traceback

try:
    # Obtain a license for the ArcGIS 3D Analyst extension
    arcpy.CheckOutExtension("3D")
    # Set environment settings
    env.workspace = "C:/output"
    # Set Local Variables
    inSurface = "extract1"
    #Execute SurfaceVolume
    result = arcpy.SurfaceVolume_3d(inSurface, "", "BELOW", "2", "1", "5")
    print result.GetMessage(0)

except arcpy.ExecuteError:
    print arcpy.GetMessages()
except:
    # Get the traceback object
    tb = sys.exc_info()[2]
    tbinfo = traceback.format_tb(tb)[0]
    # Concatenate error information into message string
    pymsg = 'PYTHON ERRORS:\nTraceback info:\n{0}\nError Info:\n{1}'\
          .format(tbinfo, str(sys.exc_info()[1]))
    msgs = 'ArcPy ERRORS:\n {0}\n'.format(arcpy.GetMessages(2))
    # Return python error messages for script tool or Python Window
    arcpy.AddError(pymsg)
    arcpy.AddError(msgs)

The script does not work at all. The error i can't seem to avoid is

'Result' object has no attribute 'GetMessage'.

SurfaceVolume_3d works when the output is a file and the result is not printed.

I just don't understand why it does not work in the Python-Script. I couldn't find any helpful information on this behavior.


texts: Displaying different ID label in Google Earth using GEarthView plugin?

texts:
When using GEarthView to show KML file in GE of a layer from QGIS, the label that shows up in GE is the Feature ID.

How can I change what is shown as the item label in Google Earth?


texts: pdf creator software to interface with Civil 3D

texts:
I've experienced a number of problems with various pdf Creator packages throughout my career, especially with different versions of AutoCAD I've used.  For example, blobby, illegible text upon export to pdf, conflicts causing the process to crash, converting north arrows to a big black blob, etc, etc.  The pdf print driver inherent to CAD is also sub par (in my experience).

I found the free pdf Creator at THIS LINK to be the most highly recommended to interface with higher level AutoCAD licences (according to the AutoDesk forum).  Feedback would be greatly appreciated regarding this software, or any other reccomendations any of you may have.  Purchasing the Adobe Pro is acceptable, however:

Is it coded to work well with AutoCAD Civil 3D?
Will the free option at the link I posted serve my needs (free is obviously preferred)?
Is there a better option than the two I've mentioned?

I'm looking for actual user experience (with CAD), not random free pdf creator options/suggestions.

Thanks in advance

CORRECTION I noticed at the AutoDESK forum that the answers were not listed by highest rank.  There are no vote-ups whatsoever on the suggestions posted there.  I assumed the answer at the top of the page was "the best" answer, while in actuality, it was just the "most recent" answer.  I should have looked more closely.


texts:
I have a group layer. Its purpose is to show feature layers over a base map. Those feature layers are logically grouped into categories and we provide the ability to add/remove whole categories through a menu.

The categorization is done using a dictionary

layerCategories <string [category name], List<IFeatureLayer>>

When I want to delete a category, I do something like the following

void deleteCategory( string categoryName){
  foreach(var layer in layerCategories[categoryName]){
    theGrouLayer.Delete(layer); // delete layer from group layer
  }
  layerCategories.Remove(categoryName); // remove layers (along with the category)
}

The code works, what I ask if it could be done better. I mean, are the layers disposed properly this way or should I do something explicitly?

Thanx in advance


texts: How to modify imposm data mapping in order to include the user name?

texts: I want to import OSM data into postgresql with imposm. How do I have to modify the import process in order to have the user name as an attribute in my tables? Do I have to change that in the data mapping file?


texts: Create time series graph with multiple rasters

texts:
I need to create a time series graph for a few random cells with known underlying attribute such as "forest" or "agriculture" based on the VCI (Vegetation Condition Index). I have been reading and looking into time series in QGIS but have not found anything that can do what I need it to do, as most only draw a time series movie of the area and usually only work with vector data, which I need neither of those 2 types of applications.

The time series is 16 years worth of data, based on 16 day NDVI composites which gives a total of about 368 Rasters. I then need to choose some random cells within areas that have known vegetation landuse classes such as a forest. Each cell selection should then have 368 total values (for each Raster). From these 368 values I would like to extrapolate them into a time series line graph format. Something similar to the below image:

So in theory I should be able to make several graphs that would have 368 values on the X-axis representative of the changes through time of the VCI of those specific cells.

Ideally it would be good if the application can either draw this graph and save the values to a .txt or .csv file but also even if it can only save the values to a .csv I can then plot this time series easily in Excel.


texts: Remove intersect/touch from polygons, to making one with JTS/GeoTools

texts:
Is there a way to remove an overlap from multiple polygons to make one final polygon?..programmatically, I am using the GeoTools or/and OpenGIS/JTS.. any ideas to put me on the right path...

example initial:

final:


texts: Intersect two same extent shapefiles, but got error 000117

texts:
I have used intersect (Analysis Tool) to intersect two polygon layers. But I got error 000117 after ran intersect tool. Same error message happened when I use clip (Analysis Tool).  The answer and two shapefiles property shown as below.

So I was wondering how to amend this?

@GetSpatial The first layer is the high-density area shapefile of Washington-Arlington-Alexandria MSA polygon shapefile. And the second layer is multiple ring buffer of MSA Virginia part roads. I want to calculate the high-density area within each buffer. So I used intersect tool. The red area is the high-density area and yellow area is road buffer.

I tried to use other spatial operations and there was no problems happened before.


texts:
We have some GIS data, which was converted from AutoCAD dwg to TAB using mapinfo, then converted to shapefiles as that's what we are required to create. Editing the data fields was done in the shape files using QGis.

Now for some of the files, the coordinate system settings created using MapInfo were wrong, so I tried to solve this by using the 'save as' option on the edited shapefiles from within QGis, setting the correct projection.

For me, all lines up as it should. However, for the one requiring the data (using Esri software) it doesn't.

I think this is due to QGis creating a .qpj file along the .prj file, and QGis knows to look at the qpj but esri software doesn't.

First of all: Could this interpretation be correct?

Second: if it is correct, how do I 'fix' this so it aligns as it should in the Esri software?


texts: How to Run Ad Hoc Database Queries in Open Jump?

texts: Which plugin allows you to run ad hoc queries from postgis in Openjump? How do you install it in Ubuntu?


texts: Extracting Table data from Raster created by iterated visibility test

texts: I would like to run a visibility test using the observer function on each of 500 random points in a region, and end up, ultimately with a table showing the visible area results (count cell for each value of 1 from each table). First I have not successfully iterated the visibility test as there are more than 16 points. I want each point to be extracted from the feature layer, its visibility analyzed and out-putted to a table with the value and count. How can I get model builder to look now at the attribute tables of each raster, and extract only the row with "value =1"?


texts: Open attribute table script in QGIS 3.6 graphical modeler

texts:
I am looking for a script which will allow me to open an attribute table of the final output of a graphical model in QGIS.

I know that you can use the following lines of code in the console but cannot find work out how to create a script from this that is accepted by the model:

layer = iface.activeLayer()
iface.showAttributeTable(layer)

I want the output of my last stage in the model to be the input to this script so that the attribute table opens when finished.


texts: Ecoregions shapefiles: subsetting, extracting and summarizing data from raster

texts:
I am attempting to subset and project an ecoregions shapefile from : EPA ECOREGIONS
my data set is found here: SHAPEFILES and RASTER
When I run this code, two major problems occur: 1) I cannot get the subset of ecoregions to plot on top of my raster file and 2) I cannot extract the values (summary information) from the subset of ecoregions from the raster. Can anyone help with this?

require (raster) raster<-raster('FG_RCP26_2050.tif')
plot (raster)

ecoregions.all<- readOGR("NA_CEC_Eco_Level2.shp",layer="NA_CEC_Eco_Level2")
head (ecoregions.all)
plot (ecoregions.all)

#this is the subset of ecoregions that I need
NA_L2<-c('5.1','5.2','5.3','8.1','8.3','8.4','8.5')

#subset the ecoregion.all to just those in the northeast
#this must be where the problem is
ne.ecoregions = coregions.all[match(toupper(NA_L2),toupper(ecoregions.all$NA_L2CODE)),]

#plot ne.ecoregions on top of raster and extract summary statistics for each ecoregion
#summary statistics =distribution of values within each ecoregion
plot (raster) plot (ne.ecoregions, add=T)

ecoregions.raster<-mask (raster, ne.ecoregions)
summary (ecoregions.raster)


texts: Given a BBOX and Zoom Range, Return All TMS Addresses

texts:
Given a bounding box using WGS84 lat long pairs (represented any way, not important to me) and a range of zoom levels, I'd like to be able to generate a CSV representing all standard z/x/y TMS addresses for a web map in a spherical Mercator projection.

For example, given the city of Seattle, with a BBOX of: (47.75, -122.14), (47.53, -122.53) and desired zoom levels 12-16, I'd like to see a delimited list something like:

z  |  x   |  y
----------------
12 | 656  | 1429
12 | 655  | 1429
...
14 | 2626 | 5724

Does anyone know of any way to do this, preferably in Python?


texts:
Elevation data isn't integrated into OSM streets directly, but both OSM street paths and elevation resources are separately available.
My objective is to find the largest hills on roads in any given area. What is the best set of tools and elevation data available to do this? Has something like this been done before? I'm coming from an R/C++ and JavaScript background.


texts: Where have all the projected national CRSs gone in QGIS 3.4.3?

texts:
I have been working with different geographic projections in QGIS 3.4.3.  I cannot find any Australian national geographic projections, such as GDA94Z55 in the new version. They were there in QGIS 3.4.2.

Any clues as to how I can access them so they are part of the default CRSs in the program?


texts: Changing CRS from Geographic to Projected on vector layer

texts:
I am trying to run a Iso-Area as Interpolation (from Point) calculation using QNEAT3 in QGIS 3.8.2.

To do so, I need to change the CRS from Geographic to Projected. I did that to my raster layer, but when I try to do the same to any of my vector layers, they just disappear.

I am trying to change it from EPSG:4326 - WGS 84, to EPSG:25832 - ETRS89 / UTM zone 32N.


texts: How to perform an attribute-weighted sum of multiple polyline features from a single vector layer?

texts:
I have a shapefile of ~8000 polyline vector features representing traffic journeys on road routes. Each has an attribute representing volume of traffic for that journey. The routes are all generated by the same routing engine so the node coordinates should match exactly where routes follow common sections of road. I am wanting to effectively add these together in a "map algebra" sort of way to give a representation of network usage. Keeping the result in vector would be preferable, but ultimately not essential.

I can see how this could be done conceptually simply by converting each feature to a separate raster of pixel values = traffic volume, then summing all the rasters but this would practically be very painful to combine 8000 high-resolution rasters. Ironically, QGIS is effectively already doing this visually for me - I am displaying the layer using combining of features by "addition" method.. but I can't see any way of getting that directly out into a much higher resolution raster format than my display:

Surely there must be a vector way to do this? (QGIS or ArcGIS) If it makes it easier to not go the shapefile route I have all the route nodes in CSV format, one line per route. The polygon "union" then "dissolve" method outlined in How to union polygons and add attribute values of combined output features looked promising, but Union does not work on line features it seems. Flicking through the ArcGIS Network Analyst "Route Analysis" pages looked potentially promising if rather OTT for this purpose and I'm guessing would require a lot of messing round to get there - I'm sure there must be a simpler way!


texts: Freezing (Avoid Re-Drawing) a MapView In a Multi-Tab/Multi-MapView WPF Application

texts:
I have a WPF application with two different MapViews in two different WPF tabs.

When I switch tabs between the first and second MapView items, the MapViews completely reload the base map and all layers.You can see it happening.

Is there a way to make the MapViews static so that they don't have to re-draw every time you switch from one MapView tab to another?


texts: Adding more information to a PostGIS table via QGIS

texts:
I have a table within PostgreSQL/PostGIS that contains all point data. I need to add another large dataset to this existing table without overwriting the table. When I go into DB Manager within QGIS and try to import the new information I keep getting an error saying that the table I'm trying to write to already exists.

How can I add more information to this table without over writing the existing information already in the table? I wish to be able to do this with QGIS because it does all of the GEOM information for me.


texts: Getting error while adding few layers to GeoPackage using PyQGIS?

texts:
I get error message when trying enter few layers to geopackage using PyQGIS 3.6.2.

This is my code:

from qgis.core import   QgsVectorFileWriter
options = QgsVectorFileWriter.SaveVectorOptions()
options.actionOnExistingFile = QgsVectorFileWriter.CreateOrOverwriteLayer
options.EditionCapability = QgsVectorFileWriter.CanAddNewLayer
options.onlySelectedFeatures = True
GpName = 'c:\\qgis_dat\\proj\\proj1'

for layer in  iface.mapCanvas().layers() :
    namee = layer.name()  # it can be 'roads' or 'buildings' ...
    options.layerName =  namee
    _writer = QgsVectorFileWriter.writeAsVectorFormat(layer, GpName, options )
    print (_writer)
    # print the error's message :
    #    (2, 'Opening of data source in update mode failed (OGR error: )')  .
    # the geoPackage and layers are not created


texts: Using IN Clause to Select By Attributes in ArcPy?

texts:
I am working in ArcMap 10.1, trying to write a script to automate the selection of a feature class. I have two input parameters:

stop1 = arcpy.GetParameterAsText(0)
stop2 = arcpy.GetParameterAsText(1)

I want to use those input paramenters with an IN statement. So if stop1 and stop2 are found in the "Address" field, the expression would look something like this:

expression = ''' "Address" IN '{0}','{1}' '''.format(stop1,stop2)
arcpy.SelectLayerByAttribute_management(layer,"NEW_SELECTION",expression)

I've manipulated the expression equation several ways, but it errors out every time.


texts: How to keep point number when converting to a line

texts:
I have a layer of point that I convert to a line with the plugin Point2One to be able to do apply a smoothing with the generalizer tool. It works well! then I convert my layer back to points using the vector tool extract vertices and when I look at the attribute table everything looks good except that all my points number became 201. The thing is I need the point numbers.

My point number is the field 6 in the image

I tried to duplicate them, add some letter before, add an other row of point number starting from 1 (mine does't), change the field type. Nothing worked I get 201 in all of these cases...


texts:
I got a GIS Data-set for swedish cities. The problem is that they have strange names such as:

TrÃ¤det

JÃ¤la

MjÃ¶nÃ¤s

Ã„lmhult

Did someone have the samoe problems? I think "Ã¶" for instance mean "ö" for example - but I am not sure if there is a system behind that.


texts: Is there a Java library to load data from a WFS server?

texts: I am writing a java client to display remote GIS data. I want to load data provided in GML by a WFS server. Do you know a java library I could use for that?


texts: Get latitude and longitude of each point forming a MultiLineString

texts:
A MultiLineString is formed by a sequence of points and each point has latitude and longitude.
Use the sql below to get the latitudes and longitudes of each point so that my MultiLineString:

SELECT
    St_AsText (ST_Transform("public"."Shapes".geom,4326)),
    "public"."Shapes"."id"
FROM
    "public"."Shapes"

I want to get the sequence of points of each line in separate rows and lat / long in separate columns example:

   ID            pt_lat            pt_lon
   1001             lat1             lon1
   1001             lat2             lon2
   1001             lat3             lon3
   1001             lat4             lon4
   1001             lat5             lon5
   1001             lat6             lon6
   1002             lat1             lon1
   1002             lat2             lon2
   1002             lat3             lon3
   1002             lat4             lon4
   1002             lat5             lon5
   1002             lat6             lon6

I know I can do that consuming one GeoJSON for example, but want to know if there is a native function of PostGIS for this.


texts:
i'm a python beginner.
I have a table with 4 fields def_score, PA_score, GPS_score and coop_score.
I would like to calculate coop_score based on the first 3 fields with the following code block, but it's not working. Where am i wrong?

defor = 'def_score'
PA = 'PA_score'
GPS = 'GPS_score'   # GPS values are [1, 1.5, 2]

def calc(coop_score):
    if PA == 0 and defor == 0:
        return GPS
    else:
        return ((PA + defor) * GPS)


texts: Displaying attribute values in IntraMaps Roam?

texts:
I am having troubles with displaying attribute values in IntraMaps Roam. Let me give you more details:

My QGIS project is based on an SQLite database
The attribute table in QGIS is displaying the attribute value (which is text) and not the unique ID's (numbers) that belong to each possible value
I created a form for one specific feature type and attribute, so it converts the unique ID's to the actual attribute values, because apparently IntraMaps Roam only displays the unique ID's linked to an attribute value. The form is also in the attachment as a screenshot.
When importing in IntraMaps roam, and select a certain feature to check the attribute values, it only displays the unique ID's. For example (see screenshot in attachment) Dieptenauwkeurigheid: 1. This should be Dieptenauwkeurigheid: tot30cm. When I click on the 'edit' button in the information table, the text is appearing and it uses the unique ID's no more.

So, how come that the form stars working ONLY when you click on the edit button?  Why doesn't it display the attribute values in the information table?  Is this somekind of bug in the codelists?


texts: OpenLayers in QGIS - Getting country map

texts:
My project consists of geocoding all Muslim mosques in Portugal (I have this data in CSV file) using QGIS.

I'm using OpenLayers to get my base map (Google Streets), but I want to get only the map of Portugal, instead of World Map. Is it possible?

I've tried to use OpenStreetMap and GoogleMaps pluggins, but I don't know how to isolate Portugal from the World Map.

I have other shapefiles of Portugal, but it doesn't show all the info as Google Maps.


texts: How to overcome missing data in Time Manager

texts: I am relatively new to QGIS, and while trying to use Time Manager, I am having the following trouble - the date field in my attribute table has some missing values, which are mentioned as '---'. When I tried to run Time Manager, it gave me an error. How can I over account for this?


texts: Postgis KNN vs ST_DWithin for nearest neighbour search with a radius

texts:
We are trying to get nearest neighbors within a user given radius. As far as I understand ST_DWithin degrades as distances are increased. The problem with using KNN is that it keeps on searching even though the radius is crossed.

For example,
let's say  we are searching for 3 petrol stations in 1 km. KNN keeps searching for 3 candidates and we had to filter out those outside the radius after KNN is done. This leads to unnecessary time consumption.

Is there a way that we can do nearest neighbors search within a user given radius which combines the benefit of both. Or is there an alternative better than these 2 ?


texts: How to create raster from points classified as ground using external Python script

texts:
I am writing a Python script that I want to run outside of ArcGIS Pro. I want it to import a LAS, classify the ground points and then create a DEM using only the ground points. The code I have written below creates a LAS dataset and classifies the ground points as expected. It is there I get stuck. I can't get it to use only the classified ground points to create the raster. I've tried using the 'Make LAS Dataset Layer' tool but I can't seem to reference this layer in the external script. It doesn't appear to be in the geodatabase. As such I am getting an error trying to make a raster from this layer.
Not sure if it's relevant but I am using PyCharm to write the script.
import geopandas as gpd
import numpy as np
import pandas as pd
import rasterio
import matplotlib.pyplot as plt
import arcpy
arcpy.env.workspace = r"C:\Tool\DSM_Processing\DSM_Processing.gdb"  # enter path and name of geodatabase workspace here

#First make las dataset from pointcloud
arcpy.management.CreateLasDataset(r"C:\Tool\DSM_Processing\Pointcloud_Data\Pointcloud.las",
                                  r"C:\Tool\DSM_Processing\Pointcloud_Data\Pointcloud.lasd", "NO_RECURSION", None,
                                  "PROJCS['British_National_Grid',GEOGCS['GCS_OSGB_1936',"
                                  "DATUM['D_OSGB_1936',SPHEROID['Airy_1830',6377563.396,299.3249646]],"
                                  "PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]],"
                                  "PROJECTION['Transverse_Mercator'],"
                                  "PARAMETER['False_Easting',400000.0],PARAMETER['False_Northing',-100000.0],"
                                  "PARAMETER['Central_Meridian',-2.0],PARAMETER['Scale_Factor',0.9996012717],"
                                  "PARAMETER['Latitude_Of_Origin',49.0],UNIT['Meter',1.0]]",
                                  "COMPUTE_STATS", "ABSOLUTE_PATHS", "NO_FILES")

##Classify the ground surface

arcpy.ddd.ClassifyLasGround(r"C:\Tool\DSM_Processing\Pointcloud_Data\Pointcloud.lasd", "AGGRESSIVE", "RECLASSIFY_GROUND", None, "COMPUTE_STATS", "DEFAULT", None, "PROCESS_EXTENT")

##make a layer from from the points classified as ground
arcpy.management.MakeLasDatasetLayer(r"C:\Tool\DSM_Processing\Pointcloud_Data\Pointcloud.lasd", "Pointcloud_LasDatasetLayer", "2", "'Last Return'", "EXCLUDE_UNFLAGGED",
                                     "EXCLUDE_SYNTHETIC", "EXCLUDE_KEYPOINT", "EXCLUDE_WITHHELD", None, "INCLUDE_OVERLAP")

#Make raster from layer
arcpy.conversion.LasDatasetToRaster(r"C:\Tool\DSM_Processing\Pointcloud_Data\Pointcloud_LasDatasetLayer", r"C:\Tool\DSM_Processing\Pointcloud_Data\Pointcloud.tif",
                                   "ELEVATION", "BINNING AVERAGE LINEAR", "FLOAT", "CELLSIZE", 0.1, 1)

I get the following error
arcgisscripting.ExecuteError: Failed to execute. Parameters are not valid.
ERROR 000732: Input LAS Dataset: Dataset C:\Tool\DSM_Processing\Pointcloud_Data\Pointcloud_LasDatasetLayer does not exist or is not supported
Failed to execute (LasDatasetToRaster).


texts: Creating Strip Map Index using ArcGIS Pro?

texts:
The Strip Map Index Features tool is listed amongst the Tools that are not available in ArcGIS Pro at version 2.2.

A workaround when needing a strip map index to be created, with an Angle field being populated so that each page can be rotated correctly, would be to use the Strip Map Index Features tool from ArcGIS Desktop but is anyone aware of a workaround for creating a strip map index on a machine when only ArcGIS Pro (with its ArcPy) is installed?


texts: Creating vector features onload in OpenLayers

texts:
Using OpenLayers and JavaScript I wish to create two vector layers and add a polygon to one of them and a circle with a centre point to the other.

My polygon will be in WKT and my circle will be a point with a radius (either as a distance or converted into approx DMS first).

I've seen this page http://dev.openlayers.org/examples/vector-formats.html but  I need to it onload as the variables will be passed from an HTML form into a Perl script which in turn writes the OpenLayers

This is as far as I've got:

Var box = new OpenLayers.Format.WKT().read("POLYGON((1.3 52.1,1.4 52.1,1.4 52,1.3 52,1.3 52.1))”);

This fails error message is: vector.div undefined

Do I need to create the layer first then add the feature?

Also can I get the style from GeoServer or do I need to set it here?
Haven’t a clue how to create the other layer.

Could someone give me and example of how to do this? By the way my base layer / coordinate system is WGS84 and so are the given coordinates.

UPDATE STYLE CODE

If I wished to draw the new vector polygon as a green box with a black border would this work?

var Green_Box ={
    fillColor: "#33CC00",
    strokeColor: "#000000",
    strokeWidth: 1

};
var vectorLayer = new OpenLayers.Layer.Vector("myPolygonLayer");
var format = new OpenLayers.Format.WKT({
    'internalProjection': map.baseLayer.projection,
    'externalProjection': new OpenLayers.Projection("EPSG:4326", {
    'style': Green_box,
    'renderers': renderer
});

Note: was playing with this at work came home found the above but now unable to test it until after Easter.


texts: Adjusting format labelling view through map scale

texts:
I am using 3.10 QGIS. I want to set the label view adjusting with map scale when zoom in zoom out.
Here is an example.

I want to make the zoom out, on the 2nd picture, label adjusting to the map scale so it should be smaller text size. How can achieve this?


texts: I am able to successfully import gpx data to Google My Maps. Most of the gpx files I have imported are in English and I am able to move them into same layer. However, some gpx files are in another language and they can't be placed on the same layer as the English ones. Are there some tricks or tools to convert the gpx file to same format as the others?


texts:
I have a postGIS database in SRS 4326 (WGS84) with several Europe layers working well, but when I load a shapefile projected to the same SRS with QGIS, the layers don't match.

To explain it better I've made two screenshots:
The pink one are the German postcodes from my postGIS database, that are matching with all the rest of my data

The second one, green polygons, are German postcodes too, but from a shapefile.
Both of them are in WGS84 (EPSG:4326, I've reprojected it with QGIS, and shows it in layer propperties). I've marked the coordinates below, to see that the it seems that the second one's coordinates are not correct.

Some Idea that how could I fix the green dataset to overlap with the rest? (pink)
The source of the green dataset is this URL: http://arnulf.us/PLZ

Thanks!


texts:
I am attempting to add a scale bar into QGIS Map Composer for a single layer showing a map of Jervis Bay, NSW Australia. I have changed both the project and layer CRS to 'MGA56 (EPSG:28356)' and have tried some other CRS types but the scale is inaccurate. Each time I insert a new scale it shows an area that would be over 3 km long as only 0.2 m or so. I need the scale to be in km and I need it to be an accurate representation of the area I have on the map.

I am completely new to QGIS and have little experience with any of this software.


texts: QGIS field calculator: Subtracting values from specific fields

texts:
I am trying to create an Isopach map of two subsurface formations. If a well has the two formations in question, then I would like a new layer with the data created so I can contour the points.

I cant figure out how to create the expression needed in the field calculator.

I have two fields, Form Alias and Top Subsurface. I need to specify what two Form Alias's subsurface values I want subtracted and create a layer with this new information. I would think it would look like this: ('CHAT''Top Subsea') - ('LNSG''Top Subsea') = CHAT_to_LNSG_ISO


texts: How to create one shapefile with categories?

texts:
In my QGIS 2.12.3 I want to create shapefiles of different lakes. I know how to do it one by one and this way, each lake gets a different color. But it should also be possible to create one (Over)-layer with all the lakes and their names below it so that I just have to click and all the lakes appear to disappear at once.

How can I create one shapefile with categories?


texts: Is there a way to turn off all labels in ArcMap without clicking on every layer?

texts: I rarely bother with labels, but sometimes it would be nice to toggle them on/off.  Something like the Pause Drawing button, except that it would only turn off labels (all labels in all layers...not just a single layer).


texts: What do I need to build a C++ plugin on windows 7?

texts:
I want to make a plugin in C++ and I need to write and build it in windows (Im better with *nix systems).

While I have skimmed the various updated C++Plugin docs I can understand if I need to build QGIS from source in windows in order to build extra plugins.

What do I need to have in order to build a plugin on windows 7? Is there any tutorial for windows users (especially windows compile-dummies!)


texts: ImageMosaic granule sorting by resolution in GeoServer

texts:
I'm having trouble ensuring that the highest resolution raster is always rendered first in my ImageMosaic.

In the folder .../svfc I have a series of overlapping geoTiff elevation rasters at different resolutions. Each raster contains only patches of data (i.e. is ~50% NoData values), though their overall bounding boxes overlap considerably.

First, I create an external imageMosaic via restAPI using:

-v -u 'admin:geoserver' -XPUT -H "Content-type: text/plain" -d "file:///srv/mal_data/external_data/processed_lidar/svfc" "http://localhost:8080/geoserver/rest/workspaces/mal/coveragestores/svfc/external.imagemosaic"

This creates the datastore and layer in GeoServer as expected, and creates the svfc.shp and svfc.properties files within .../svfc

Next, I create an indexer.properties file within .../svfc containing the following:

PropertyCollectors=ResolutionExtractorSPI(resolution)
Schema=*the_geom:Polygon,location:String,resolution:String
Caching=false
AbsolutePath=false

Finally, I manually set the Granule Sorting (WFS like syntax) field on the svfc layer configuration page in GeoServer to: resolution A

When I examine the layer in QGIS via WMS the highest resolution granules are not rendered first (i.e. 2m res granules are covering 1m res granules). At the bottom of the image below you can see the smaller pixels of the 1m res granule peeping out from beneath the 2m res granule.

Mosaic index file

Contents of directory /srv/mal_data/external_data/processed_lidar/svfc:

indexer.properties
oxfordshire_lidar-composite-dtm_dtm-2m_svfc.tif
oxfordshire_lidar-composite-dtm_dtm-2m_svfc.tif.aux.xml
oxfordshire_lidar-composite-dtm_dtm-50cm_svfc.tif
oxfordshire_lidar-composite-dtm_dtm-50cm_svfc.tif.aux.xml
oxfordshire_merged-products_dtm-1m_svfc.tif
oxfordshire_merged-products_dtm-1m_svfc.tif.aux.xml
sample_image.dat
svfc.dbf
svfc.fix
svfc.prj
svfc.properties
svfc.qix
svfc.shp
svfc.shx
wberks_lidar-composite-dtm_dtm-2m_svfc_20200602_032452.tif.aux.xml
wberks_merged-products_dtm-1m_svfc.tif.aux.xml
wsussex_lidar-composite-dtm_dtm-2m_svfc.tif
wsussex_lidar-composite-dtm_dtm-2m_svfc.tif.aux.xml
wsussex_merged-products_dtm-1m_svfc.tif
wsussex_merged-products_dtm-1m_svfc.tif.aux.xml