texts: Leaflet only one placemark visible at exact location, whereas is more than 1 there

texts:
I have a problem, because sometimes on my GeoJSON layer plotted into the Leaflet map more than 1 placemarks has been assigned to one exact location (by postcode).

The problem is, that only one of them is visible, whereas there is at least 3-4 another ones.

For example in Google Earth the problem is solved like this: when you hover one placemark, then you can see another ones hidden under the same location, as per in the picture below:

In leaflet this situation not occur and I see only 1 placemark only instead.

I was searching some stuff regarding to this problem and found a following examples:

http://jawj.github.io/OverlappingMarkerSpiderfier-Leaflet/demo.html

https://github.com/jawj/OverlappingMarkerSpiderfier-Leaflet

http://ivansanchez.github.io/Leaflet.Marker.Stack/demos/color_ramps.html

https://stackoverflow.com/questions/50813139/leaflet-map-with-multiple-marker-on-exact-same-position

but they don't really solve this problem

In turn, the links below:

https://stackoverflow.com/questions/22168558/multiple-marker-on-exact-same-position-on-a-leaflet-map/22173475

they are focused on the markercluster, which I don't need.

My situation looks like this:

where as you can see (after the onclick marker trigger) you see 2 jobs on the same location, whereas only 1 appear (always the same in addition).

Could you give me some hint, what should I tweak in the code, to make it pretty much the same, as it's in Google Earth?


texts: Can GDAL with "Kml Super Overlay" output a georeferenced geoTiff from a kmz?

texts:
I've got a kmz containing an unrefrenced jpg and a boundary box in kml.

ogrinfo is happy to report (filenames redacted):

INFO: Open of `Sim[...].kmz'
  using driver `LIBKML' successful.
1: Sim[...]

but, gdalinfo reports:

  ERROR 4: `Sim[...].kmz' not recognised as a supported file format.

I've run these info commands with the generic gdal/ogr tools on 14.04 and a specially compiled 1.11.1 using this answer.

Is there a way, using some approximation of these tools (barring manual unzipping and georeferencing of the jpg), to export a geoTiff containing appropriately referenced raster data from this data?


texts:
i have buildings data from osm and ward boundaries. i need to find which building lies in which ward and add that information in the building data. I prefer use osm output format.

INPUT:

  <way id='242474181' action='modify' timestamp='2013-10-17T11:42:37Z' uid='1700096' user='GautamPratik' visible='true' version='1' changeset='18401787'>
    <nd ref='2499343546' />
    <nd ref='2499343543' />
    <nd ref='2499343539' />
    <nd ref='2499343541' />
    <tag k='building' v='yes' />
  </way>

OUTPUT:

  <way id='242474181' action='modify' timestamp='2013-10-17T11:42:37Z' uid='1700096' user='GautamPratik' visible='true' version='1' changeset='18401787'>
    <nd ref='2499343546' />
    <nd ref='2499343543' />
    <nd ref='2499343539' />
    <nd ref='2499343541' />
    <tag k='building' v='yes' />
    <tag k='district' v='kathmandu' />
    <tag k='vdc' v='kathmandu metropolitan city' />
    <tag k='ward' v='4' />
  </way>


texts: Creating and re-using map template in ArcMap?

texts:
I am using ArcMap 10.1 and working on a project which requires 5 different analyses on the same data. I created five map documents (.mxd) and completed my analyses. Now, the problem is that I want to export my results in map format but I want them to be uniform in display i.e. having a same map template. Is there a way to save one of the maps as a template and then use it for creating other 4 maps?

I tried saving one document in 'Map Templates' folder in ArcGis folder but a message box appears which says that I don't have permission to save in that location and I must contact administrator for that, even though I am the only user and administrator of the computer. I don't want to manually add all map elements for 5 different maps, thats too tedious!


texts: Fixing raster with wrong Eastings using QGIS and GDAL?

texts:
I get an error when trying to warp a tif file in QGIS:

  ERROR 1: Too many points (35721 out of 35721) failed to transform,
  unable to compute output bounds.

I am using QGIS 1.7.0 Wroclaw with gdal 1.8. I am trying to reproject from ETRS89 to HD72 EOV. There was an error in the data, with the UTM zone added to the eastings. That is, the x coordinates read 34327213.050 when they should read 327213.050 in UTM zone 34N. I have tried changing the .tfw file but this did not help. Then I translated the tif to a .vrt and edited the vrt, removing the rogue 34, and fixing the extent. This also did not help - I still get the same error when I try to warp. In layer properties, this raster now has the same origin, extent and SRS as other similar rasters, but only this one gives the error.

How do I work this out?

A similar question (Setting projection of tif file in QGIS?) was asked 4 months ago, but it was not the exact same problem, and the OP did not explain how the problem was resolved.


texts: Error when using r.reclass in QGIS - GRASS

texts:
I am trying to run the r.reclass tool but I got the following error:

  Traceback (most recent call last):

  File "C:/OSGEO4~1/apps/qgis/./python/plugins\processing\algs\grass7\Grass7Algorithm.py", line 431, in processAlgorithm
  Grass7Utils.executeGrass(self.commands, feedback, self.outputCommands)

  File "C:/OSGEO4~1/apps/qgis/./python/plugins\processing\algs\grass7\Grass7Utils.py", line 389, in executeGrass
  for line in iter(proc.stdout.readline, ''):

  File "C:\OSGEO4~1\apps\Python37\lib\encodings\cp1253.py", line 23, in decode
  return codecs.charmap_decode(input,self.errors,decoding_table)[0]
  UnicodeDecodeError: 'charmap' codec can't decode byte 0x9a in position 25: character maps to

I changed the path in Windows but still I get the above error.

Do you know something about it and how can I fix it?

I am using the last version of QGIS 'A Coruna' 3.10.2.


texts: Network analysis, fastest path in QGIS 3.2

texts: I used this plugin and got results in the cost field, but I do not know in which units the cost is displayed (€, $ ...). Whether it's money or something else?


texts: How can I convert graphic text to attribute field within polygon that encloses it

texts:
I was able to convert graphic polygons to GDB. Similarly, I can convert the graphic text to annotation within same GDB. But I need to have that text as an attribute or field within the polygon table. How can I do this?

The graphic text was drawn within the polygons.


texts: Merge/Combine two adjacent contour lines with ArcPy?

texts:
I am a Python beginner; I found a script in Append or merge consecutive line segments (contours) as one single line of the same elevation value

I am trying to use it to merge 2 shapefiles (1 ft -contours):
N1_Terrain_SurfaceContour ( 0 -12) ft elevetions and
N1_Terrain_SurfaceContour (0-11)

These are located at my workspace:
C:\Users\efa8\Documents\AA_NRMD_Bvrd\CONTOUR PJ\Contour.mdb

Each one has the same fields: OBJECTID, Shape*, Contour and Shape_Lenght. I used "Contour" instead "ELEV_FT" in the script (as you suggested).

I got the following ERROR (Parameters are not valid):

Traceback (most recent call last):
  File "C:\Users\efa8\Documents\AA_NRMD_Bvrd\CONTOUR PJ\Contourlines merge BREVARD.py", line 35, in <module>
    arcpy.Merge_management(contourList, "contour_merge.shp")
  File "C:\Program Files (x86)\ArcGIS\Desktop10.1\arcpy\arcpy\management.py", line 3788, in Merge
    raise e
ExecuteError: Failed to execute. Parameters are not valid.
ERROR 000735: Input Datasets: Value is required
Failed to execute (Merge).

This is the script that I adapted:

Contour_Merge_And_Dissolve.py - CONTOUR PJ BREVARD COUNTY

# Import modules
import arcpy
arcpy.env.workspace = "C:\Users\efa8\Documents\AA_NRMD_Bvrd\CONTOUR PJ\Contour.mdb"
from arcpy import env

# Define functions
def min_max(featureClass, fieldName):
    listItems = []
    rows = arcpy.SearchCursor(featureClass)
    for row in rows:
        listItems.append(row.getValue(fieldName))
    del rows
    listItems.sort()
    listMin = listItems[0]
    listMax = listItems[-1]
    return listMin, listMax

# Arguments
contourList = arcpy.GetParameterAsText(0) #input shapefiles, ElevationContours
outFolder = arcpy.GetParameterAsText(1) #output folder

# set workspace and overwrite option
arcpy.env.workspace = outFolder
arcpy.env.overwriteOutput = True
arcpy.env.outputZFlag = "Disabled"
arcpy.env.outputMFlag = "Disabled"

# local variables
dissolveFields = ["Contour"] # change field name here

# Process: Merge
arcpy.AddMessage("Merging contours...")
arcpy.Merge_management(contourList, "contour_merge.shp")

# Get Min/Max values of Contour field
ft_min, ft_max = min_max("contour_merge.shp", "Contour") # change field name here
ft_range = ft_max - ft_min

arcpy.AddMessage("Contour MIN: " + str(ft_min))
arcpy.AddMessage("Contour MAX: " + str(ft_max))

# Set progressor
arcpy.SetProgressor("step", "Breaking apart contours...", 0, ft_range, 1)
arcpy.AddMessage("Breaking apart contours...")

# Break apart and dissolve for each Contour value
mergeList = []
counter = ft_min
while counter <= ft_max:
    arcpy.SetProgressorLabel("Breaking out " + str(counter) + "ft contours...")
    arcpy.AddMessage("Attempting " + str(counter) + "ft contours...")
    counter_str = str(counter)
    counter_str = counter_str.replace("-", "neg")
    thisLyr = "elev_" + counter_str + "_lyr"
    thisDis = "elev_" + counter_str + "_dissolve.shp"
    where_clause = '"Contour" = ' + str(counter) # change field name here
    # Make Layer
    arcpy.management.MakeFeatureLayer("contour_merge.shp", thisLyr)
    # Select by attributes
    arcpy.SelectLayerByAttribute_management(thisLyr, "NEW_SELECTION", where_clause)
    # Dissolve Layer to new shape
    arcpy.Dissolve_management(thisLyr, thisDis, dissolveFields, "", "SINGLE_PART", "DISSOLVE_LINES")
    # Add to mergeList
    mergeList.append(thisDis)
    counter += 1
    arcpy.SetProgressorPosition()

arcpy.ResetProgressor()

# Merge contours
arcpy.AddMessage("Merging dissolved layers...")
arcpy.Merge_management(mergeList, "contour_dissolve.shp")

# Multipart to singlepart
arcpy.AddMessage("Converting multipart to singlepart...")
arcpy.MultipartToSinglepart_management("contour_dissolve.shp", "contour_singlepart.shp")

# Clean up
for shape in mergeList:
    arcpy.Delete_management(shape)


texts: List group of line segments that are not connected to any other lines

texts:
In the image below, the system of lines includes 6 segments. This connectivity gives the appearance of an "island" as the internal segments are connected to each other and are not to the outside world.

Given that the geometry is correct and I'm not looking to alter it, in a large line network, how can I identify those groups of lines or "sub networks" that are not connected to the main network?


texts: ST_Contains(geometry A, geometry B) doesn't return polygons of B sharing external border with A

texts:
I have a large shapefile with 400k census sections in a table called census_poly. I want to return all the sections from census_poly contained in a district. Visually, the external border of districts matches perfectly with the borders of the census sections.

Still when I run this query

SELECT "census_poly".geom
   FROM "census_poly"
   JOIN "district_poly"
   ON ST_Contains("district_poly".geom, "census_poly".geom)

I get this

which excludes all the bordering census districts (FALSE NEGATIVE) such as these

I wonder why they are excluded and how to tweak my query to include them.

I could use the centroids of the census sections but that would not indicate whether a section is truly over two districts (TRUE NEGATIVE). Should I define a small buffer around my districts?


texts: CRS of QGIS CAD tool and measure line tool

texts:
Does any one know what CRS the QGIS digitizing toolbox (CAD tool) and the measure line tool use?

Because I couldn't find any way to manually set the CRS for the tools in QGIS 3, my assumption is that it will use the CRS set for the QGIS project. But I am seeing inconsistent values of the two tools when setting the project CRS to EPSG 3857 (Pseudo-Mercator). While the CAD tool shows a distance of 5m, the measure line tool only measures a distance of 3.3m.


texts: Geoserver: how to filter a layer related to another layer?

texts:
I have the following problem:

I have 2 vector layers from PostGIS datastores: "firepoints" and "threats". The corresponding tables "firepoints" and "threats" have IDs and other attributes. They are related via a foreign key "id_firepoint" in "threats" table. One threat can have multiple firepoints, but one firepoint can have only one threat.

What I need to do is to filter WMS layers in Openlayers using CQL filters: first I want to filter the "firepoints" layer using some attribute ("satellite", for example) and then apply the derived IDs to filter "threats" layer using the "id_firepoint" attribute.

The problem is that I don't know how to derive the list of IDs of filtered "firepoints".

I have tried using the Cross-layer filtering with "querylayer" module, but I don't think that I understood how to use it correctly.

This is the filter that I try to apply to the "GSS:threats" layer:

  id_firepoint IN queryCollection('GSS:firepoints','id','satellite =
  AQUA')

Is there any way to filter a layer related to another layer using dynamically changing attributes?


texts: How can I convert PostGIS geography to geometry and split polygons that cross the antimeridian

texts: I have a PostGIS 2.2 table that stores polygons as geography type.  The table contains several polygons that cross the antimeridian.  I would like to split these polygons at the antimeridian and convert the table's topology from geography to geometry.


texts: Why aren't Mapbox vector tiles rendered in the browser?

texts:
Would it be technically possible to have Mapbox's compact vector tile format sent to the client (instead of PNG tiles), then rendered inside the browser with JS code?

Second question: when using TM2/Mapbox Studio, why are the tiles rendered remotely? Couldn't Mapnik be directly used locally to render the tiles? (still using remote data sources like Mapbox Streets)


texts: What angular units does the QGIS raster calculator use: Degrees, Radians, Gradians?

texts:
I am working with QGIs and am converting a slope raster with values in degrees to a slope raster with values in percent using the raster calculator and the following formula:

(arctan(slope_raster))*100

My original slope raster in degrees has no values outside of the 0-90 degree range which is expected but the percent slope raster it generates gives me a low of -62000 which is pretty silly.

I guess my question is what angular units does the QGIS raster calculator use? Degrees? Radians? Gradians?


texts: How to convert this geometry location format to a qgis usable format?

texts: I am a beginner with regard to dealing with gis information, so please bear with me. I would be obliged if someone could help. I have data on asssembly constituencies of Delhi, with these constituencies being represented as polygons. But I am unfamiliar with the format in which it's geometry is given. Can someone tell me how to convert this to a format which can help me use the data on qgis or google fusion table.


texts: How to convert from UTM Zone 36N Mercator to WGS 84 degrees of lat and long?

texts:
I have acquired many coordinates from a points I extracted from polygons. The points are in WGS 84 UTM Zone 36N, and I want to use them on Google Maps. I found out that they use different systems (I am no GIS expert :S) "that uses a Mercator Projection appropriate for that zone, with Units of Meters", and I want to convert them to units of degrees of Latitude and Longitude.

Example of points I have:

3390003.50135    3390003.50135


texts: Cannot open TIF files using TauDEM in ArcMap:

texts:
I am trying to use TauDEM in Arcmap.
When I try to start the "Pit Remove" script it does not accept my .tif Input (see picture).
Same with every other script from the tbx.
If I press OK a empty error messagebox occurs.

Any idea ?


texts: Converting Matplotlib contour objects to Shapely objects

texts:
I would like to convert a Matplotlib contour collection to a Shapely geometry collection, in order to calculate the area contained within closed contour lines.

it seems that all the information is there in the Matplotlib contour collection to create Shapely polygon collection, based on these stackoverflow answers on extracting data from matplotlib contours:
https://stackoverflow.com/questions/5666056/matplotlib-extracting-data-from-contour-lines and https://stackoverflow.com/questions/17051131/matplotlib-extracting-values-from-contour-lines/17056513#17056513

Has someone already done this?


texts: ArcSDE - AutoIncrementing Object ID in table stored in ORACLE when doing INSERT or UPDATE or DELETE using FME SDE Writer

texts:
I have data in an Oracle table (11 million records) and would like to insert 2000
records approximately using the FME ArcSDE Writer to maintain ObjectID put on the Oracle table by ArcSDE via ArcGIS. I have been using FME ArcSDE writer to write the insert records to the Oracle table. I have encountered a problem of Object ID not auto incrementing to show new sequential object ids for the records inserted. I would like object id to be maintained when doing INSERT, UPDATE and DELETE operations via FME

I am not sure why this is happening and how to resolve.

Should I be using FME Oracle Writer?


texts:
I've recently started working with spatial data in R. I am currently working on a project that uses data in both raster and netcdf formats, and I'd be happy to convert the latter to the former, however I have certain doubts about the procedure.

First I tried the very simple approach of just opening the netcdf using stack(), so basically just:

library(raster)
filename <- "C:/.../X.nc"
pol_st <- stack(filename)

This seems to work fine, I can plot the rasters etc. However, the resulting rasters do not have the proper extent:

> pol_st
class       : RasterStack
dimensions  : 576, 992, 571392, 4475  (nrow, ncol, ncell, nlayers)
resolution  : 1, 1  (x, y)
extent      : 0.5, 992.5, 0.5, 576.5  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
names       : X0, X1, X2, X3, X4, X5, X6, X7, X8, X9, X10, X11, X12, X13, X14, ...

Now, I know the actual extent should be somewhere around 73, 135, 18, 54. I suppose it wouldn't be a problem to change the extent of the raster (right?), but:

Can I be sure that the projection is correct? (Do I understand correctly that there is no projection in case of netcdf, just an array of longitudes, latitudes and values?)
I still need to get the exact extent from the netcdf data.

I was looking for some answers and found threads such as this and that, which made me think that actually the solution might be more complicated.

I've tried to recreate the methods shown in those threads, but have difficulties in creating rasters of lat and lon, that is, for:

library(raster)
filename <- "C:/.../X.nc"
vname <- "lat"
pol_st <- stack(filename, varname = vname)

I get:

Error in .varName(nc, varname) :
  varname: lat does not exist in the file. Select one from:
NO2

I've tried different names for the variable too, but the "Select one from: NO2" makes it look like the data is missing lon and lat?

Unfortunately I am unable to share the data, but if I run:

library(ncdf4)
polnc <- nc_open(filename)
polnc

I get:

File C:/.../X.nc (NC_FORMAT_NETCDF4_CLASSIC):
 3 variables (excluding dimension variables):
    float latitude[lat]
    float longitude[lon]
    float NO2[lon,lat,time]

 3 dimensions:
    lat  Size:576
    lon  Size:992
    time  Size:4475


texts: Reading LiDAR attributes from a .mat file?

texts:
I am working on some LiDAR dataset. I am used to importing the LiDAR data from .las files. I generally use liblas library in Python to import the same.

Recently, I found the dataset in .mat format which is actually a structure in MATLAB. What I am able to see in that file is only a single cell which contains all the x,y,z,intensity fields as attributes.

So, to be clear, my query is how to convert this .mat file to .las file so that I can directly read the attributes using liblas?


texts: Using What3Words addresses for 911 addressing in US?

texts:
In the US 911 districts spend a lot of money creating rural addresses for their emergency dispatch systems.

Over the past few years I've worked with several districts implementing tools that assign addresses to "structures" based on a distance along a linear feature (e.g. a road).  With the adoption of 911, rural addresses based on postal route and box numbers have been converted to a street and house (or structure) number.  This can be complicated.

I've been told by these agencies that they are following procedures prescribed by NENA.

I searched for mention of What3Words in the final report from NENA's National Address Database Summit, but found nothing. (Perhaps being sponsored by the Department of Transportation steers them toward road oriented addressing systems.)

At first glance it appears W3W might be easier to maintain than systems relying on linear features.

Has anyone used What3Words for 911 rural addressing in the US?

If so, please post a link to a whitepaper describing how it was done.


texts: Calculating field with hour from date time field with ArcPy?

texts:
This expression works in the field calculator   to  give me the hour from the datetime field called dtrepor

arcpy.time.ParseDateTimeString(!dtrepor!).strftime('%H')

But how do I call it in python?
I have tried the following:

hour = arcpy.time.ParseDateTimeString(!dtrepor!).strftime('%H')
arcpy.CalculateField_management(saved_Incidents, "mil_time", hour , "PYTHON")

But it gives me a syntax error?


texts: Does the Postgis UpdateGeometrySRID function update the column structure as well as the data inside it?

texts:
I had to manually insert  a geometry from a shapefile which had SRID 0, and my column had SRID 4326, so I had to use the UpdateGeometrySRID function to update my column to SRID 0. I inserted the geometry and it worked fine, but now I need to transform everything back to SRID 4326, the column as well as the data inside it. Does UpdateGeometrySRID do that?

Additional question: how can I  transform the geometry's SRID before inserting so that I don't have to modify my column?


texts: GeoJSON coordinates transformation to Web Mercator coordinates

texts: I'm not sure how OpenLayers is converting the raw coordinates from a GeoJSON's geometry.coordinates property into a Feature's extent property. I can't really seem to find a lot of documentation on it either, does anyone understand what's going on to make this transformation?


texts: Polygons get lost when doing Join attributes by location

texts:
I am trying to do the operation "Join attributes by location" in QGIS. I have a point layer and a few polygon layers, all of them in shapefile format. I want to sum the values of the attributes of the points inside each polygon, and add that sum as a new field in the output polygon layer. The problem is that in the output polygon layer, several of the polygons from the input layer have gotten lost in the process. Why could this be?

This problem happens with some polygon layers, but not with others. And for one polygon layer, still, the operation worked fine on one computer, but not on three others.


texts: Designing Update Cursor that updates rows found between two values?

texts:
I am new to python, and have not successfully coded an Update Cursor.  I am unable to understand the documentation.

I have a geodatabase table of wells that i have built using python.  I now want to loop through each row and update a field.  the table contains a field called 'MD' which is type double and contains a depth value.  I need to look at a second field called 'FM_Name' for two text strings that indicate the TOP and BOTTOM of a certain perforation zone.  So, the cursor would loop though the rows which are sorted in ascending order, and update 'FM_Name' field with the word "Perf" in between the TOP and BOTTOM for each well.

This table has over 50K records, representing roughly 300 wells, so i really dont want to do this by hand.


texts: How to simplify line geometry and preserving topology in Oracle?

texts: I have some polylines in an Oracle database which are topologically connected and I want to simplify these lines based on a tolerance directly in Oracle. I need something like st_simplifypreservetopology in PostGIS to use in Oracle.


texts:
I've been trying to figure this out, but have come up with nothing. I have this code:

....
var lat = 7486473;
var lon = 4193332;
var lonlat = new OpenLayers.LonLat(lon, lat);
....

It uses floats for the lat/lon values. My issue is I don't know how to convert regular latitude/longitude coordinates into these numbers. I've searched for the answer, but I haven't been able to make sense of anything since I'm pretty new to GIS as this isn't my field.

The location given represents Moscow. 55.7500° N, 37.6167° E
Moscow, Coordinates


texts: Transformation functions for EPSG:3395 projection vs. EPSG:3857

texts:
The transformation functions for the Web Mercator projection EPSG:3857 for latitude and longitude in radians to projected coordinates in meters are defined as

R = 6378137
x = R * longitude
y = R * ln(tan(pi/4 + latitude/2))

I'm unable to find corresponding functions for EPSG:3395, the so-called World Mercator projection.

Running a WMS that supports both, x values appear to be identical in both projections. There is however a deviation in y values that increases with increasing distance from the equator.

So how would the function for calculating y have to be written?


texts: Do you know how it's possible to convert a layer symbology or an *.Lyr file to an *.style file?


texts: How can I put OSM data in a geography view in PostGIS?

texts:
I loaded the *.osm file to postgis with osm2pgsql. The data are stored in lat/long, so the same time I expect them to should be available in geography view, but they are not (only present in geometry column). Is there any way the osm tables to be visible in geography view?
As far as I know I cannot define other unit in function like ST_buffer than the original table is stored in? So for example if the data are in lat/long the function St_buffer will return the distance in degree. Instead of that I would like to be able to execute queries with meter parameters even if the data are stored in lat/long. To do that actually I need the osm tables to be visible in the geography column, Am I right? If so ... how to?

Thanks a lot for an answer,
Gregory


texts: Making layer visible without having it be snapped to in ArcMap?

texts:
I have one feature created called 'parcel' that shows the borders of houses and another feature called 'conduit'. I need to attach a point feature to the conduit preferably where the parcel boundaries meet.

The problem I run into is that when I go to place the point feature it always seems to want to attach to the parcel field.

I don't want to turn the field off as I want to make sure the point is placed at the correct spot but it can get annoying when I am placing 100's of points and have to undo a couple because they stopped getting attached to the conduit line feature.


texts: Iterating over layers and export them as PNG images with PyQGIS in a standalone script?

texts:
I would like to export specific layers as a png from a QGIS .qgs map. From the table of contents the script would first export the "boundary, climits, and Div1_Irrig_1956_0" layers as one png. The standalone script would then iterate to export "boundary, climits, and Div1_Irrig_1976_1" layers as the next png and so on. I am working with the script below to begin with but am only getting one png with all of the layers exported.

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

import sys, os
import glob

qgs = QgsApplication(sys.argv, True)
qgis_prefix = "C:\\OSGeo4W\\apps\\qgis"
QgsApplication.setPrefixPath(qgis_prefix, True)
qgs.initQgis()

#layers = glob.glob((configData['destination_folder'])+"\\*.shp")
layers = glob.glob(r"E:\IrrigatedLands\FC_qgis\*.shp")

for layer in layers:
    print layer
    irrig = QgsVectorLayer(layer, "testlayer_shp", "ogr")
    print irrig.isValid()
    layerset = []
    QgsMapLayerRegistry.instance().addMapLayer(irrig)
    layerset.append(irrig.id())

    # create image
    imageType = "png"
    pathToFile = "C:\\Users\\James\\Desktop\\"
    name = "render"
    img = QImage(QSize(800, 600), QImage.Format_ARGB32_Premultiplied)

    # set image's background color
    color = QColor(255, 255, 255)
    img.fill(color.rgb())

    # create painter
    p = QPainter()
    p.begin(img)
    p.setRenderHint(QPainter.Antialiasing)

    render = QgsMapRenderer()

    # set layer set
    layer_set = [irrig.id()]  # add ID of every layer
    print layer_set
    render.setLayerSet(layer_set)

    # set extent
    rect = QgsRectangle(render.fullExtent())
    rect.scale(1.1)
    render.setExtent(rect)

    # set output size
    render.setOutputSize(img.size(), img.logicalDpiX())

    # do the rendering
    render.render(p)

    p.end()

    img.save(pathToFile + name + "." + imageType ,imageType)


texts: Convert .dwg file to .shp file using ArcGIS Desktop?

texts:
I have a map in CAD 2002 format with 20 layers all with different features.

I need the procedure and steps on how to convert these layers to shape files so as to be able to work on it in ArcGIS Desktop.


texts: Does GPS usually use more than 4 satellites for a 3D fix

texts:
In reading around it seems clear that 4 satellites are required for a 3D fix, and while it is often stated a GPS chipset will track more than 4 satellites, literature often implies it will only use 4.

Is this true? And if so, is there a reason not to use more?


texts: gdal.jar cannot find gdalalljni on Ubuntu

texts:
I'm trying org.gdal library on Ubuntu.

<dependency>
  <groupId>org.gdal</groupId>
  <artifactId>gdal</artifactId>
  <version>2.3.3</version>
</dependency>

And got runtime error:
Native library load failed.
java.lang.UnsatisfiedLinkError: no gdalalljni in java.library.path

1) Before this I installed gdal:

sudo apt-get install gdal-bin
sudo apt-get install libgdal-dev

There are no any gdalalljni files installed.

2) export LD_LIBRARY_PATH=/usr/lib
which contains files libgdal.a, libgdal.so, libgdal.so.20 libgdal.so.20.3.2 and not contains any gdalalljni

3) export GDAL_DATA=/usr/share/gdal/2.2

4) set VM options -Djava.library.path="/usr/lib/"

I tried some advices such as http://geoexamples.blogspot.com/2012/05/running-gdal-java.html
Here is about 1.9 release. I not found swig/java directory in my installation of 2.2.

How to use last releases of gdal in java on Ubuntu?


texts: Using parameterAsSink in QGIS processing script?

texts:
Note I have edited this question to make it more direct and hopefully get a response.

I am attempting to write a series of processing algorithms to automate and extend a workflow in QGIS 3.4.  I have been using the processing script template and following a number of examples including.

pyqgis - post-processing of layer output of QgsProcessingAlgorithm
Defining fields and adding features to output layer using QGIS
Python processing script with clip in QGIS 3.2

I also reference the porting_processing.dox by nyalldawson and the workshop by NathanW.

My problem is ( I think) how to properly use parameterAsSink as described by nyalldawson.  Unfortunately none of the listed resources, some of which are unanswered, quite answer my question.  There seem to be many variations of how to output results from a processing script.  I am sure I am missing something very simple but it has me beat.  My best guess is I am not properly assigning the results the output sink.

I have tested the code and the file is generated correctly if I send the result to a file it works.  What I can't do is get the memory layer to work - the same as if you are usign the inbuild processing scripts

The following code is a modification of the script template attempting something simple to get started.  It runs and as above the result is valid - I just can't get it to the result to the QGIS interface.     At present I am just implementing the existing processing difference algorithm but I will extend this when I can get it working properly.

# -*- coding: utf-8 -*-

from PyQt5.QtCore import QCoreApplication
from qgis.core import (QgsProcessing,
                       QgsFeatureSink,
                       QgsProcessingException,
                       QgsProcessingAlgorithm,
                       QgsProcessingParameterFeatureSource,
                       QgsProcessingParameterFeatureSink)
import processing

class ceaEz(QgsProcessingAlgorithm):
    INPUT = 'INPUT'
    OVERLAY= "OVERLAY"
    OUTPUT = 'OUTPUT'

    def tr(self, string):

        return QCoreApplication.translate('Processing', string)

    def createInstance(self):

        return type(self)()

    def name(self):

        return self.tr('ceaez')

    def displayName(self):

        return self.tr('CEA minus EZ')

    def group(self):

        return self.tr('Sample plan algorithms')

    def groupId(self):

        return 'sampleplanning'

    def shortHelpString(self):

        return self.tr("Example algorithm short description")

    def initAlgorithm(self, config=None):
        self.addParameter(
            QgsProcessingParameterFeatureSource(
                self.INPUT,
                self.tr('CEA'),
                [QgsProcessing.TypeVectorAnyGeometry]
            )
        )
        self.addParameter(
            QgsProcessingParameterFeatureSource(
                self.OVERLAY,
                self.tr('EZ'),
                [QgsProcessing.TypeVectorAnyGeometry]
            )
        )
        self.addParameter(
            QgsProcessingParameterFeatureSink(
                self.OUTPUT,
                self.tr('OUTPUT'),
                #QgsProcessing.TypeVectorAnyGeometry
            )
        )

    def processAlgorithm(self, parameters, context, feedback):
        diff = processing.runAndLoadResults("native:difference", {
                 'INPUT': parameters['INPUT'],
                 'OVERLAY': parameters['OVERLAY'],
                 'OUTPUT': 'memory:'
        }, context = context, feedback = feedback)['OUTPUT']

        return {self.OUTPUT: diff}

Can anyone show me how to get the correct output?


texts:
I'm trying to access a WMS endpoint through Openlayers to display the data in a slippy map.

Here's my code:

var basemap = new ol.source.TileWMS( " Layer Name 1", "http://localhost:8080/geoserver/oswgc-india/wms", {layers: 'oswgc-india:india_grouplayer'} );

         var map = new ol.Map({
                layers: [
                        basemap
                        ],
                target: 'map',
                view: new ol.View({
                        center: [312406.33, 6150500.21],
                        zoom: 5
                })
          });

In Geoserver, I have created a grouped layer called india_grouplayer in the workspace oswgc-india. However, my index.html with the embedded JavaScript just returns a blank page. I suspect that I'm missing configuration options in Geoserver, as I'm not exactly very familiar with the software.


texts: QGIS/PostGIS: geom as Point in -2147483648

texts:
Based on QGIS 3.4 and PostGIS 2.5:

I imported a shapefile in my PostGIS database using shp2pgsql.

When I connect to my PostGIS database, I see that 2 tables were imported. One with the hint "geom as Point in -2147483648" and one with the hint "geom as Point in 0". When I add the one with -2147483648 to the map, I see no features and the attribute table is empty. The other table (with "geom as Point in 0") works fine (I can see the features and the records in the attribute table)

Also when I create a new table from scratch (with a geometry field) using the database manager in QGIS, 2 tables are created with the same result. Why?

What is the table/layer with -2147483648 for? Why is it created and can I delete it somehow or is this not recommended? In pgAdmin, I only see one table.

When importing a shapefile with the database manager in QGIS, this table is not created.

Coordinate system: EPSG 4326

Windows 7, 64 bit

Here is an image how it looks when I connect to PostGIS in QGIS (ignore the white space):


texts: How do duplicate script on Google Earth Engine code editor

texts:
I'm following https://docs.google.com/document/d/1bRCEB8ybRUp2D6vFc7mjqDMG2bpmCd2iTPQIEern0mk/edit#, its about classification and regression.

I followed the regression part then I'll do the classification part.

To keep things organized I prefer to do the classification part in another file. But a lot of code can be reused. So I would like to duplicate script, how to?

I can copy paste code from one file to another. Problem is with the imports, at least I don't know how to copy paste them to another file.


texts:
I have a polygon file buffer_stands.shp with 20 features and I have a line file OSMroads.geojson with 9 features.

I expected the script to loop through the buffers and print out for each buffer loop the the point count of each of the OSMroads. I know it doesn't sound useful, but the script is a simplification of a bigger script.

The outer loop loops all the way through. But the inner loop only loops once. What is causing that and how can I fix it?

import ogr

bufferfn = ogr.Open('buffer_stands.shp')
buffer_lyr = bufferfn.GetLayer()

roadfn  = ogr.Open('OSMroads.geojson')
road_lyr = roadfn.GetLayer()

for buffer_feat in buffer_lyr:
    buffer_geom = buffer_feat.GetGeometryRef()
    print 'buffer'

    for road_feat in road_lyr:
        road_geom = road_feat.GetGeometryRef()
        print road_geom.GetPointCount()

The output looks like this

>>>
buffer
10
62
17
46
10
28
97
78
121
buffer
buffer
buffer
buffer
buffer
buffer
buffer
buffer
buffer
buffer
buffer
buffer
buffer
buffer
buffer
buffer
buffer
buffer
buffer

It shows the inner loop is running only once.


texts: Batch Contour Process Broken in New Update?

texts: When I try to run the batch processing tool for contours, it does not give me the option to save the contours as a file. This stops me from running the tool all together. Has anyone else encountered this bug in 3.12.2?


texts: How to convert shapefiles into ascii grid .asc in QGIS?

texts:
I downloaded multiple zips that cover a small island from the Rhode Island Geographic Information System webpage, specifically the geographic files that make up Block Island.

Here's the source: http://www.edc.uri.edu/rigis/data/download/2011USGS/muni.aspx?ID=22

I'm trying to convert these files into .asc. From what little information I managed to find online I must make a the vector into a raster and make it a  point/grid image and then save it as a .asc

I can get the point style raster image but I'm unable to save/export as the.asc.

Does anyone know how to do it? I also have geotiff and grass gis if it helps the process.


texts: Calculate percentage coverage datapoints

texts:
I have about 8000 datapoints of giraffe observations and would like to calculate the surface of the area they spend 90% of their time in.

Anyone who can help me with solving this question?

I use QGIS 2.14.3 Essen.


texts:
I want to find the area of a set of shapes in QGIS. I can do this step(Vector>Geometery Tools>Export add Geometry Columns or alternatively add feature $area) but I need the output to be in square meters or alternatively square kilometres.

I'm struggling to ascertain what projection I should be using, (I've tried a range, WGS 84 EPSG:4326, WGS 84 /UT< zone 59S EPSG:32759 & NZGD2000 / New Zealand Transverse Mercator 200 EPSG:2193  but at no stage have I got values that are fit for purpose i.e. in square meters or kilometers). The area I'm interested in is on the SE coast of New Zealand or http://gpsinformation.net/utm-zones.gif 59 G on this map.

The other possibility is i'm using/have used the correct the projection but have not set up QGIS correctly to measure area. My process is to set the project CRS > resave the project with the new CRS > add layers > change layers to the selected crs.


texts: Custom sprite and glyphs not showing in Maputnik

texts:
I'm completely new to vector tiles and mapbox. In order to test different icons, fills, etc. I found Maputnik would be the best option. Working with regular outlines, fills and other primitives is just fine but when it came to implementing custom icons and text labels I became confused.

Using spritezero-cli I created a custom sprite: https://raw.githubusercontent.com/yaras-phoenix/test/master/sprite/test_sprite

And according to this tutorial I created a simple onefont-fontstack: https://developer.tomtom.com/maps-sdk-web/tutorials-advanced/creating-custom-glyphs

I placed my links in Maputnik's Style Settings and expected to fill "Pattern" option to be activated and to see my font in text "Font" drop-down list but none of that happened.

Sprite's .json and .png files are fetched with code 200 but Maputnik seems to expect one more file I don't have. I see "GET https://raw.githubusercontent.com/yaras-phoenix/test/master/glyphs.json 404 (Not Found)" in Chrome DevTools. So I guess I need some extra json file.

Could anyone point me in some direction? I'm quite sure I do something wrong.


texts:
A bit of a strange question perhaps, but let me give you a breif explanation of the background before my actual questions:

Atomic force microscopy (AFM) is a method that, in short (and to my limited knowledge) lets researchers scan areas at the micro- and nanoscale. It works by "scanning" an area using a probe of sorts. More that is difficult for me to explain, as I have no real understanding of it. What I do know, and what triggered my curiosity was that the result is in fact a "grid" of "height" values (a matrix of say 512x512 values describing the height of the probe at that point).

I then thought: Well, apart from the scale, this is in fact a digital elevation model! And, this means that if I can manage to create a DEM file as understood by GIS tools I could apply GIS analysis to it!

As it is, my significat other works at a lab that has a AFM-machine, and is using it in one of her projects. I got some scan files from her, and have managed, using Python (struct and numpy), to parse these binary files and what I have now is a numpy array of size 512x512 filled with int16-values.

What I'm planning on next, and what I need some help with is the "mapping to a proper DEM"-part. I have some knowledge about DEMS, but when it comes to actual generation of them I am quite new.

What I'm thinking is that I have to georeference my data in some way, and for this I need a custom (planar) coordinate system. I envision my coordinate system would use micro- or nano-meters as units. Then it's just a matter of finding the size of the area scanned with the AFM (this I believe is somewhere in the binary file, assume this known).

update: I also have several scans at different resolutions, but of the same area. For example i have this information about two scans:

larger image:

Scan Size: 51443.5 nm
X Offset: 0 nm
Y Offset: 0 nm

smaller (detail) image:

Scan Size: 5907.44 nm
X Offset: 8776.47 nm
Y Offset: 1486.78 nm

What I', thinking is that my custom coordinate system should have an origin in 0,0 and for the larger image I cam assign pixel 0,0 the coordinate value of (0,0) and pixel 512,512 the coordinate value (51443.5, 51443.5) (Guess you get the picture for the other points needed).

Then, the larger image would map pixel (0,0) to  (8776.47, 1486.78) and (512,512) to (8776.47 + 5907.44, 1486.78 + 5907.44)

1st Question is then: How do I create a proj4 def for such a coordinate system? I.e: how do I assign these "real world coordinates" to my custom coordinate system (or, if I follow whubers suggestion and using a local coordinate system and lying about units (i.e. treating my nanometers as kilometers)

Then I have to transfer my numpy 2-dimensional Array to a Georeferenced DEM file format. I was thinking of using GDAL (or, rather the Python bindings).

2nd Question is then: How do I create a georeferenced DEM from "arbitrary" data such as mine? Preferably in Python and using open source libraries.

The rest should then be fairly easy going, just a matter of using the right analysis tools. The problem is, this task is driven by my own curiosity, so I'm not quite sure what I actually should do with a nanoscale DEM. This begs the

3rd question: What to do with a nanoscale DEM? What kind of analysis can be done, what are appropriate tools for DEM analysis and finally: is it feasible to make a map with hillshade and contour lines from this data?:)

I welcome all suggestions and pointers, but keep in mind that I am looking for free alternativies, as this is a strictly hobby-based project with no budget or funding (and I do not have access to any lisenced GIS-applications). In addition, I do know that Bruker, the company that sells these AFM machines, ships some software, but using that would not be any fun.


texts:
I have a layer with symbology that includes geometry generators. It would appear that because of these the symbol is not visible in the layers panel.

Is there a way to show the up most symbol (green square) in the panel?

Following the answer by NettaB

As suggested I tried locking the marker but the after clicking okay the symbol flashes up for a short period then disappears.

Geometry generator expressions

For the line:

make_line(make_point($x, $y), make_point("xcoord", "ycoord"))

For the marker:

 make_point("xcoord","ycoord")


texts:
I found a query in the web which create catchments by pgrouting.
I tried the query several times and modified a bit. But it dosen´t work.
It alwys shows me the same error.

    DO $$
DECLARE
    c integer;
BEGIN
    DROP TABLE IF EXISTS tmp;
    CREATE TABLE tmp (
        agg_cost integer,
        geom geometry(MULTIPOLYGON, 4326)
    );

    FOR c IN SELECT agg_cost FROM temp1 GROUP BY agg_cost HAVING COUNT(1);
    3 ORDER BY 1 DESC LOOP
        RAISE INFO '%', c;
        INSERT INTO tmp (agg_cost, geom)
        SELECT
            c,
            ST_Multi(ST_SetSRID(pgr_pointsAsPolygon(
                'SELECT
                        temp1.id::integer,
                        ST_X(temp1.the_geom)::float AS x,
                        ST_Y(temp1.the_geom)::float AS y
                FROM temp1
                WHERE agg_cost = ' || c
            ), 4326));
    END LOOP;
END$$;

the query is from the page:
https://www.zimmi.cz/posts/2017/routing-with-pgrouting-catchment-area-calculation/
and I got the error:

  ERROR:  missing "LOOP" at end of SQL expression
  LINE 11: ...ELECT agg_cost FROM temp1 GROUP BY agg_cost HAVING COUNT(1);

Does someone know what I need to change?


texts: Is it possible to automatically create categories with mapserver?

texts:
I have a mapfile connecting to a postgres/postgis database. To create a thematic map I can set a class item and a number of classes, for example:

CLASSITEM  "aantal"
    CLASS
        NAME "0"
        EXPRESSION ([aantal] = 0)
        STYLE
            COLOR 116 255 117
            OUTLINECOLOR 178 178 178
        END
    END
    CLASS
        NAME "100-250"
        EXPRESSION ([aantal] >= 1 AND [aantal] <= 2)
        STYLE
            COLOR 0 255 0
            OUTLINECOLOR 178 178 178
        END
    END
    CLASS
        NAME "250-500"
        EXPRESSION ([aantal] >= 3 AND [aantal] <= 10)
        STYLE
            COLOR 0 210 0
            OUTLINECOLOR 178 178 178
        END
    END

However I have to do this 25 times with different attributes and different values. So I would like to know whether it is possible to do this automatically.

Similar to the option in QGIS where you can set the number of classes and mode.


texts: How to synchronize layer orders in QGIS "layers()" and in the Layers Panels GUI

texts:
I'm trying to change the order of layers in qgis3 (python). But the code I can come up with only changes the order of layers in the list iface.mapCanvas().layers(), and not the displayed layer order in the lower left Layers Panel in the GUI.

My code below, e.g. swaps the orders of two arbitrary loaded layers. The print() statement shows that orders in layers() are indeed swapped with each run of the code.

lyr1 = iface.mapCanvas().layers()[0]
lyr2 = iface.mapCanvas().layers()[1]
treeRoot = QgsProject.instance().layerTreeRoot()
treeRoot.setHasCustomLayerOrder(True)
order = treeRoot.customLayerOrder()
order.append(order.pop(0))
treeRoot.setCustomLayerOrder( order )
iface.layerTreeCanvasBridge().setCanvasLayers()

print (iface.mapCanvas().layers()[0].name())
print (iface.mapCanvas().layers()[1].name())

I expected this to change the layer orders in QGIS GUI as well. However, the Layers Panel stays the same after running the scripts multiple times.

How can I synchronize the Layers Panel to reflect the order changes in the underlying layers()?


texts: ArcGIS converting list of enumeration in single field of s-57 data?

texts:
Im storing s-57 data in an MDB which is being symolised via a layer. Some attributes in s-57 data can be a list of enumerations and I want to find a way to show the actual values the enumeration represents, please see example below.

My feature has a column called colour, this can be represented as a single number (i.e. 1 = white) or it can be represented by a list (ie. 2,3,2 = black, red, black).

I have included the enumeration table in the MDB which has values similar to those shown below.

VALUE   COLOUR

1   white

2   black

3   red

where a single colour is listed in the field a simple join to the enumeration table works fine and I can set it up so the users see the full attribute values instead of the enumeration but when a list of enumerations I still need to convert it to full attribute values, is this possible?

so in summary if I have a colour field that contains '2,3,2' using the enumeration table above I would like to be able to get the layer to display black, red, black instead of 2,3,2 in the attribute table, is this possible?


texts: How to cast data types from GeoJSON feed (aka: how to sync when everything is a string)?

texts:
Just getting started with the CartoDB Way and have run up against first problem. Looks like I'm not alone.

We have a GeoJSON feed coming from a Fulcrum data collection app. We want to see current results of that collection on a map, so I have set up a GeoJSON feed from Fulcrum as a starting point.

Over in CartoDB, I create a table and provide the GeoJSON URL.

The data appears...all as strings, except for the_geom (thankfully).

I want to create a sort of "set and forget" transfom process so that
the column types are always converted properly and made available to CartoDB as their correct data types.

What is the next step?

Thanks indeed!


texts: How to make a cloud mask from a daily MOD09 image?

texts:
I'm currently working with MODIS data and i've been having some trouble with several aspects of the steps I need to take.

For example, I need to make a cloud mask from a daily MOD09GQK image, and i'm confused on how to do it. I have the QC_250m layer in a tif, but I don't understand the pixel values. I know I should convert, say 4096 to bits and get: 100000010010110, where the 3rd and 4th bits would represent the 00 'clear' flag.

But my main doubt is why do I get the same value where the image is cloudy and where the image is clear?

Help, please?

Thanks a lot


texts:
I have WMS map, after applying CQL filter to layer how can I zoom map to see filter result?

UPDATE:

Here is my code what I'm trying to do

var layerArray = map.layers;
var someLayer;

for (var i = 0; i < layerArray.length; i++) {
    if (layerArray[i].name == "myLayer"){
        someLayer = layerArray[i];
        break;
        }
}

someLayer.mergeNewParams({
    "CQL_FILTER": celIDCondition /* (CID =1351 and LAC = 401) */
});

// instead of somelayer.maxExtend I want to extent map acording to CQL_FILTER result
map.zoomToExtent(someLayer.maxExtent);


texts: Use R to publish GeoTIFF as WMS in GeoServer

texts:
I’ve been trying to automate the publishing of GeoTIFF rasters as a Web Map Service (WMS) to GeoServer using R.
So far I’ve only been able to establish a new datastore from inside R using this:

system(paste('curl -u admin:geoserver -v -XPOST -H "Content-type: application/xml" -d "<coverageStore><name>myraster</name><workspace>myworkspace</workspace><enabled>true</enabled><type>GeoTIFF</type><url> /usr/share/geoserver/data/data/myworkspace/myraster.tif</url></coverageStore>" "http:// my_mv_ip/geoserver/rest/workspaces/myworkspace/coveragestores"', sep=""))

Similarly, you can also do this using the geosapi package in R. (https://github.com/eblondel/geosapi/wiki).
But progressing to the next step of adding and publishing a GeoTIFF as a Web Map Service (WMS) has been problematic and raster management is currently not supported with geosapi.

How do I publish GeoTIFFs as WMS to GeoServer in R?


texts: Decode encoded polylines to usable format for QGIS

texts:
Beginner alert.

I was given a set of encoded polylines that I need to import into a QGIS project. I understand what encoded polylines are but now I am stuck with how to decode them into some usable formats (I assume csv files?) in order to import into QGIS.

I did some research but have only found resources around encoding the data but not the other way around.

https://developers.google.com/maps/documentation/utilities/polylineutility

^seems to do what I want but I need to download the polylines in a format that can be imported into QGIS.

Here's an example:
"_{reFxmjVpAeBn@{@w@aAc@m@yAsBJuAC}BIeAO_Ak@sBCQOg@g@aAk@}@qAwAs@s@}I{JgCqCsC_D}DeEyDgEsA}AwB{BcCwCoGyGqBwByEwFgHqHgDqDkIaJyFuGeEsEeEcEcGkGkEaFyDyDgBqBcAgAm@y@w@gAuCcCmFeHk@i@eAgAaDaDcC_CoBaCeBkCwAiCkCuFoPy]eAkCy@eC{@eDq@kDk@cEwAyQiBsVkAsN}Egn@q@sJ_@iGWgFg@mMcAuWQcCg@yFkBkRc@{EiAmKuBiT]sDWiEy@kPg@_Gg@_E{AuOkAgKo@gFe@oDq@oEe@qBmAeEWkAOaAKoAC}ADiBPeBFc@\wA~@oC@sAh@CVgBJoA@c@EWBkAAcBGgB[yC}AuNu@iHG}B@uBJcBJmAjEeZbAuGz@uFh@}DV]H[RkAfAqGJe@l@mBl@eAv@w@n@]TIj@Ml@In@?Df@n@Fh@@|H|AzDj@zDl@l@F^KhAZhCp@pAXnCr@bBh@nBbAfCdB|@z@vCfDPTjAvAhArAjBjBlBdBlEnC\TCzAtEpCtDtBrEdCdExB~DpBnAj@z@XdAPtBLr@CH?pAQZKfBs@t@e@x@s@|@cAt@mAh@gA@mAlAmEFELSpAyD|A_GfAeC|B}HT_At@gDjCaJZq@dB{F~AeFmCuAgCyAv@mCd@}AfA}C~@oBr@iA|BcDn@gA~@iBTo@~@qCh@gCXoB??P{AP_BU_@e@w@d@v@nArBl@v@ClBlBtAdATnAXLCFAXiALc@Ng@X}A^eCvB{T\qDZsBSeATmBr@mE?@ASVkAv@iDlBqHvA_FhAcDv@qBtBaFbIyQpBsDrAoBdAiAp@k@r@g@v@a@lAi@dBk@dG{AfD}@v@YzAq@|BsAfA{@[|A_BvAiBv@oAAiB@A@cAf@_BzBaJrM}j@x@eCz@mBbAgBlAeBhAmArAmAhCkB|MeJnD_CjByAxE}DtG{FlNsL~WaUjMmK|GyFjR}OdI_HhYcV~SoQfCwBnGeF~I}HbHcGnDyCDmCbJyHhA}@~FeF~GwFnNmLzCkCtJkJ|MeNDeDjGkGxGmGtAyApHcHrBiBlAsAdF_GRyZhI{MlJcPBeC~PaYdPaXjOqVbLcRdXmc@dDqFhLgRtGmKj@aAlB}Cz@oAR]tAeC~EqHtWwb@fRc[bDgGzB}Df@{@fIuLxGsKzBwCtCyCbB{AfAw@bC}AlGgDbEsBdIgEl@]fNkH~VcN~f@wWzP_Jr@c@zIeFAo@xEwCbHmEdGyDnIsF~MyIrJmGrUaOlSmMzDgC~k@e@bf@{Z~n@ca@|QqLtRcMdp@ab@naAkn@fBeADaBfDwANBJ?nAWhCm@^Gf@@XHZPp@f@xDfCpF}NpAcDpA_DjF}MtGyPLY]UyDkCoFqDaAC_@z@aApBoAzC]W??\Vt@iBd@gAdAuBNa@r@gBz@h@~B~AdBhAhChBZTLYnCeHtAmDnHgRjHcRtCqH~FcOvBuFdEyKpLoZvL_[v@mBjH}QXaAfAkCx@sBXaAlEwKpA}CfDgIx@mB^y@FE@Ap@uAfCaGrIsSA_CxAeDnGyOlFiMdBeECcGhG_O??AaCnDrC~DvCrNtKz_@~YfUfQtFhE~DD|@p@\Pn@RbATjBBpABzANv@ThAh@|@n@@^ZXf@p@NfA@t@KXSV_@LWA[KU]Ie@?WJk@Le@hA}CBK@g@bAwBlCkGzFoN~KqYpOk@zF_OzAaDfB_DpCmEzJwMfPqTT}XjLoOvJwMh@q@vK}NjEwFbDsErKiOvWi_@d\{d@V_@r@cA~@qAb@u@xB{CrIyLhIoLxKsOnG_JjB}Bp@u@fC_CxCyBnC{A|D}AfFuBrOeGG}Bfh@mShc@}PdIeDhG{BbE{AjD_AzDw@pFq@zJiA|@MxAKbNcBfJaAtOkBbKiA~TwCzFo@nC]bKaAhGm@~Hi@d@AbPaA~AOpCe@pEkAdA]zEkBj@UfB}@|Aw@xBs@|Ac@xJ_Cf@Dr@GrBa@fBYdBQvDQv@MlA]pG{BrAc@xCu@~Cm@zCg@|@WSqAOeA_@oDMqDT}`@JaC^qBZcAZw@xD}HvAsCpBiE~BkGJYg@c@eCaBO^i@tAc@dAo@~A"

Which lead my other question: most of the examples I found online were definitely not as long as the above. is there something wrong with it?


texts: How to restrict search area by city in cloudmade?

texts:
I use cloudmade geocoding api for search some geographic entity in a certain city, but i don't know how to restrict search area by city. I try to use structured query, for example:

http://beta.geocoding.cloudmade.com/v3/api key/api/geo.location.search.2?q=[city=my_city_name]

but, response, what i retrieved, doesn't looked like what i need.

I know about bbox, but i don't know how to determine bbox for a certain city.

I will welcome any suggestions.


texts: Calculation between Raw GPS data to latitude,longitude points?

texts: I am working on the GPS related project. The objective is to understand the calculation process of the raw GPS data (server side calculation). It means " How did we get the latitude and the longitude points from GPS raw data (the converion between the Raw GPS data to latitude and longitude points).I have been searching for the material related to this. But still can't find one.


texts:
I have a lat / long and I'm trying to convert that point to a screen point (the pixel location on the current map). Current version of arcgis-js-api: "3.24.0". Here is my implementation:

var point = new Point(track.geometry.coordinates, map.spatialReference); // 4326 spatial reference
return ScreenUtils.toScreenPoint(map.extent, map.width, map.height, point);

sometimes the screen point will display the correct location, then if you pan around the world it no longer has the correct pixel location. Or if I have two points across the IDL then one is in the correct location and the other is not.


texts: Getting buildings that house specific amenities using R?

texts:
I am new using OSM and I am little bit lost here. I managed to get all restaurants in Winterthur, a city in Switzerland using amenity=restaurant. I used the R-package osmdata as follows:

q1 <- getbb ("winterthur, switzerland") %>%
      opq () %>%
      add_osm_feature("amenity", "restaurant")

I extract the locations as follows:

out <- osmdata_sf(q1)
out$osm_points

Is there an easy way to find the corresponding buildings (e.g. osm ids), belonging to the amenities.


texts: Access individual feature properties when running l.esri.find or l.esri.query

texts:
Here is what my site looks like: https://brandonlprice.github.io/landrecords/. How can I access the feature properties for a layer when running a find or a query on a layer? I have the the feature.properties.length being displayed in a .innerhtml to show the total count of features returned. I also would like to be able to display the individual attributes for each feature in a layer in a .innerhtml when they are clicked on in the map. The onEachFeature would be one way although I do not know how to link that to the query.run or find.run. How can I do that? Are there any other way to access the individual attributes. My query snippet are below.
Query:
function tmqueryFunction() {
     latlng = markerLayer.getLatLng();
     tmbufferResult.clearLayers();
document.getElementById("tractcounts").innerHTML = '';
    if (document.getElementById("tm1").checked == true){
       tmquery.nearby(latlng, distance);
       tmquery.run(function (error, featureCollection, response) {
  if (error) {
    console.log(error);
    return;
  }
  tmbufferResult.addData(featureCollection).addTo(map)
  document.getElementById('tractcounts').innerHTML = 'Tract Map Records: ' +
  featureCollection.features.length;
  });
  map.fitBounds(L.geoJSON(featureCollection).getBounds());
});
}
}


texts: Using the UpdateCursor to change the values in fields added with arcpy.Addfield_management

texts:
I am working with some ArcGIS data and trying to add fields using AddField_management.

>>> arcpy.AddField_management("huc4","Class","TEXT")

I want to populate this new field using the UpdateCursor and two other fields that already exist.

>>> cursor = arcpy.da.UpdateCursor("huc4",["SHAPE_Area","SHAPE_Leng","Class"])
>>>
>>> for row in cursor:
...   if row[0]/row[1] > 40000:
...         row[2] = "Low"
...         cursor.updateRow([row])
...
Runtime error
Traceback (most recent call last):
  File "<string>", line 4, in <module>
TypeError: sequence size must match size of the row
>>>

I don't really understand this error, but it's happened with every attempt I've made to use the UpdateCursor to populate my new fields.


texts: OpenLayers 3 display WFS EPSG:4326 on EPSG:3857 map

texts:
Is it possible in OpenLayers 3 to project an EPSG:4326 WFS vector layer on to say a MapQuest EPSG:3857 map?

I've tried the following code (modifying this example OpenLayers 3 WFS example, but the EPSG:4326 layer coordinates get interpreted as EPSG:3857 and so I get a tiny drawing situated around 0 long, 0 lat.

There's obviously a reprojection needed somewhere, but where?


texts: Unable to connect to data source attendance "ODBC:attendance"

texts:
I'm having an issue connecting to my MS Access DB (.md file) under the  ODBC  USER DSN configuration in Windows.

Setting up my DSN configuration under the System DSN failed, rather USER DSN works(which to me is somewhat strange). I am able to list the layers(tables) within the .md file pointed to by my DSN configuration.

Just pointing out: My DSN is rightly configured for the right for a 64-bit Windows virtual machine.

The issue here is this:  I'm trying to create the server in my pgAdmin or psql session using the  SQL returned by running:  ogr_fdw_info -s "ODBC:attendance" -l USERINFO.

But doing that raises a connection error that reads like so:


texts: How to find the coordinates of the intersection points between two SpatialLines objects (or SpatialLines and SpatialPolygon edges)

texts:
I am try to find the co-ordinates of points where SpatialLines cross the edge of a SpatialPolygon. I have searched and searched, but can only find results such as over which tell me which lines cross but not where, or gIntersection which cut the lines down to sections within the polygon.

Here is a subset of my data:

# Required packages
library(rgeos)
library(sp)
library(maptools)
library(raster)
library(grDevices)

# Create convex hull
# Points to create hull
Hull_pts <- structure(list(x = c(38.87584, 38.89215, 38.87062, 38.86157, 38.72808, 38.22315, 38.12702, 38.05936, 37.95169, 37.96915, 38.14758, 38.22325, 38.34001, 38.6394, 38.6447),
                           y = c(4.17092, 4.05521, 3.93639, 3.9168, 3.85601, 3.89487, 4.17066, 4.38951, 4.98351, 4.99706, 5.12187, 5.17153, 5.16003, 4.77422, 4.76607)),
                      .Names = c("x", "y"),
                      row.names = c(1574L, 1540L, 1490L, 1482L, 1457L, 1473L, 294L, 1718L, 1131L, 974L, 2838L, 2101L, 111L, 1914L, 1909L),
                      class = "data.frame")

# Create SpatialLines of Hull
Hull_lines <- list()
for (i in 1:length(Hull_pts$x)-1) {
  Hull_lines[i] <- Lines(list(
    Line(rbind(Hull_pts[i,1:2], Hull_pts[i+1,1:2])))
    , ID=i)
  Hull_lines[15] <- Lines(list(
    Line(rbind(Hull_pts[15,1:2], Hull_pts[1,1:2])))
    , ID="15")
}
# Make it a SpatialLines object
Hull_spLines <- SpatialLines(Hull_lines, proj4string = CRS("+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0"))

# Create SpatialPolygon of Hull (as alternative)
Hull_poly <- SpatialPolygons(list(
  Polygons(list(
    Polygon(Hull_pts))
    , "ID"))
  , proj4string=CRS("+proj=longlat +datum=WGS84"))

# Outer points
Site <- c(seq(1,10,1))
x <- c(38.21467,37.22799, 38.22852, 39.42621, 37.85457, 38.10600, 37.84077, 38.02475, 38.07228, 37.40716)
y <- c(6.327209, 5.691638, 6.237393, 5.382760, 6.042394, 6.078254, 6.007846, 5.973651, 6.044092, 5.746368)
Sites <- data.frame(Site, x, y)
xy <- Sites[2:3]
Sites2 <- SpatialPointsDataFrame(coords=xy, data=Sites)
projection(Sites2) <- CRS('+proj=longlat +datum=WGS84')
crs(Sites2) <- "+proj=longlat +datum=WGS84"

# Create lines
Radial_lines <- list()
for (i in 1:length(Sites$x)) {
  Radial_lines[i] <- Lines(list(
    Line(rbind(c(38.35419, 4.483533), Sites[i,2:3])))
    , ID=Sites[i,1])
# Make SpatialLines object
Radial_spLines <- SpatialLines(Radial_lines, proj4string = CRS("+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0"))

# Plot
plot(Hull_spLines, xlim=c(37,40.5), ylim=c(2.5,6.6), axes=T)
points(Sites2$x, Sites2$y, pch=20, col="blue")
plot(Centre, add=T, col="red", pch="+", cex=1.5)
plot(Radial_spLines, add=T)

What I need to find are the co-ordinates of the points where the radial lines cross the polygon, as achieved here: http://uk.mathworks.com/help/map/ref/polyxpoly.html. I know this is possible in ArcGIS, and it doesn't seem like it should be that difficult.

Does anyone have any ideas?


texts: Adding multiple vector files in a folder as an input in QGIS processing modeler using iterator tool

texts: I have shapefiles in a folder that I intend to use as an input in QGIS Processing Modeler. I would like to iterate within the folder to add each of the shapefiles as an input. It did it in ArcGIS using the 'Iterator" tool but I prefer to do it on QGIS using similar tools like 'Iterator'.


texts: Error while running SAGA algorithms in QGIS 3.10.14-A Coruña

texts:
I've been trying to run several hydrological algorithms (Strahler, Fill Sinks etc) from SAGA in QGIS ltr and I keep running into the same error:
The following layers were not correctly generated.C:\Users\HANNES~1.BYS\AppData\Local\Temp/processing_54346125960a43788d22bdb05d351ac6/a33e5a9db9f7418e8b1c75f5b08eac12/STRAHLER.sdatC:\Users\HANNES~1.BYS\AppData\Local\Temp/processing_54346125960a43788d22bdb05d351ac6/ffe2e0caaf2c4a5a9f12648fa77f4364/STRAHLER.sdatYou can check the 'Log Messages Panel' in QGIS main window to find more information about the execution of the algorithm.
I've reinstalled QGIS via OSGeo4W and stand alone installer and I've also reinstalled SAGA.
Here is a picture of my workflow:


texts: Extracting Z data from point layer using ArcGIS Desktop?

texts:
I have created a point layer with ~150 points. I have a DEM beneath it, which is where I have taken the X and Y values from.

Does anybody know how I can retrieve the Z values, without using the info button?

I am using ArcGIS Desktop 10. I have licences for the spatial analyst etc tools.


texts: Automating Buffer Wizard in ArcGIS Desktop by working in Shell

texts:
I am using ArcGIS Desktop 10.x.

I am trying to automate the Buffer Wizard since I have several large sets of the data that I need to buffer at multiple distances. Doing it by hand will take a long time and be quite tedious, but I've already compared buffer times: even arcpy scripting with the buffer tool takes longer than using the Buffer Wizard by hand.

I do not have experience working in the Shell, so I am wondering whether this is feasible?


texts: Merging CAD and .mdb table into a feature class

texts:
I have cad file (.dwg) which contains the geometry and the attributes are in the tables of access (.mdb) database. The common FIELD in both of them is "HANDLE".

How do I convert these into a geodatabase taking the geometry from the CAD and the attributes from the .mdb?

I have tried first by converting the cad into feature class and then JOIN it with the mdb tables, but the fields are getting repeated.
How to update or copy the attributes from the tables into the geodatabase?


texts: Looking for advice regarding solar potential applications

texts:
I have the task of looking into whether it is possible to build a solar potential application which would inform the user whether a roof is suitable for solar panels and a brief report about the roof´s potential output.

I have seen quite a few of these applications and they are pretty user friendly and easy to use. Being a GIS Analyst with experience in vector/raster data analysis and presentation using opensource desktop/webgis tools ( and ESRI products of course!) I am not a newcomer to this kind of project...but I must admit, the consequences of getting the methodology wrong would be quite embarassing.

Does anyone have any "real" experience in this area and can report or supply advice regarding the methodology and processes involved. Is this something that can be realistically done by someone without any experience in the solar world...or should this be left to the universities?

I look forward to hearing from any comments.


texts: Importing kml google earth polygons to QGIS

texts:
I've created some polygons in google earth, saved it to kml and then I loaded these files into QGIS. The problem seems to be the projection, these polygons are slightly out of the place. I've tried everything I know:
- enabling "on the fly" projection
- saving the kml to shape file with other projections
- changing my project's projection and all the layers to WGS84 (EPSG:4326)

Nothing seems to work, they always appear in the wrong place.


texts: GeoServer WMS Cascading how to define a custom parameter

texts:
I am attempting to configure our internal GeoServer to cascade an external WMS server.  The external server has some custom properties which need to be passed when making a call to that WMS.
http://externalwms.com/?VERSION=1.3.0&REQUEST=GetMap&LAYERS=ENC&CRS=EPSG:4326&BBOX=53.8,8.4,54.1,9.2&CSBOOL=384&CSVALUE=,,,,,3&WIDTH=1280&HEIGHT=800&FORMAT=image/png
as you can see, the custom properties are
&CSBOOL=384
&CSVALUE=,,,,,3

I am unsure how to configure our WMS service with these parameters so the remote service is called with them.
I tried adding them to the service definition / GetCapabilities without luck.  it appears to ignore the additional parameters.

I then tried to configure in the layers but there is no clear place to add custom parameters.
How can I achieve this?


texts: Survey123 Python Module Feature Report Names

texts:
I am using Survey123 version 3.9 and the ArcGIS API for Python. I am successfully downloading reports, but I seem to have no control over the report file name. The test file name I am using is "Feature Report Title", but the file name of the report ouptut is "SurveyNameOIDXXXYYYYMMDDxxxxxx". When I don't specify the URL of my organization's ArcGISOnline account for GIS(), I see "Feature Report Title" in the survey's Report --> Recent tasks pane on Survey123 for a failed report download, but the report failed to download, and I get this message, "An error occurred when uploading result files to Portal for ArcGIS. 498. Invalid token.". Can I control the file name of the generate_report() outputs? I read this, https://community.esri.com/t5/arcgis-survey123-questions/survey123-report-feature-how-to-change-the-word-doc-file-name/m-p/829348#M17993, which seems to indicate that in early 2019, controlling file names was not an option.
Here is my code:
username = 'myname_organization'
password = 'password'
template = 1
utc_offset = '-07:00'

gis = GIS("https://organization.maps.arcgis.com", username=username, password=password)
sm = arcgis.apps.survey123.SurveyManager(gis)

survey_id = 'some unique identifier'
where_clause = 'objectid=someOID'
report_title = "Feature Report Title"
templates = survey.report_templates
survey.generate_report(templates[template], where_clause, utc_offset, report_title)

for x in survey.reports:
    if survey_id in x.description:
        data_item = gis.content.get(x.id)
        data_item.download(save_path=r'DownloadFolderPath')


texts: Auto-increment ID in shapefile table when adding new feature

texts:
In QGIS I have a shapefile containing about 2000 features where I added a unique ID using the $id expression . I want to continue adding features to the same table, using the attribute feature pop-up after finishing the geometry creation.

My question is: is there an expression or function that auto-increment on the id numbers already in the id column? E.g. if you have id's from 1-2000, then QGIS automatically suggesting id number 2001 for the next feature added?


texts: Using ArcPy to ListFields inside SearchCursor?

texts:
The current code I have simply lists ALL fields in the feature class and seems to ignore my logic of if row[0] is not None:
Why does it do this?

fields = arcpy.ListFields('\\\\Bvweb104draco\\e\\03_MapServices\\SDE_CONNECTIONS\\bvcolumbusdev.sde\\bvcolumbusdev.sde.Telecom\\bvcolumbusdev.sde.AREAS')
featureClass = '\\\\Bvweb104draco\\e\\03_MapServices\\SDE_CONNECTIONS\\bvcolumbusdev.sde\\bvcolumbusdev.sde.Telecom\\bvcolumbusdev.sde.AREAS'

with arcpy.da.SearchCursor(featureClass, str(fields)) as cursor:
    for field in fields:
        if row[0] is not None:
            print(field.name)

I want the cursor to loop through each field and if it stumbles upon a row that has a value other than NULL, print that field name.


texts:
These are 3 Isochrones created with the OSM-Tool. They are marking catchment areas for footwalking.
Since im pretty new to QGIS I have no idea how to unite the isochrones into 1 Isochrone. Any idea?


texts: Default UI Theme is still set to Dark Mode (QGIS 3.8)

texts: I've noticed that the default for version 3.8 appears to be dark mode. I've set it Night Mapping and back to Default and it's still dark mode (I restarted both times). Any suggestions for setting ui to a white back ground (Blend of Gray is not cutting it for me either)?


texts: Template for Feature Edit Popup in GeoExplorer

texts:
I need to customize GeoExplorer (OpenGeoSuite). I have to create an input mask for a point layer in editing (create and modify tools).
When I insert a point, GeoExplorer shows a popup with attributes as in DB (Postgis) in a grid form. I'd like to let users choose between few values, using combobox, propose defaults for some attributes,  and so on.
I've found the Freemarker Template for identify (GetFeatureInfo Template) and it's ok but I need something similar in editing.

I'm using GeoExplorer in Geonode


texts: Is it possible to change the column order of a table in postgresql/postgis 2.0? If so, how is it done?


texts: Stream Gradient - two types - using NHD high resolution data

texts:
I need to calculate slope, or stream gradient, two different ways:

Main channel slope from the pour point to the head/source point along the vector stream network
a local slope consisting of 100 meters up from the pour point

I am using NHD High Resolution data (need to use this dataset as some other calculations are dependent on it so I don't want to go the direction of deriving my own streams from DEM data) and have just shy of 500 "pour points" that are snapped to the NHD High Resolution data. I fully understand the intent and concepts but am struggling with how to compose the methods to automate this for all the sites within the confines of ArcGIS (version 9.3.1). I do have 10m DEM data for the whole study area for getting the elevations.

Obstacles I'm wrestling with mentally:

How to deal with identifying all the stream segments that make up the "main channel" so that I can sum their lengths. In some cases near the head of the stream there may be a first-order tributary, but then the main channel also has the same stream order (I've calculated both Strahler and Shreve stream orders using a tool called RivEx so the stream segments are attributed with these as well as their from-nodes and to-nodes via the RivEx tool, and I have the from-nodes and to-nodes). There doesn't seem to be a common attribute in the NHD flowline attributes (e.g., GNIS ID) to essentially dissolve common, or at least to somehow select all of them. Many of the segments simply don't have GNIS ID's.
How to measure just 100 m up from a pour point and then to generate a new "from node" to represent that 100 m upstream location to obtain Elevation data and thus calculate the local stream slope.

I've been working some long hours so my brain is probably making this harder than it is...

and I have come across these two other posts on GIS StackExchange but have been struggling with how to adapt them to my needs:

Workflow for determining stream gradient?

Assigning points to their respective Reaches using Stream Gradient?

It seems like all the pieces are there ("m-aware vector streams"), etc. but I'm starting to wonder if that's not the case without, again, a common identifier to identify which upstream segments should be considered part of the downstream "main channel".

Here is an example:

With the shown example, the darker red dots are the nodes, and the bright red dots are the stream sampling points. To calculate a main channel slope, if I write an automation program I need some sort of attribute to cue off of to be able to identify which of the shown upstream segments are part of the "main channel". This case is an instance where the segment that has the sampling points (bright red) and the upstream segement that branches to the left have the same GNIS ID so I could use that. But the NHD data have many instances where smaller streams of this nature simply don't have GNIS ID's at all. And the from/to node info or even the NHD Flow table which I have simply tells you that both upstream segments flow into the sampling point segment. With this example as a demonstration of other instances that lack a GNIS ID, I need a way to automate the decision of WHICH upstream segment to use for calculating the main channel slope in a scripted program. What else could be used to do this?

Another example of where the GNIS ID could not be used as a way to hydrologically identify all of the segments that make up the main channel is where an upstream section has one GNIS ID and unique GNIS Name (e.g., "stream X", and then it takes on a different name (e.g., stream "y") from a confluence with another trib. A point on stream "y" hydrologically still receives main channel flow from stream x but the GNIS ID can't be used to connect the two.


texts: OpenLayers WMS disappears when including OpenStreetMap

texts:
Setup:

ArcSDE 9.3.2
Oracle 11g
GeoServer 2.1.4
OpenLayers (latest build)
OpenStreetMap
InternetExplorer 7.0/ Google Chrome

EDIT: I have updated information just below the screen shots

Problem:

I am trying to take some feature classes from our arcsde and overlay them with OSM, originally i was trying to keep them in their native projection of BC_Albers (EPSG:42102) but i couldn't get them to overlay properly. So i consigned my self to having to do a reprojection of my data from 42102 to 4326.

I did the reprojection this morning and then modified my openlayers test page to view the data. I can view my featureset by itself just fine, I can view the OSM layer by itself just fine but when I try and combine the two all I get is OSM.

<title></title>

<form id="form1" runat="server">

location

</form>

it appears it may have something to do with the way that openlayers is changing the EPSG:

without OSM on map:
    http://wldev1.canfor.ca:8080/geoserver/WIM/wms?LAYERS=WIM%3AFOREST_WIM.WEB_COPI_TRAPPER&FORMAT=image%2Fgif&SERVICE=WMS&VERSION=1.1.1&REQUEST=GetMap&STYLES=&SRS=EPSG%3A4326&BBOX=-180,-90,0,90&WIDTH=256&HEIGHT=256

WITH OSM on map:
    http://wldev1.canfor.ca:8080/geoserver/WIM/wms?LAYERS=WIM%3AFOREST_WIM.WEB_COPI_TRAPPER&FORMAT=image%2Fgif&TRANSPARENT=true&SERVICE=WMS&VERSION=1.1.1&REQUEST=GetMap&STYLES=&SRS=EPSG%3A900913&BBOX=-180,-90,10018574.17,10018664.17&WIDTH=256&HEIGHT=256


texts:
I have a geopandas dataframe that contains bus boarding data. I'm trying to display a subset of the data on a plot using a Contextily basemap but the basemap does not appear.
I created the geometry attribute this way.
geometry = [Point(xy) for xy in zip(df['longitude'], df['latitude'])]
crs = {'init':'epsg:3857'}

gdf = gpd.GeoDataFrame(df, crs=crs, geometry=geometry)

My data looks like this.
stop name   |   latitude    |   longitude   |  geometry
----------------------------------------------------------------------
stop a      |  36.09541     |  -79.5464     | POINT(-79.5464,36.09541)
unknown     |  0.0000       |  0.00000      | POINT(0.0000,0.0000)
unknown     |  36.14894     |  -80.5465     | POINT(-80.5465,36.14894)

Some of the data has unknown stop names and I'm trying to plot those that have a lat/long. It is possible that some data do not have coordinate locations.
unknown_transfer_loc = transfer_stops[(transfer_stops['stopName'] == 'unknown') & (transfer_stops['latitude'] != 0)]

ax = unknown_transfer_loc.plot(figsize=(12,7))
cx.add_basemap(ax, zoom=10, source=cx.providers.CartoDB.Voyager)

Using the code above, I can successfully plot the points but the basemap is blank, see image. How can I get the basemap to appear?


texts: 2D export renderer creating anomalies using ArcScene

texts:
A lot have been improved in the 3D Analyst 10 including the ability to export large scenes as large 2D images. However, some problems still exist one of which is:
When I export a scene that contains a very large TIN the resulting image appears to have some triangles shaded darker which causes a "spotty" appearance. It's almost as if some "faces - triangles" were shaded incorrectly, or reversed?
My current graphics card is NVIDIA Quadro 4000.  (Used to have the NVIDIA Quadro FX 1700 which caused similar issues with the 9.3 version of 3D Analyst) I should also mention that it does not matter what type of image i export, compressed or not the result it more less the same.
Does anyone have any ideas as to how to rectify this?
It seems that this is only happening when the TIN is exaggerated, in this case 5x
If I drape imagery over the terrain and turn the IN off, there are no issues. This problem pertains only to the display of a TIN.
This is an ArcScene screen capture:

Exported Image (jpeg)


texts: How to load h2 database from geonetwork in netbeans?

texts:
I'm using Geonetwork 2.10.3. I have all needed files and i want to connect to h2 database in NetBeans, so I could add new language. Howerver, I can't connect to h2. I have followed these steps:
1. In Services tab, Right-Click on Databases and click New Connection;
2. In new window, in in padded menu, i have chosen New Driver and added .jar file from lib folder named h2-1.3.152.jar, and then Next;
3. I have typed user and pass in fields, url fieled has been filled with jdbc:h2:, and then Next.

After that, it throws me exception:
 Cannot establish a connection to jdbc:h2:  using org.h2.Driver (IO Exception: "java.io.FileNotFoundException: C:\Program Files\NetBeans 8.0.2\ .lock.db (The filename, directory name, or volume label syntax is incorrect)"; "C:\Program Files\NetBeans 8.0.2\ .lock.db" [90031-152]).
How can i solve this?


texts: Land Use classification and canopy cover using LT8 imagery and ground truth points in QGIS

texts:
I am trying to run a land use classification with QGIS on LANDSAT8 imagery to classify the various land uses within an area of northern Vietnam.

I want to distinguish between primary forest, secondary forest, sub canopy cultivation, monoculture crops, built up areas and riparian habitat.

I also want to classify the imagery for canopy cover (%).

I have recorded over 100 ground truth points in situ for land use type and canopy cover percentage but I have no knowledge on how to import this information into QGIS so that it can train my classification in order to create spectral signatures.

The only way I can think at the moment would be to import the UTMs with IDs then zoom to each one and create a spectral signature manually, but this will take forever.

I fear it may also create some error.

I am new to QGIS.


texts: How to set the maxBound for the map in OL3

texts:
I need to show the Map only for the World

like https://stackoverflow.com/questions/16558231/leafletjs-setmaxbounds-wont-let-zoom-out-to-see-whole-world

How to set this in OL3,

tried by giving extent to the view but not working

view : new ol.View({
                    center : center.getCoordinates(),
                    zoom : this.mapOptions.zoom,
                    extent : [-180,-180,180,180]
                })


texts: How to join a specific cell in excel to the corresponding polygons in QGIS/MapInfo

texts:
I would like to be able to create an excel database that has hyperlinks into a QGIS or MapInfo table so that if I clicked within the cell (which is linked to a particular polygon) I would be taken directly to the polygon on MapInfo/QGIS (I am saying either as I currently work on both platforms). The Site ID  in excel would match the cell to the polygon.

Can this be done?


texts: Make custom function using dynamic variable

texts:
In the function editor, I am trying to make this function :
@qgsfunction(args='auto', group='Custom')
def size_reduction( x, feature, parent):
    s = @map_scale/1000/x
    y = if( s<1, s, 1+ln(s) )
    return y

But I get the error shown below.
Can anyone help me ?
An error occurred during execution of following code:
from qgis.core import *
from qgis.gui import *

@qgsfunction(args='auto', group='Custom')
def size_reduction( x, feature, parent):
    s = @map_scale/1000/x
    y = if( s    return y

  File "", line 6
    s = @map_scale/1000/x
        ^
SyntaxError: invalid syntax

Python version:
3.8.5 (default, Jan 27 2021, 15:41:15)
[GCC 9.3.0]
QGIS version:
3.16.4-Hannover 'Hannover', exported


texts: Labeling Poles in Sequential Order on Data Driven Pages Along Curvy Route?

texts:
Is there a way to automate the process of labeling poles #'s on data driven pages along a curvy route in sequential order?

The many bends and turns of the route is throwing off the ordering. the snippet is a visual representation of our problem.


texts: Adding multiple icons to a single csv layer

texts:
I am loading the NumberIcon1 image for every point added through the csv file, is it possible to loop through the images so that way they go in order like NumberIcon 2 for the second point in the csv file and so on.

    layer1 = new CSVLayer("2.5_week.csv", {
      copyright: "USGS.gov"
  });
  var marker[] = new PictureMarkerSymbol("resources/markers/red/NumberIcon1.png", 20, 20);
  var renderer = new SimpleRenderer(marker);
  layer1.setRenderer(renderer);
  var template = new InfoTemplate("${type}", "${place}");
  layer1.setInfoTemplate(template);
  map.addLayer(layer1);


texts: Run-time pointer NULL error in C with OGR

texts:
I want to find out the geometry column and then want to export the geometry to WKB format from a shapefile. I wrote a code in C to find the geometry column. But a run-time error occurred. Can you please explain?

#include <stdio.h>
#include <stdlib.h>
#include <ogr_api.h>
#include <cpl_conv.h>
#include <cpl_port.h>
#include <gdal.h>
#include <string.h>
#include <gdal_utils.h>
#include <cpl_progress.h>
#include <cpl_minixml.h>
#include <ogr_core.h>
int main()
{
    char path1[200],filename1[]="/1.shp";
    GDALDatasetH inShapeDS;
    OGRLayerH inShapeLayer;
    OGRFeatureH inShapeFeature;

    GDALAllRegister();

    getcwd(path1, sizeof(path1));
    strcat(path1,filename1);

    inShapeDS = GDALOpenEx(path1, GDAL_OF_VECTOR, NULL, NULL, NULL );
    if( inShapeDS == NULL )
    {
        printf( "Open failed.\n" );
        exit( 1 );
    }

    inShapeLayer = GDALDatasetGetLayerByName(inShapeDS, path1);

    printf("number of layers:%d\n",GDALDatasetGetLayerCount(inShapeDS));

    printf("geom col name:%s",OGR_L_GetGeometryColumn(inShapeLayer));
    GDALClose(inShapeDS);

    return 0;
}

Error:

number of layers:1
ERROR 10: Pointer 'hLayer' is NULL in 'OGR_L_GetGeometryColumn'.

geom col name:(null)


texts: Add GeoJSON features from an API call in OpenLayers 4, layer not displayed

texts:
I'm trying to create a little demo to show the data from our API into a map but I cannot see my layer

Here my code, is just basic for the moment, we just need to display the data:

const projection = new ol.proj.Projection({
  code: 'EPSG:27700',
  getPointResolution: (r) => {
    return r
  },
  units: 'm'
})

var view = new ol.View({
  projection: projection,
  center: [383000, 398000],
  extent: [0, 0, 700000, 1250000],
  resolution: 28,
  resolutions: [140, 70, 28, 14, 7, 2.8, 1.4, 0.7],
})

var maplayer = new ol.layer.Tile({
  extent: [0, 0, 700000, 1250000],
  preload: Infinity,
  source: new ol.source.TileWMS({
    url: 'https://www.buchananmapping.co.uk/mapcache/?',
    params: {
      'LAYERS': 'bccache2',
      'TILED': true
    },
    ratio: 1,
    projection: 'EPSG:27700',
    tileGrid: new ol.tilegrid.TileGrid({
      resolutions: [140, 70, 28, 14, 7, 2.8, 1.4, 0.7],
      extent: [0, 0, 700000, 1250000],
      tileSize: [128, 128]
    })
  })
})

var appSource = new ol.source.Vector({
  projection: 'EPSG:27700',
  format: new ol.format.GeoJSON()
})

var appLayer = new ol.layer.Vector({
  source: appSource,
  style: new ol.style.Style({
    fill: new ol.style.Fill({
      color: 'rgba(255, 255, 255, 0.1)'
    }),
    stroke: new ol.style.Stroke({
      color: 'rgb(85, 172, 238)',
      width: 2
    })
  })
})

var map = new ol.Map({
  target: 'map',
  layers: [
    appLayer,
    maplayer
  ],
  view: view
})

  $.ajax({
    type: 'GET',
    url: 'https://www.traffwebdev.uk/tro/v1/restrictions?epsg=&x=383666+&y=398197&area=&radius=500&output=',
    dataType: 'json',
    success: function (data) {
      // reads and converts GeoJSon to Feature Object
      var features = (new ol.format.GeoJSON()).readFeatures(data.geometry)
      appSource.addFeatures(features)
    }
  })

Both my map and GeoJSON use EPSG:27700 but I cannot see any features in the map, I can see that the call to the API is made.

Sorry I don't have a demo but if you copy the whole code is working, there is no other code a part from that.


texts: Why should Google Earth's satellite imagery differ from what is available via the Openlayers plugin on QGIS?

texts:
I suspect that there's likely to be a relatively straightforward explanation of why, in some instances, satellite imagery on Google Earth may differ from that which is available from Google Maps on QGIS, but this oddity has puzzled me for quite sometime.

I'm involved in researching the built environment and would like a clear grasp of any limitations.  I'll illustrate this point in two screen-dump images of the same building.  The first one, immediately below, was taken from Google Earth software, and the other, beneath this, from QGIS using the Openlayers Plugin

 .

In this particular case the satellite imagery that I sourced from Google Earth seems to be a lot more recent as it shows a building that has been extended and is now considerably larger from what is shown via the OpenLayers Plugin. I'm guessing but perhaps this issue is due to how Google updates its maps. Can I use something other than the OpenLayers Plugin in QGIS, which may potentially be able to access more recent Google satellite imagery?


texts: What are the advantages of using Rails for a GIS/mapping application?

texts:
I've heard some good things about Rails, but I don't know much about it yet. There are only a few questions about Rails on GIS.SE (and there's not even a Rails tag).

Is it worthwhile to learn about Rails, from a GIS professional's perspective? (eg, what advantages does Rails give when building a website containing GIS functionality, compared to not using Rails?)
Are there any examples of good GIS/mapping websites built on Rails?

(I don't have any specific requirements in mind yet - I just want to know whether to invest time in learning Rails.)


texts: ArcGIS JavaScript Web Application: Workflow for Tokens?

texts:
This may appear to be a duplicate but I can't find any online source which is able to help me. Essentially, I want to protect a web Service from unauthorized access but still need to allow a web application to call the Service.

So, using ArcGIS Manager, I made a Service as "secure". After that, I launched: http://mygisserver/arcgis/tokens page, logged in and selected http referrer as default, entered the url of the web page hosting the web application like http://mywebserver/mygisapp/ and then generated a Token for 1 year.

Then, in the web application, I have http://mygisserver/arcgis/mysecureservice/mapserver?token=xxxxxxx for a ArcGISDynamicMapServiceLayer load parameter.

But running the web application I see an "Invalid Token" 498 error.

What could be happening? Is this even the right workflow or do I have to use some Proxy.

I don't want any login challenge for the users.

Update I am finding that the token is indeed working and the invalid token error I am getting are for the Legends and TOC widgets--somehow I need to pass the token to them. Hmmm.

Update 2 I am now trying to use the proxy.ashx. So I have defined the secure service's url in proxy.config, granted a new Arcgis Server user (say, proxy_culand) access to the locked service in Arcgis Manager, entered those credentials along with the url in the proxy.config file. After these, when I go to

http://mywebserver/proxy/proxy.ashx?http://mygisserver/arcgis/rest/services/CULAND/CU_LAND_v2a/MapServer then I still prompted for a password. Why is that?


texts: How to add attribute in shapefile guidedly

texts:
I have a point layer (31 by 31 grid = 961 points) created by fishnet tool > feature vertices to point tool (but with duplicate points). NOx Concentration data are in an Excel file.  I want to populate an NOx Concentration field in the point layer using the corresponding values in the Excel file automatically.  How can this be done?
Excel File

Map

I am writing this question in connection with my previous post.


texts: How does ArcGIS modify a TIF file while georeferencing?

texts:
I have some images, . TIF files that I am georeferencing using ESRI's software ArcGIS. I would like to know how does the program modify the initial .TIF file when one uses the Update Georeferencing tool after finding adequate Control Points.
Is there a way to revert such a georeferenced TIF file to its initial (ungeoreferenced) state by modifying it?


texts: ArcGIS IOS load TPK from Shared Container AGSLocalTiledLayer

texts:
I have 2 apps that are trying to share the same TPK offline base map.  I have App Groups enabled and have the TPK file stored in the shared container but it does not seem to work when referencing the file with AGSLocalTiledLayer.

The error message I see says "Cache at path '(null)' not found"

Does anyone know if this approach is possible?


texts: Cannot install qgis-mapserver, on Deb8.3,' Dependencies issue'

texts:
I am trying to install qgis-mapserver.

I started at Qgis org Installers website - Latest Release   2.12.x Lyon,
with a fresh Debian 8.3 jessie install.

The repos's are as instructed, which gives me:
qgis-mapserver 2.4.0+jessie1, that depends on:
qgis-providers (=2.4.0+jessie1) but 1:2.12.3+13jessie is all I can install.

Do I need to install 2.8.x Wien, to avoid the issues??

Also, I've read:
"The package is called "qgis-server" (no "map" part anymore) now.",
by Michal Mackiewicz, in this two year old Q&A: Cannot Install QGIS Mapserver On Ubuntu 12.04

No recent notes at the Qgis website, or from other Debian users,
on this issue?


texts: Why don't the Clip Raster by Mask show proper band values/interpolated values in QGIS?

texts:
I'm trying to clip a raster interpolated data (using IDW) with a polygon in QGIS 3.10. It worked with the first file (as shown in the linked here [desired Clipped (masked) layer] and this would be the result I would want.
For the rest of the interpolated layers that I have (which all went through similar steps to the first layer mentioned), I used the exact same methods to Clip Raster by Mask Layer with the desired polygon. However, the resulting layer would not display any values. I have included the screenshots of what I did in order as below:

I've tried to reupload my datasets and redo the interpolations, but nothing has worked so far. I have checked the projections as well and all were in the desired projection. It just seems odd to me that these steps worked well before except for now. Not sure if it is a bug on the QGIS software or if something else went wrong.


texts: list all authorized services for a given user in ArcGIS 10.1 API FLEX

texts:
I'm trying to get the list of all authorized services for a given user in ArcGIS 10.1 API for FLEX.
I can't find a way to automatically list all of the services granted to my account.

I can fix the services in the configuration file but I need to get them automatically in case of any authorizations modifications made by the administrator.

Thank you in advance.
N.


texts: Tabular display of features not working as expected in QGIS?

texts:
I have just installed QGIS 2.18 and was dismayed to find that "open attribute table" no longer displays data in a tabular format.  Instead I get list of 'filter' values with a check box and the rest of the pane taken up with the details of the selected value as a form.

I have just looked at the 2.18 Docs and this is clearly not what is meant to happen.

Any hints on how to get back the old behaviour?


texts: ArcMap 10.5 graticule incorrectly displays positive long tick labels

texts:
In ArcMap 10.5 I created a map of my study site and added a graticule. I realized that the longitude tick labels are missing the "-" sign e.g., displayed label 79°19'W should actually be -79°19'W, since it's western hemisphere.
How can I correct this in the graticule options?


texts:
I cannot find correct usage for using custom mapbox tilesets with folium.

From the documentation, "Folium also supports Cloudmade and Mapbox custom tilesets- simply pass your key to the API_key keyword" (http://python-visualization.github.io/folium/quickstart.html).

The docs provide this example...

folium.Map(location=[45.5236, -122.6750],
           tiles='Mapbox',
           API_key='your.API.key')

So I have tried...

#Create Map
my_map = folium.Map(location=[45.372, -121.6972],
                    zoom_start=12,
                    tiles='Mapbox',
                    API_key = 'pk.xxxxxxxxxxxxxxxxxxxxxxxxxxx') #mapbox map id

#Create Map
my_map = folium.Map(location=[45.372, -121.6972],
                    zoom_start=12,
                    tiles='Mapbox',
                    API_key = 'pk.xxxxxxxxxxxxxxxxxxxxxxxxxxx')#mb api key

#Create Map
 my_map = folium.Map(location=[45.372, -121.6972],
                    zoom_start=12,
                    tiles='user.2fb7fc73',
                    API_key = 'pk.xxxxxxxxxxxxxxxxxxxxxxx',
                    attr='Mapbox Data Attribution')

I have not been able to find any more example usage or explanation in the documents. If anyone has had success please advise?