texts: Python code to compare two fields in different feature layers?

texts:
I have two shapefiles which encompass 3897 features each. Corresponding shapefiles (village and watershed that flows into it) have the same ID (ws_id in one, OBJECTID in the other). Many watersheds overlap. Because of errors in the watershed layer, I want to compare the areas of each watershed with its corresponding village (settelment_area for villages and Shape_Area for the watersheds). If the watershed has a smaller area than the village it flows into, it must be deleted and the corresponding village's data added to that of the village into which the next watershed it flows in (data includes urban area and population).

I've tried doing this using cursors and temporary tables, but my python skills are not great and I seem to loose a lot of the related data along the way (ending up with empty lists).

using arcgis desktop 10.2 with advanced license

Code

#the aim of this code is to get rid of erroneous data by creating a new serviceshed_v0 shapefile which would contain only serviceshed which are larger than the settlement they run into.
#it could then be extended to have different minium size requirements for servicesheds

import arcpy
arcpy.CheckOutExtension("spatial")
arcpy.env.overwriteOutput = True

workdir = r'C:\Users\xxx\beneficiarylayers.gdb'
arcpy.env.workspace = workdir
fc1 = workdir + r'\servicesheds_v0'
fc2 = workdir + r'\miyun_settlements'
fc1table = arcpy.CreateTable_management(workdir, 'fc1table')
fc2table = arcpy.CreateTable_management(workdir, 'fc2table')
template = workdir + r'\template'
arcpy.CreateTable_management(workdir, 'servicesheds_v1', template)
newservicesheds = workdir + r'\servicesheds_v1'

rows = arcpy.SearchCursor(fc1)
for row in rows:
    arcpy.Append_management('Shape_Area', fc1table)

del row, rows

rows = arcpy.SearchCursor(fc2)
for row in rows:
    arcpy.Append_management('settlement_area', fc2table)

del row, rows

rows = arcpy.SearchCursor(fc1table)
for row in rows:
    if fc1table > fc2table:
        arcpy.Append_management(row, newservicesheds)

del row, rows

that's a first attempt. Problem is I realise tables are not the way to go because the data relating each tuple to a polygon is lost. i'd like to get something like this to work (below) but i'm not sure if it's allowed in python synthax

arcpy.CreateTable_management(workdir, 'servicesheds_v1')
newservicesheds = workdir + r'\servicesheds_v1'

cursor1=arcpy.da.SearchCursor(fc1, "Shape_Area")
cursor2=arcpy.da.SearchCursor(fc2, "settlement_area")

for row in cursor1:
    if 'Shape_Area' in cursor1 > 'settlement_area' in cursor2:
        arcpy.Append_management(row, newservicesheds)

FINAL CODE
here's a method which allows a rapid creation of a list giving all the polygons that fulfil the function of being too small.

import arcpy
arcpy.CheckOutExtension("spatial")
arcpy.env.overwriteOutput = True

workdir = r'C:\Users\xx\beneficiarylayers.gdb'
arcpy.env.workspace = workdir
fc1 = workdir + r'\servicesheds_v0'
fc2 = workdir + r'\miyun_settlements'

#set up cursors
cursor1 = arcpy.da.SearchCursor(fc1, ["ws_id", "Shape_Area"])
cursor2 = arcpy.da.SearchCursor(fc2, ["OBJECTID", "settlement_area"])
wrongsheds = []
#make a dictionary and store values from watershed table
serviceshed_area = {}
for row in cursor1:
  serviceshed_areas[row[0]] = row[1]

#loop through other table
for row in cursor2:
  if row[1] > serviceshed_areas[row[0]]:
    # if serviceshed_areas < village_Areas add to list and print
    wrongsheds.append(row[0])
    print "serviceshed {} is wrong".format(row[0])


texts: Can I use a cursor to select two rows and compare attributes

texts: Can I use the cursor search function to look at the first two rows in a table and select them so I can compare attributes and possibly merge them or change attributes? I am using ArcGIS 10.0


texts: Points In Attribute Table Will Not Display In Data View

texts:
I am using ArcMap 10.2.2 and opening a shapefile of points. My issues is this: the points do not appear in the data view. In other words, I cannot view the points in the display where they usually appear. I opened the attribute table and there are indeed several hundred points listed there; they just don't appear in the data view as they normally would. I have also tried "zooming to layer" just to be sure that I'm not off in la la land while the data is sitting somewhere else. This is does not help either.

***NOTE: I was given this shapefile from another source, and it did not come with a spatial reference. Hence, when I add this shapefile to the Table of Contents, it alerts me to the fact that there is no spatial reference. Perhaps this has something to do with it?

Also, the person who gave me this shapefile only provided me with a .shp, .shx and .dbf. Perhaps this has something to do with it as well...


texts:
I'm working on a project where we need our tiles rendered with tilt perspective.
I have installed and running openstreetmap server and it's working, but I need this projection to work:

  +proj=tpers +h=500000.000 +lat_0=19 +lon_0=19 +x_0=0 +y_0=0 +azi=30
  +tilt=40 +R=6371000 +units=m

if I insert this in style.xml I get this error from rendered
Unknown projection string, using web mercator as never the less

how can I get the tiles rendered with this projecting string or any other
perspective view


texts:
I want to view a .osm file that was generated from an Overpass API query.

I tried the steps given in QGIS tutorial on OSM wiki. The tutorial says

  Start QGIS Desktop
  Go to: Vector (menu) → OpenStreetMap → Import Topology from XML.
  Select the .osm file, and it'll suggest a name for the SpatiaLite database that will be created from this data.

But when I opened the QGIS, the vector menu does not give drop-down submenu.

Am I missing something in the steps or the QGIS tool on my machine is not installed correctly? Any suggestions?

This is the most related question I could find. But it was not very helpful when it comes to my problem.


texts: Topography - represent result in WKT as polygons with rings or multipolygon without rings?

texts:
I'm new to GIS vector, especially for GIS tools like ArcGIS, QGIS,...
My question is, for this image below with poylgon 1 and polygon 2, typically how a WKT (Well-known text) is represented in GIS tools?
I can see two possible cases, please speak from your experience what case should be used normally by GIS tools (e.g: ArcGIS, QGIS) when I have multiple polygons nested:

By using polygon with interior ring. Is that correct, WKT only supports 1 nested polygon as interior ring?
POLYGON ((35 10, 45 45, 15 40, 10 20, 35 10), (20 30, 35 35, 30 20, 20 30))

By using multipolygon containing polygons without any rings:
MULTIPOLYGON ( ((30 20, 45 40, 10 40, 30 20)), ((15 5, 40 10, 10 20, 5 10, 15 5)) )


texts: Geoserver SLD string functions for GetFeatureInfo requests

texts:
I have several Geoserver WMS layers styled with SLD.  I also run GetFeatureInfo requests against those layers.  My understanding is that Geoserver checks the layer's SLD rules even for GetFeatureInfo requests.

I'd like to be able to do string functions (for example StrReplace) on the attribute names to control how they are look when they are returned by GetFeatureInfo requests.  For example if the column name is 'CONAME', I'd like the ability replace it with "County" using an SLD rule.

I could easily do this as well using the Freemarker templates or a Postgres view, but the we'd like to be able to frequently change these through the Geoserver web admin interface.

I can't find a good example of the complete syntax I need to make a rule like that and everything I've tried either doesn't validate or breaks the request.


texts: How to retrieve latitude/longitude from Tiff with Geotools

texts:
I am new to Geotools and I have to work on project which includes a set of tiff image maps. In that I need to get latitude/longitude for a given point.

Ex: When mouse click occurs I want to get the lat/long for the clicked point.

I tried following code snippet

frame.getMapPane().addMouseListener(
            new MapMouseAdapter(){
                public void onMouseClicked(MapMouseEvent ev) {

                    //DirectPosition2D pos = ev.getMapPosition();
                    //System.out.println(pos.x + " " + pos.y);
                    DirectPosition2D pos = ev.getWorldPos();
                    System.out.println(pos.x + " " + pos.y);

                }
            }

    );

Are these simply the X/Y coordinates for a given point?

I appreciate your help, thanks.


texts: Batch intersect and multiprocessing with arcpy

texts:
I admit what I am working on is way over my head, especially as I've only been dabbling in python for a few months. I am attempting to write a script that allows me to use Intersect analysis on about sixty raster images for each of their sides. To make this more efficient, I would like to use multiprocessing. I have read many tutorials on this, but I have decided to follow the framework that ESRI used in their blog post as it appears to be the least confusing.

What I have written so far is only for two processes, but I will eventually be using a computer that will allow six.

import os
import re
import multiprocessing
import arcpy

def intersect_1():
    inFeatures = ["file1", "file2"]
    intersectOutput = ["intersect_result"]
    clusterTolerance = 1.5
    arcpy.Intersect_analysis(inFeatures, intersectOutput, "", clusterTolerance, "")

def intersect_2():
    inFeatures = ["file2", "file3"]
    intersectOutput = ["intersect_result2"]
    clusterTolerance = 1.5
    arcpy.Intersect_analysis(inFeatures, intersectOutput, "", clusterTolerance, "")

def main():
    workspace = 'c:/Users/Hannah/Desktop/wrs2_d'
    arcpy.env.workspace = workspace

    pool = multiprocessing.Pool()
    pool.map(intersect_1,intersect_2)

    pool.close()
    pool.join()

if __name__ == '__main__':
    main()

These are the errors I continue to get, but my grasp of basic python is so shaky I cannot understand why. Can anyone point me in the right direction?

Traceback (most recent call last):
  File "C:/Users/rs1/Desktop/possol_1.py", line 31, in <module>
main()
  File "C:/Users/rs1/Desktop/possol_1.py", line 25, in main
pool.map(intersect_1,intersect_2)
  File "C:\Python27\ArcGIS10.3\lib\multiprocessing\pool.py", line 251, in map
return self.map_async(func, iterable, chunksize).get()
  File "C:\Python27\ArcGIS10.3\lib\multiprocessing\pool.py", line 304, in map_async
iterable = list(iterable)
TypeError: 'function' object is not iterable


texts:
In QGIS 3.4 MADEIRA, I get no active "data defined override"; I click on the sybol right toio "Size", that I want to modify with a function. The words "data defined override" are pale grey and italic, and when I click on "description " I get this

I suppose I must activate it somewhere, but I was unsuccessfull to find where.


texts: IDL - Cropping GeoTiff with IDL (no Envi)

texts:
How can I crop a GeoTiff (EPSG:4326) in IDL by using pixel coordinates, and then saving the cropped files with the correct georeferencing?

I assume I have to re-set the information in the MODELTIEPOINTTAG in the GEOTIFF-Structure?
If I understand that tag correctly, it contains (in my case) the coordinates of the Upper-Left corner in Degrees positive in East and North direction:

MODELTIEPOINTTAG={0, 0, 0, 90.08, 23.63, 0}

So basically I can just crop the image by pixel coordinates, and then convert the pixel-coordinates of the upper-left crop-point into Lat-Lon coordinates (in the original image frame) and write them into the MODELTIEPOINTTAG, but leave every other tag the same, correct?

Do I take the corner of the for that, or the center (i.e. add 0.5 to the pixel coordinates)?

How can I accomplish this in IDL? Are there already methods available for that?


texts: Making polygon for every point in set using voronois?

texts:
I'm trying to make polygons from points, were each point will be represented by one polygon. I tried using voronois, but the results are incomplete. While I more or less understand where this uncompletness comes from, I don't know how can I solve such problem.

I also have another, bigger polygon that is an area containing all the points. In the end, I want all the small voronois to cover 100% of this big polygon area.
In other words:

I want all the empty white spaces filled with voronois/polygons and I want each voronoi/polygon to contain just one point.

My code:

from scipy.spatial import Voronoi
import numpy as np
import shapely
import geopandas as gpd

#pointset is geodataframe with points, region is geodataframe with big
#polygon I want to fill with voronois

x = pointset.geometry.x.values
y = pointset.geometry.y.values
coords = np.vstack((x, y)).T
vor = Voronoi(coords)

lines = [shapely.geometry.LineString(vor.vertices[line]) for line in
    vor.ridge_vertices if -1 not in line]
polys = shapely.ops.polygonize(lines)
voronois = gpd.GeoDataFrame(geometry=gpd.GeoSeries(polys))

voronois = gpd.overlay(voronois, region)

So basicly I just use Voronoi and then intersect the result with my region polygon. Any advice about how to make result complete?


texts:
In reading about DBSCAN clustering it mentions the SEQUOIA 2000 Benchmark database was used for benchmarking.  Is this data set available anywhere?  Or, is there another set like it?

So far my searching hasn't turned up anything.


texts: Authenticated WMS not accessible through ArcCatalog?

texts:
I have a Geoserver WMS which has basic authentication for several layers that have data security restrictions against them.

After putting in the username and password in ArcCatalog's add WMS server facility, I am noticing these layers are not showing (layers that do not have any restrictions against them do show).
Are there any known issues with connecting to authenticated WMS services? Do I need to configure anything?

I am using ArcCatalog 10.2.


texts: How to migrate a specific workspace, store and layer in Geoserver

texts:
Geoserver is running on a development and production environment. I would like to migrate just one workspace, store, layer, style, etc from development to production. How can I do this?

I know I can copy the entire data-dir to production and restart geoserver; but than all workspaces, stores, etc are migrated. The point is I want just one workspace, store, etc.

In the development environment I have several test workspace, store, and layers, which I do not want on production.

In the data-dir I see a directory gwc-layers, with a lot of xml-files. But how do I know which xml is which layer?


texts: Adding built-in tools to python add-in toolbar

texts:
I want to make an ArcGIS Desktop add-in toolbar that has a regular ESRI built-in tool on it (like the buffer tool, or select tool, etc), instead of a custom tool. Is that possible or are we talking ArcObjects at that point?

This would be simple if I was looking to add tools for me only. But that's not the case. I'm looking to create a toolbar with standard tools that anyone in the organization could add to ArcMap via an add-in. We have a situation where sometimes the ArcMap profiles get reset in Citrix, and it would be nice when that happens, if the user could simply launch an add-in and get the standard tools back fairly easy. So using my above example, I want to make a toolbar with the buffer tool and the select tool, and we'll call it "Toolbar A". And anyone in the organization could load "Toolbar A" and get those tools. Is this possible using the ESRI python add-in solution?


texts: Getting current center of view in OpenLayers 3?

texts:
I would like the get the center (coordinates) of the current view. I need it for the geoserver print module.

I've found the code var center= ol.extent.getCenter({}); but it unfortunately doesn't work. Furthermore I've found a code in the API (getview.getcenter()), but it only shows the initial center and it is not changing when I explore on the map.

This is my code, it works fine except the center params... It gives back tha same values everytime.

var configs = {

    "layout":"A4 portrait",
    "srs":"EPSG:3857",
    "units":"meter",
    "dpi":300,
    "outputFilename": "valami",
    "layers":[
     {baseURL: 'xy',
            singleTile: false,
            type: 'OSM',
            maxExtent: [-20037508.3392,-20037508.3392,20037508.3392,20037508.3392],
            tileSize: [256, 256],
            extension: 'png',
            opacity: 1,
            resolutions: [156543.0339, 78271.51695, 39135.758475, 19567.8792375, 9783.93961875, 4891.969809375, 2445.9849046875, 1222.99245234375, 611.496226171875, 305.7481130859375, 152.87405654296876, 76.43702827148438, 38.21851413574219, 19.109257067871095, 9.554628533935547, 4.777314266967774, 2.388657133483887, 1.1943285667419434, 0.5971642833709717]
        },
{"baseURL":"xy",
        "opacity":1,
        "singleTile":true,
        "customParams":{},
        "type":"WMS",
        "layers":["mission:corine_ps"],
        "format":"image/png",
        "styles":[],
        "overview":true},
{"baseURL":"xy",
        "opacity":1,
        "singleTile":true,
        "customParams":{},
        "type":"WMS",
        "layers":["mission:utak_ps"],
        "format":"image/png",
        "styles":[],
        "overview":true},
            ],
    "pages":[
        {
        "center":[map.getView().getCenter()[1],map.getView().getCenter()[0]],
        "mapTitle":"",
        "comment":" ",
        "scale": 1000000,
        "rotation":0,
        "icon":""}],
    "legends": [
        {
            "classes": [
                {
                    "icons": [
                        ""
                    ],
                    "name": "",
                    "iconBeforeName": true
                }
            ],
            "name": "utak"
        }
    ]
    };

 var printUrl= 'xy?spec='+JSON.stringify(configs);


texts: Divide polygon into parts of certain size?

texts:
I have a shapefile polygon which is 13 acres. I want to divide it for 3, 4 and 6 acres polygons.

Is there any plugin or something which I can divide the original polygon?

I use Qgis 2.18.5

I added a photo how I want to divide


texts:
I have a lot of DWG files like basemap, and water and wastewater network.

How do I import these files into my QGIS project?


texts: The GetLegendGraphic query does not return legends for the specified bbox value

texts:
I have a GetLegendGraphic request like:
http://localhost:8080/geoserver/ows?service=WMS&request=GetLegendGraphic&width=19&height=19&layer=Workspace:Layer&rule=&bbox=8126200.9087728467,9581252.8025055025,9577653.389397664,10414397.006185452&srcwidth=2743&srcheight=1574&srs=EPSG:3857&legend_options=hideEmptyRules:true&format=image/png.

Layer is based on PostGIS. GeoServer version 2.19.
I cannot get the legend in this BBOX, although I get the correct image when I request an image by GetMap. The GetLegendGraphic and GetMap queries are based on the same principle of filtering features, as I understand it, but it is not clear why then they give different results.
There are 4 screen shots. BBOX in GetMapCorrect gives me correct legend in GetLegendGraphicCorrect query. BBOX in GetMapIncorrect gives me error legend in GetLegendGraphicIncorrect.
GetLegendGraphicCorrect:
Request
http://localhost:8080/geoserver/ows?service=WMS&request=GetLegendGraphic&format=image/png&width=19&height=19&layer=Workspace:Layer&rule=&bbox=7830826.700968759,9960757.157649865,9252323.14903124,10776706.38235013&srcwidth=2743&srcheight=1574&srs=EPSG:3857&legend_options=hideEmptyRules:true

Response

GetLegendGraphicIncorrect:
Request
http://localhost:8080/geoserver/ows?service=WMS&request=GetLegendGraphic&width=19&height=19&layer=Workspace:Layer&rule=&bbox=8126200.9087728467,9581252.8025055025,9577653.389397664,10414397.006185452&srcwidth=2743&srcheight=1574&srs=EPSG:3857&legend_options=hideEmptyRules:true&format=image/png

Response

GetMapCorrect:
Request
http://localhost:8080/geoserver/ows?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&FORMAT=image%2Fpng&TRANSPARENT=true&LAYERS=Workspace:Layer&CRS=EPSG%3A3857&STYLES=&WIDTH=2743&HEIGHT=1574&BBOX=7830826.700968759,9960757.157649865,9252323.14903124,10776706.38235013

Response

GetMapIncorrect:
Request
http://localhost:8080/geoserver/ows?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&FORMAT=image%2Fpng&TRANSPARENT=true&LAYERS=Workspace:Layer&CRS=EPSG%3A3857&STYLES=&WIDTH=2743&HEIGHT=1574&BBOX=8126200.9087728467,9581252.8025055025,9577653.389397664,10414397.006185452

Response


texts:
I can't download QGIS 2.8.2 via googlePlay. Device is not compatible. Device is Nexus 7 2012 with android 5.1.1
Do you have same issue? On old android 4.4 it was working fine.


texts:
I have been facing an issue to capture the delta updates on SDE.

We have written triggers on A and D tables to capture the delta updates(to capture the version name based on state id's), however I couldn't able to get the very first delete of any feature from a new version that has been created on Default.

After debugging, I have observed that, after the very first delete of a feature, the State Id in versions table is not being updated however the State Id in SDE_States is updated. So, the state id in versions table and the state id in SDE_States table is different.

I would like to say that the same query is working fine from second deletes and so on.... Having issue with getting the version name for only the first delete.

For your reference, please do find the query below.

Select @ObjectId=SDE_DELETES_ROW_ID, @SdeStateId=DELETED_AT, @SdeCurrentStateId = SDE_STATE_ID FROM inserted;

SELECT @Version= Name FROM sde.SDE_versions v,sde.SDE_states s INNER LOOP JOIN
    sde.SDE_state_lineages l ON l.lineage_name = s.lineage_name WHERE s.state_id = v.state_id and lineage_id = @SdeStateId and s.state_id = @SdeStateId

IF @version IS NULL
SELECT @Version= Name FROM sde.SDE_versions v,sde.SDE_states s INNER LOOP JOIN
    sde.SDE_state_lineages l ON l.lineage_name = s.lineage_name WHERE s.state_id = v.state_id and lineage_id = @SdeCurrentStateId


texts:
I'm trying to use the ArcMap Standardize Address tool (http://desktop.arcgis.com/en/arcmap/latest/tools/geocoding-toolbox/standardize-addresses.htm). I'd like it to change "West" to "W" in the PreDir field as it shows in the example, but my results are actually the opposite. Addresses with either "West" or "W" both get "West" in the PreDir field.

Why is this happening and how can I get it to function as the tool describes?


texts: How can the flow accumulation in a hydrological model/stream network be maintained across tiles of a large area of interest?

texts:
My study area is too large to be processed in its entireity by any hydrological modelling software I have found yet. I will have to break it up in to smaller areas, but there will be flow across these tiles. My outcome needs to include the total flow for each stream.    Are there any tips or best practices for approaching this?

EDIT: My  study area is 3740 square km and the DEM is 1 metre resolution (from lidar). It can be split into two major watersheds of roughly equal size. I'd assumed these would be too large to process on their own as I have trouble even loading them into most software.


texts: Why does ST_Area(g.geom) + ST_Area(p.geom) - ST_Area(ST_Union(g.geom,p.geom)) produce negative result when using ST_Overlaps(g.geom, p.geom) = true?

texts:
I have two set of shape data. Let's say file gov.shp contains list of post codes and their boundaries. File post.shp contains list of government area name and their boundaries. I want to find out the the union (may not be the correct word here) between two files (i.e. post code 1, 2 and 3 overlaps with government area a; post code 1 also overlaps with government area b). Bear in mind that government area is bigger than post code area. So I have a query like below:

select  g.id, g.name, p.id, p.code
from    gov g
inner join
        poa p
on      ST_Overlaps(g.geom, p.geom) = true

However if I apply this where clause, to my surprise the count(*) actually returns results. I am confused by why total area of g + p - union_area of(g,p) can be negative...

where   ST_Area(g.geom) + ST_Area(p.geom) - ST_Area(ST_Union(g.geom, p.geom)) < 0

Am I using the wrong function (ST_Overlaps)? Should I be using ST_Intersecs/ST_Equals/other functions?

I am new to this area so I am not too sure what the differences are in terms of solving my particular problem.

I imported shape file as geometry and MULTIPOLYGON object and I tried different results in both PostGIS and SQL Server 2008 R2 and both produced slightly different results. The above sql query syntax is in PostGIS though.

--SQL Server 2008
-- STOverlaps
-- DIFF < 0 = 27
-- DIFF > 0 = 1646

--SQL Server 2008
-- STEquals()
-- 42

--PostGIS
-- DIFF < 0 = 730
-- DIFF > 0 = 5312

--PostGIS
-- STEquals()
-- 0

SQL Server 2008 syntax

select  g.id, g.name, p.id, p.code
from    gov g
inner join
        poa p
on    g.geom.MakeValid().STOverlaps(p.geom.MakeValid()) = 1
--for some reason I had to use *.geom.MakeValid()
where g.geom.MakeValid().STArea() + p.geom.MakeValid().STArea() - g.geom.MakeValid().STUnion(p.geom.MakeValid()).STArea() < 0

I used Shape2SQL to import into SQL Server 2008, Planer Geometry, didn't check SRID, create spatial index.

And I used PostGIS Shapefile Import/Export Manager which came with PostGIS 2.0 to import into PostGIS, and all default settings.


texts: Uploading shapefile onto Garmin GPS, and displaying it over basemap?

texts:
I have a roads shapefile layer in ArcMap and I want to be able to view it on my Garmin64st. Ideally, it would be overlaid on a basemap so I would always be able to see it. It seems like that means I need to convert my .shp file to .img.

I'm finding various programs that claim to help with some step of this, but some of them are really expensive, and none of the free ones I've tried have gotten me to .img so far. I did download an ArcGIS Toolbox which converts .shp to .gpx, but I haven't been able to find a way to convert that to .img for free yet.

Am I missing something?

Does anyone know a simple and free way to convert a .shp to be viewed in a Garmin GPS like I've described?

If I absolutely must purchase a program, do you know of one that works?


texts: Source for cloud cover map layer over UK

texts:
I'd like to add a layer to display historical data on cloud cover in a Leaflet project I'm currently working on. This is specific to a given day in 2012 for a city in England, and ideally I'd like to have 1-hour updates. This way users could watch the changes in cloud cover - it shouldn't take much work to put in a couple of buttons to skip back and forth one hour.

I seek suitable data sources/layers, preferably with an open license.


texts: One label per polygon in Mapzen's Tangram JS (MVT)

texts:
I have problem with Mapzen's TangramJS labelling. I am using MVT tiles and labels as an attributes to polygons that are multiplying in large scales. Below an example:

I would like to put one label per feature. Here is an example solving code but only for TopoJSON and GeoJSON tiles - not working for MVT.

Do you have any idea how to resolve this problem for Mapbox Vector Tiles?

I have found mapzen's code for GeoJSON i TopoJSON:

key: "preprocessFeatures",
value: function(e) {
if (this.config.generate_label_centroids) {
        var t = [],
            n = {
                label_placement: !0
            };
        e.forEach(function(e) {
            if (null != e.geometry) {
                var i = void 0,
                    a = void 0;
                switch (e.geometry.type) {
                    case "Polygon":
                        a = s(i = e.geometry.coordinates, e.properties, n), t.push(a);
                        break;
                    case "MultiPolygon":
                        i = e.geometry.coordinates;
                        for (var r = -1 / 0, o = 0, l = 0; l < i.length; l++) {
                            var u = f.default.polygonArea(i[l]);
                            u > r && (r = u, o = l)
                        }
                        a = s(i[o], e.properties, n), t.push(a)
                }
            }
        }), t = t.filter(function(e) {
            return e
        }), Array.prototype.push.apply(e, t)
    }
    return e

What should I do to implement it to MVT?


texts:
Update for QGIS 3.2 of previous queries raised for earlier versions
Symbols above labels in map QGIS

Labels always appear on the top label layer, as previously described.  I want to label multiple layers, but don't want every label to appear in the top layer.  I use some labels for background context, but don't want them to take priority over the main layers of maps.  I've looked at some workarounds suggested, but none are really working for me.  My best solution at the moment is to use transparency & other formatting to lessen impact of some labels.

Below is an example, where I am using place names used for context, but with the point labels removed. I can format the point labels differently and make them take preference to differentiate, but ideally I'd like these to appear under my point layer. The points are the key element of the map, not the place name and should be taking priority - both the marker and the label.

I can also tweak individually label positions to fix specific issues, but not ideal as I need to use the map at different scales and don't want to manually have to adjust for each case.

Have I missed anything in 3.2 which allows more control over labels/layers, or any other workarounds that allow better control of labels or allow other layers to appear above labels?

I am also playing around with using multiple maps in print composer, but ideally I'd like to find a way to do on one map.


texts:
I am trying to download building footprints from OpenStreetMap so I can import them into ESRI's CityEngine but I am only able to download road files (as routes), as seen below.

I am trying to download a section of Seoul (Gangnam) to display the area in CityEngine.

I am new to downloading OSM data.


texts:
I was able to download and install the Kartograph.py library. When I call it from the command line I am getting:

kartograph is not an external command

Usually from what I understand this is solved by setting the Path on Windows 7. I have tried everything I can think of and still get the same error.

For example if my path is:

c:\python\katrograpy I have tried to set it to c:\python\katrograpy, c:\python\katrograpy\build, c:\python\katrograpy\bin with no luck.

I am not a python programmer. I just need to use this to convert SHP files to SVG.

I am calling it this way. in c:\python\katrography I saved my mapConfig.json file and call it like this: kartograph mapConfig.json


texts: Debugging RuntimeError: Object: Error in executing tool from arcpy.RasterDomain_3d()?

texts:
I´m always getting an Error with my Python-Code. Even though it is almost the same as on the ArcGIS Help Website.

The Error occurs in line 19:

  Traceback (most recent call last):

  File "C:_geodata\check.py", line 19, in

  arcpy.RasterDomain_3d(rasterList,out_file,geometry)

  File "c:\program files (x86)\arcgis\desktop10.1\arcpy\arcpy\ddd.py", line 998, in RasterDomain

  raise e
  RuntimeError: Object: Error in executing tool

  Failed to execute (IterateRasterDomain).

My Code looks like this:

import arcpy
from arcpy import env

arcpy.CheckOutExtension("3D")

env.workspace = arcpy.GetParameterAsText(0)
rasterList = arcpy.ListRasters()

if rasterList:
    for raster in rasterList:
        geometry = "POLYGON"
        out_file = "domain_" + raster[:-4] + ".shp"
        arcpy.RasterDomain_3d(rasterList,out_file,geometry)
    print "Finished"
else:
    print "No Data available"


texts: How to use Scale tool when Georeferencing Raster in ArcMap?

texts:
I'm georeferencing a bunch of maps in ArcMap 10.0.  The Shift and Rotate commands are pretty easy to use, no problems there.  But whenever I use the Scale command, the image completely changes position.  Is there a way to put an anchor down (similar to the rotate tool)?  I feel like I must be missing something but haven't figured out what it is yet.  I don't want to rubbersheet them into place, so I'm trying not to use control points until I have it in position, but first I need to make sure the scale bar on the map matches the distance it's supposed to.

Edit: Asking a question does wonders in finding an answer, even when nobody has answered yet.  I found the Rescale_management tool and will use that to set the scale properly first, then position.  If the length of the scale bar is 1256ft and it should be 1200ft, I'll use a factor of .9554 in the X and the Y directions.  Is there a better way?


texts: Is it possible to have a customized zoom level?

texts:
I am setting up my own Map Tile Server according to this guide and it is working well now:
http://switch2osm.org/serving-tiles/building-a-tile-server-from-packages/

The map currently has 19 zoom levels and it is initialized in the OpenLayers object "map" as below code shows:

map = new OpenLayers.Map ("map", {
            controls:[
                new OpenLayers.Control.Navigation(),
                new OpenLayers.Control.PanZoomBar(),
                new OpenLayers.Control.Permalink(),
                new OpenLayers.Control.ScaleLine({geodesic: true}),
                new OpenLayers.Control.Permalink('permalink'),
                new OpenLayers.Control.MousePosition(),
                new OpenLayers.Control.Attribution()],
            maxExtent: new OpenLayers.Bounds(-20037508.34,-20037508.34,20037508.34,20037508.34),
            maxResolution: 156543.0339,
            numZoomLevels: 19,
            units: 'm',
            projection: new OpenLayers.Projection("EPSG:900913"),
            displayProjection: new OpenLayers.Projection("EPSG:4326")
        } );

I am wondering if I want to zoom beyond the maximum zoom level 19, say zoom into level 23 or even better, have my own zoom level reconfigured, how to do that? This is because I want to display an indoor room map on my current map server. The indoor room map is small it needs zoom level higher than 19.

Please point me to some resources to achieve this task  as well.


texts: How do Wikipedia maps add a coordinate point to a base map?

texts:
On wikipedia pages for localities, there is often a map on the right side of the page in the introduction, and the map is tagged with the point of the locality. For example, see below:

https://en.wikipedia.org/wiki/Buenos_Aires

How do I achieve something similar in a website (use one map image, and use code to plot a coordinate point and labelling on the map)?


texts: What is difference between hosted and non-hosted feature layer in ArcGIS Online?

texts:
When I add a slpk to ArcGIS Online there is an option if I want to host it as a hosted feature layer but this applies to 2D feature layers as well.

What is the difference between a hosted feature layer and a non-hosted feature layer?


texts:
I'm using QGIS Version 3.16.4. I thought one can do heatmaps by going to the properties of a layer and selecting it under symbolization. But as you can see, it's not available in the dropdown menu. Does anyone know what's the reason for that or how I can activate the heatmap option? The Vector-Layer only consists of polygons which have 3 attributes, so nothing fancy there. Maybe there's a plugin I can use? The two which I found need a timestamp unfortunately, which is not part of my layer

I know you can do something similar by using graduated in the drop down menu, it's explained in this thread but I still want to know why I can't do it with heatmaps.


texts:
In the image below, the green points are associated with radar wave amplitude data. They're generated from RasterToPoint, so they're at the centroid of each pixel and have only the value of the amplitude at that point. The red points are a associated with surface motion at that location, so they have values for displacement, acceleration, coherence, etc.

What I want to do is associate each red point (motion) with its nearest green neighbor (amplitude). The resulting dataset will have the attributes of both datasets - except that for green points not near a red point (this will be most of them) the data for the motion attributes will simply be empty.

I suppose the better way to do it is to snap the datset with the smaller number of points to the datset with the greater number of points.

I am using ArcGIS 10.2 for Desktop.


texts: How to get l.control layers or copyright into panel / div

texts:
I am tring to get l.controls like layers, copyright or draw into panel/div like:

...HERE LAYERS...

...HERE COPYRIGHTS...

So my question is, how to define new position or function that will get layers into panel/div?


texts: Adding some background information to a layer

texts:
I have 2 Python Qgis plugins. The first plugin is creating a QgsVectorLayer via the ogr dataprovider. The second plugin is offering some services on the features in the layer.

I need to pass some additional information (a single string) which is not in the file that ogr reads to the second plugin. What's the best /easiest way to attach it to the layer so the second plugin can easily extract it?


texts:
I am using QGIS 2.2 Valmiera to try to find the average slope of a water course. I've drawn in the water course as a polyline over a large number of contour lines (the water course is long with a very large catchment area).  What I've tried to do so far is create a tin file. I've used the Raster > interpolation function using the heights attribute of the contour line shapefile to make the tin and each time I try it gives me a mini-dump error.  I've tried to simplify the geometry (from 392515 to 66801 vertices) but I still just get the mini-dump error.

Is there another way to get the slope of this polyline without having to make a tin file or am I doing something wrong in trying to make the tin?


texts: Formula for maximum coverage of a geostationary satellite based on earth station elevation

texts:
I am trying to create a JavaScript function to represent the maximum coverage of a geo satellite for different antenna elevations (0 degrees, 5 and 10 degrees) with 3 leaflet polylines.

The result is similar to http://reference.wolfram.com/language/ref/GeoVisibleRegion.html (see "The visible region from geostationary orbit", the contour in the example is for an elevation of 0 degree).

Any hint is welcome, I only found partial formulas on google.. it is not very clear


texts: Constructing a function to pass to lidR::grid_metrics

texts:
I am attempting to write an R function, as part of a package, that takes a lidar tile (LAS object) and a vector of cut-points defining vertical layers, and returns a RasterStack where each layer has a layer point count as the cell value. The function makes use of the grid_metrics function in the lidR package. Here is a simplified version of the code to illustrate:

# Main package function
get_layer_counts <- function(las, breaks, res = 10) {
  f <- function(Z) layer_count_metrics(Z, breaks)
  stats <- lidR::grid_metrics(las, func = f(Z), res = res)
  lidR:::as.raster.lasmetrics(stats)
}

# Metrics function to count points within height intervals
layer_count_metrics <- function(Z, breaks) {
  zcat <- findInterval(Z, breaks, left.open = TRUE)
  as.list(table(zcat))
}

Unfortunately this gives an error: "Error in eval(call) : object 'f' not found" which comes from the lidR:::lasaggregate function.

It seems (I think) that unless the function passed to grid_metrics is defined in the global environment it will not be found. Is that correct and, if so, can anyone suggest a work-around?


texts: Zonal statistics (mode, variety) in PostGIS

texts:
I would need to get different statistics than are available in ST_SummaryStats.
In particular I need mode (majority) and variety (number of different values) for a discrete raster and according polygonal zones. This calculations are available for example in zonal statistics tool in QGIS, but I am looking for solution with PostGIS query.
Any clue?


texts: Letters change when saving tab into shp in QGIS

texts:
I am working in QGIS 2.18. When I open and then save a .tab file as .shp file I get different symbols in the letters of the attribute table. I don't know what coding option I should use in order to see the original letters (Spanish alphabet symbols: ñ, accents, etc...).

By default the "save vector layer as..." dialog window chooses UTF-8.

I already tried some other options but I get unreadable text, so I rather ask before I do something that hurts the file.

In fact, when I open the .tab file in MapInfo I can see it correctly. Is it something related to the system settings of the program (QGIS/MapInfo)?


texts: Determine start/end points for line slope calculation

texts:
This is a follow up question based off of this question.

I have a line feature class representing storm drains and a point feature class representing where storm drains are connected.  The point feature class has an elevation for each point allowing me to determine the slope of the storm drain as the difference in elevation from start to end divided by the length of the storm drain.  However, most of my points act as starting points for one storm drain and end points for another at the same time.  I do know that the start point should always be the point with a higher elevation to facilitate water flow.

How can I calculate the slope for each line segment and designate which points are the start/end points for each line?

Would it be better to calculate new lines based off of my points layer that would include a slope field or is there a way to automate this calculation between the two layers?

I have an ArcMap Advanced license but no 3D or Spatial analyst extensions.  I am open to Python and QGIS solutions.  I just don;t want to have to enter these values in manually as I know there must be a way to automate it.


texts: How to create header files for multiple rasters in R

texts:
I have nearly 1000 generic binary images, i.e., headerless, band sequential and binary (16-bit signed integer).  I have been succesfully able to load these with package 'raster' only when I create an associated .hdr file with the same name as the image.

Hence I would like to create an associated .hdr file for each image, with the following information;

ENVI
description = {
  AVHRR_GIMMS 8km global NDVI x1000}
samples = 2160
lines   = 4320
bands   = 1
header offset = 0
file type = ENVI Standard
data type = 2
interleave = bsq
sensor type = AVHRR
byte order = 1

I am drawing a blank on how to accomplish this.  I don't have any programming skills outside of R.

Any suggestions? Is there a more efficient way to accomplish this (e.g., specifying a single header for all files)?

This question is similar to this (presently) unanswered question, except I suspect it may be more straightforward in R than Arcgis.


texts: ImportError: No module named qgis.core from PyQGIS?

texts:
Trying to load QGIS --

I have a working copy of Python: 2.7.9 --  It looks like QGIS is not finding PyQGIS.

I am aware that this is a problem with the PATH but do not see how to fix it.

 > My OS: Debian GNU/Linux 7 (wheezy)

I got a warning when I issued the "apt-get update" command

  W: GPG error: http://qgis.org wheezy Release: The following signatures couldn't be verified because the public key is not available: NO_PUBKEY 3FF5FFCAD71472C4

When I launch the application ("~> qgis") I get two warnings

  Warning: loading of qgis translation failed [/usr/share/qgis/i18n//qgis_en_US]

  Warning: loading of qt translation failed [/usr/share/qt4/translations/qt_en_US]

and a window with the following info:

  Couldn't load PyQGIS.

  Python support will be disabled.
  Traceback (most recent call last):
    File "", line 1, in

  ImportError: No module named qgis.core

  Python version:
  2.7.3 (default, Mar 13 2014, 11:26:58)
  [GCC 4.7.2]

  QGIS version:
  2.10.0-Pisa 'Pisa', exported

  Python path:
  ['/usr/share/qgis/python', u'/home/cjohnson/.qgis2/python', u'/home/cjohnson/.qgis2/python/plugins', '/usr/share/qgis/python/plugins', '/usr/local/anaconda/lib/python2.7/site-packages/pyxnat-0.9.5.3-py2.7.egg', '/usr/local/anaconda/lib/python2.7/site-packages/VTK-6.0-py2.7.egg', '/usr/local/anaconda/lib/python2.7/site-packages/distribute-0.6.28-py2.7.egg', '/usr/local/anaconda/lib/python2.7/site-packages/MySQL_python-1.2.4b4-py2.7-linux-x86_64.egg', '/usr/local/anaconda/lib/python2.7', '/usr/local/anaconda/lib/python2.7/plat-linux2', '/usr/local/anaconda/lib/python2.7/lib-tk', '/usr/local/anaconda/lib/python2.7/lib-old', '/usr/local/anaconda/lib/python2.7/lib-dynload', '/usr/local/anaconda/lib/python2.7/site-packages', '/usr/local/anaconda/lib/python2.7/site-packages/PIL', '/usr/local/anaconda/lib/python2.7/site-packages/Sphinx-1.2.3-py2.7.egg', '/usr/local/anaconda/lib/python2.7/site-packages/runipy-0.1.1-py2.7.egg']


texts: Converting Points to Line using ArcGIS Desktop?

texts:
I have a shapefile that contains points (xy coordinates) that display several different routes for a transit agency.  In the attribute field I'm given information on what route each point corresponds to and also the sequential numbering for each point in each route (starting back at 1 for every new route).

How do I convert these points into lines separated by each of the routes?

I used the point to line tool but it created extra lines connecting each of the routes together.  I even tried specifying the Line Field and Sort Field to the fields with the route numbers/sequential numbering and got an error.


texts: Exporting Contours from QGIS into AutoCAD

texts:
I'm trying to export contour lines from QGIS into standard AutoCAD but I'm not having much luck.

It does not seem to matter what type of file I select as the output during the contour creating (eg. DXF or CSV) all I get is a folder containing 3 files (an SHX,a compiled shape source file and a DBF file all titled as contour. I am then unable to open these files in AutoCAD.

If I then export the contour lines from QGIS as a DXF that file will open in AutoCAD but the lines appear to be missing their Z property.

Is there a simple way to export contours from QGIS that retains their xyz information?


texts: Why can't QGIS recognize the right Layer CRS?

texts:
I have a vector layer (shapefile), that I projected with ArcMap 10.4 into UTM33 (EPSG: 25833). If I load this layer in ArcMap, 25833 ist always recognized as the CS of the layer, as expected.

But when I load the same layer into QGIS 2.14, it can't recognize 25833. Instead it tells me the layer is allegedly 3006 (Sweref), which should be wrong.

What is the reason behind that problem and how can I prevent it? Happens quite often....

Here is the *.prj file content:
 PROJCS["ETRS89_UTM_zone_33N",GEOGCS["GCS_ETRS_1989",DATUM["D_ETRS_1989",SPHEROID["GRS_1980",6378137,298.257222101]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",15],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]]


texts: Creating raster object in R (raster package) by browsing directory

texts:
I am working with raster data (.tif) in R using the Raster Package, and I know how to create a raster object through
r <- raster(filename)

I am working with models with some long filenames, and so having to type out the filename is slowing me down. e.g.
filename <- "wc2.1_2.5m_bioc_CNRM-CM6-1_ssp245_2041-2060.tif"

Is there a way to create a raster object using my windows browser, or a simpler way to create raster objects, when I have multiple filenames with unwieldy names?


texts: Arcade Expression to Intersect Features and add Attribute to Popup

texts:
New to Arcade. I have two features in an AGOL map, "SS Parcels" and "SS Zoning". I would like to configure the popup for "SS Parcels" to also display the zoning code from "SS Zoning".

It seems like this is possible from this article:
Overlapping Features in Pop-Ups

However I can't seem to make my code work at all. It looks like the function "Intersects" actually returns a Boolean, maybe this was changed since the article was written. If I instead use the function "Intersection", which returns a geometry, my code looks like this:

var intersectLayer = Intersection(Buffer($feature, -10, 'feet'),FeatureSetByName($map, 'SS Zoning'))
for (var f in intersectLayer){
return f.ZONING
}


texts: Using OpenLayers to collect data by CrowdSourcing?

texts:
I'm currently making my first steps in OpenLayers, trying to build a page showing a map where people can add data, which is then saved back in our database.

(Background: Mapping of the distribution of indigenous vegetables. We have some data, but the "crowd" should be able to add their information, too).
Similar to this example (but in my case polygons should get added): http://ediblecities.org/

Is this possible with OpenLayers?

And this is maybe a too open question, but how? (a reference to a simple example would already be helpful).
It does not seem possible to somehow use the tools offered in the toolbar (OpenLayers.Control.EditingToolbar(vlayer)) for this.

A more concrete follow-up question as to HOW to collect data via Openlayers, Geoserver and PostGIS can be found here: How to collect data via WFS and Popup


texts: Creating VRT or index of large amount of aerial images to create mosaic in QGIS?

texts:
I have about 1600 aerial images that I want to create a mosaic of. I tried doing about 400 at a time in a VRT, it's extremely slow (been building polygons for about 15 minutes now + the 20 minutes it took to build the VRT). Is there a more efficent way to do this, like by creating an index somehow?

To be honest I don't know exactly what that entails, but I noticed at my job we keep large amounts of tiles indexed in an excel file, you load the file into arcmap and bam you've got a snappy large scale mosaic. Can I do something to the same effect in QGIS?


texts:
I dowloaded the climate data from this site http://catalogue.ceda.ac.uk/uuid/3f8944800cc48e1cbc29a5ee12d8542d and they have the extension .dat .
I'd like to import the information in QGIS as raster data.


texts: Mapserver and Tilecache configuration issue

texts:
Currently I am using mapnik2 with tilecache for map generation. As there is shift in some of the tiles generated,I am trying to configure mapserver with tilecache for caching tiles. I have installed mapserver 6.2.1 and tilecache 2.11. In the tilecache.cfg file , I configure as below to access map data:

[Layername]
type=MapServerLayer
bbox=255923,904698,1050731,1437039
data_extent=255923,904698,1050731,1437039
layers=layer1,
mapfile=/path to mapfile/
srs=EPSG:32643

When I check if tiles are generated by using the url:

http://tilecache.org/docs/README.html

an xml file is generated with error being,

<?xml version='1.0' encoding="ISO-8859-1" standalone="no" ?>
<!DOCTYPE ServiceExceptionReport SYSTEM "http://schemas.opengis.net/wms/1.1.1/exception_1_1_1.dtd">
<ServiceExceptionReport version="1.1.1">
<ServiceException code="LayerNotDefined">
msWMSLoadGetMapParams(): WMS server error. Invalid layer(s) given in the LAYERS parameter. A layer might be disabled for this request. Check wms/ows_enable_request settings.
</ServiceException>
</ServiceExceptionReport>
Cache-Control: max-age=10, must-revalidate
Content-type: text/plain

An error occurred: msWMSLoadGetMapParams(): WMS server error. Invalid layer(s) given in the LAYERS parameter. A layer might be disabled for this request. Check wms/ows_enable_request settings.
  File "/var/www/tilecache/TileCache/Service.py", line 343, in cgiHandler
    format, image = service.dispatchRequest( params, path_info, req_method, host )
  File "/var/www/tilecache/TileCache/Service.py", line 208, in dispatchRequest
    return self.renderTile(tile, params.has_key('FORCE'))
  File "/var/www/tilecache/TileCache/Service.py", line 138, in renderTile
    data = layer.render(tile, force=force)
  File "/var/www/tilecache/TileCache/Layer.py", line 444, in render
    return self.renderTile(tile)
  File "/var/www/tilecache/TileCache/Layers/MapServer.py", line 50, in renderTile
    wms.loadOWSParameters(req)
  File "/usr/lib/python2.7/dist-packages/mapscript.py", line 1770, in loadOWSParameters
    def loadOWSParameters(self, *args): return _mapscript.mapObj_loadOWSParameters(self, *args)

Please help me in solving this issue.Thanks in advance.


texts: PostGIS - when I add a polygon, delete overlapping areas in other layers

texts:
I'm starting to work with PostGIS and I'd like to solve the following problem:

In my database, I have multiple land use class layers, for example: "residential", "industrial", "commercial".
The polygons in the layers don't overlap. Some areas are not covered by any layer, which means that the area is undeveloped land.

Now, let's say that I want to simulate a change of the landuse: a new "residential" area polygon is created instead of a part of an industrial area and surrounding undeveloped land (see picture).

I know that to insert the new polygon to the "residential" table, I must write a query like this:

INSERT INTO residential (name, the_geom)
   VALUES ('new development', ST_GeomFromText('POLYGON((-71.1776585052917 42.3902909739571,
   -71.1776820268866 42.3903701743239,
   -71.1776063012595 42.3903825660754,
   -71.1775826583081 42.3903033653531,
   -71.1776585052917 42.3902909739571))'))

Now, my question:
When I insert the new polygon, how do I make sure that all overlapping areas (whole polygons or parts of polygon) in all other layers are deleted?


texts: Most optimized way of extracting data from planet.osm

texts:
I want to cut extracts from planet.osm.bz2 .
I have found this command

bzcat downloaded.osm.bz2 | osmosis  --read-xml enableDateParsing=no file=-  --bounding-box top=49.5138 left=10.9351 bottom=49.3866 right=11.201 --write-xml file=- | bzip2 > extracted.osm.bz2

Only concern I have the disc space. planet.osm uncompressed is 400GB!

How much space would use this command from WIKI?


texts: I cannot able to popup the marker from the geojson in the leaflet

texts:
var json1 = {
  "type": "FeatureCollection",
  "features": [
    {
      "type": "Feature",
      "properties": {},
      "geometry": {
        "type": "Point",
        "coordinates": [
          83.94299268722534,
          28.218707821027166
        ]
      }
    }
  ]
}
var geolayer = L.geoJSON(json1).addTo(map);

geolayer.bindPopup(features.properties.name).openPopup()


texts:
What are the various formats/storage technologies which fall under the name Esri Geodatabase?
What are the most important differences between them?
and (briefly) what are the mainstream means to utilize them?

There are conversations all over the place about personal and file and enterprise and (?) geodatabases. These are each different beasts which need distinct handling, but there is a lot of confusion in the answers as to which means apply to which formats of geodatabases.

update: I should add that don't feel any single answer has to address the whole range of possibilities. It would be okay to say "the two single user gdb formats are personal and file gdb, they are appropriate for xxx, have these limits xxx, and the primary differences between them are xxx" etc.


texts: Trilateration using 3 x,y points and 3 distances, all in cm

texts:
I want to find out the position (x,y) of a given point, knowing the distance from 3 x,y points.

Below you can view the x,y coordinates and their given distances in cm.

B1x = 795
B1y = 332
B1d = 245.735630302

B2x = 795
B2y = 569
B2d = 327.946642001

B3x = 516
B3y = 569
B3d = 220.635445928

I have had a look at this question, but this is for earth, and I get the error:

z = math.sqrt(pow(DistA,2) - pow(x,2) - pow(y,2))
ValueError: math domain error

FIXED the above error with abs to allow negative numbers but getting the wrong coordinates.

I want my algorithm to calculate in a 2D environment, with the given values above.

Here is the map:

You can see the green X, is where I want to locate, from the three given beacons (x,y) locations.


texts:
I´m creating a map from a WMS with time dimension.

The user should have the opportunity to change the time, so that other features are displayed on the map.

I have a start time and the user should fill an input-field with his desired time.

A string will be created containing the time in the right format.

With that time another layer should be generated and just be laid over the existing layer.

function CreateLayer(myOtherTime)
        {

            var newlayer =  new ol.layer.Tile({source: new ol.source.TileWMS({
              url: 'http://webgis-kogge.auf.uni-rostock.de:8080/geoserver/sm/wms?service=WMS&version=1.1.0&request=GetMap&layers=sm:twit_corp_live_test&time='+myOtherTime

              })

            });

            map.addLayer(newlayer);
        }

Problem: I get the error message: map.addLayer is not a function.

The map is created in the following:

function CreateMap(center, myTime){

          new ol.Map({
            layers: [
              new ol.layer.Tile({source: new ol.source.OSM()}),
              new ol.layer.Tile({source: new ol.source.TileWMS({
              url: 'http://webgis-kogge.auf.uni-rostock.de:8080/geoserver/sm/wms?service=WMS&version=1.1.0&request=GetMap&layers=sm:twit_corp_live_test&time='+myTime
              })
              })
            ],
            view: new ol.View({
              center: center,
              zoom: 11.2
            }),
            target: 'map'
          });


texts:
Trying to use integrate a simple leaflet map within Bootstrap grid / column.  While I am trying to follow tutorials and templates, for some reason the
 id="map"
does not seem to be nesting within the column. The column only seems to be 1px in height.  In fact, when looking at Firefox Inspector, the full extent of the map window doesn't even appear to be registering as part of the body (apologies for any dubious terminology - very much starting out and probably running before I can walk).

My simple code is as follows:

Simple Mapping Page

Some witty lead in text...

Some test text

I also apply the following CSS:

#map {
width: 100%;
position: absolute;
height: 600px;
margin:0;
padding:0;
border: 1px solid black;
border-radius: 8px;
}

While the second column appears to the right of the map, there is no padding etc.


texts: OpenLayers2 - setting opacity to a webservice layer (if statement)

texts:
I've got an openlayers application which calls the layers from MSSQL via a web service.
The problem is that I want to set an opacity level for one of the layers (If layer name = X, set opacity level to Y).

How do I do that ?


texts:
I know is a silly question but I thought I should ask since you guys know more about postgresql. Sorry for being analytical, but I think its the only way to make my question clear.

I'm using postgresql 9.1 /postGIS 2.0 /openlayers/geoserver and apache

I've managed to get a feature's fid when clicked from a vector layer. I edited it using javascript and kept only the number, for example 64.

Then I set this as a value of a hidden html form. Using ajax I m trying to pass this value to a PHP value and perform a query. Openlayers, javascript and ajax work fine but php never gets the value. Here is my code

            function selected_feature(event){
var sf = event.feature.fid;
    var sfs = sf.split(".");
    //that's the number out of fid
    var sff = sfs[1];
    //feed it to the form
    document.getElementById("blahhh").value=sff;
    //alert it-works
    alert(document.getElementById("blahhh").value);
    //call ajax
    loadXMLDoc();
    }

    <!--the form-->
    <form id="hifo" name="hifo" method="post">
    <input type="hidden" id="blahhh" name="blahhh" />
    </form>

    //outside openlayers code
            function loadXMLDoc(){
        var xmlhttp2;
        alert("inside");//works
        var t2lFname=document.getElementById("blahhh").value;
        var parb="blahhh="+t2lFname;

        if (window.XMLHttpRequest)
          {// code for IE7+, Firefox, Chrome, Opera, Safari
          xmlhttp2=new XMLHttpRequest();
          }
        else
          {// code for IE6, IE5
          xmlhttp2=new ActiveXObject("Microsoft.XMLHTTP");
          }
        xmlhttp2.open("POST","map.php",true);
        xmlhttp2.setRequestHeader("Content-type", "application/x-www-form-urlencoded");
        xmlhttp2.send(parb);

        xmlhttp2.onreadystatechange=function()
          {
          if (xmlhttp2.readyState==4 && xmlhttp2.status==200)
            {
        alert("suces");//works
            }
          }

//and then php
$blah = $_POST['parb'];
echo $blah;
 $queryd='
  SELECT
   pins.p_name
  FROM
    pins
    WHERE
    pins.p_id='.$blah;
 //fetching to fill array
 $resultd=pg_query($conn, $queryd);
 $nord=pg_num_rows($resultd);
 $d=0;
  while($arrayd=pg_fetch_array($resultd)){
      $name[$d]=$arrayd['p_name'];
    $d++;
        echo $arrayd['p_name'];
                echo $arrayd[$d];
        echo $d;
    }

So basically use the fid as an id, via ajax to query the db. All this code is in the same file (named map.php).Problem is php's echos not working and query not working at all. I ve searched the web, I red tutorials and still no luck. Everything works except php. What should I do? Please advise

Thanks

PS I also tried OpenLayers.Request, found here and still no luck


texts: Adding different animated widget at same map in cartobuilder

texts:
I have 4 layers in Carto and I want to make animated maps for each of them but it seems I can only make one.

Is there any way that I can make more than one?


texts: Django client install postgresql/postgis on Ubuntu 12.04

texts:
Does anyone know what apt packages one should install for a postgres(8.4)/postgis(1.5) django (1.4) client on ubuntu 12.04?

All recipes seem to assume that I want to install the full postgres/postgis database server on the machine, but I have a external database server so I only want the packages that can connect to that.

For a full server I use:

apt-get install binutils gdal-bin libproj-dev postgresql-8.4-postgis postgresql-server-dev-8.4 python-psycopg2 python-setuptools


texts: Create image using raster photo layer with vector shapefile overlaid (landuse)

texts:
I am a beginner with QGIS.

I have a raster aerial photo and wish to overlay a vector shapefile (I have already created) showing polygons of land use.

Can this be done in QGIS, or do I need to use photoshop?


texts: Using pgRouting for 'multiple pedestrians' on an existing table of 'students'?

texts:
I'm trying to use the pgRouting sample for 'multiple pedestrians going to the same destination' to determine the 1) length and 2)geometry of students to their schools.

Essentially I'm trying to write functions that use pgRouting the way ST_Distance and ST_MakeLine accept geometry from students and schools which we use frequently.

Here is the logic we use to get students 'straight line' distance to a particular school (works for multiple schools too):

select

stu.studentid
, sch.schnum
, ST_Distance(stu.geom, sch.geom)

from dpsdata."Geo_Active_Student" as stu

, dpsdata."Schools_Current" as sch

where stu.studentid in ('750909', '818767', '768952')

and sch.schnum = '215'

The way ST_Distance accepts two geometries, and the where clause specifies which geometries to pass in (schnum, studentid's) is the same way I'm trying to get the pgRouting tools to work in a function.

I can get the 'single pedestrian' samples to work pointing to my existing student and school tables:

--one student to one school
SELECT * FROM pgr_dijkstra('
    SELECT gid AS id,
         source,
         target,
         length AS cost
        FROM osm.osm_ways',

     ---student
    (SELECT id FROM osm.osm_ways_vertices_pgr ORDER BY the_geom <->
     (select geom from dpsdata."Geo_Active_Student" as stu where stu.studentid = '750909') limit 1 )::integer,
---school
 (SELECT id FROM osm.osm_ways_vertices_pgr ORDER BY the_geom <->
 (select geom from dpsdata."Schools_Current" as sch where sch.schnum = '215') limit 1 )

   , directed := false
     );

However, I can't think of the logic that would generate the same in the 'multiple pedestrians' sample using the same method of passing in multiple geometry instances for each student:

(note the use of the in statement to select which students I want to run  'multiple pedestrian' against)

    --many students to one school
    SELECT * FROM pgr_dijkstra('
        SELECT gid AS id,
             source,
             target,
             length_m AS cost
            FROM osm.osm_ways',

--instead of passing in an array of ID's, get them from the source tables
        --ARRAY[84234, 38905, 107851], 149952

--student

        (SELECT id FROM osm.osm_ways_vertices_pgr ORDER BY the_geom <->
         (select geom from dpsdata."Geo_Active_Student" as stu where stu.studentid in ( '750939', '814767', '768922')) limit 1)

    --school
         , (SELECT id FROM osm.osm_ways_vertices_pgr ORDER BY the_geom <->
     (select geom from dpsdata."Schools_Current" as sch where sch.schnum = '215') limit 1 )

        , directed := false
        );

The error returned is: ERROR:  more than one row returned by a subquery used as an expression

Is there either an error in logic or SQL that I have missed?


texts: How to get latitude and longitude of Point plotted on the map using mapbox

texts:
How can I retrieve the longitude and latitude of a point plotted on map? I am testing one website in which I need to verify the longitude and latitude of a point plotted on the map.

Map has been developed using Mapbox.


texts: Display OSGB (EPSG: 27700) ArcGIS tile MapServer in OpenLayers 3

texts:
I am trying to add this ArcGIS OS Open Background map service in OpenLayers 3.7.0.

http://tiles.arcgis.com/tiles/qHLhLQrcvEnxjtPr/arcgis/rest/services/OS_Open_Background_2/MapServer

I set the projection and the resolutions for this service in my OL setup like so:

proj4.defs('EPSG:27700', '+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs');

ol.proj.get("EPSG:27700").setExtent([0, 0, 700000, 1300000]);

var layerOSGBtiles = new ol.layer.Tile({
    source: new ol.source.TileArcGISRest({
        url: 'http://tiles.arcgis.com/tiles/qHLhLQrcvEnxjtPr/arcgis/rest/services/OS_Open_Background_2/MapServer',
        projection: 'EPSG:27700'
    })
});

var view = new ol.View({
    projection: 'EPSG:27700',
    center: [300000, 500000],
    extent: [0, 0, 700000, 1300000],
    resolutions: [132291.9312505292, 66145.9656252646, 26458.386250105836, 19843.789687579378, 13229.193125052918, 6614.596562526459, 2645.8386250105837, 1322.9193125052918, 661.4596562526459, 264.5838625010584, 132.2919312505292, 66.1459656252646, 26.458386250105836, 19.843789687579378, 13.229193125052918, 6.614596562526459, 2.6458386250105836, 1.3229193125052918, 0.6614596562526459],
    zoom: 6
});

map = new ol.Map({
    target: 'map',
    controls: [],
    layers: [layerOSGBtiles],
    view: view
});

So far I have been unable to get any love from TileArcGISRest. Any help will be welcome.

Here is a JSFiddle which shows the right coordinates and a control geometry for the outline of the British Isles.


texts:
I am facing some difficulties configuring properly my geoext tree. I'll appreciate if anybody can help me.

When I try to open the "images" folder, the node stays rounding, trying to load something forever.

Here is my code:

var DNI2008_dec = new OpenLayers.Layer.WMS("2008_DNI_dec",
        "http://localhost:8080/geoserver/wms", {
            layers: "personal:2008_DNI_dec",
            transparent: true,
            format: format,
            tiled: true,
            tilesOrigin : map.maxExtent.left + ',' + map.maxExtent.bottom
        }, {
            buffer: 0,
            displayOutsideMaxExtent: true,
            ratio: 1,
            opacity: 1,
            visibility: false
        }
    );

map.addLayers([DNI2008_dec]);

    var LayerTree = new Ext.tree.TreePanel({
    title: 'Map Layers',
    renderTo: 'layers',
    loader: new Ext.tree.TreeLoader({
                    applyLoader: false
                }),
    root: {
        nodeType: "async",
        text: "Layers",
        leaf: false,
        expanded: true,
        children: [
            new Ext.tree.TreeNode({
                nodeType: "folder",
                text: "images",
                leaf: false,
                children: [
                       new GeoExt.tree.LayerNode({
                           text: "test",
                           nodetype: "gx_layer",
                           layerStore: mapPanel.layers,
                           layer: "2008_DNI_dec",
                           leaf: true
                       })
                ]
            })
        ]
    }
});

(CORRECTING A MIS-EDIT)


texts: how to convert shapefile in gpx route not gpx waypoints

texts:
When I do "save as gpx", it creates a gpx file with waypoints only, although I force route and track.
  Besides other non-gpx schema columns all relavant gpx columns exist based on a attribute table merge of a 'gpx to point shapefile' layer and 'custom point shapefile' layer.


texts:
I have a map where I have a KMLLayer that defines a geographic boundary and a CSVLayer that defines points within that boundary. I do not want any infoWindow information popping up for the KMLLayer and want to take some actions on the click event on the CSVLayer.

In building this page, I have found that I can no longer apply a click event to the CSVLayer once I have added the KMLLayer to the map.

Seems like it would not be too difficult... but I have played with multiple instances of this and keep coming back to the same place.

var map, kml, csv;
        require(["esri/map", "esri/layers/KMLLayer", "esri/layers/CSVLayer", "esri/renderers/SimpleRenderer", "esri/InfoTemplate", "esri/urlUtils", "esri/symbols/PictureMarkerSymbol", "esri/geometry/Point", "esri/SpatialReference", "esri/geometry/webMercatorUtils", "esri/geometry/screenUtils", "esri/geometry/ScreenPoint", "dojo/dom", "dojo/_base/array", "dojo/domReady!"], function(Map, KMLLayer, CSVLayer, SimpleRenderer, InfoTemplate, urlUtils, PictureMarkerSymbol, Point, SpatialReference, webMercatorUtils, screenUtils, ScreenPoint, dom, array) {
            map = new Map("map", {
                basemap : "hybrid",
                center : [-83.71, 36.05],
                zoom : 8
            });

            map.on("load", function() {
                map.infoWindow.set('anchor', 'bottom-right');
                map.infoWindow.resize(625, 425);
            });
            map.infoWindow.on("hide", centerMap);

            urlUtils.addProxyRule({
                proxyUrl : "/proxy/proxy.php",
                urlPrefix : "www.srh.noaa.gov"
            });

            var template = new InfoTemplate("${name}", "<iframe src=${fcstlink} class="fcst"></iframe>");
            var csvUrl = "http://www.srh.noaa.gov/images/rtimages/mrx/js_local/riverpoints.csv";
            csv = new CSVLayer(csvUrl, {
                copyright : "weather.gov",
                infoTemplate : template
            });
            var logo = new PictureMarkerSymbol("http://www.srh.noaa.gov/images/mrx/surf.png", 20, 20);
            var renderer = new SimpleRenderer(logo);
            csv.setRenderer(renderer);
            map.addLayer(csv);
            csv.on("click", mapFcstAction);

            var kmlUrl = "http://www.srh.noaa.gov/images/rtimages/mrx/kml/RiversLakes.kml";
            kml = new KMLLayer(kmlUrl);
            map.addLayer(kml);

            kml.on('load', function(lyr) {
                var layers = lyr.layer.getLayers();
                array.forEach(layers, function(l) {
                    l.setInfoTemplate(null);
                });
            });

            function centerMap() {
                alert("CENTER");
                map.graphics.clear();
                map.centerAt([-83.71, 36.05]);
            }

            function mapFcstAction(evt) {
                alert("FCST");
                map.infoWindow.show(evt.screenPoint);
                var screenX = evt.screenPoint.x;
                var screenY = evt.screenPoint.y;
                var newScreenX = screenX + 225;
                var newScreenY = screenY + 250;
                var newScreenCenter = screenUtils.toMapPoint(map.extent, 800, 550, new ScreenPoint(newScreenX, newScreenY));
                map.centerAt(newScreenCenter);
            }
        });


texts: OpenLayers transform EPSG 27700 coordinate to Google 900913

texts:
Community,

I'm currelty capturing my mouse click event and trying to transform the coordinates into Long Lat. I need to get at the Long Lat as strings so I can build a url and throw it at Google Street view. Current transform code only prints the BNG coords.

lonlat = map.getLonLatFromPixel(e.xy);

brit = new OpenLayers.Projection("EPSG:27700");
google = new OpenLayers.Projection("EPSG:900913");

lonlatclone = lonlat.clone()
lonlatclone.transform(brit, google)
//alert("You clicked near " + lonlatclone.lat + "N," + +lonlatclone.lon + "E");

Thanks!


texts: QGIS not saving layer in Lo 19 / EPSG:2048

texts:
I am reprojecting a shapefile from lat/long to Hartebeeshoek Lo19 using 'set layer crs', however when saving the shapefile ('export as') the output remains in lat/long?

Set CRS:

Save the reprojected shapefile:


texts: How to attach GP script tools, standard tools, and Python script to a toolbar?

texts:
I need to build a custom toolbar that can launch the following items:

A Python script tool - shouldn't be a problem, as you can already
attach these to a toolbar
Out of the box functions such as the Editor Toolbar maybe or a standard GP tool - again, should be do-able as you can already attach these to a toolbar
A Python script that launches a PyQT/wxPython/Tk GUI. This is the tricky one I think. Now, this GUI does not need to interact with ArcMap. My other custom tools will do all of that, this GUI is for the user to setup parameters that will get written to a XML config file. It will also do some arcpy calls in the background and slurp some data and values out of a geodatabase and push them to the same XML config file.

So the question is - how to do this as simply as possible? It appears I cannot attach a Python script (that is not encased in a tool) to a toolbar in ArcMap. I'm thinking I'm going to have to go the Add-In route and create a simple Toolbar with buttons to kick off each of the items listed above. VBA is not an option. Thoughts?


texts: Seeking Leaflet-compatible tile layers for Palestine?

texts:
Where can I find a Leaflet-compatible tileset that correctly displays Palestine with its native borders and Arabic-language placenames?

I've searched high and low and I've found all kinds of tilesets, but the only ones I've found display the Occupied Territories with post-occupation toponyms.

If the tileset were otherwise pretty similar to the standard OpenStreetMap tileset, that would be great.


texts: Identifying Primary Key field of Layer in ArcGIS Server using ArcGIS API for JavaScript?

texts:
I need to get the primary key column of the Particular layer added in the Map.
Usually 'FID' be found for a Shapefile Layer Source and OBJECTID  found for a SDE source layer.

Is there any way to Identify which one is the primary key columns so that I can access those values.

I need to get this using ArcGIS JS API Version.3.3 and I am using ArcGIS Server 10.


texts: I am trying to symbolize data as  bar charts in ArcMap. I have a feature class with two sets of data, and I want the bar charts next to each other on the map. The problem I am having is that the map will only display the bars if both sets of data are not Null. Is there a way to display all the bars with data? If I put 0 for the Nulls, then a flat bar appears, and I don't want that either.


texts: Merge a polygon from a layer with polygons from another layer

texts:
I have some subbasins as polygons. I know that an area is totally isolated so that it builds a separate subbasin (Polygon - dark purple - Fig.2).

I wanted to merge them in a way that the boundary of the isolated area (polygon) cuts the underlaying polygons.


texts: Setting output of spatial ETL tool in ModelBuilder?

texts:
I am using the data interoperability extension and building spatial ETL tools, but I can't seem to specify the output location of the ETL tool so it's impossible to use this dynamically.

Has anyone figured out how to do this?


texts:
I have pandas dataframew ith two fields: ID and geometry (geometry is already in WKT). I want to change this dataframe into geodataframe so I can save it as shapefile.
I have tried to change it into geodataframe in the followong way:

covex = gdf(covex, crs="EPSG:4326", geometry='geometry')
covex.to_file("covex.shp")

that didn't work and I got this error message:

TypeError: 'GeoDataFrame' object is not callable

covex = gdf(covex, crs="EPSG:4326", geometry=covex.iloc[:,1])
covex.to_file("covex.shp")

but I got the same error messege.
What do I do wrong? I thought I should pass the geoemtry column there but seems like it should be somethng else. how can I change this into geodataframe and then to shapefile?


texts: Plotting graphs showing precipitation of area on map?

texts: I have precipitation data of an area and I need to create a bar graph or something on the map of the area showing the precipitation. Is there a way I can do it?? I am a beginner so I really don't have an idea about it?


texts: Auto Generate kml files along the roads for more than 8 points

texts:
What I am trying to do: Auto Generate kml files for bus routes.

Input: Bus stop list for each route with name and lat/long (every list has more than >30 bus stops.)

Output: kml file of the route passing through the bus stops and along the roads.

Research so far:

Found this tool from geocodezip.com to plot routes with more than 8 points.
Found this excellent tool to write kml files from manually clicked points. It works beautifully for less than 8 way-points.

Possible direction: Ideally I want to combine the functionality of both the the above tools and have a tool that writes kml code for routes without restriction on no. of waypoints.


texts: QGIS reading Sentinel 2 image .zip gives error?

texts:
I have a very primitive question.
The sentinel-2 .zip tiles downloaded from Copernicus Open access Hub are being recognized and opened with STEP SNAP toolbox (obviously)! However this is not necessarily true with QGIS 3.4. When I try to open the same .zip file with QGIS I receive an error : Raster Layer Provider is not valid

I wonder if it is something wrong with the GDAL attached to QGIS?


texts: Why does QgsFeature has no function to retrieve the corresponding layer?

texts:
Every QgsFeature should have one corresponding parent layer, right? But there is no function to retrieve this layer?

Do I overlook something here?

@qgsfunction(args=0, group='Custom')
def getLayer(feature,parent):
    layer = feature.layer()
    return layer

The Expression function sample doesn't make much sence, but it shows my intention.


texts:
I am using JTS in my application, and I want to transform the geometry from longitude and latitude to pixel.

And the algorithm for this translation has been providered by an interface:

public double[] fromLngLatToPixel(double lat, double lng);

Now I have to read the coordinates for the source geometry ,transform the coordinates and then create a new geometry.

However, I found that JTS has a transfrom method like this:

public static Geometry transform(Geometry geom, MathTransform transform);

But I do not know to to implement the MathTransform according to fromLngLatToPixel.

Any suggestion?


texts: Shapely adds more vertices to a polygon for smother snapping

texts:
I'm using shapely to snap together some polygons. The problem is that snapping in shapley only supports vertices to vertices snapping (i.e. a vertices from a polygon can only ever snapped to the vertices of another polygon) and not like other snapping tools where it is possible to snap to the whole geometry. This leads to undesirable effects as seen below.
First the polygons unsnapped:

and then snapped:

Since snapping to other geometric features (like for example in postgis) in shapely is not possible I thought if it is possible to introduce simply more vertices to the polygon so that snapping has more snapping targets. Is that possible and if so how?
my code currently looks like this:
input_gdf = input_gdf.copy()
for index, row in input_gdf.iterrows():
    # first we select polygons that are close to our
    # row but not the geometry of the row itself
    temp_gdf = input_gdf[input_gdf.index != index]
    # this selects everything close by within buffer distance
    temp_gdf = temp_gdf[temp_gdf.geometry.intersects(row['geometry'].buffer(0.0001))]#
    # the closeby polygons get put together into one big multipolygon
    other_polys = MultiPolygon(temp_gdf['geometry'].to_list())

    # then we snap the geometry of row to the just created multipolygon
    # the snaped geometry overwrites the exisiting geometry in input_gdf
    input_gdf.loc[index,['geometry']]= snap(row['geometry'],other_polys,0.00003)

as you can see I iterate over the Geodataframe and select for every polygon the close polygons and then snap the current polygon to the vertices of the other polygons.


texts:
Direct comparison between software may not be the best thing to do but as I can't arrive at a conclusion thought of asking here.

Is there any gvSIG user in this GIS community?

Looking at the screenshots gvSIG seems to have lot of features. I particularly liked the single window interface when compared to GRASS. So I wanted to know how gvSIG would suit for image processing and GIS operations in general.

How are the map composing features in gvSIG?

I would also like to know which fares better amongst gvSIG or GRASS for image classification task - particularly digitizing training sites which is improving in GRASS these days(even then working with QGIS is easy for the purpose), ability to handle large rasters and vectors in the process, ease of map composing etc.


texts:
I am attempting to write a small bash script that will mass process a large number of raster .tif files. It should do the following:

Clip to a smaller extent
Reproject to UTM

This is what I have, so far:

ODIR="utm19"

if [ ! -d "$ODIR" ]; then
  mkdir $ODIR
fi

for i in `ls *.tif`;do
  gdalwarp -s_srs EPSG:4326 -t_srs EPSG:32718 \
  -te -80 -55 -15 -60 -te_srs EPSG:4326  \
  -r bilinear -wm 6048 -multi -srcnodata "-9999" \
  -srcnodata "-9999"  $i $ODIR/${i}.tif
done

However, I am a bit puzzled by the different opinions existing at:
http://www.gdal.org/gdalwarp.html
vs.
§4.5 in http://www.geos.ed.ac.uk/~smudd/TopoTutorials/html/tutorial_raster_conversion.html

So, what is the SRS units we should use for the -te option in gdalwarp


texts: OpenLayers 4 select a feature on click in a list

texts:
I have a project where I load a GeoJson map and when I click on it I launch a function, till here everything is Ok, I used the same code in the example.

https://openlayers.org/en/latest/examples/select-features.html

I just add a few extra bit to get the name of the ward

  var areasel = [];

  var select = null;  // ref to currently selected interaction
  var selectSingleClick = new ol.interaction.Select();

  var selectClick = new ol.interaction.Select({
    condition: ol.events.condition.click
  });

  // select interaction working on "pointermove"
  var selectPointerMove = new ol.interaction.Select({
    condition: ol.events.condition.pointerMove
  });

  var selectAltClick = new ol.interaction.Select({
    condition: function(mapBrowserEvent) {
      return ol.events.condition.click(mapBrowserEvent) &&
          ol.events.condition.altKeyOnly(mapBrowserEvent);
    }
  });

  var selectElement = document.getElementById('type');

  var changeInteraction = function() {

    if (select !== null) {
      map.removeInteraction(select);
    }

    var value = selectElement.value;

    if (value == 'singleclick') {
      select = selectSingleClick;
    } else if (value == 'click') {
      select = selectClick;
    } else if (value == 'pointermove') {
      select = selectPointerMove;
    } else if (value == 'altclick') {
      select = selectAltClick;
    } else {
      select = null;
    }

    if (select !== null) {
      map.addInteraction(select);
      select.on('select', function(e) {
        document.getElementById('status').innerHTML = '&nbsp;' +
            e.target.getFeatures().getLength() +
            ' selected features (last operation selected ' + e.selected.length +
            ' and deselected ' + e.deselected.length + ' features)';

        if (e.selected.length === 0 && e.deselected.length > 0) {
          for (x = 0; x < e.deselected.length; x++) {
            var int = areasel.indexOf(e.deselected[x].getProperties().Name);
            if(int != -1) {
             areasel.splice(int);
            }
          }
        }

        if (e.selected.length > 0 ) {

        for (i = 0; i < e.selected.length; i++) {
          if (e.deselected.length !== 0 ) {
            for (x = 0; x < e.deselected.length; x++) {
              var int = areasel.indexOf(e.deselected[x].getProperties().Name);
              if(int != -1) {
               areasel.splice(int);
              }
            }
          }
          areasel.push(e.selected[i].getProperties().Name);
        }
      }
          $('#wardsInput').val(areasel);

         if ($('#wardsInput').val() !== "") {
          update();
        } else {
          mainView();
        }

      });
    }
  };
  selectElement.onchange = changeInteraction;
  changeInteraction();

in the website I also create a list(it's a table) with the wards name displayed in the map, I manage to get the name on click and launch the same function I launch when I click on the wards in the map.

$('.tdclick').on('click', function() {
  $('#wardsInput').val($(this).find('td').first().text());
  update();
});

How can I select the ward in the map when I click on the name of that wards in the table?

Here a screenshot of how it looks


texts: Different popup for every polygon (leaflet)

texts:
I'm making a map (you can see here)  with different markers and layers to group these markers. I added polygons in another layer (Poligono layer), and I want to know if exist a way to add different information on each polygons that after can be show on a popup.

I can add different information on each marker using this

        L.marker ([21.017697, -89.641353], {icon:hospitalIcon}).bindPopup('Info info info').addTo(hospitales),

Exist a way to add information for each polygon using something like .bindPopup for the polygons I am using?

L.polygon ([
        [20.6894, -88.2017],
        [20.572997, -88.730264],
        [20.3457, -88.1677],
    ]).addTo(polygon),

Something like the image. When you click, show information, I need to show different information for each polygon.

and you can see the complete code


texts: Topology Checker - Polyline layer, How to detect overshoot/undershoot errors?

texts:
I have created a river network as polyline shapefile for training purposes. No I wanted to integrate undershoots/overshoots that needs to be detected and corrected from the participants.

I have used the snapping option to generate the polyline shapefile. I used as rule within the Topology Checker Plugin: "must not have dangles". The QGIS User Guide is saying (QGIS User Guide) that this rule is "show the overshoots in the line layer". Actually all end- and strating points are detected as error, even if they share a vertex and are connected to each other.
I think dangles are per definition something else. Have a look to the definition of GIS-Lounge, dangles are: "Dangles or dangling nodes are lines that are not connected but should be.  With dangling nodes, gaps occur in the linework where the two lines should be connected.  Dangling nodes also occur when a digitized polygon doesn’t connect back to itself, leaving a gap where the two end nodes should have connected, creating what is called an open polygon."

I am right, that Topology Checker Plugin is not able to detect undershoots/overshoots? If yes, is there any other possibility how to easily to it within the QGIS environment?


texts: A software architecture for viewing/editing roads and routes

texts:
I have a routing software module which uses PostGIS/pgRouting to calculate routes.

The road network is generated from an .osm file with osm2pgrouting tool.

I need to enhance the functionality of the software to meet the following requirements:

The user should be able to edit the roads network in a graphical editor (add/edit/delete roundabouts, lanes, road signs, highway tags, etc.);
The edits should be immediately updated in the geometry columns - basically, if there is a new road in the OSM map, we have a new line in the respective PostGIS column. pgRouting will use the updated roads network);
The user should be able to view the underlying OSM area and the edited roads network as two layers in the graphical editor.

Note:

I tend to think, that JOSM is a very convenient tool to edit roads, but it can not save_to/query a PostGIS-enabled database.

On the other side, there is QGIS, OpenJUMP and uDig, but the support for OSM editing is far less mature than in JOSM.

I thought about having a bash script which runs osm2pgrouting each time the .osm file is updated. This way I can have the .osm file edited via JOSM and the roads which are inside this .osm stored in a PostGIS-enabled database and viewed with QGIS (the most recent roads network would be accessible by QGIS).

But:

I would have to force the user to use two graphical editors;
and QGIS does not render OSM data in a very attractive way.

How can I meet the requirements 1, 2 and 3 in a more elegant way?


texts:
In my project I needed to convert TIFF file to SRTM data and came up with following warning/error.

gdal_translate -of SRTMHGT Terrain.tif Terrain.hgt

Output:

    Warning 1: The source projection coordinate system is PROJCS["RGF93 / CC48",GEOGCS["RGF93",DATUM["Reseau_Geodesique_Francais_1993",SPHEROID["GRS 1980",637813
    7,298.257222101,AUTHORITY["EPSG","7019"]],TOWGS84[0,0,0,0,0,0,0],AUTHORITY["EPSG","6171"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.017
    4532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4171"]],PROJECTION["Lambert_Conformal_Conic_2SP"],PARAMETER["standard_parallel_1",47.25],PARAMETER[
    "standard_parallel_2",48.75],PARAMETER["latitude_of_origin",48],PARAMETER["central_meridian",3],PARAMETER["false_easting",1700000],PARAMETER["false_northing"
    ,7200000],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["X",EAST],AXIS["Y",NORTH],AUTHORITY["EPSG","3948"]]. Only WGS 84 is supported.
    The SRTMHGT driver will generate a file as if the source was WGS 84 projection coordinate system.

    Warning 1: The corner coordinates of the source are not properly aligned on plain latitude/longitude boundaries.

    ERROR 1: Image dimensions should be 1201x1201, 3601x3601 or 1801x3601.

Here is the gdalinfo of the Terrain.tif file.

Driver: GTiff/GeoTIFF
Files: Terrain.tif
Size is 375, 375
Coordinate System is:
PROJCS["RGF93 / CC48",
    GEOGCS["RGF93",
        DATUM["Reseau_Geodesique_Francais_1993",
            SPHEROID["GRS 1980",6378137,298.257222101,
                AUTHORITY["EPSG","7019"]],
            TOWGS84[0,0,0,0,0,0,0],
            AUTHORITY["EPSG","6171"]],
        PRIMEM["Greenwich",0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.0174532925199433,
            AUTHORITY["EPSG","9122"]],
        AUTHORITY["EPSG","4171"]],
    PROJECTION["Lambert_Conformal_Conic_2SP"],
    PARAMETER["standard_parallel_1",47.25],
    PARAMETER["standard_parallel_2",48.75],
    PARAMETER["latitude_of_origin",48],
    PARAMETER["central_meridian",3],
    PARAMETER["false_easting",1700000],
    PARAMETER["false_northing",7200000],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AXIS["X",EAST],
    AXIS["Y",NORTH],
    AUTHORITY["EPSG","3948"]]
Origin = (1350983.998099997639656,7223238.974300000816584)
Pixel Size = (2.000000000000000,-2.000000000000000)
Metadata:
  AREA_OR_POINT=Point
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  ( 1350983.998, 7223238.974) (  1d41'25.74"W, 48d 6'48.59"N)
Lower Left  ( 1350983.998, 7222488.974) (  1d41'23.54"W, 48d 6'24.35"N)
Upper Right ( 1351733.998, 7223238.974) (  1d40'49.55"W, 48d 6'50.06"N)
Lower Right ( 1351733.998, 7222488.974) (  1d40'47.35"W, 48d 6'25.82"N)
Center      ( 1351358.998, 7222863.974) (  1d41' 6.54"W, 48d 6'37.21"N)
Band 1 Block=375x5 Type=Float32, ColorInterp=Gray
  NoData Value=0


texts: ArcGIS (ArcObjects) Windows 7 - 64 bit -support

texts:
I have been developing with ArcGIS Cross Platform C++ SDK (with Qt 4.8.2) on 32 bit Windows XP SP3 (development machine). I used to use MapTool Control within ArcObjects on the development machine. Moreover I compile my application with Visual Studio 2010 in both release and debug mode. I installed Engine Runtime on Windows 7 (both 32 bit and 64 bit) on the test machines, however my application simply crashes on both machines whenever

new QAxCtl(AoPROGID_MapControl, &mapctl, Map Control");

is called in my application. However, my application works on another machine running Windows XP SP3 which makes me believe that the cause of the problem is Windows 7 (or something I forgot about Windows 7). I am not familiar with ArcGIS and even I am not sure if I installed everything necessary for the test machines (I just installed just engine runtime and use a valid license). I wonder if using QAxCtl for controls is possible on Windows 7?


texts: How to position custom button controls on Google Maps API?

texts:
I've this simple map built on Google Maps APi v3 ....

    <title>Simple Map</title>
    <meta name="viewport" content="initial-scale=1.0">
    <meta charset="utf-8">

The map result is the follow

In my code the styles for the two buttons, marker and trash, are set in a way that shoud be aside the "Satellite" button in this manner

Why my style settings are not kept?

Suggestions / examples?


texts: Shadow impact of building on certain area

texts:
I have done some research but was not able to find a solution to my task!

What do I want:

Pretty simple - I have a building and I want to know how its shadow is influencing a certain site at different days of a year.

Questions asked:

--> How long (hours) is the building shading the site on a specific day?

--> When starts/ends (daytime) the impact by the buildings shadow on the site?

Whats the best way to do this? Can it be done using QGIS/GRASS? Or is ESRIs City Engine the only feasible solution?


texts: Is qualityMosaic the same as MVC (max value composite)?

texts:
I want to make an annual Landsat 5 image for area A in Google Earth Engine using MVC. The growing season in this area is from April to August.

I get data from April to October.
I masked the clouds.
I calculated the NDVI value.
I used qualityMosaic.
I extracted the NDVI band and named it as the Annual Data.

The problem is, the Annual NDVI chart I got is very strange.

From literature, and according to my knowledge, the NDVI peak in area A is in June and July. I assume that I made some mistakes in my data process, but I don't know where.


texts:
I have a shapefile that contains a road network. A while ago I could easily see the whole of my country.

There is no filter and every and the table contains all rows.

When displaying the file in Alteryx, it shows the complete layer.

What could be wrong with QGIS?


texts: Replacing multiples words using field calculator in attribute table python

texts:
I need to replace multiple words with one word in a field at one time in one script.

Example in my LettersField I need to replace all the "Aa" and all the "Bb" and all the "Cc" with "Rr".

I am aware of the replace fxn and I cannot use the find and replace tool - which would be MUCH easier but I have to create a python script.

How can I do this - what would the code block be and the 'LettersField =' be?


texts: Displaying properties from GeoJSON file in choropleth map using Leaflet.js?

texts:
I'm using the Leaflet.js library to generate a choropleth map of zip codes colored by some property values stored within an externally loaded GeoJSON file. I'm following a tutorial here but have run into some trouble, since the tutorial shows how to access the data through assigning it to a variable within the code. I'd like to proceed without explicitly assigning a variable (i.e. using the AJAX function I currently have in place to load the data).

I've been able to display the raw data (viewable here, source code available upon "View Source" in your browser) on the map, but have been unsuccessful in displaying the properties I'm interested in showing. I suspect I need to add some lines within the function that is loading the data, but can't understand how to do that.


texts: Creating a 'brick wall' style for polygons

texts:
I am trying to create a brick wall style for a vector polygon in QGIS 2 and 3. I can render a svg drawing of a brick wall with SVG fill (image below). But I do not want to use SVG since the style is part of a free geological map, where the SVG file will complicate the distribution and rendering of the data for the public. Saving the SVG through a default QML forces relative paths to the SVG file in the QML file and the SVG file must be in the QGIS system SVG folder.

Anyone created a polygon vector style like:

That does not need a external SVG file for rendering.


texts:
I'm trying to read ESRI personal geodabase.mdb through OGR. I have my driver "PGeo" but when I write

  ogrinfo INFLOR090812.mdb

I have this bug:

ERROR 1: Unable to initialize ODBC connection to DSN for DRIVER=Microsoft Access
 Driver (*.mdb);DBQ=INFLOR090812.mdb,
[Microsoft][ODBC Microsoft Access Driver] Could not find file '(unknown)'.
ERROR 1: Unable to initialize ODBC connection to DSN for DRIVER=Microsoft Access
 Driver (*.mdb);DBQ=INFLOR090812.mdb,
[Microsoft][ODBC Microsoft Access Driver] Could not find file '(unknown)'.

Someone could connect to a database.mdb to read the shapefiles in MapServer?