texts: Migrating to the use of databases over shapefiles in QGIS

texts:
I mostly use ArcGIS at work and am digging into QGIS at home and for my business. At work we are phasing out the use of shapefiles and moving towards geodatabases. So I am thinking I should look into this at home and for the business. having very limited knowledge of this, this is what I assume I need to do.

Install Postgresql on my Ubuntu Server
Install PostGIS on my Ubuntu Server
Connect to PostGIS from my Windows Machine that currently has QGIS 2.0 installed.

Once this is done I have to make the connection to the server, then figure out how to turn some of my existing shapefile data into a database.

So here is my question:

Does the above process sound correct and are there any step by step resources for the beginner to set implement this process?

Be gentle I have no formal training in GIS, its all on the fly.


texts: How to drape a simple jpeg image over dem of some area?

texts: How can we drape a simple jpeg image which doesnot have any geo referencing information over the digital elevation model in order to visualize it in 3D in QGIS?


texts: Creating Rasyer object from .adf file using ArcPy gives ERROR 000732

texts:
I am trying to use Raster class in Arcpy, to create a raster object from a .adf file

My folder contains all the related .adf files but I am unsure how can I "point" the correct one to create the raster object

for example I tried

import arcpy
from arcpy.sa import *
out_raster = Raster("C:/folder/hdr.adf")

But also I tried with w001000.adf and w001001.adf
and the result I get is an error 007032 that my file is not recognised


texts: GRASS: v.drape, v.split and v.to.db on spearfish sample data

texts:
I am new to GRASS and started training on the spearfish sample data set in GRASS 6.4.3. I selected a bicycle tour and transformed it into a polyline. My next goal is to read out the slope and elevation data from the elevation.dem in order to create a 3D polyline. I want to use the v.drape function, but although the input polyline had an attribute table containing rows for all the street segments that form part of the bicycle tour, my output 3D-polyline doesn't have an attribute table anymore.

v.drape --overwrite input=tour1poly@foo rast=elevation.dem@PERMANENT output=tour1poly3D

Warning: Database connection not defined for layer <0>

First of all: Why is that?

Going on, I split the 3D-Polyline (using v.split along the vertices) into over 300 segments, so I can calculate the slope of each segment with v.to.db. But, of course, I first of all need an attribute table for all the 300+ segments to write their slope into.

Is there a function to help me write the results (coordinates, elevation, length etc. of each line or segment) from v.drape and v.split into an attribute table? So v.to.db can add the slope of each segment to said attribute table?

I know there is v.db.addtable, but how can I access the necessary information to create the proper attribute table with over 300 rows that v.to.db needs?


texts: Ceating cross-sectional styled visualisation of multiple layers in map

texts:
Is it possible to create a map using multiple layers to visual them in a cross-sectional manner in QGIS, ArcGIS or R?
I have attached two examples to illustrate what I want to achieve.


texts: QGIS: No easy way to change .csv columns to integers?

texts:
I am trying to make a heat map of every town in my state. After joining my data (a .csv) to the shapefile, the column I need to map is not visible on the graduated color option under Symbology.

I have read many responses to similar questions here, but none of them have worked for me. Changing the .csv to a .csvt or .txt or other file type does not work.

I also saw posts that it was recommended to convert using the Field Calculator ... but the message says, "This layer does not support adding new provider fields. You can only add virtual fields."

I somehow need to get the data file to be read as integers, not text, if I am correct. However, 12 columns are available to be mapped on my .csv file ... just not the joined column that I need. So I am quite confused.

FWIW, under Symbology, the correct column DOES appear when I tried to map by "categorized," but that is not what I need to do, I need it "graduated."


texts: Some polygons failing to fill with color

texts:
I have several Spatial Polygon Data Frames I am adding successively to a plot in R 4.0.2. The polygons represent US States, National Forests, National Parks, and wildfire perimeters. While it probably doesn't make sense to provide links to all the data for a fully reproduced working example, the following code is a simplified version of the plot I made:
# equal area projection
proj_coords <- "+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=23 +lon_0=-96 +x_0=0 +y_0=0 +datum=NAD83 +units=m +no_defs +ellps=GRS80 +towgs84=0,0,0"

# load libraries
library(rgdal)
library(rgeos)

# Make a bbox of the total extent
e <- as(raster::extent(-1400000,-850000,2300000,2750000), "SpatialPolygons")
proj4string(e) <- proj_coords

# Make a smaller box to add to the plot for Stack overflow example
e2 <- as(raster::extent(-1300000,-900000,2400000,2650000), "SpatialPolygons")
proj4string(e2) <- proj_coords

# Plot
lt.green <- adjustcolor("lightgreen", alpha.f = 0.5)
plot(e,axes=T,xlim=c(-1150000,-1100000),ylim=c(2350000,2700000),lwd=2,border="grey",xlab="m",ylab="m")
plot(e2,add=T,col=lt.green,border="grey")

My question:
In my real plot, e2 are National Forest boundaries, some of which extend beyond the extent of the plot axes. When I try to add them to the plot with
lt.green <- adjustcolor("lightgreen", alpha.f = 0.5)
plot(e2,add=T,col=lt.green,border="grey")

all of the National Forest boundaries plot, but some of them are not filled with color. This is true both when I print the graphic in the console and when I try to export it as a png.
Below is the actual png. The polygons with empty, thin gray borders are the National Forests that are failing to fill with color, while the green ones are correctly filled.


texts: Is there a FME Spatial Filter test for "centroid is completely within"?

texts: Using SAFE FME I am trying to select all polygons from featureclass A whose centerpoints are "completely within" polyon B. Being used to ArcGIS I would have thought this would be available as a test with the Spatial Filter but there are only limited options here. Does anyone know how I could create this filter?


texts: Why is ArcMap deleting data from a field when I add a new one?

texts:
My copy of ArcMap 10 is acting pretty strangely. I have a shapefile and all I've done is added a new field as a FLOAT and for some reason all data in another one of my fields is being deleted! Does anyone know why this is happening?

Thanks,
Alan!


texts: Shared edge between two geojson polygons

texts:
I am interested in getting the coordinates of a shared edge (if exists) between two geojson polygons. I am (currently) using turfjs and geojson-equality libraries

I tried this but I consider it quite inefficient

function sharedEdge(poly1, poly2) {
    for (p1=0;p1<poly1.geometry.coordinates[0].length-1;p1++)
        for (p2=0;p2<poly2.geometry.coordinates[0].length-1;p2++) {
            var point1 = turf.point(poly1.geometry.coordinates[0][p1])
            var point2 = turf.point(poly2.geometry.coordinates[0][p2])
            var point3 = turf.point(poly1.geometry.coordinates[0][p1+1])
            var point4 = turf.point(poly2.geometry.coordinates[0][p2+1])

            if ((equality.compare(point1,point2) && equality.compare(point3,point4))
                || (equality.compare(point1,point4) && equality.compare(point3,point2)))
            {
                    return [point1,point3]
            }
        }
    return false
}

Is there a shorter/more efficient way to do this?


texts: Getting user drawn line returned to tool in python script tool?

texts:
I'm trying to have a user draw a line on the map and then use that line as the input to the fence diagram tool (3D). Code is below and I've tried it several different ways but just can't make sense of it. I can get the tools to allow the user to draw the line, but once the tool runs, it fails.
This is a script tool running in ArcGIS Pro 2.1.

Using this help as guidance:

#variables defined above but excluded.
#get input from user
#what to call the cross section
crossSectionName = arcpy.GetParameterAsText(0)
if crossSectionName == '#' or not crossSectionName:
    crossSectionName = "CrossSection1" # provide a default value if unspecified

#What disstance to search from the line created by the user
Distance = arcpy.GetParameterAsText(2)
if Distance == '#' or not Distance:
    Distance = "50 Meters" # provide a default value if unspecified

#get user drawn line
def getParameterInfo():
        param1 = arcpy.Parameter(
            displayName="Cross Section Line",
            name="in_feature_set",
            datatype="GPFeatureRecordSetLayer",
            parameterType="Required",
            direction="Input")

        param1.value = crossSectionLineLyr

#check out extension
arcpy.CheckOutExtension("3D")

#create fence
arcpy.FenceDiagram_3d(getParameterInfo(), [fenceSurf1, fenceSurf2, fenceSurf3, fenceSurf4], (fenceout + "\\" + crossSectionName))

#check in extension
arcpy.CheckInExtension("3D")


texts: C/C++ function to return minimum set of points required to represent a geometric shape

texts:
I am looking for a method from Grass, GEOS, or similar that will delete those points in a geometry which do not affect its shape.

Example 1

+----+----+
p0   p1   p2

p1 can be deleted without affecting the shape of the line

Example 2

   p0   p1
   +----+
   |     \
p4 +      \
   |       \
p3 +--------+ p2

p4 can be deleted without affecting the shape of the polygon

I am sure this has been asked/answered and there is probably a name for the output point set (akin to convex hull) but I am not sure what to search for.

It would be relatively straightforward (no pun intended) to write my own function but obviously silly if there is something out there already.


texts: How to remove items from Layer List in ArcGIS Operations Dashboard?

texts:
After making a web-map I used some data that I don't want to be displayed in the final Dashboard (as I made tile cache layers for them).  I want to use their statistics but for them not to be included in the layer list on the Dashboard.
I have used the remove from Legend - so it doesn't appear in the Dashboard but I don't want it in the layer list.
Is there a way to have these layers only in the background of the Operations Dashboard? And remove them from the Layer List? (see picture - I want these removed)

Tom


texts: Creating buffer doughnuts in ArcGIS Desktop?

texts:
I have tried Multiple Buffer in ArcMap and I get this image. What I want however is the interval between buffers. I.e. 0-10, 11-20, 21-30 meters.

This tool just uses the point as origin. Hopefully I made my question clear enough.


texts: How can I find EPSG when it's not obvious from the .PRJ file?

texts:
Using Geoserver, I'm trying to get the correct projection for a map layer that is provided by the State of Idaho. I've tried looking up the EPSG code from http://prj2epsg.org/search. However, the search results throws an error.

Here is the location of the data - https://fishandgame.idaho.gov/ifwis/portal/opendata/controlled-hunt-areas-2015-elk

Here are the results of the .prj file

PROJCS["NAD_1983_Transverse_Mercator",GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Transverse_Mercator"],PARAMETER["False_Easting",2500000.0],PARAMETER["False_Northing",1200000.0],PARAMETER["Central_Meridian",-114.0],PARAMETER["Scale_Factor",0.9996],PARAMETER["Latitude_Of_Origin",42.0],UNIT["Meter",1.0]]

I found this page that should help, but I don't know what it all means.

https://www.idwr.idaho.gov/GIS/IDTM/


texts: QGIS Server WMS not displaying correctly in Google Earth Pro

texts: So I'm serving up WMS using QGIS Server in Ubuntu 20 using Apache. It works fine when being consumed by QGIS or even online on OneGeology, but when I use Google Earth Pro there seems to be a problem with boundaries or projection. Any ideas?


texts: More Efficient Spatial join in Python without QGIS, ArcGIS, PostGIS, etc

texts:
I'm attempting to do a spatial join much like the example here: Is there a python option to "join attributes by location"?. However, that approach seems really inefficient / slow. Even running this with a modest 250 points takes almost 2 minutes and it fails entirely on shapefiles with > 1,000 points. Is there a better approach? I'd like to do this entirely in Python without using ArcGIS, QGIS, etc.

I'd also be interested to know if it's possible to SUM attributes (i.e. population) of all the points that fall within a polygon and join that quantity to the polygon shapefile.

Here is the code I'm trying to convert. I get an error on line 9:

poly['properties']['score'] += point['properties']['score']

which says:

  TypeError: unsupported operand type(s) for +=: 'NoneType' and 'float'.

If I replace the "+=" with "=" it runs fine but that doesn't sum the fields. I've also tried making these as integers but that fails as well.

with fiona.open(poly_shp, 'r') as n:
  with fiona.open(point_shp,'r') as s:
    outSchema = {'geometry': 'Polygon','properties':{'region':'str','score':'float'}}
    with fiona.open (out_shp, 'w', 'ESRI Shapefile', outSchema, crs) as output:
        for point in s:
            for poly in n:
                if shape(point['geometry']).within(shape(poly['geometry'])):
                    poly['properties']['score']) += point['properties']['score'])
                    output.write({
                        'properties':{
                            'region':poly['properties']['NAME'],
                            'score':poly['properties']['score']},
                        'geometry':poly['geometry']})


texts: Looking for PDF map downloading capability in OpenGeo GXP

texts:
I am looking for a plugin to OpenGeo GXP in order to download maps in PDF format.

Could you help me?


texts: How to check a polygonal area for diagonal lines?

texts:
I currently have an algorithm which turns a polygon drawn on a map into a contiguous set of 100m^2 chunks. Sometimes the algorithm fails and I end up with chunks missing and a diagonal line remains part of my polygon.
I'm working on fixing the algorithm but in parallel, I'm looking for a way to check that my new polygon bordering my collection of chunks contain no diagonal lines.

Before applying the algorithm:

After: (in a case where it worked)


texts:
I have a raster representing rivers/drainage basins within a terrain:

I have a raster DEM of the terrain:

In the basin raster, everything in black is 9999.  White is nodata.

The DEM is obviously a range of values, from 1-500.

I want to combine these two rasters, in either of two ways:
The drainage raster has its values added to those of the DEM.  So if the course of the drainage basin in the DEM is originally a range of values between 2-300, the new range would be between 10199-10299.

Or:
The drainage raster has its values replace those of the DEM.  So if the course of the drainage basin in the DEM is originally a range of values between 2-300, it would now be a static 9999.

This seems like it should be straightforward, but I can't come up with a method in raster calc to give me the desired output.  The closest I've come is an output raster which adds the two values together, but all of the surrounding DEM data is missing.  i.e., the output looks identical to the basin raster, but with the combined value of DEM + basin raster.

I found the following solution from about a year ago:
Updating/replacing raster values in one raster values from another in QGIS?

However, I can't replicate this behaviour in QGIS 3.4.  When attempting to point the vectorization output at the existing DEM, it does not update, my ultimate output is the same DEM i began with.

My goal is to be able to run r.drain on the resulting map and get paths that do not follow the drainage basins, but instead find the next most efficient path through the terrain.

If someone knows of a better method to achieve this result, I'd welcome that as well.  My DEM is based on vector contours/breaklines/point data, and I considered trying to re-interpolate the DEM based on the original vector drainage map, but wasn't able to come up with a way to do that either.


texts: QGIS: moving one side of polygon parallel

texts:
I've got a project with different polygon shapes.
Now I´m trying to move only one side of a polygon parallelly and about an exactly defined distance. The result has to be as exactly as possible (a few centimeters).
Is there an add-on which offers me this function?


texts: Why do I get coordinates in QGIS3 as nnnnnn; nnnnnn instead of Lat/Lon when I upload a shapefile supposedly projected in WGS84

texts:
I received a couple of shapefiles made in France supposed to be projected in WGS84. When I upload them to QGIS3, it recognizes them as WGS84. However, their coordinates are displayed in a strange format with 6 digits. My shapefiles projected in WGS84 display well with the usual Lat/Lon degrees format.
Anybody has any ideas of what might be going on?
The CRS layer data in the shapefiles I received is exactly the same as that of my own:
WGS 84 WKT GEODCRS["WGS 84", DATUM["World Geodetic System 1984", ELLIPSOID["WGS 84",6378137,298.257223563, LENGTHUNIT["metre",1]]], PRIMEM["Greenwich",0, ANGLEUNIT["degree",0.0174532925199433]], CS[ellipsoidal,2], AXIS["geodetic latitude (Lat)",north, ORDER[1], ANGLEUNIT["degree",0.0174532925199433]], AXIS["geodetic longitude (Lon)",east, ORDER[2], ANGLEUNIT["degree",0.0174532925199433]], AREA["World"], BBOX[-90,-180,90,180], ID["EPSG",4326]] Proj4 +proj=longlat +datum=WGS84 +no_defs Extensão -180.00, -90.00, 180.00, 90.00
I am baffled and stuck...


texts:
I'm using a dictionary created by a search cursor from one feature class to update another feature class with an update cursor. I'm getting an invalid syntax error for the cursor.updateRow(row) and..have no idea why.

There must be something blatant I'm missing, yet I am simply not seeing it. Can someone point out why this code snippet would be providing an invalid syntax error for cursor.updateRow(row)? The format, from what I see, matches all the examples I've found.

    with arcpy.da.UpdateCursor(testingSHP, ['SERVICEMXL','STREETADDR']) as cursor:
        for row in cursor:
            mxLoc = row[0]
            oldAddr = row[1]
            if mxLoc in svcDict:
                print("Looking at row {0} with an address of {1}.".format(mxLoc, row[1]))
                row[1] = svcDict(mxLoc)
                print("Address is now {0}.".format(row[1])
                cursor.updateRow(row)
    del cursor


texts: Finding intersections for interstates using Overpass API?

texts:
I found this query, which works for streets:

way[highway][name="Mexico Road"];node(w)->.n1;
way[highway][name="Bryan Road"];node(w)->.n2;
node.n1.n2;
out meta;

I assumed it would work for interstates but when I plug in "I 55" and "I 270" for example, it gives me nothing. The problem appears to be related to exit ramps being their own structures in the data model.

Does anyone know an equivalent that works for interstates or just all roads?


texts:
I would like to add a raster wms layer to Mapbox GL but in similar way as OpenLayers single tiles example
Mapbox example https://docs.mapbox.com/mapbox-gl-js/example/wms/ use the 256px tile's size declared in code:
map.on('load', function () {
map.addSource('wms-test-source', {
'type': 'raster',
// use the tiles option to specify a WMS tile source URL
// https://docs.mapbox.com/mapbox-gl-js/style-spec/sources/
'tiles': [
'https://img.nj.gov/imagerywms/Natural2015?bbox={bbox-epsg-3857}&format=image/png&service=WMS&version=1.1.1&request=GetMap&srs=EPSG:3857&transparent=TRUE&width=256&height=256&layers=Natural2015'
],
'tileSize': 256
});
map.addLayer(
{
'id': 'wms-test-layer',
'type': 'raster',
'source': 'wms-test-source',
'paint': {}
},
'aeroway-line'
);
});

On the other hand, OpenLayers set the single tile size automatically depending on container's size (full screen in my case):
new OpenLayers.Layer.WMS(
    "Single Tile",
    "http://vmap0.tiles.osgeo.org/wms/vmap0",
    {layers: "basic"},
    {singleTile: true, ratio: 1}
)

Is possible to make something similar using Mapbox?


texts: geom_line() does not plot points from a shapefile in its correct order

texts:
I have been trying to load and plot a few shapefile in R using ggplot() and am facing an issue with a line vector shapefile. It seems that the readOGR() command is importing the line vectors ordering the points of the lines in the x-axis order so that the plot does not correspond to its correct shape. The "curvy" part of the lines are plotted wrongly, especially those in the right-most part of the plot.

The first figure attached illustrates how that shapefile is, plotted in QGIS. The second the output by loading the .shp with readOGR() and plotting it with ggplot(). I do so by running

library(rgdal)
library(ggplot2)
areas <- readOGR('.', 'areas')
ggplot(fortify(areas)) + geom_line(aes(long, lat, group = group))

Any thoughts on how to fix it? Here the source shapefile.


texts: QGIS server installation: Do I need to create a qgis user?

texts:
I am following the instructions in the official docs to install QGIS Server (3.8) on my Ubuntu server (18.04).

Throughout the docs, a /home/qgis/ folder is mentioned.

12.1.3. Follow Along: HTTP Server Configuration

(qgis.demo.conf)

  # default QGIS project   SetEnv QGIS_PROJECT_FILE

  /home/qgis/projects/world.qgs

  # QGIS_AUTH_DB_DIR_PATH must lead to a directory writeable by the
  Server's FCGI process user   FcgidInitialEnv QGIS_AUTH_DB_DIR_PATH
  "/home/qgis/qgisserverdb/"   FcgidInitialEnv QGIS_AUTH_PASSWORD_FILE
  "/home/qgis/qgisserverdb/qgis-auth.db"

  # See
  https://docs.qgis.org/testing/en/docs/user_manual/working_with_vector/supported_data.html#pg-service-file
  SetEnv PGSERVICEFILE /home/qgis/.pg_service.conf   FcgidInitialEnv
  PGPASSFILE "/home/qgis/.pgpass"

  # default QGIS project   SetEnv QGIS_PROJECT_FILE
  /home/qgis/projects/world.qgs

and then again (same paragraph)

  mkdir /home/qgis/qgisserverdb chown www-data:www-data
  /home/qgis/qgisserverdb

However, I don't have a qgis user in my system nor the installation has created one for me.

So the question is:

Do I need to create a qgis user or should I change all the references to /home/qgis to /home/myuser?

I cannot seem to find anything useful about this not even in the Advanced configuration page.


texts:
How can I save the results from a MapServer file as an image? I would like to be able to write a png file as a result of a query. My query which I can display in browser is the following:

http://mysite.com/mapserv.cgi?map=../../mapfiles/demo.map&layer=land&mode=map&FORMAT=image/png

I have tried the following (which does not seem to work):

wget http://mysite.com/mapserv.cgi?map=../../mapfiles/demo.map&layer=land&mode=map&FORMAT=image/png

and also

wget http://mysite.com/mapserv.cgi?map=../../mapfiles/demo.map&layer=land&mode=map&FORMAT=image/png >> myfile.png

I have seen the shp2img function but that appears to be used just for testing mapfiles. In this case I would like to use a SQL filter with my mapfile, which I don't think I can do using shp2img.

Included first wget attempt that I tried.


texts: In earlier version we could install the PostGIS using the Stack Builder but now it is not showing the option of spatial extension. Can anyone help me install it?


texts: QGIS Simple fill symbology for vector layer not exporting to pdf

texts:
I'm trying to use simple fill i.e. Diagonal X, on vector layers, the layers look fine on the composer when zoomed out (for example to a scale of 1:6,000) however when I export to pdf the distance between the lines decreases to the extent it is difficult to see the pattern. I have been using the line pattern fill to create my own Diagonal X symbology which with a large enough distance between the lines is fine.

There are a number of us using a wide range of symbols daily.

Is a quick way of fixing this for all symbols possible?

I've tried changing the units to pixels, mm or map units which I've seen mentioned a few times and this is not making a difference


texts:
Coming from an ESRI background I'm used to ESRI geodatabases (.mdb or .gdb).
In those databases, I could create feature classes, tables and relationships between them. Every time I loaded a feature class from one of my geodatabases into a new ArcMap project, relationships were automatically loaded as well: in example, if I loaded a point file which had two related tables, the INFO feature showed me the information coming from the point feature AND from the two related tables.

I'm trying to reproduce this in QGIS. I created a Spatialite database, with the same point feature, added the two tables and created the relationships between them. I then loaded my point file in QGIS, and with the INFO tool I can only retrieve the point info, NOT the related tables too. So I loaded the tables to my QGIS project and manually recreate the relationships.

So my question is: is this the only way to do this? Every time I start a new QGIS project do I have to load the point, the two tables and recreate the relationships?

EDIT: I'm talking about one-to-many relationships, every point has two or more related rows in the related tables


texts: Exporting Multiple Raster Data with custom Cell Size using ArcPy?

texts:
I Have 172 Raster, and need to change the cellsize. So I use Export Data from TOC

and change the cellsize manually (example I want to change it to "3")

I tried using this python script to create input and output folder (because batch output in resample had to be filled one by one and I just want to let all value as default including file name except for cellsize)

import sys, os, arcpy

InRaster = sys.argv[1]
OutRaster = sys.argv[2]

InFolder = os.path.dirname(InRaster)
arcpy.env.workspace = InFolder

for Ras in arcpy.ListRasters():
    arcpy.AddMessage("Processing " + Ras)
    #arcpy.RasterToOtherFormat_conversion(InFolder + "\\" + Ras , OutRaster ,"TIFF")
    arcpy.CopyRaster_management(InFolder + "\\" + Ras , OutRaster+ "\\" + Ras ,"DEFAULTS","","","","","")
    arcpy.AddMessage("Done")

Setting Properties :

When I try to add Input folder and Output Folder it has different submit

Result in error :

Is there any way to do it automatically or fixed for my python code to process multiple raster at once?


texts: Execution times when executing v.net.steiner is too high

texts:
I have a question related to execution times when running v.net.steiner.
When I try to connect 200 points by running v.net.steiner, it doesn't take more than 2 - 5 minutes. However, when I try to connect 2000 points, execution time rises up to 10-15 hours, is there any alternative way to decrease these times?
I am using Windows 10, I have 16 GB of memory, Intel Core i7 vPro. GRASS version is 7.8.2, and executing the code from the Python console of QGIS with the following code:
# Compute steiner(v.net.steiner)
alg_params = {
    '-g' : False,
    'GRASS_MIN_AREA_PARAMETER' : 0.0001,
    'GRASS_OUTPUT_TYPE_PARAMETER' : 0,
    'GRASS_REGION_PARAMETER' : None,
    'GRASS_SNAP_TOLERANCE_PARAMETER' : -1,
    'GRASS_VECTOR_DSCO' : '',
    'GRASS_VECTOR_EXPORT_NOCAT' : False,
    'GRASS_VECTOR_LCO' : '',
    'acolumn' : '',
    'arc_type' : [0,1],
    'input' : outputs['ExplodeLines']['OUTPUT'],
    'npoints' : -1,
    'output' : pSteiner_temp,
    'points' : lRAhouses,
    'terminal_cats' : '1-100000',
    'threshold' : 50
}

Find attached a picture of the Steiner I am trying to compute. In this particular case, I am trying to connect 4252 points, since this number is by far too high, I've divided the set of points in 4 sets of ~1000 points each one, but even though the running time is too high, is there any alternative way to shorter the execution time?


texts: How to add the Python Help and Syntax Box to ArcGIS Pro?

texts:
ArcGIS Desktop features a "Help and Syntax Window" in the Python prompt:

ArcGIS Pro has removed this feature, which was useful as it both listed and explained parameters:

Is there any way to add a similar window to ArcGIS Pro? At the moment the help() function is the only option I am aware of, which requires more typing and clutters the transcript.


texts: Vincenty's formulae to calculate in MySQL and PHP

texts:
I was reading Haversine formula on wikipedia and at the end of article its state that "More accurate methods that consider the Earth's ellipticity are given by Vincenty's formula and the other formulas in the geographical distance article.". Now It made me think that Vincenty's formula has better scope when calculating distance between 2 pair of latitude and longitude (Distance between two points ).

So i need to calculate minimum distance between 2 points in FEET

My questions are

1) Is Vincenty Formula the better option to calculate distance between 2 points than Haversine formula when generating latitude and longitude from Android Device

2) If Vincenty Formula is the better option then how to calculate distance between 2 points in FEET in MySQL

3) For example:Latitude 28.5312119 Longitude 77.2569046 and Latitude 28.5312138 Longitude 77.2569046 . these locations is generated from same location (with few feet difference).As you can see longitude is same but latitude is bit different in last decimals places only.

When i calculate these two location with Haversine formula i got 0 meters difference,When i know difference is few meters.So how can i achieve results in Feet with Vincenty Formula to get actual difference between these two locations. Please Also Read


texts: Transferring annotations from ArcMap to MAPublisher?

texts:
I'm collaborating with a designer who's using Avenza MAPublisher, while I'm using ArcMap Basic. We're still figuring out where the handover between GIS and Illustrator will be.

Can I create labels or annotations using ArcMap and Maplex, and pass these to the designer to refine in Illustrator? The designer will need access to the individual letters and words in the annotations, in order to adjust each letter's size, placement, etc.

I've tried exporting in *.AI format, and saving the labels as annotations in a file geodatabase. But the designer isn't able to import either of these into MAPublisher/Illustrator. In the case of the AI file, the layer has lost its georeferencing and needs to be placed manually, which isn't suitable for this high-precision work. The "add data" function in MAPublisher simply doesn't see the annotation featureclass in a file geodatabase.

Does the designer need to generate the labels himself using the Label Pro extension to MAPublisher?


texts: Create square buffers around points in QGIS with Python?

texts:
I have x, y coordinates in lat/long and I need to create 5x5 degree square cells around them, with the lat/long coordinates being the centroids.

My first option is to create a buffer around the centroids with 1 segment and a distance of 1/2(5°) * sqrt(2) (must multiply by sqrt(2)  bc the tool uses centroid to the corner of the square as the buffer distance, as opposed to centroid to edge), which results in sideways squares around my points, then rotate each feature by 45 degrees. I'd prefer not to do this as the distance isn't as precise and I do not know how to rotate each individual buffer feature quickly.

My second option, which seems much simpler, is to create a buffer around the centroids with the distance I need ((1/2)*5°) and then use something like ArcMap's Feature to Envelope tool.

I see that someone has the same question here and an answer was provided, but I have no idea how to do it programmatically.


texts:
Am new to gis and I hope to find help here. I am working on a research that requires that I overlay conditions of drainage systems and quality of healthcare in a coastal state to get areas which have potential for increase in the incidence of malaria. Each of these two attributes of areas in the state have more than one value. For instance, under quality of health care, I have data on:presence of healthcare service, drug provision and attitude of healthcare personnel.  Under conditions of drainage system, data on: presence of open drainage, effectiveness of waste management. The data is derived from answers from respondents in a household survey. The data is in numeric form

Am using arcgis 9.2

To my question, can I use frequency analysis in gis to summarize the data and go on to overlay the two attributes


texts: Changing number of classes in Symbology for ArcMap?

texts:
The goal is to style a polygon layer using graduated colors, something like

-9999 < value < 0 --> some color
0 =< value < 0.1 --> some other color
0.1 =< value 0.2 --> some other color

However, for some reason the number of classes is grayed out.

and when I click on classify:

classes is still grayed out.

Any thoughts?

ArcMap 10.5, the datalayer is using a join of a shapefile and csv table.

Number of classes is 7 and colors are custom (cannot be generated by colorbar).


texts: How to show hours and minutes in the time slider on CartoDB?

texts:
CartoDB automatically display only days and months on my custom map. However, I would like to show the following format (included in my data): dd-mm hh-mm.

Here is my code:

<script>
function main() {
  var map = new L.Map('map', {
    zoomControl: true,
    center: [0,0],
    zoom: 7
  });

  L.tileLayer('URL', {
  }).addTo(map);

      cartodb.createLayer(map, {
      type: "torque",
      order: 1,
      options: {
        table_name: "table_name",
        user_name: "",
       tile_style: $('#cartocss_template').html()
        }
    }).done(function(layer) {
      map.addLayer(layer);
    });
  }
  window.onload = main;
</script>


texts: Is there any swamps/marshes symbol using Quantum GIS?

texts: I am having hard time to find nature symbols for polygons. Any tips to find more symbol (any plugins?)


texts: Set the spacing above group in a Legend - PyQGIS

texts:
I am in the process of creating a script which needs to add a caption to an already existing layout.
So I created a PyQGIS code in order to do it :
legend = QgsLayoutItemLegend(layout)
layout.addLayoutItem(legend)

legend.setLinkedMap(map) # pass a QgsLayoutItemMap object

legend.setLegendFilterByMapEnabled(True)
#legend.adjustBoxSize()
legend.setResizeToContents(False)

legend.attemptMove(QgsLayoutPoint(3.5, 12, QgsUnitTypes.LayoutMillimeters))
legend.attemptResize(QgsLayoutSize(5.5, 2, QgsUnitTypes.LayoutCentimeters))

legend.setLineSpacing(15)
legend.setSymbolHeight(6)

#legend.setSpacingAboveGroups(15)

The legend is displayed well but I would like there to be space before the different groups of the legend.
I saw that this was possible in the layout editor (screenshot) but after studying the documentation for QgsLayoutItemLegend, I can't find a "setSpacingAboveGroups()" method that would allow this value to be modified. (setSpacingAboveGroups does not exist, this is an example)

How can I change the distance between groups in the legend using PyQGIS?


texts: Self terminating Python script in ArcMap?

texts:
I have tried sys.exit, but it doesn't seem to work in ArcMap 10.

Anyone know how to do it?


texts:
I'm trying to use Folium but am having difficulty adding multiple points to a MarkerCluster in one fell swoop.  Documentation is... sparse, so I'm not sure how to proceed.

This recreates the problem in an IPython notebook:

import folium print folium.__version__
# __version__ == 0.2.1

Joburg=[-26.2041, 28.0473]
locations=[[-26.15048, 28.10404],  [-28.67343, 32.1051],  [-25.97249, 28.12367]]

# Adding the markers one-by-one works fine:
map1 = folium.Map(location=Joburg,zoom_start=9)
my_marker_cluster1 = folium.MarkerCluster()
for location in locations:
    folium.Marker(location=location).add_to(my_marker_cluster1)
map1.add_children(my_marker_cluster1)

# Adding them together fails (i.e., the map renders, but with no markers on it)
map2 = folium.Map(location=Joburg,zoom_start=9)
my_marker_cluster2 = folium.MarkerCluster(locations)
map2.add_children(my_marker_cluster2)


texts: QGIS opening file - very slow on PostGIS database with many layers

texts:
When opening a QGIS file with many layers (20 or so) it takes +5 minutes to load. All the layers come from a single PostGIS database.

When looking into the server logs, I see that there are many requests like:

LOG:  SELECT typname,typtype,typelem,typlen FROM pg_type WHERE oid=1700
LOG:  SELECT pg_catalog.format_type(atttypid,atttypmod) FROM pg_attribute WHERE attrelid=34757 AND attnum=3
LOG:  SELECT description FROM pg_description WHERE objoid=34757 AND objsubid=3

Since there are thousands (+3000) of these queries, I suppose the problem lies there...

Any idea?


texts: Reclassify raster values of 9999 to NoData, leaving all other values alone, in ArcGIS Desktop?

texts:
I have a raster in ArcGIS which has a number of cells with the values 9999, showing that the data could not be generated for that cell. I want to convert all of these 9999 values to NoData, so that then I can do statistics on the dataset without getting crazy results.

How should I do this?

I have tried to use the reclassify tool, just adding one reclassification from 9999 -> No Data, but it seems to change all of the other values as well. Is there a way to do a reclassification which only does the changes you specify, and leaves all other values alone?

In case it matters, my dataset is a TIFF, with floating point values in it.


texts:
I can't filter out small areas (e.g. 200 sq.m) from the segmented NDVI map.

Is it possible to do that with the GEE means?
The script is here: https://code.earthengine.google.com/?scriptPath=users/sharovgennady/FieldData:Zones


texts: Determining which field is returned as NEAR_FID by ArcGIS Near tool?

texts:
I have a point layer and a polygon layer (both containing hundreds of records) and I want to know how far those points are from a nearest polygon feature,

I've used the Near tool and it returns two columns, NEAR_FID (identifies which polygon feature is the nearest) and NEAR_DIST (distance between a point and the nearest polygon feature). However, I would like to use a different field as the identifying field to show which polygon feature is the nearest (my own ID field instead of ArcGIS's automatic FID field).

Point layer:
FID ID
0   2
1   7
2   10
3   213
4   476

Does anyone know if this is possible, and if yes, how?


texts: Is it possible to get coordinate system name from SQL Server geometry type?

texts: I need to get the name of the coordinate system of a geometry type field in SQL Server 2008 R2 (for example NAD_1983_UTM_Zone_14N or GCS_WGS_1984 or other). Is this possible?


texts: Openlayers coordinate grid with text lables

texts:
Need to display the coordinates in the middle of the cell of grid. I do this grid, but can not to display coordinates. At zoom this coordinates must change and display on the grid (example write on leaflet).
My example in Openlayers
https://jsfiddle.net/andrii555/setmtbsq/
Grid is write on this example, but in leaflet
http://nakarte.tk/#m=7/48.84664/26.58142&l=O/Ng


texts: How can one solve the "Handle Unavailable Layers" error on a csv file in QGIS?

texts:
I've been working on this project for over a year and I have never had any substantial problem. Two days ago I was working on the project and everything worked fine. Today, when I try to load the project, I get the following:

The layer in question depicts points on my map. If I click on "Keep Unavailable Layers" and hover my mouse over the layer in question in the layers pane, I get the following message: Unavailable layer! Layer data source could not be found. Click to select a new data source.

I click and select the exact same data source, i.e. I "remind" QGIS where the .csv file is located, but then almost nothing seems to happen, apart from the layer's symbol changing from three dots to a cute cinderblock. The points on my map do not come back.

I am certain that neither the .qgs project file nor the .csv have been moved. Nor have I updated my operating system.
The problem surfaced with QGIS 3.10 running on MacOS 11.1 (Big Sur), but I have tried with QGIS 3.16 too and I get the same error. The source file is readable by QGIS, because I can create a new layer with the same csv file and it works.
Why is this happening? How can I get my layer back?
Edit: If I try to copy the style (Styles > Copy Style > All styles categories) of the cinderblock layer to a newly created layer with the same source file (Styles > Paste style > All styles), I get the following error:

But if I copy the style before the old layer becomes a cinderblock (i.e. before clicking "click here to select a new data source"), I do get back my old points...

...but I don't get back the callout lines (the lines between the labels and the points). Before the problem, it looked like this (from a PDF file):

If I try to zoom in on this region in QGIS...

...I get this:

Meanwhile the callouts on the new layer are defined like this (it used to say 0 mm of minimum length, I changed it to 3mm to no avail, nothing changes).

The content of the non-working Datasource (the one in the first screenshot) is file:../../../Articles/[Some folder name]/[some file name].csv?type=csv&delimiter=%5Ct;&detectTypes=yes&xField=longitude&yField=latitude&crs=EPSG:4326&spatialIndex=no&subsetIndex=no&watchFile=no The content of the Datasource in the working layer is: file:///[Some path]/[Some file name].csv?type=csv&maxFields=10000&detectTypes=yes&xField=longitude&yField=latitude&crs=EPSG:4326&spatialIndex=no&subsetIndex=no&watchFile=no
The delimiters appear different in these two versions (%5Ct; in the first, and not mentioned in the second so presumably commas). But when I try to paste the working Datasource path in the problematic layer, the problem is not resolved.


texts: Projecting polygon when calculating its area using R?

texts:
I have a broad range of polygons (in vector format and in SpatialPolygonsDataFrame class), some going east-west, others going south-north, while even others span polar regions or very large ranges as mid-to-high/low-to-mid latitudes/longitudes. In addition, as a single polygon can consist of several merged sub-polygons, they do not necessarily have to be aggregated or continuous, but can instead be chopped-up across large ranges. I have until now calculated the area without projecting the data using gArea {rgeos} as follows
gArea(SpatialPolygons(spdf@polygons)). Since I am not interested in visualization of any maps but only comparing polygon area and range sizes, I reckoned that the unit (square degrees vs. square meters) did not matter... but now I am not so sure.

So my questions are;

(i) Calculating polygon areas, would I need to project my data using an e.g. equal area projection?

If yes,

(ii) Since I do not care about any map visualizations and my polygons can span very broad ranges making it hard to find a single best projection, can I use the same projection for all polygons e.g. Lambert azimuthal equal-area projection when calculating the area?


texts:
I have labeled a map.

The problem is that the text is outside of the map in some places.

Is there a way to fix this?


texts: Understanding Slope Line and Flow Direction in ArcGIS Desktop?

texts:
I would like to create a Flow Direction (Polyline - not raster) which needs to follow slope direction (as arrow in picture I attached)

I'm working in mike urban (SWMM Model) and the software assign catchmente slope in according with flow direction. Specifically it need as inpunt a DEM Raster and Flow direction line. The polygon in this img rapresent the sub-catchment and in order to address the flow in each specific subcatchment-nodes I need to assign slope.

Can you suggest me how I can work in ArcMap to get them?


texts: Converting projected coordinates to google maps pixels

texts:
I have a single polygon in a local projection system.  I would like to convert that polygon into an image that can be used as an overlay on google maps.  I understand that it may be better to use existing projects such as Mapnik, but I am still interested in the process of turning coordinates into pixels.

My scripts current flow is as follows:

Convert polygons to Google Spherical Mercator Projection.
Divide the coordinates by the map resolution (metres per pixel) for the zoom level in question.
Transform the coordinates of the polygon to be relative to the origin of the image.
Render the coordinates to an image file.

Will such a process result in a cartographically correct rendering of the data?  Am I leaving anything out?


texts: Exporting training samples object statistics in eCognition

texts: I have training samples collected for my four sites. I need to run CART for feature selection off eCognition. The only way I know is to recollect samples using manual classification which is not desired. Is there any way which I can directly convert samples to object statistics with classes and features? or at least convert samples to classification as it is easier to export from classification to object statistics. Please help.


texts: Map orientation in print composer in QGIS

texts: I need to create a layout in which the vertices of the map in the layout coincide with the vertices of a rectangle of a grid that is oriented differently from the map on which it is overlayed: what can be the correct way to do it?


texts:
I am trying to get field values from one table and use them to perform calculations on a row in another fc table.  Does anyone have any ideas on how to make this run faster? I've tried a dictionary but am having problems when i try to reassign a value in the dictionary.

import arcpy
import numpy as np

arcpy.env.workspace = "C:\Users\michellep\Documents\project.gdb"

print "Arcpy & Numpy imported successfully and workspace set"
print "Start Processing"
start = time.time()
#define variables
group = "Test_Group"
fields = ['ObjectID', 'R', 'A', 'Total']
pc = "PC_Table"

cursor = arcpy.da.UpdateCursor(group, fields)
    #row = cursor.next()
for row in cursor:
    R = row[0]
    A = row[1]
    arcpy.CalculateField_management(pc, "R", R, "PYTHON", "")
    arcpy.CalculateField_management(pc, "A", A, "PYTHON", "")
    arcpy.CalculateField_management(pc, "WVal", "!R! / !A! * (!B! / !A!)**(!R!-1), "PYTHON", "")
    arcpy.CalculateField_management(pc, "Hours", "!WVal! * !K!", "PYTHON", "")
    arr = arcpy.da.TableToNumPyArray(pc, "Hours")
    arr.dtype = np.dtype(float)
    row[2] = np.sum(arr)
    cursor.updateRow(row)
    del row

del cursor
end = time.time()
processTime = end - start
print "Processing Completed - {:.1f} seconds".format(processTime)


texts: ArcGIS spatial query returned wrong count and ID list

texts:
I am currently encountering a very perplexing issue:

I have an ArcGIS feature layer of polygons
I am trying to query if any polygons in the feature contains a specific point
When I call the API to ask for the whole feature (returnGeometry: false, where: 1=1, geometryType: esriGeometryPoint, spatialRel: esriSpatialRelWithin) it works fine
When I call the API to ask for feature counts (returnCountOnly: true) or feature ID only (returnIdsOnly: true) it works for most cases; but there are some places where the query returns 0 counts, no feature IDs but 1 feature.
When I take note of the specific feature's Object ID and query it directly (objectIds: {objectId}), all 3 returned results are consistent

Is returned feature count and ID list unreliable like this? Do I really have to ask for whole features just to get the correct result? It seems to me it shouldn't behave this way, and I cannot find any mentions of others encountering similar issue.


texts: Set map service in idle mode in ArcGIS Server 10.1

texts:
I am using ESRI suit 10.1 (ArcGIS Desktop, ArcGIS Server, ArcGIS JSAPI 3.3/3.4) at my end and having few performance questions regarding ArcGIS Map service.

-We have three simple map services published on our server.

-Each map services contains query layer and with the help of query layer, we have created few xy events source layer in each map service.

-Integrated those map services in ArcGIS Java script application.

-While querying on the map we are facing load at our sql server database (I think the request coming from client side as we used query layer)

-Interesting thing is that in down time we got request only from arcgis server machine (On
DB server machine, checked the machine details and host name it’s our arcgis server machine)

So is it possible to set map service in idle mode when not in used? I checked the pooled and non-pooled option but in 10.1, all map services are pooled.

Going through System Design Strategies

Any help regarding this will be great !!

EDIT

According to @Alex Tereshenkov answer I can set the minimum number of instance to be zero so

-Is it take time to UP the instance?

-How many no. of instances to be set in Maximum number of instances (default is two)?


texts:
I do not understand the command line when I try to work with gdal_translate in QGIS. I found this at the forum: Normalize raster values to 0-255

And  tried this:

"If not, from within QGIS choose Raster -> Conversion - > Translate. Choose your files and options. At the bottom of the window the full gdal_translate command will be displayed. Click edit and add the -scale option."

When I type the following command into the line I get an error message "wrong parameters"??? Maby the commands are wrong.

My current grid file has a range from 1 to 9. I need to have a range from 1 to 5 (5 classes).


texts: warping and outputting coordinates of GeoTiff data: advice and error checking

texts:
So I have a massive TIF file (projected in Lambert azimuthal equal-area), and I'm trying to get lat/lon coordinates for each point of the grid.

My first step is to warp the projection to the WGS84 datum using:

>gdalwarp giganticfile.tif warpedcroppedfile.tif -t_srs "+proj=longlat +ellps=WGS84" -te -120 25 -60 55

Then I make that tif file into the following xyz format file:

>gdal_translate -of XYZ warpedcroppedfile.tif latlongvalues.dat

I can read a manual page as well as anyone else, but I'm new to this, so my two questions are:

does anyone have suggestions on how do do this better? I don't have ArcGIS access.
How do I error check myself.

I plotted the warped/cropped file, and it looked 'reasonable'(i.e., I could see a map of the US), but I'm doing some numerical analysis with the outputs, and I want to make sure that my data are good. So how do I error check myself? Manually go through the math for a few points?


texts:
I have been trying to use brew to install QGIS 3.0 dev on my Mac running 10.13 (High Sierra) without any success. I have reported the details on GitHub as Cmake errors installing on mac os 10.13.1 but have not had any response.  If anyone has managed this it would be great know if there were anything special that was needed.

I am using this command at present:

brew install --no-sandbox qgis/qgisdev/qgis3-dev

which is referenced at Install QGIS3 (QGIS-dev) on MacOSX as working recently.

It is frustrating that brew removes the build directory including the output from the qgis compile so one can't examine the entrails :(

Alternately is there any other way of building QGIS on a Mac?


texts:
Looking to  find the exact type of map that the following is :

It looks like a hill-shaded topographic map, most likely DEM, but I'm not sure. I need to know for certain so i can contact the map provider and order this exact style, for use in my projects.

EDIT: After some research, I am wondering if anyone knows were to purchase this exact type of map. I have looked into the opencyclemap, but that map is not this 3-d nature, and my client is very keen on getting this 3-D effect. Once again, the "must be available offline" restriction is the biggest problem for me here.

I have looked into a program called Google Maps Terrain Downloader, and I am wondering if anyone has had any experience with using it. The program downloads "small tiles", so I`m concerned that it may not be able to help me with rendering most of the country.


texts: Our company tracks multifamily development and is relying on opensource satellite imagery to track physical changes at development sites over time. Are there any open source providers (aside from Google Earth, Landsat - which is 30m) that we can leverage internally?


texts:
I have been trying to find out how to make the SAGA and GRASS kernel density tools work in QGIS 2.8. I have had a lot of trouble trying to get them to work as they produce useless results seen below.

GRASS' v.kernel gives me a solid black square with one upper value and 'nan' (the red box is the region). I have searched for reasons this could have occurred and found this. It says the likely cause of error is my region's number of columns and rows. In my case, I have tried numerous variations in the number of columns and rows and have not found any change to the results.

In SAGA's KDE, I get this weird 4 square result with only two values for the black and white squares. Each of these squares extend to probably 1000x the size of the point data and I'm unsure as to why. Also, is there a way to perform the analysis with no weights?

I need one of these tools to work because the heatmap plugin can not be used in processing modeler and  I need to automate maps for people who do not know how to use QGIS/GIS in general.


texts: Opengeosuite blank and postgresql errors

texts:
I had the unfortunate idea to upgrade opengeosuite from 4.0 to 4.1.
I followed the docs found on opengeosuite site but something has gone wrong.

I uninstalled opengeosuite and all connected packages and trying a brand new install but nothing seems to work anyway.

The dashboard page works but I can't access the Geoserver admin page and the Geoexplorer demo page. All I get is a blank page (400 Bad request).
I solved this problem by installing the single opengeosuite packages.
I don't understand, I followed the upgrade instructions but something didn't work as expected.

I had installed phppgadmin on the same server but now I can't login anymore. I get a login failed error. I solved the phppgadmin by running the following command:

setsebool -P httpd_can_network_connect_db 1

The command must be executed at each reboot of the server.

I'm not sure if the 2 problems are linked together. The postgresql problem is probably related to a failure in initdb when I reinstalled opengeosuite (yum install opengeo) since there were databases in the data folder.
But I have no clue about the geoserver/geoexplorer problem.

Any suggestion is welcome.


texts: Interpolation outside the contour of a set of points

texts:
I have a field contour in a shape layer. I have another layer compound by a set of points which are the locations where soils samples have been taken within the field. I have in the attribute table of this second shape layer the Phosphorus content of the soil in those points.

As I have said, the contour of the set of points is smaller than the field contour, being the contour of the set of points enclosed within the field contour.

I have interpolated this set of points by creating a TIN in order to have a value of Phosphorus content for any point in the field. But through this procedure, the area between the contour of the set of points and the contour of the field is excluded from the TIN.

I have also tried the Kriging method (this method allows to enter the field contour shape as the Processing Extent in the Environments section, being the referred problem solved) but through this interpolation method a new problem arises: the sample points are lying over classes which value ranges don’t include the real values of the sample points. Then, the original data is in some way adulterated.

I guess the IDW method isn’t good for this purpose: there isn’t any point in the field contour, and thus the value will go down towards the field contour in an unreal way due to this fact.

Any idea or solution for this issue (interpolating the area excluded from the TIN)?


texts: Dynamic IPs for polling and results requests using PyWPS?

texts:
I'm using PyWPS 3.2.4 that I containerized into Docker Containers to run into Mesos Marathon.

The external IP's and ports of my WPS Servers are therefore dynamic for being assigned by the Marathon framework.

I am able to send a GetCapabilities or Execute requests to these distributed WPS Servers in the Cluster since I know the 'public' IPs and ports. However the XML documents sent back to the client by the PyWPS servers contain 'polling' and 'results' IPs refering to the values set in 'pywps.cfg' (and not the IP and port used by the client to communicate with the server).

For example:

[wps]

...
serveraddress=http://172.16.10.40
...

[server]
...
outputUrl=http://172.16.10.40/wps/wpsoutputs/
...

Thus, the problem arise once it gets to poll for the end of the execution.
How can I configure PyWPS to force it to provides polling and results URLs on the same public URLs that the client request?


texts:
I am new to using ArcMap and I'm trying to get a simple greyscale/monochrome heightmap.

How would I go about doing this?

Right now, the only thing I really managed to do is to find the topo maps

I've read the instructions on how to remove the labels/overlay, but it mentions turning off the layer with the labels, which I don't actually get, when I save the map.


texts: From .shp to .kml to .shp Without Losing Attributes

texts:
I've converted a shapefile with attributes to a .kml for a client.  The .kml and attributes were then edited within google earth.  I received the updated .kml and now I need to get it back to .shp.  The problem is: all of my attributes are now in html format, which makes the popup info box in google earth.  When I try to convert back to .shp my attributes get lumped into a (google created) "Description" field.  And ALL of my attributes and information is lumped into a single cell.  Is there some way to extract the data from the "Description" field.  Even if I do that, is there an easy way to get rid of the html format and make it usable?

This is what I get using the "Fusion Table" method I found on some forums.  You can see the html in the left, it continues on further than it shows.

I am using ArcGIS Desktop 10.2.


texts:
from my Plugin I've loaded a rasterlayer using the following code:

fileName = "/path/to/raster/file.tif"
fileInfo = QFileInfo(fileName)
baseName = fileInfo.baseName()
rlayer = QgsRasterLayer(fileName, baseName)
if rlayer.isValid():
    QgsMapLayerRegistry.instance().addMapLayer(rlayer)
    ...

Now I wonder how I can move the freshly loaded layer to the back, that it lies behind all previously loaded layers.

Thanks and regards


texts: How to reproject from Michigan Georef to WGS84 in QGIS?

texts: I am trying to reproject a michigan state land shapefile which is the NAD83/ Michigan Oblique Mercator EPSG:3078 to the WPG84 EPSG:4326 projection. However when i try to do this my map gets squished down. I have tried the commands in "Reprojecting from MGRm to WGS84 using open-source tools" with no luck. I need to get these coordinates into a decimal degrees format. Thanks


texts:
I'm new to GIS. I have a shapefile (shp) that I need to add lat/lng coordinates too.

I used QGIS to open the file and selected:

Processing->ToolBox->QGIS geoalgorithms->Vector geometry tools->Extract nodes to add nodes to the file. I set the coordinate system to WG84G, when I go to add the $x coordinates I get output of 585292.76.

I believe it's the coordinate system not being correct but not sure which one I should be selecting.


texts: I have polygon with an attribute which has numbers populated in it. I would like to create random points for each of the polygon depending on the attribute. How could this be done ??


texts: Moving Esri grid format and info folder inside working directory using ArcPy?

texts:
Using the below code, I want to extract a new raster from an original raster using value field and I want to automate the process for 70 more rasters. However, for testing purpose, I'm using only 3 rasters.

Out of 3 rasters inside the working directory, this code works only for one raster which has both " ESRI grid" and "info" folder inside the working directory. I was not able to keep the "info" folders of other two raster inside the working directory as the "info" folder of 1st raster is already in the working directory.

How can I keep the "info" folders of other two rasters inside the working directory, so that code will work and I can automate the process?

Code:

import arcpy, os
from arcpy import env
from arcpy.sa import *

#To overwrite output
arcpy.env.overwriteOutput = True

#Set environment settings
env.workspace = "C:/Subhasis/Test/Neshanic_Python"

outws="C:/Subhasis/Test/Neshanic_Python/extract"

#checkout ArcGIS spatial analyst extension license
arcpy.CheckOutExtension("Spatial")

inraster = arcpy.ListRasters("*", "GRID")

for i in inraster:
    flds = ("VALUE", "COUNT")
    dct = {row[0]:row[1] for row in arcpy.da.SearchCursor(i, flds)}
    sumcnt = sum(dct.values())
    dct1 = {k:v for (k,v) in dct.items() if k >= 15}
    sumcnt1 = sum(dct1.values())
    percentage=(float(sumcnt1)/float(sumcnt))
    print percentage
    newraster = ExtractByAttributes(str(i), "VALUE>=15")
    outname=os.path.join(outws,str(i))
    newraster.save(outname)


texts: Is there a reliable source for Electronic Navigation Charts, where they get updated regularly? I am able to find for American waters on the NOAA website, but not for other regions.


texts: How to add a label based on multiple fields and brackets in GeoServer SLD

texts:
Is it possible to create a label which is comprised of multiple fields, in GeoServer using SLD?
For example, on the Points page of the SLD Cookbook, the point with default label section shows how to create a label based on the point attribute name:
<TextSymbolizer>
  <Label>
    <ogc:PropertyName>name</ogc:PropertyName>
  </Label>
  <Fill>
    <CssParameter name="fill">#000000</CssParameter>
  </Fill>
</TextSymbolizer>

How can I create a label based on two fields, which contains brackets around one of the fields? For example,

name (population)
Borfin (157860)

Testing shows that it's possible to add a new <ogc:PropertyName> element for the second attribute, and put the brackets on new lines:
<TextSymbolizer>
  <Label>
    <ogc:PropertyName>name</ogc:PropertyName>
    (
    <ogc:PropertyName>pop</ogc:PropertyName>
    )
  </Label>
  <Fill>
    <CssParameter name="fill">#000000</CssParameter>
  </Fill>
</TextSymbolizer>

but this results in a label containing a line break. How can I generate a label with all of the elements on a single line, as in Borfin (157860)?


texts: Creating shapefile without ArcGIS Desktop?

texts:
I have to create a shapefile but I do not have access to ArcGIS Desktop, only to ArcExplorer, which will not allow me to create a shapefile.

What other way can I do this?

I have GRASS but am not familiar with it.


texts: Mapping territory boundaries based on county lines and street?

texts:
I don't have any experience with GIS and mapping (other than drawing polygons in Google Earth). I need to create custom territory boundaries on a map, with county lines and street/highways as boundaries. These territories that I want to create will be next to each other, so I want to prevent empty spaces between them also.

Are there free tools that will allow me to do this? I attempted to do this using Google Earth's polygon tool, but tracing the county lines along with major highways over a large area is a massive task, and the result is not clean (i.e. my traces don't "snap" to the streets perfectly).


texts:
I'm using Qgis 2.18.17 with Postgres 9.6.

I have a feature line (below) and I want to create a feature AREA table VIEW with the following SQL using DB manager of Qgis.

select   distinct  row_number()  OVER () as id2 , (st_dump(st_polygonize(a.shape))).geom as shape
FROM   aaa_teste_qgis a.

The result table is:

When I try to create the VIEW and put on screen I receive the following msg.

  "Camada não é válida: A camada dbname='bdd' host=xx.xx.x.x port=5432
  user='TEST' password='xxxxxxxx' sslmode=allow key='id2' srid=4674
  type=Polygon table="TEST"."aaaaaaa" (shape) sql= não é uma camada
  válida e não pode ser adicionada ao mapa"

Then a note that ID2 column has the same number one.

Could you give a help?


texts:
I've got two layers with houses. One has housing groups (MultiPolygon), which has a mouseover (hover function) to show, which houeses belong together. The other layer contains single houses. On mouseover, a tooltip should show its address.

Tooltip works: https://jsfiddle.net/48n4sk4q/2/
Hover works: https://jsfiddle.net/48n4sk4q/3/
Combination doesn't (try fiddle's revision /1/ insead of /2/ or /3/).

It think, the problem is that both layers to be topmost?
Is there a work-around?


texts:
I have some problem with the representation of the coordinates.
I need to convert the coordinates from "EPSG: 32633" to "EPSG: 4326", and then create geometries in PostGIS.
For the conversion (from java) I use GeoTools.

double xc = 586073.93799999997;
double yc = 4644834.9100000001;
Coordinate c = new Coordinate(xc, yc);
Point p = gf.createPoint(c);
CoordinateReferenceSystem utmCrs = CRS.decode("EPSG:32633");
CoordinateReferenceSystem tCrs = CRS.decode("EPSG:4326");
MathTransform mathTransform = CRS.findMathTransform(utmCrs, tCrs, false);
Point p1 = (Point) JTS.transform(p, mathTransform);
System.out.println(p1.getCoordinate());

with the output:

  (40.271915834122524, 17.83405351767443, NaN)

However, this representation is wrong for PostGIS.

create table coordinate(
id int  primary key,
tipo varchar(15)
);
SELECT AddGeometryColumn('coordinate','luogo','4326','POINT',2);
INSERT INTO coordinate(id, tipo, luogo) VALUES(1,'calcare', ST_GeomFromText('POINT(40.271915834122524, 17.83405351767443)', 4326));

returns

  ERROR:  parse error - invalid geometry HINT:
  "POINT(40.271915834122524, " <-- parse error at position 26 within
  geometr

how can I solve my problem?

...
the wkt of GeoTools is the same as that contained in spatial_ref_sys:

crs = CRS.decode("EPSG:4326");
wkt = crs.toWKT();
System.out.println("\n\nwkt for EPSG:4326");

wkt for EPSG:4326
GEOGCS["WGS 84",
  DATUM["World Geodetic System 1984",
    SPHEROID["WGS 84", 6378137.0, 298.257223563, AUTHORITY["EPSG","7030"]],
    AUTHORITY["EPSG","6326"]],
  PRIMEM["Greenwich", 0.0, AUTHORITY["EPSG","8901"]],
  UNIT["degree", 0.017453292519943295],
  AXIS["Geodetic latitude", NORTH],
  AXIS["Geodetic longitude", EAST],
  AUTHORITY["EPSG","4326"]]


texts: GEE 'featureCollection.union()' keeps failing

texts:
I need to perform a union on two separate feature collections (one of species ranges, one of protected areas), which are unavoidably large and complex. My understanding is that this is done using featureCollection.union(). It repeatedly fails. See reproducible example for the Atlantic puffin:
https://code.earthengine.google.com/?scriptPath=users%2Frasenior%2Fpublic%3Apuffin_problem
// Puffin range
var puffinRange = ee.FeatureCollection('users/rasenior/temp/puffin');
Map.addLayer(puffinRange,{color:'blue'},'original');

// Are there different kinds of range polygon?
var rangeTypes = puffinRange.aggregate_array('seasonal').distinct();
print(rangeTypes); // Yes: breeding & non-breeding

// Union of range
var puffinUnion = ee.Feature(puffinRange.union().first());

// Range bounding box
var bbox = puffinUnion.bounds();

// Import PA polygons
var wdpaColl = ee.FeatureCollection('WCMC/WDPA/current/polygons')
  // Filter bounds to species range
  .filterBounds(bbox.geometry())
  // Filter status
  .filter(ee.Filter.inList('STATUS', ['Designated', 'Proposed']))
  // Exlude UNESCO Man and Biosphere reserves because these include
  // buffer and transition areas which do not meet the IUCN definition of protected areas.
  .filter(ee.Filter.neq('DESIG_ENG', 'UNESCO-MAB Biosphere Reserve'))
  // Exclude 100% marine reserves
  .filter(ee.Filter.neq('MARINE', '2'));

// Union of range
var wdpa = ee.Feature(wdpaColl.union().first());

// Range overlap with all PAs
var rangePA = puffinUnion.intersection(wdpa).area().floor();

// Export the results
Export.table.toDrive({
  collection: ee.FeatureCollection(ee.Feature(null, {'rangePA':rangePA})),
  description: 'puffinPA',
  fileFormat: 'CSV'
});

The union of the species' range cannot be plotted, and exporting the task results in Execution failed; out of memory.. Is there anything I can do to make this run? I've tried simplifying the features by 1 km (with cylindrical equal area projection), but so far that has not worked either. Perhaps I have to work with image collections instead?
EDIT: The ultimate goal is to calculate the area of overlap between species ranges and Protected Areas. To do this accurately, both feature collections need to be flattened, plus intersection() only works on a feature.


texts: How to replace VectorWriter with QgsVectorFileWriter?

texts:
Within QGIS 2.8.1 I found within the Processing Toolbox a script called "Points on crossing lines" which should exactly as it is named. Unfortunately it needs processing.core.VectorWriter which I guess has disappeared in this latest installment of QGIS. I tried to replace the writer with QgsVectorWriter but get the following error message:

"Algorithm [Unnamed algorithm] starting...
Index lines...
Find crossing points...
arguments did not match any overloaded call:
QgsVectorFileWriter(QString,QString, QgsFields, QGis.WkbType,
QgsCoordinateReferenceSystem, QString driverName="ESRI Shapefile",
QStringList datasourceOptions=QStringList(), QStringList
layerOptions=QStringList(), QString newFilename=None,
QgsVectorFileWriter.SymbologyExport
symbologyExport=QgsVectorFileWriter.NoSymbology): argument 3 has unexpected
type 'list' QgsVectorFileWriter(QgsVectorFileWriter): argument 1 has
unexpected type 'unicode' See log for more details"

Can you point me at what I am doing wrong?

##Lines=vector
##Results=output vector

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

cutLayer = processing.getObject(Lines)
cutPrder = cutLayer.dataProvider()
n = cutLayer.featureCount()
l = 0

# build spatial index of lines

index = QgsSpatialIndex()
geom_ix = {}
progress.setText("Index lines...")

for feat in processing.features(cutLayer):
progress.setPercentage(int(100*l/n))
l+=1

index.insertFeature(feat)
geom_ix[feat.id()] = feat.geometry().asWkb()

# find points on crossing lines
progress.setText("Find crossing points...")

l = 0
i = 0
ptindex = QgsSpatialIndex()
pt_ix = {}
secgeom = QgsGeometry()
featgeom = QgsGeometry()
resfeat = QgsFeature()

for feat in processing.features(cutLayer):
progress.setPercentage(int(100*l/n))
l+=1

near = index.intersects(feat.geometry().boundingBox())

for f in [x for x in near if x != feat.id()]:   # exclude self

    featgeom = feat.geometry()
    secgeom.fromWkb(geom_ix[f])

    if featgeom.crosses(secgeom):
        crosspts = feat.geometry().intersection(secgeom).asGeometryCollection()

        for pt in crosspts:
            i += 1

            # index point
            resfeat.setGeometry(pt)
            resfeat.setFeatureId(i)
            ptindex.insertFeature(resfeat)
            pt_ix[i] = pt.asPoint()

feat = QgsFeature()
fields = [QgsField("nodeid", QVariant.Int)]
writer = QgsVectorFileWriter(Results,None, fields, QGis.WKBPoint, 27700)

progress.setText("test")

# only save unique points
progress.setText("Save unique points...")
n = len(pt_ix)
featgeom = QgsGeometry()

while len(pt_ix) != 0:
progress.setPercentage(int(100*(n-len(pt_ix))/n))

i = pt_ix.keys()[0]

# write point

attrs = [i]
feat.setGeometry(featgeom.fromPoint(pt_ix[i]))
feat.setAttributes(attrs)
writer.addFeature(feat)

# delete close points
near = ptindex.intersects(buffRect(pt_ix[i], buff))

for pt in near:
    feat.setFeatureId(pt)
    feat.setGeometry(featgeom.fromPoint(pt_ix[pt]))
    deleted = ptindex.deleteFeature(feat)
    del pt_ix[pt]

del writer


texts:
It could be user error, but for some reason, the reclass tool is just hanging and hanging...

QGIS 1.8
Compiled against GDAL/OGR 1.9.2
Ubuntu 12.04.1

Running the tool from the Sextante toolbox
Input is a 10' cellsize raster ~10mi x 10mi in area
grass region cellsize = 10
Reclass file is in order from what I can gather from the help

When run:
- There is no process operating in the system monitor
- Dialogue hangs and cursor shows progress spinny thing...

Not sure what else to say?

Does the GRASS tool run from this dialogue on a .tif? Or does it have to be run from the GRASS tools specifically?

Thanks!


texts: How do I query OpenStreetMap for a list of US counties?

texts:
I would like to query OpenStreetMap (remotely) for a list of US counties (and their boundary vectors if possible), in a scriptable automated fashion.

Preferred format is CSV.

How do I do so? I am new to querying OpenStreetMap


texts: How to perform Visibility Analysis using QGIS?

texts:
How can I calculate Visibility using QGIS?

I have not found any suitable plugins for QGIS to do this.


texts:
I am displaying and labelling a polygon Layer using a GeoJSON file and the field nom_comm. I would like to label the polygons between two scales i.e. two zooms.
(So to speak, labels are to be invisible when scale of the map is too high or too little.)

How can I do that?

Sample of HTML:

var styleDS = new OpenLayers.Style({
                        'fillOpacity': 0,
                        'strokeColor': "red",
                        'strokeWidth': 0.8,
                        'fontFamily': "Verdana",
                        'fontWeight': "bold",
                        'fontSize': "6.5px",
                        'labelOutlineColor': "white",
                        'labelOutlineWidth': 1,
                        'label': "${getLabel}"/* options */
                                        },
                        {
                        context: { getLabel:function(feature) //affichage nom commune en fonction echelle
                                    { if (map.getZoom() > 9 && map.getZoom() < 12.5)
                                    { console.log('affiché'); var label = feature.attributes['nom_comm'];
                                    return label; }
                                    else { console.log('caché'); return ''; }
                                    }
                                }
                            });

Then I applied the style like this:

var commune = new OpenLayers.Layer.Vector(
        "commune", {
        strategies: [new OpenLayers.Strategy.Fixed()],
        protocol: new OpenLayers.Protocol.HTTP({
            url: "communes.geojson",
            format: new OpenLayers.Format.GeoJSON()
        }),
        //autre maniere de creer un style
        styleMap: styleDS,
        }
);
  map.addLayer(commune);

And I don't have a link to the map because it belongs to a domain.


texts: Filtering an ImageCollection returns no elements (Earth Engine)

texts:
I'm trying to replicate this python script to calculate NDVI mean for a feature collection (https://stackoverflow.com/questions/47633088/get-results-in-an-earth-engine-python-script). It seems like the code calculates a result, but a cant see result.
This is the code:

import datetime
import ee
ee.Initialize()

#Feature collection
fc = ee.FeatureCollection("ft:1x290vohty0Wgdn5jL3RlpzryK7dfOPtG6yY213e0");
fc_filtered = fc.filter(ee.Filter.equals('NUM_DECS', 1))
#Image collection
Sentinel_collection = (ee.ImageCollection('COPERNICUS/S2')
    .filterBounds(fc_filtered)
    .filterDate(ee.Date('2017-01-01'),ee.Date('2017-08-01')))

def GetSeries(feature):
  def NDVIcalc(img):
    red = img.select('B4')
    nir = img.select('B8')
    ndvi = nir.subtract(red).divide(nir.add(red)).rename(['NDVI'])
    return (feature
            .set(ndvi.reduceRegion(ee.Reducer.mean(), feature.geometry(), 10))
            .set('date', img.date().format("YYYYMMdd")))

  series = Sentinel_collection.map(NDVIcalc)
  list = series.reduceColumns(ee.Reducer.toList(2), ['date', 'NDVI']).get('list')
  return feature.set(ee.Dictionary(ee.List(list).flatten()))

result = fc_filtered.map(GetSeries)
print(result.getInfo())

Give me this result:

{'columns': {}, 'type': 'FeatureCollection', 'features': [], 'properties': {'name': 'kmltest', 'DocID': '1x290vohty0Wgdn5jL3RlpzryK7dfOPtG6yY213e0'}}

So it does not looks like something is coming out of it?

I havent yet tried the next where the following is included:

# Get all possible dates.
dates = ee.List(Sentinel_collection.map(function(img) {
      return ee.Feature(null, {'date': img.date().format("YYYYMMdd") })
}).aggregate_array('date'))

# Make a default value for every date.
header = ee.Feature(null, ee.Dictionary(dates, ee.List.repeat(-1, dates.size())))
output = header.merge(result)
ee.batch.Export.table.toDrive(...)

Any suggestions for what a do wrong? Will prefer to have the result in the feature collection or exported as a list.


texts: How to export a training shapefile as an XML spectral signature file via QGIS?

texts:
I recently used the answer from Using Orfeo Toolbox object based classification in QGIS tutorial/workflow? to perform an object base segmentation on a SPOT image I acquired a while ago. I tried both of the solutions given, they both work mostly but I have an issue with the second half of the first one.

What it says, in short, is that after using the "Exact large scale meanshift segmentation" tool from the Orfeo Toolbox via Qgis, you need to create and edit training areas manually using a shapefile. The thing being that the said shape is needed afterwards as a signature XML file containing mean values and variance for each of the training areas. I did think about using the semi automatic classification tool, creating the ROIs and converting them into an XML signature file directly. The issue being here that any manual edit on the ROI shapefile makes QGIS insta-crash 8 times out of 10, and when I'm actually able to do it, it crashes again when I try to export the shapefile as a spectral signature via the inbuilt conversion tool.

Is there any work-around that might help creating the XML file from a shapefile?
Is semi automatic classification the best way to create it?
If not, what should I use?


texts: width is bounded when applying SLD in Geoserver?

texts:
I applied a SLD style to color a road network in Geoserver (see the picture, before the styling, all lines are in the same color). Before the styling, I can query the WMS layer with arbitrary width and height values. But after the styling, I cannot get the right picture if I increase width to certain number.

For example, consider the following query to WMS

http://"hostname"/geoserver/Workspace/wms?service=WMS&version=1.1.0&request=GetMap&layers=Lines&styles=&bbox=-91.44184599959314,36.99574399968581,-87.52404399962063,42.508301999999986&width=410&height=512&srs=EPSG:4269&format=image%2Fpng

Before applying the color style, I can set the width value to any number;

After applying the color style, I can only set width no larger than 410. I do not understand why the styling matters.

Any ideas?


texts: Converting JSON data to GeoJSON data from an Arduino GPS tracker

texts:
I am working on a school project that tracks the location of someone via an Arduino GPS tracker. The issue I am running into is trying to convert JSON data into GeoJSON data. When the console is opened the error presented is invalid GeoJSON object. I have seen posts here describing the process to convert the data, but due to my rudimentary understanding of coding I am not sure how to implement this conversion into code. The code is based on this: https://github.com/perliedman/leaflet-realtime

index.js

var map = L.map('map');
    realtime = L.realtime({
        url: 'https://data.sparkfun.com/output/3J0bK7jKYvs6OL6V5zZE.json',
        jsonp: 'callback',
        cache: true,
        crossOrigin: true,
        type: 'GET',
    },
    {
        interval: 10 * 1000
    }).addTo(map);

L.tileLayer('http://{s}.tile.osm.org/{z}/{x}/{y}.png', {
    attribution: '&copy; OpenStreetMap contributors'
}).addTo(map);

realtime.on('update', function() {
    map.fitBounds(realtime.getBounds(), {maxZoom: 3});
});

The url https://data.sparkfun.com/output/3J0bK7jKYvs6OL6V5zZE.json describes the provided JSON data from my Arduino, but in order for the code to work it seems the data needs to be in GeoJSON format.


texts: Labeling nodes with letters but removing last letter using QGIS

texts:
I want to label nodes of my polygon with upper letters. I got to the point where I am able to label them with the right style. What I still need to do is to dismiss the final letter because it overlays the first letter, see the following screenshot where you can see that Q (the last label) overlays A (the first label):

Here is my expression :
CASE WHEN @geometry_part_num <27 THEN upper(char( @geometry_part_num-1+97)) WHEN @geometry_part_num <677 THEN concat(upper(char( to_int( (@geometry_part_num-1)/26 + 96))) , upper(char( to_int((@geometry_part_num-1)%26 + 96)))) END


texts: GeoFence is securing GeoServer quite nicely, but I can't find how to secure the GeoFence REST API. I've been able to access the REST API without any authentication.


texts: Creating probability of occurrence map in ArcMap Pro

texts:
I am trying to figure out how to create a map representing probability of some phenomena occurrence. I have a binary raster (1- phenomena exists; 0 - doesn't exist) and several raster that explain this phenomena (temperature, precipitation, slope and etc.).
Could you advise on steps I need to take to have this map?
I tried to use different regression tools in ArcMap Pro, but seems that they don't suit my situation?


texts: GeoServer Layer from PostGIS not updating

texts:
I created a WMS layer from a PostGIS table with a fraction of total polygons in GeoServer. Everything displayed fine.

Then, I added the rest of the polygons to the table in postgis. For some reason its not updating the layer. I tried restarting geoserver and its still not updating. How do I get geoserver to update the layer?

I use Open layers WMS layer and OpenLayers.Protocol.WFS.fromWMSLayer(polygons) to get polygon data.

-N

UPDATE: I tried to click Compute from data/Compute from bounding box links to get the new MBR. Unfortunately, it did not update, so here's as much info as I can think to provide.

I have a set of polygons scattered around the US. They are from a pair of SHP files in WGS_1984. I started testing the system with a small set of 100 polygons and now that the system works, I wanted to add the the other polygons to the data. I used shp2sql to create a new table in postgis. Then I used a INSERT INTO polygons (SELECT from ...) to add them to the table that geoserver uses. Now there is about 50k polygons in the table.

I tried updating the MBR as answers below suggest. Here are the current settings of the layer in geoserver:

Native SRS: EPSG:4326
Declared SRS: EPSG:900913

When I hit Compute from data for Native bounding box it gives me:

(Min x,Min y, Max x, Max y)

-10Mil, 5Mil, -10mil, 5Mil

and for the Declared:

-95.xxx,40.xxx,-95.xxx,43.xxx

For whatever reason, they are now updating correctly. The MBR should cover most of the US now and not just that little portion in the middle.


texts: I've been tasked with generating some samples from a larger raster file to use as inputs to a computer model. Since we already have rasters prepared, the hope was that this would just be a matter of opening the file and selecting x by y pixels from the raster using something like the "Select Features" tool. However, none of the selection features are enabled in the menu (see below). How can I select an x by y portion of the raster (ideally in the UI) so that it can be exported to a new raster file?


texts: OSMNx Extracting blocks from point and save coordinates of corners of each block

texts:
I am playing with OSMNx and I would like to extract street blocks,

I have this code which should display the blocks

import osmnx as ox
import geopandas as gpd
import shapely
from shapely.geometry import Polygon
import matplotlib.pyplot as plt
import numpy as np

point = (50.090988, 14.425114)

streets_graph = ox.graph_from_point(point, distance=500, network_type='drive')
streets_graph = ox.project_graph(streets_graph)

streets = ox.save_load.graph_to_gdfs(streets_graph, nodes=False, edges=True,
                                     node_geometry=False, fill_edge_geometry=True)

point = streets.unary_union.centroid

polygons = shapely.ops.polygonize(intersections.geometry)
polygons = gpd.GeoSeries(polygons)
target = polygons.loc[polygons.contains(point)]

This is how the area near the point looks like on Google maps:
This is what I get:

Then I tried to isolate the intersections which I would then use to create the polygons I need

intersections = ox.clean_intersections(streets_graph, tolerance=15, dead_ends=False)

gdf = gpd.GeoDataFrame(geometry=intersections)
gdf.crs = streets_graph.graph['crs']
gdf_projected = ox.project_gdf(gdf, to_latlong=True)
polygons.plot()
plt.show()

And I got something more realistic:

And my question is:

How can I connect those intersections in order to get these polygons as shown on the first image from google maps AND (this is the most important part for me: How to get coordinates of each polygon edge?

Or at least, how can I fix the drawings that the 2nd picture would resemble first picture?


texts: Are there ArcIMS instances still around or phased out?

texts: My understanding is there are few companies still using ArcIMS and I am wondering if it is still alive or phased out?


texts: Cropping GeoTIFF with GeoTools using integer indexes

texts:
I know how to crop a raster by providing coordinates in a CRS:

Envelope envelope2D = new ReferencedEnvelope(envelopeStartX, envelopeEndX, envelopeStartY, envelopeEndY, obj.targetCRS);

GridCoverage2D tile = (GridCoverage2D) ops.crop(obj.coverage, envelope2D);

How do I crop by providing raster indices (ranging from 0..height,
0..width) instead of coordinates in a CRS?


texts:
I am trying to carry out a FLAASH atmospheric correction on radiometricaly corrected GeoEye-1 spatial image (working within ENVI).
I am getting the wrong results from my FLAASH correction and believe that this is being caused by the Scale factor I have been using or the radiometrically corrected image being in the wrong units for the FLAASH correction to work.

I am using the Scale factor of 100 as the image is a 16 byte pixel value (also used 1 and 10 as recommend in ENVI help but had no luck).

As far as I can read the units in the GeoEye image are in mW/(cm2-srt-nm-dn) at the radometrically corrected level. However to carry out the FLAASH correction it should be in µW/(cm2 * sr * µm) and I do not know how to convert the one to the other.

Could you give me some advice about either of these options (converting the units or using a separate scale factor for each band)?


texts: Projection Transformation from ETRS89 LAEA to WGS1984?

texts:
I am having difficulties to transform a Shapefile in ArcGIS 10 from ETRS 1989 LAEA to WGS 1984 N33.

I have 2 files:

in WGS1984 N33, with a grid of 10x10km,
in ETRS 1989 LAEA, with small 100x100m polygones looking somewhat gridlike.

If I transform the ETRS-Shapefile (on the fly or using the transformation tool in Data Management), it remains skewed and does not align with the WGS1984 file (see picture with green the WGS file and lilac the transformed ETRS-file, where the horizontal lines are not identical)

Properties of the ETRS1989 before transformation:

Projected Coordinate System:    ETRS_1989_LAEA
Projection: Lambert_Azimuthal_Equal_Area
False_Easting:  4321000.00000000
False_Northing: 3210000.00000000
Central_Meridian:   10.00000000
Latitude_Of_Origin: 52.00000000
Linear Unit:    Meter

Geographic Coordinate System:   GCS_ETRS_1989
Datum:  D_ETRS_1989
Prime Meridian:     Greenwich
Angular Unit:   Degree

Properties of the WGS-File

Projected Coordinate System:    WGS_1984_UTM_Zone_33N
Projection: Transverse_Mercator
False_Easting:  500000.00000000
False_Northing: 0.00000000
Central_Meridian:   15.00000000
Scale_Factor:   0.99960000
Latitude_Of_Origin: 0.00000000
Linear Unit:    Meter

Geographic Coordinate System:   GCS_WGS_1984
Datum:  D_WGS_1984
Prime Meridian:     Greenwich
Angular Unit:   Degree


texts:
I have two data sets, one is a polygon data set (Norwegian counties) and the other is points (roadkill incidence data). I want to do a spatial join where the points dataset inherits attributes (name/number etc) from the county they are within.

I am setting up a batch job which will be supplied with regularly updated data. The tool is FME 2015.0.

I have tried using the SpatialFilter transformer to do this using my point data as candidate data and the county data set as filter data. I have tried Contains and within filters with no luck. I only the county polygons from the passed output. Swithching around does not help much either. Still polygons.

Is there an other way to do this using FME?


texts:
I'm trying to create a multiband raster from several component rasters, which have the same alignment. The first step is to identify which tiles completely overlap (i.e., have the same extent and upper left coordinates).
From the docs, it would seem that r1.rast=r2.rast would be a quick way to identify these tiles, i.e. are the bounding boxes the same.
But the = operator is identifying tiles in completely different parts of the world. Moreover, the ~= operator (which should be the same as = but with lower precision) is returning False.
SELECT r1.rast~=r2.rast, r1.rast=r2.rast
FROM raster1 AS r1, raster2 AS r2
WHERE g1.rid=66585 AND g2.rid=170411

returns False, True for the row in question. How is that possible?


texts: Creating and serving pdf files from ArcGIS server 10.0

texts:
I am tasked with sending map data from the Flex client to an MXD then to serve that as a PDF to the client and/or a server side location.  Creating the PDF isn't an issue, however, I am having trouble finding a resource to let you now how to serve it to the client.  Does anyone have any experience of this?

Thanks in advance.


texts: Mapshaper fails when converting files to .topojson due to file size?

texts:
I have a .gejson file that's about 378 meg (378,776k)
and I'm trying to convert it to topojson.

I was successful with files that sized to about 262 meg
but files greater than that failed leaving the error:

Error: "toString()" failed
at Buffer.toString (buff.js:480:11)
etc....

I'm using the command line format:
ex)

mapshaper - i bigfile.geojson -o little.json format=topojson

and I'm starting to think that this is a buffering issue for memory that needs to be changed/tweeked.

What exactly do I  need  for files >~262 meg?


texts: ArcGIS InfoTemplate setTitle setContent is not a function

texts:
I'm using ArcGIS JS API 3.24 and I got an error saying that

Uncaught TypeError: cellInfoTemplate.setTitle is not a function

But I'm sure that I have already included esri/InfoTemplate in my require.
Here is more:

 require([
  "esri/map", "esri/layers/FeatureLayer", "esri/InfoTemplate",
  "esri/tasks/GeometryService", "esri/layers/ArcGISDynamicMapServiceLayer",
  "esri/layers/TableDataSource", "esri/layers/QueryDataSource", "esri/layers/LayerDataSource", "esri/layers/DynamicLayerInfo",
  "esri/symbols/SimpleFillSymbol", "esri/renderers/SimpleRenderer", "esri/layers/LayerDrawingOptions",
  "esri/symbols/PictureMarkerSymbol", "esri/symbols/SimpleLineSymbol",
  "esri/symbols/SimpleMarkerSymbol", "esri/renderers/UniqueValueRenderer",
  "esri/config",
  "esri/Color", "dojo/domReady!"
], function (
  Map, Layers, FeatureLayer, InfoTemplate,
  GeometryService, ArcGISDynamicMapServiceLayer,
  TableDataSource, QueryDataSource, LayerDataSource, DynamicLayerInfo,
  SimpleFillSymbol, SimpleRenderer, LayerDrawingOptions,
  PictureMarkerSymbol, SimpleLineSymbol, SimpleMarkerSymbol, UniqueValueRenderer,
  esriConfig,
  Color
) {
    var cellInfoTemplate = new InfoTemplate();
    cellInfoTemplate.setTitle("Details");
    cellInfoTemplate.setContent(getCellHTML());

    var featureLayer = new esri.layers.FeatureLayer("myUrl", {
        "opacity": 0.6,
    });

    featureLayer.setInfoTemplates(infoTempDict);

    function getCellHTML(){ ... }
    ... other codes ...
 }

And if I avoid using .setTitle and .setContent with

var cellInfoTemplate = new InfoTemplate("Details",getCellHTML());

It will return me

Uncaught TypeError: featureLayer.setInfoTemplates is not a function

So can anyone tell me how to solve the error?


texts:
When i select polygon with select tool- the whole polygon colored:

i seek a way to select the shape only with outline, as it possible in ArcGIS. I know i can change the opasitiy in "Map Tool" but this isn't a solution for me.

i red QGIS 2.2 Identify tool highlight settings but didn't found an answer.

Is it possible with QGIS?


texts: Reclassifying NDWI using Google Earth Engine?

texts:
I'm trying to reclassify a NDWI raster into 1 for values >0.4 and 0 for <0.4. I understand that its possible to use .remap in case of discrete values, but not sure how to reclassify NDWI raster using a threshold. Ideally, I need a function which will reclassify each image in a collection.

Attached is the earth engine code I used (and struck)

// Import MODIS 8 day reflectance collection.
var collection = ee.ImageCollection('MODIS/006/MYD09A1')
                   .filterBounds(roi);

//roi refers to an area of interest
//start data and end date
var iniDate = ee.Date.fromYMD(2016,1,1);
var endDate = ee.Date.fromYMD(2017,12,31);
print('MODIS', collection);

// MODIS bands
var modisBands =
['sur_refl_b03','sur_refl_b04','sur_refl_b01','sur_refl_b02',
'sur_refl_b06', 'sur_refl_b07'];
var lsBands = ['blue','green','red','nir','swir1','swir2'];

// helper function to extract the QA bits
function getQABits(image, start, end, newName) {
// Compute the bits we need to extract.
var pattern = 0;
for (var i = start; i <= end; i++) {
pattern += Math.pow(2, i);
}

// Return a single band image of the extracted QA bits, giving the band
// a new name.
return image.select([0], [newName])
            .bitwiseAnd(pattern)
            .rightShift(start);
}

// A function to mask out cloudy pixels.
function maskQuality(image) {
// Select the QA band.
var QA = image.select('StateQA');
// Get the internal_cloud_algorithm_flag bit.
var internalQuality = getQABits(QA,8, 13, 'internal_quality_flag');
// Return an image masking out cloudy areas.
return image.updateMask(internalQuality.eq(0));
}

// create cloud free composite
var noCloud = collection.filterDate(iniDate,endDate)
                        .map(maskQuality)
                        .select(modisBands,lsBands);

// Function to calculate and add an NDWI band
var addNDWI = function(image) {
return image.addBands(image.normalizedDifference(['green', 'swir1']));
};

// Add NDWI band to image collection
var noCloud = noCloud.map(addNDWI);

// create NDWI image collection
var NDWI = noCloud.select(['nd']);

// Reclassify each image in NDWI collection and calculate yearly composite


texts: Resizing shapefile features proportionally in QGIS

texts:
In QGIS (1.8.0-Lisboa) on Windows (OSGeo4W Install), how would you resize a feature proportionally?

The only way I've found to do this is by hand using either the Freehand Editing (0.2.6) plugin or the Reshape Features button. Neither of these provides a means to uniformly scale (as in resize) the selected feature up or down.

For example, I want to increase the features size by 130% so that you retain the shape of the feature but expand the area of the shape itself proportionally in all directions.


texts: what is the fastest way to plot a list of Canadian Addresses?

texts: I have a list of roughly 500 addresses that I'd like to plot on a map either as a point layer or as a kml/placemarks. The addresses are in form 123 QGis Street, Province, Canada, A1A 1A1


texts: Error saving polygon with hibernate spatial and postgis

texts:
I'm trying to persist a com.vividsolutions.jts.geom.Polygon in my PostgreSQL (with PostGIS enabled) database, using Hibernate and Hibernate-Spatial plugin.

Caused by: org.postgresql.util.PSQLException: One can not infer the SQL type to use for an instance of org.postgis.PGgeometry

It fired the below exception. What am I doing wrong?

My pom.xml

<dependency>
            <groupId>org.hibernate</groupId>
            <artifactId>hibernate-spatial</artifactId>
            <version>${version.hibernate.spatial}</version> //version 4.3
            <exclusions>
                <exclusion>
                    <groupId>postgresql</groupId>
                    <artifactId>postgresql</artifactId>
                </exclusion>
                <exclusion>
                    <groupId>org.postgis</groupId>
                    <artifactId>postgis-jdbc</artifactId>
                </exclusion>
                <exclusion>
                    <groupId>org.hibernate</groupId>
                    <artifactId>hibernate-core</artifactId>
                </exclusion>
            </exclusions>
        </dependency>

        <dependency>
            <groupId>com.vividsolutions</groupId>
            <artifactId>jts</artifactId>
            <version>${version.vividsolutions.jts}</version> //version 1.13
        </dependency>

        <dependency>
            <groupId>org.postgis</groupId>
            <artifactId>postgis-jdbc</artifactId>
            <version>${version.jdbc.postgis.version}</version> //version 2.1.3
        </dependency>

persistence.xml

<persistence-unit name="Test-PU">
      <provider>org.hibernate.ejb.HibernatePersistence</provider>
        <jta-data-source>java:jboss/datasources/testDS</jta-data-source>

        <properties>
            <property name="hibernate.dialect"
             value="org.hibernate.spatial.dialect.postgis.PostgisDialect"/>
        <property name="hibernate.hbm2ddl.auto" value="update" />

            <!-- Debug -->
            <property name="hibernate.show_sql" value="true" />
            <property name="hibernate.format_sql" value="true" />
            <property name="hibernate.generate_statistics" value="true" />

...

My module in Wildfly

modules/org/postgresql/main

module.xml
postgis-jdbc-2.1.3.jar
postgresql-9.3-1101.jdbc4.jar

standalone.xml

<datasource jndi-name="java:jboss/datasources/testDS" pool-name="testDS" enabled="true"
 use-java-context="true">
<connection-url>jdbc:postgresql://localhost:5432/testDB</connection-url>
  <driver>postgresql-driver</driver>
    <security>
       <user-name>postgres</user-name>
       <password>postgres</password>
    </security>
 </datasource>

<driver name="postgresql-driver" module="org.postgresql">
  <driver-class>org.postgresql.Driver</driver-class>
</driver>

Mapping

@Column(name = "perimeter", columnDefinition = "Geometry")
@Type(type = "org.hibernate.spatial.GeometryType")
private Polygon perimeter;